00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #include "EvtGenBase/EvtPatches.hh"
00023 #include <stdlib.h>
00024 #include "EvtGenBase/EvtParticle.hh"
00025 #include "EvtGenBase/EvtDecayTable.hh"
00026 #include "EvtGenBase/EvtPDL.hh"
00027 #include "EvtGenModels/EvtPycont.hh"
00028 #include "EvtGenModels/EvtPythia.hh"
00029 #include "EvtGenBase/EvtId.hh"
00030 #include "EvtGenBase/EvtReport.hh"
00031 #include "EvtGenModels/EvtGlobalSet.hh"
00032 #include <string>
00033 #include <iostream>
00034 using std::endl;
00035
00036
00037 extern "C" struct {
00038 double maximum;
00039 } cbbeam_;
00040
00041
00042 extern "C" struct {
00043 int myisr;
00044 } isrflag_;
00045
00046 extern "C" struct
00047 {
00048 int mint[400];
00049 double vint[400];
00050 }pyint1_;
00051
00052 extern "C" struct
00053 {
00054 int mstp[200];
00055 double parp[200];
00056 int msti[200];
00057 double pari[200];
00058 }pypars_;
00059
00060 extern "C" {
00061 extern void pystat_(int &);
00062 }
00063
00064 extern "C" struct
00065 {
00066 int dc[18];
00067 } decaych_;
00068
00069
00070 extern "C" {
00071 extern void initpythia_(int *);
00072 }
00073
00074
00075 EvtPycont::~EvtPycont() {
00076
00077
00078 }
00079
00080 void EvtPycont::getName(std::string& model_name)
00081 {
00082 model_name="PYCONT";
00083 }
00084
00085 EvtDecayBase* EvtPycont::clone()
00086 {
00087 return new EvtPycont;
00088 }
00089
00090 void EvtPycont::init()
00091 {
00092
00093 if ( getNArg() != 12 && getNArg() != 0 && getNArg() != 1 ) {
00094 report(ERROR,"EvtGen") << "EvtPYCONT expects "
00095 << " 12 arguments (d u s c b t e nu_e mu nu_mu tau nu_tau) but found: "
00096 << getNArg() <<endl;
00097
00098 }
00099 checkNArg(0,1,12);
00100 isrflag_.myisr=1;
00101 for( int i=0; i<18; i++) decaych_.dc[i]=0;
00102 if ( getNArg() == 12 ) {
00103 decaych_.dc[0]=(int)getArg(0);
00104 decaych_.dc[1]=(int)getArg(1);
00105 decaych_.dc[2]=(int)getArg(2);
00106 decaych_.dc[3]=(int)getArg(3);
00107 decaych_.dc[4]=(int)getArg(4);
00108 decaych_.dc[5]=(int)getArg(5);
00109 decaych_.dc[10]=(int)getArg(6);
00110 decaych_.dc[11]=(int)getArg(7);
00111 decaych_.dc[12]=(int)getArg(8);
00112 decaych_.dc[13]=(int)getArg(9);
00113 decaych_.dc[14]=(int)getArg(10);
00114 decaych_.dc[15]=(int)getArg(11);
00115 } else if(getNArg()==1){
00116 isrflag_.myisr=(int)getArg(0);
00117 decaych_.dc[0]=1;
00118 decaych_.dc[1]=1;
00119 decaych_.dc[2]=1;
00120 decaych_.dc[3]=1;
00121 }
00122
00123 }
00124
00125 void EvtPycont::initProbMax()
00126 {
00127 noProbMax();
00128 }
00129
00130 void EvtPycont::decay( EvtParticle *p)
00131 {
00132
00133 double energy=p->mass();
00134 cbbeam_.maximum = energy;
00135
00136
00137
00138 if(energy==0){std::cout<<"Particle "<<EvtPDL::name(p->getId())<<" has zero mass"<<std::endl;abort();}
00139 EvtPythia::pythiaInit(0);
00140 int mydummy=0;
00141 initpythia_(&mydummy);
00142
00143 EvtVector4R p4[100];
00144
00145
00146 int i,more,nson;
00147 int ndaugjs;
00148 int kf[100];
00149 EvtId id[100];
00150 int type[MAX_DAUG];
00151
00152 double px[100],py[100],pz[100],e[100];
00153
00154 if ( p->getNDaug() != 0 ) { return;}
00155 do{
00156 EvtPythia::pythiacont(&energy,&ndaugjs,kf,px,py,pz,e);
00157
00158 double toteng=0;
00159
00160 for(i=0;i<ndaugjs;i++)
00161 {
00162
00163 id[i]=EvtPDL::evtIdFromStdHep(kf[i]);
00164
00165 type[i]=EvtPDL::getSpinType(id[i]);
00166
00167
00168
00169
00170
00171
00172 toteng += e[i];
00173 if (px[i]*px[i]+py[i]*py[i]+pz[i]*pz[i]>=e[i]*e[i])
00174 e[i]=sqrt(px[i]*px[i]+py[i]*py[i]+pz[i]*pz[i])+0.0000000000001;
00175
00176 p4[i].set(e[i],px[i],py[i],pz[i]);
00177
00178 }
00179
00180
00181 int channel=EvtDecayTable::inChannelList(p->getId(),ndaugjs,id);
00182
00183 more=((channel!=-1)&&(channel!=p->getChannel()) );
00184
00185
00186 }while(more);
00187
00188 p->makeDaughters(ndaugjs,id);
00189
00190 for(i=0;i<ndaugjs;i++)
00191 p->getDaug(i)->init( id[i], p4[i] );
00192
00193 if(energy<2.0) { EvtGlobalSet::ConExcPythia=0;}else{EvtGlobalSet::ConExcPythia=1;}
00194 return ;
00195 }
00196