00001
00002
00003
00004
00005
00006
00007
00008 #include "BesTwogam/BesTwogam.h"
00009
00010
00011 #include "HepMC/GenEvent.h"
00012
00013
00014
00015 #include "HepMC/HEPEVT_Wrapper.h"
00016 #include "HepMC/IO_HEPEVT.h"
00017
00018 #include "GaudiKernel/MsgStream.h"
00019 #include "GaudiKernel/ISvcLocator.h"
00020 #include "GaudiKernel/AlgFactory.h"
00021 #include "GaudiKernel/DataSvc.h"
00022 #include "GaudiKernel/SmartDataPtr.h"
00023
00024 #include "GeneratorObject/McGenEvent.h"
00025
00026 #include "cfortran/cfortran.h"
00027
00028 #include <stdlib.h>
00029
00030 typedef struct{
00031 int ISEED;
00032 }ISEEDC_DEF;
00033 #define ISEEDC COMMON_BLOCK(ISEEDC_DEF,iseedc)
00034 COMMON_BLOCK_DEF(ISEEDC_DEF,ISEEDC);
00035
00036 typedef struct{
00037 int i_cardName;
00038 char cardName[255];
00039 }FNAMES_DEF;
00040 #define FNAMES COMMON_BLOCK(FNAMES_DEF,fnames)
00041 COMMON_BLOCK_DEF(FNAMES_DEF,FNAMES);
00042
00043 typedef struct{
00044 int NEVTOT,NTUNW,NSMALL,IPROC1,LUNOUT,IDEBUG,LUN,NTRY;
00045 }INPARC_DEF;
00046 #define INPARC COMMON_BLOCK(INPARC_DEF,inparc)
00047 COMMON_BLOCK_DEF(INPARC_DEF,INPARC);
00048 typedef struct{
00049 float RLUMI;
00050 float CST;
00051 }LUMINO_DEF;
00052 #define LUMINO COMMON_BLOCK( LUMINO_DEF,lumino)
00053 COMMON_BLOCK_DEF(LUMINO_DEF,LUMINO);
00054 typedef struct{
00055 int NEVENT,NEVCUT,NEVUNW;
00056 }EVENUM_DEF;
00057 #define EVENUM COMMON_BLOCK(EVENUM_DEF,evenum)
00058 COMMON_BLOCK_DEF(EVENUM_DEF,EVENUM);
00059 typedef struct{
00060 double MAXESW,SUMW[5],MEANW[5],SUMW2[5],MAXW;
00061 double SUMWC[5],MEANWC[5],SUMW2C[5];
00062 }UNWEIT_DEF;
00063 #define UNWEIT COMMON_BLOCK(UNWEIT_DEF,unweit)
00064 COMMON_BLOCK_DEF(UNWEIT_DEF,UNWEIT);
00065 typedef struct{
00066 double WT,SUMWT,SUMWTC;
00067 }WEIGHTED_DEF;
00068 #define WEIGHTED COMMON_BLOCK(WEIGHTED_DEF,weighted)
00069 COMMON_BLOCK_DEF(WEIGHTED_DEF,WEIGHTED);
00070 typedef struct{
00071 int GOODEV;
00072 }CFLAGS_DEF;
00073 #define CFLAGS COMMON_BLOCK(CFLAGS_DEF,cflags)
00074 COMMON_BLOCK_DEF(CFLAGS_DEF,CFLAGS);
00075
00076 typedef struct{
00077 float OM1MNE, OM1MXE, OM2MNE, OM2MXE,
00078 T1MINE, T1MAXE, T2MINE, T2MAXE,
00079 Q12MNE, Q12MXE, Q22MNE, Q22MXE,
00080 W2MINE, W2MAXE, TKMINE, TKMAXE, EBEAME, MAXESW;
00081 }EXCUTS_DEF;
00082 #define EXCUTS COMMON_BLOCK(EXCUTS_DEF, excuts)
00083 COMMON_BLOCK_DEF(EXCUTS_DEF,EXCUTS);
00084
00085 typedef struct{
00086 int N;
00087 int K[5][4000];
00088 float P[5][4000],V[5][4000];
00089 }LUJETS_DEF;
00090 #define LUJETS COMMON_BLOCK(LUJETS_DEF,lujets)
00091 COMMON_BLOCK_DEF(LUJETS_DEF,LUJETS);
00092
00093 PROTOCCALLSFSUB0(PRINTEVT,printevt)
00094 #define PRINTEVT() CCALLSFSUB0(PRINTEVT,printevt)
00095
00096 PROTOCCALLSFSUB0(INIRUN,inirun)
00097 #define INIRUN() CCALLSFSUB0(INIRUN,inirun)
00098
00099 PROTOCCALLSFSUB0(FINISH,finish)
00100 #define FINISH() CCALLSFSUB0(FINISH,finish)
00101
00102 PROTOCCALLSFSUB0(GENEVT,genevt)
00103 #define GENEVT() CCALLSFSUB0(GENEVT,genevt)
00104
00105 PROTOCCALLSFSUB0(GEN1EVT,gen1evt)
00106 #define GEN1EVT() CCALLSFSUB0(GEN1EVT,gen1evt)
00107
00108 PROTOCCALLSFSUB0( RDCUTS , rdcuts)
00109 #define RDCUTS() CCALLSFSUB0( RDCUTS,rdcuts)
00110
00111 PROTOCCALLSFSUB0( SETDEF , setdef)
00112 #define SETDEF() CCALLSFSUB0( SETDEF,setdef)
00113
00114 PROTOCCALLSFSUB0( QCDSTA , qcdsta)
00115 #define QCDSTA() CCALLSFSUB0( QCDSTA,qcdsta)
00116
00117 PROTOCCALLSFSUB0( DOCUTS , docuts)
00118 #define DOCUTS() CCALLSFSUB0( DOCUTS,docuts)
00119
00120 PROTOCCALLSFSUB0( ESTMAXCS , estmaxcs)
00121 #define ESTMAXCS() CCALLSFSUB0( ESTMAXCS,estmaxcs)
00122
00123 PROTOCCALLSFSUB0( PRTCUT , prtcut)
00124 #define PRTCUT() CCALLSFSUB0( PRTCUT,prtcut)
00125
00126 PROTOCCALLSFSUB0( USRINI , usrini)
00127 #define USRINI() CCALLSFSUB0( USRINI,usrini)
00128
00129 PROTOCCALLSFSUB0( USRHIS , usrhis)
00130 #define USRHIS() CCALLSFSUB0( USRHIS,usrhis)
00131
00132 PROTOCCALLSFSUB0( USREND , usrend)
00133 #define USREND() CCALLSFSUB0( USREND,usrend)
00134
00135 PROTOCCALLSFSUB1(LUHEPC, luhepc, INT)
00136 #define LUHEPC(ICONV) CCALLSFSUB1(LUHEPC, luhepc, INT, ICONV)
00137
00138 PROTOCCALLSFSUB0(HEPEVT_CLEAN, hepevt_clean)
00139 #define HEPEVT_CLEAN() CCALLSFSUB0(HEPEVT_CLEAN, hepevt_clean)
00140
00141 PROTOCCALLSFSUB1(LULIST, lulist, INT)
00142 #define LULIST(ICONV) CCALLSFSUB1(LULIST, lulist, INT, ICONV)
00143
00144 PROTOCCALLSFSUB0(HEPEVT_PRINT, hepevt_print)
00145 #define HEPEVT_PRINT() CCALLSFSUB0(HEPEVT_PRINT, hepevt_print)
00146
00147 extern "C" {
00148 extern void pygive_(const char *cnfgstr,int length);
00149 }
00150
00151
00152 BesTwogam::BesTwogam(const string& name, ISvcLocator* pSvcLocator):Algorithm( name, pSvcLocator )
00153 {
00154 declareProperty("InitialSeed",m_seed=1001);
00155 declareProperty("CMEnergy", m_cmEnergy = 3.097);
00156 declareProperty("MinimumW2",m_m2min=0.02);
00157 declareProperty("Unweighted",m_unw=1);
00158 declareProperty("CardFileName", m_fcard="BesTwogam.conf");
00159
00160 m_pypars.clear();
00161 declareProperty("setPythiaPars", m_pypars);
00162 }
00163
00164 StatusCode BesTwogam::initialize(){
00165 MsgStream log(messageService(), name());
00166 log << MSG::WARNING << "BesTwogam initialize" << endreq;
00167
00168 if(m_fcard.size()>254){
00169 log << MSG::ERROR << "Too big name of BesTwogam card file" << endreq;
00170 return StatusCode::FAILURE;
00171 }
00172
00173 strcpy(FNAMES.cardName,m_fcard.c_str());
00174 FNAMES.i_cardName=m_fcard.size();
00175
00176 RDCUTS();
00177 ISEEDC.ISEED=m_seed;
00178 EXCUTS.EBEAME=0.5*m_cmEnergy;
00179 EXCUTS.W2MINE=m_m2min;
00180 INPARC.NTUNW=m_unw;
00181
00182 SETDEF();
00183 QCDSTA();
00184 DOCUTS();
00185 ESTMAXCS();
00186 PRTCUT();
00187
00188 m_numberEvent=0;
00189
00190 for(int i=0;i<m_pypars.size();i++){
00191 pygive_(m_pypars[i].c_str(),strlen(m_pypars[i].c_str()));
00192 }
00193 return StatusCode::SUCCESS;
00194 }
00195
00196 StatusCode BesTwogam::execute()
00197 {
00198 MsgStream log(messageService(), name());
00199 log << MSG::INFO << "BesTwogam executing" << endreq;
00200 HepMC::HEPEVT_Wrapper::set_max_number_entries(2000);
00201 HepMC::HEPEVT_Wrapper::set_sizeof_real(8);
00202 HepMC::IO_HEPEVT HepEvtIO;
00203 HEPEVT_CLEAN();
00204 GEN1EVT();
00205
00206 if(CFLAGS.GOODEV!=1){
00207 log << MSG::ERROR<<" BesTwogam: fail to generate good event"<<endl;
00208 return StatusCode::FAILURE;
00209 }
00210 m_numberEvent++;
00211 if( log.level() < MSG::INFO )LULIST(1);
00212 LUHEPC(1);
00213
00214 HepMC::GenEvent* evt = HepEvtIO.read_next_event();
00215 evt->set_event_number(m_numberEvent);
00216 evt->set_signal_process_id(1);
00217
00218 SmartDataPtr<McGenEventCol> anMcCol(eventSvc(), "/Event/Gen");
00219 if (anMcCol!=0) {
00220
00221 MsgStream log(messageService(), name());
00222 log << MSG::INFO << "Add McGenEvent to existing collection" << endreq;
00223 McGenEvent* mcEvent = new McGenEvent(evt);
00224 anMcCol->push_back(mcEvent);
00225 } else {
00226
00227 McGenEventCol *mcColl = new McGenEventCol;
00228 McGenEvent* mcEvent = new McGenEvent(evt);
00229 mcColl->push_back(mcEvent);
00230 StatusCode sc = eventSvc()->registerObject("/Event/Gen",mcColl);
00231 if (sc != StatusCode::SUCCESS) {
00232 log << MSG::ERROR << "Could not register McGenEvent" << endreq;
00233 delete mcColl;
00234 delete evt;
00235 delete mcEvent;
00236 return StatusCode::FAILURE;
00237 }
00238 }
00239 return StatusCode::SUCCESS;
00240 }
00241
00242 StatusCode BesTwogam::finalize()
00243 {
00244 MsgStream log(messageService(), name());
00245 log << MSG::INFO << "BesTwogam finalized" << endreq;
00246 FINISH();
00247 return StatusCode::SUCCESS;
00248 }