00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #include "EvtGenBase/EvtPatches.hh"
00022 #include <stdlib.h>
00023 #include <stdio.h>
00024 #include <string.h>
00025 #include "EvtGenBase/EvtParticle.hh"
00026 #include "EvtGenBase/EvtDecayTable.hh"
00027 #include "EvtGenBase/EvtPDL.hh"
00028 #include "EvtGenModels/EvtJscont.hh"
00029 #include "EvtGenModels/EvtJetSet.hh"
00030 #include "EvtGenBase/EvtId.hh"
00031 #include "EvtGenBase/EvtReport.hh"
00032 #include <string>
00033 #include <strstream>
00034 using std::strstream;
00035
00036 extern "C" {
00037 extern void continuum_(double *,int *,int *,int *,
00038 double *,double *,double *,double *);
00039 }
00040
00041 extern "C" {
00042 extern void lugive_(const char *cnfgstr,int length);
00043 }
00044
00045
00046 EvtJscont::~EvtJscont() {}
00047
00048 void EvtJscont::getName(std::string& model_name){
00049
00050 model_name="JSCONT";
00051
00052 }
00053
00054 EvtDecayBase* EvtJscont::clone(){
00055
00056 return new EvtJscont;
00057
00058 }
00059
00060 void EvtJscont::init(){
00061
00062
00063
00064 checkNArg(1,2);
00065
00066 }
00067
00068
00069 void EvtJscont::initProbMax(){
00070
00071 noProbMax();
00072
00073 }
00074
00075
00076 void EvtJscont::decay( EvtParticle *p){
00077
00078 EvtJetSet::jetSetInit();
00079 static int first=1;
00080
00081 if (first){
00082 first=0;
00083
00084 float val=0.6;
00085 if ( getNArg()>1) {
00086 val=getArg(1);
00087 }
00088 char vak[20];
00089 sprintf(vak,"PARJ(13)=%f",val);
00090 std::string temp(vak);
00091 lugive_(temp.c_str(),strlen(temp.c_str()));
00092 }
00093 EvtVector4R p4[100];
00094
00095 double energy=p->mass();
00096
00097 int flavor;
00098
00099 int i,more;
00100 int ndaugjs;
00101 int kf[100];
00102 EvtId id[100];
00103 int type[MAX_DAUG];
00104
00105 flavor=(int)getArg(0);
00106
00107 double px[100],py[100],pz[100],e[100];
00108
00109 if ( p->getNDaug() != 0 ) { return;}
00110 do{
00111
00112 continuum_(&energy,&flavor,&ndaugjs,kf,px,py,pz,e);
00113
00114 for(i=0;i<ndaugjs;i++){
00115
00116 id[i]=EvtPDL::evtIdFromStdHep(kf[i]);
00117
00118 type[i]=EvtPDL::getSpinType(id[i]);
00119
00120
00121
00122
00123
00124 if (px[i]*px[i]+py[i]*py[i]+pz[i]*pz[i]>=e[i]*e[i]){
00125
00126 e[i]=sqrt(px[i]*px[i]+py[i]*py[i]+pz[i]*pz[i])+0.0000000000001;
00127
00128 }
00129
00130 p4[i].set(e[i],px[i],py[i],pz[i]);
00131
00132 }
00133
00134 int channel=EvtDecayTable::inChannelList(p->getId(),ndaugjs,id);
00135
00136 more=((channel!=-1)&&(channel!=p->getChannel()));
00137
00138
00139 }while(more);
00140
00141 p->makeDaughters(ndaugjs,id);
00142
00143 for(i=0;i<ndaugjs;i++){
00144 p->getDaug(i)->init( id[i], p4[i] );
00145 }
00146 return ;
00147 }
00148
00149
00150
00151