CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Py8GunBase.cc
Go to the documentation of this file.
3 
4 using namespace Pythia8;
5 
7 
8 namespace gen {
9 
10 Py8GunBase::Py8GunBase( edm::ParameterSet const& ps )
11  : Py8InterfaceBase(ps)
12 {
13 
15  ps.getUntrackedParameter<double>("filterEfficiency", -1.) );
17  GenRunInfoProduct::XSec(ps.getUntrackedParameter<double>("crossSection", -1.)) );
19  GenRunInfoProduct::XSec(ps.getUntrackedParameter<double>("crossSectionNLO", -1.)) );
20 
21 
22  // PGun specs
23  //
24  edm::ParameterSet pgun_params =
25  ps.getParameter<edm::ParameterSet>("PGunParameters");
26 
27  // although there's the method ParameterSet::empty(),
28  // it looks like it's NOT even necessary to check if it is,
29  // before trying to extract parameters - if it is empty,
30  // the default values seem to be taken
31  //
32  fPartIDs = pgun_params.getParameter< std::vector<int> >("ParticleID");
33  fMinPhi = pgun_params.getParameter<double>("MinPhi"); // ,-3.14159265358979323846);
34  fMaxPhi = pgun_params.getParameter<double>("MaxPhi"); // , 3.14159265358979323846);
35 
36 }
37 
38 // specific to Py8GunHad !!!
39 //
41 {
42 
43  // NO MATTER what was this setting below, override it before init
44  // - this is essencial for the PGun mode
45 
46  // Key requirement: switch off ProcessLevel, and thereby also PartonLevel.
47  fMasterGen->readString("ProcessLevel:all = off");
48  fMasterGen->readString("ProcessLevel::resonanceDecays=on");
49  fMasterGen->init();
50 
51  // init decayer
52  fDecayer->readString("ProcessLevel:all = off"); // Same trick!
53  fDecayer->readString("ProcessLevel::resonanceDecays=on");
54  fDecayer->init();
55 
56  return true;
57 
58 }
59 
61 {
62 
63  Event* pythiaEvent = &(fMasterGen->event);
64 
65  int NPartsBeforeDecays = pythiaEvent->size()-1; // do NOT count the very 1st "system" particle
66  // in Pythia8::Event record; it does NOT even
67  // get translated by the HepMCInterface to the
68  // HepMC::GenEvent record !!!
69  int NPartsAfterDecays = event().get()->particles_size();
70  int NewBarcode = NPartsAfterDecays;
71 
72  for ( int ipart=NPartsAfterDecays; ipart>NPartsBeforeDecays; ipart-- )
73  {
74 
75  HepMC::GenParticle* part = event().get()->barcode_to_particle( ipart );
76 
77  if ( part->status() == 1 )
78  {
79  fDecayer->event.reset();
80  Particle py8part( part->pdg_id(), 93, 0, 0, 0, 0, 0, 0,
81  part->momentum().x(),
82  part->momentum().y(),
83  part->momentum().z(),
84  part->momentum().t(),
85  part->generated_mass() );
86  HepMC::GenVertex* ProdVtx = part->production_vertex();
87  py8part.vProd( ProdVtx->position().x(), ProdVtx->position().y(),
88  ProdVtx->position().z(), ProdVtx->position().t() );
89  py8part.tau( (fDecayer->particleData).tau0( part->pdg_id() ) );
90  fDecayer->event.append( py8part );
91  int nentries = fDecayer->event.size();
92  if ( !fDecayer->event[nentries-1].mayDecay() ) continue;
93  fDecayer->next();
94  int nentries1 = fDecayer->event.size();
95  if ( nentries1 <= nentries ) continue; //same number of particles, no decays...
96 
97  part->set_status(2);
98 
99  Particle& py8daughter = fDecayer->event[nentries]; // the 1st daughter
100  HepMC::GenVertex* DecVtx = new HepMC::GenVertex( HepMC::FourVector(py8daughter.xProd(),
101  py8daughter.yProd(),
102  py8daughter.zProd(),
103  py8daughter.tProd()) );
104 
105  DecVtx->add_particle_in( part ); // this will cleanup end_vertex if exists, replace with the new one
106  // I presume (vtx) barcode will be given automatically
107 
108  HepMC::FourVector pmom( py8daughter.px(), py8daughter.py(), py8daughter.pz(), py8daughter.e() );
109 
110  HepMC::GenParticle* daughter =
111  new HepMC::GenParticle( pmom, py8daughter.id(), 1 );
112 
113  NewBarcode++;
114  daughter->suggest_barcode( NewBarcode );
115  DecVtx->add_particle_out( daughter );
116 
117  for ( int ipart1=nentries+1; ipart1<nentries1; ipart1++ )
118  {
119  py8daughter = fDecayer->event[ipart1];
120  HepMC::FourVector pmomN( py8daughter.px(), py8daughter.py(), py8daughter.pz(), py8daughter.e() );
121  HepMC::GenParticle* daughterN =
122  new HepMC::GenParticle( pmomN, py8daughter.id(), 1 );
123  NewBarcode++;
124  daughterN->suggest_barcode( NewBarcode );
125  DecVtx->add_particle_out( daughterN );
126  }
127 
128  event().get()->add_vertex( DecVtx );
129 
130  }
131  }
132 
133  return true;
134 
135 }
136 
138 {
139 
140 
141 // FIXME !!!
142 
143  //******** Verbosity ********
144 
145  if (maxEventsToPrint > 0 &&
147  {
150  {
151  fMasterGen->info.list(std::cout);
152  fMasterGen->event.list(std::cout);
153  }
154 
156  {
157  std::cout << "Event process = "
158  << fMasterGen->info.code() << "\n"
159  << "----------------------" << std::endl;
160  event()->print();
161  }
162  }
163 
164  return;
165 }
166 
168 {
169 
170  fMasterGen->statistics();
171 
172  double xsec = fMasterGen->info.sigmaGen(); // cross section in mb
173  xsec *= 1.0e9; // translate to pb (CMS/Gen "convention" as of May 2009)
174  runInfo().setInternalXSec(xsec);
175  return;
176 
177 }
178 
179 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
std::auto_ptr< Pythia8::Pythia > fMasterGen
void setFilterEfficiency(double effic)
double fMinPhi
Definition: Py8GunBase.h:75
bool initializeForInternalPartons()
Definition: Py8GunBase.cc:40
void setInternalXSec(const XSec &xsec)
virtual bool residualDecay()
Definition: Py8GunBase.cc:60
void finalizeEvent()
Definition: Py8GunBase.cc:137
static const std::vector< std::string > p8SharedResources
Definition: Py8GunBase.h:85
unsigned int pythiaPylistVerbosity
std::vector< int > fPartIDs
Definition: Py8GunBase.h:74
std::auto_ptr< HepMC::GenEvent > & event()
Definition: Py8GunBase.h:69
part
Definition: HCALResponse.h:20
unsigned int maxEventsToPrint
double fMaxPhi
Definition: Py8GunBase.h:76
std::auto_ptr< Pythia8::Pythia > fDecayer
tuple cout
Definition: gather_cfg.py:121
static const std::string kPythia8
void setExternalXSecNLO(const XSec &xsec)
GenRunInfoProduct & runInfo()
Definition: Py8GunBase.h:68
void setExternalXSecLO(const XSec &xsec)