![]() |
HepMC Reference DocumentationHepMC |
00001 00002 // Matt.Dobbs@Cern.CH, October 2002 00003 // example of generating events with Herwig using HepMC/HerwigWrapper.h 00004 // Events are read into the HepMC event record from the FORTRAN HEPEVT 00005 // common block using the IO_HERWIG strategy. 00016 00017 #include <iostream> 00018 #include "HepMC/HerwigWrapper.h" 00019 #include "HepMC/IO_HERWIG.h" 00020 #include "HepMC/IO_GenEvent.h" 00021 #include "HepMC/GenEvent.h" 00022 #include "HepMC/HEPEVT_Wrapper.h" 00023 #include "HerwigHelper.h" 00024 00025 int main() { 00026 // 00027 //........................................HEPEVT 00028 // Herwig 6.4 uses HEPEVT with 4000 entries and 8-byte floating point 00029 // numbers. We need to explicitly pass this information to the 00030 // HEPEVT_Wrapper. 00031 // 00032 HepMC::HEPEVT_Wrapper::set_max_number_entries(4000); 00033 HepMC::HEPEVT_Wrapper::set_sizeof_real(8); 00034 // 00035 //.......................................INITIALIZATIONS 00036 00037 hwproc.PBEAM1 = 7000.; // energy of beam1 00038 hwproc.PBEAM2 = 7000.; // energy of beam2 00039 // 1610 = gg->H--> WW, 1706 = qq-->ttbar, 2510 = ttH -> ttWW 00040 hwproc.IPROC = 1706; // qq -> ttbar production 00041 hwproc.MAXEV = 100; // number of events 00042 // tell it what the beam particles are: 00043 for ( unsigned int i = 0; i < 8; ++i ) { 00044 hwbmch.PART1[i] = (i < 1) ? 'P' : ' '; 00045 hwbmch.PART2[i] = (i < 1) ? 'P' : ' '; 00046 } 00047 hwigin(); // INITIALISE OTHER COMMON BLOCKS 00048 hwevnt.MAXPR = 1; // number of events to print 00049 hwuinc(); // compute parameter-dependent constants 00050 hweini(); // initialise elementary process 00051 00052 //........................................HepMC INITIALIZATIONS 00053 // 00054 // Instantiate an IO strategy for reading from HEPEVT. 00055 HepMC::IO_HERWIG hepevtio; 00056 // Instantiate an IO strategy to write the data to file 00057 HepMC::IO_GenEvent ascii_io("example_MyHerwig.dat",std::ios::out); 00058 // 00059 //........................................EVENT LOOP 00060 for ( int i = 1; i <= hwproc.MAXEV; i++ ) { 00061 if ( i%50==1 ) std::cout << "Processing Event Number " 00062 << i << std::endl; 00063 // initialise event 00064 hwuine(); 00065 // generate hard subprocess 00066 hwepro(); 00067 // generate parton cascades 00068 hwbgen(); 00069 // do heavy object decays 00070 hwdhob(); 00071 // do cluster formation 00072 hwcfor(); 00073 // do cluster decays 00074 hwcdec(); 00075 // do unstable particle decays 00076 hwdhad(); 00077 // do heavy flavour hadron decays 00078 hwdhvy(); 00079 // add soft underlying event if needed 00080 hwmevt(); 00081 // finish event 00082 hwufne(); 00083 HepMC::GenEvent* evt = hepevtio.read_next_event(); 00084 // define the units (Herwig uses GeV and mm) 00085 evt->use_units(HepMC::Units::GEV, HepMC::Units::MM); 00086 // set cross section information 00087 evt->set_cross_section( getHerwigCrossSection(i) ); 00088 // add some information to the event 00089 evt->set_event_number(i); 00090 evt->set_signal_process_id(20); 00091 if (i<=hwevnt.MAXPR) { 00092 std::cout << "\n\n This is the FIXED version of HEPEVT as " 00093 << "coded in IO_HERWIG " << std::endl; 00094 HepMC::HEPEVT_Wrapper::print_hepevt(); 00095 evt->print(); 00096 } 00097 // write the event to the ascii file 00098 ascii_io << evt; 00099 00100 // we also need to delete the created event from memory 00101 delete evt; 00102 } 00103 //........................................TERMINATION 00104 hwefin(); 00105 00106 return 0; 00107 }