CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Static Public Member Functions | Private Attributes | Static Private Attributes
gen::EvtGenInterface Class Reference

#include <EvtGenInterface.h>

Inheritance diagram for gen::EvtGenInterface:
gen::EvtGenInterfaceBase

Public Member Functions

void addToHepMC (HepMC::GenParticle *partHep, EvtId idEvt, HepMC::GenEvent *theEvent, bool del_daug)
 
HepMC::GenEvent * decay (HepMC::GenEvent *)
 
 EvtGenInterface (const edm::ParameterSet &)
 
void go_through_daughters (EvtParticle *part)
 
void init ()
 
const std::vector< int > & operatesOnParticles ()
 
void setRandomEngine (CLHEP::HepRandomEngine *v)
 
void update_candlist (int theIndex, HepMC::GenParticle *thePart)
 
 ~EvtGenInterface ()
 
- Public Member Functions inherited from gen::EvtGenInterfaceBase
 EvtGenInterfaceBase ()
 
virtual void SetPhotosDecayRandomEngine (CLHEP::HepRandomEngine *decayRandomEngine)
 
virtual ~EvtGenInterfaceBase ()
 

Static Public Member Functions

static double flat ()
 

Private Attributes

std::string decay_table_s
 
std::vector< EvtId > forced_Evt
 
std::vector< int > forced_Hep
 
std::vector< std::string > forced_names
 
int index [10]
 
HepMC::GenParticle * listp [10]
 
EvtGen * m_EvtGen
 
Pythia6Service * m_Py6Service
 
int nevent
 
int nforced
 
int nlist
 
int npartial
 
int nPythia
 
int ntotal
 
std::string pdt_s
 
std::map< int, float > polarizations
 
std::vector< int > polarize_ids
 
std::vector< double > polarize_pol
 
myEvtRandomEngine * the_engine
 
bool useDefault
 
bool usePythia
 
std::string user_decay_s
 

Static Private Attributes

static CLHEP::HepRandomEngine * fRandomEngine
 

Additional Inherited Members

- Protected Attributes inherited from gen::EvtGenInterfaceBase
std::vector< int > m_PDGs
 

Detailed Description

Definition at line 46 of file EvtGenInterface.h.

Constructor & Destructor Documentation

EvtGenInterface::EvtGenInterface ( const edm::ParameterSet &  pset)

Definition at line 48 of file EvtGenInterface.cc.

References gather_cfg::cout, edm::hlt::Exception, edm::ParameterSet::exists(), edm::FileInPath::fullPath(), edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), nevent, and groupFilesInBlocks::ntotal.

49 {
50 
51  ntotal = 0;
52  nevent = 0;
53  std::cout << " EvtGenProducer starting ... " << std::endl;
54 
55  the_engine = new myEvtRandomEngine(nullptr);
56 
57  // Get data from parameter set
58  edm::FileInPath decay_table = pset.getParameter<edm::FileInPath>("decay_table");
59  edm::FileInPath pdt = pset.getParameter<edm::FileInPath>("particle_property_file");
60  useDefault = pset.getUntrackedParameter<bool>("use_default_decay",true);
61  usePythia = pset.getUntrackedParameter<bool>("use_internal_pythia",true);
62  polarize_ids = pset.getUntrackedParameter<std::vector<int> >("particles_to_polarize",
63  std::vector<int>());
64  polarize_pol = pset.getUntrackedParameter<std::vector<double> >("particle_polarizations",
65  std::vector<double>());
66  if (polarize_ids.size() != polarize_pol.size()) {
67  throw cms::Exception("Configuration")
68  << "EvtGenProducer requires that the particles_to_polarize and particle_polarization\n"
69  "vectors be the same size. Please fix this in your configuration.";
70  }
71  for (unsigned int ndx = 0; ndx < polarize_ids.size(); ndx++) {
72  if (polarize_pol[ndx] < -1. || polarize_pol[ndx] > 1.) {
73  throw cms::Exception("Configuration")
74  << "EvtGenProducer error: particle polarizations must be in the range -1 < P < 1";
75  }
76  polarizations.insert(std::pair<int, float>(polarize_ids[ndx], polarize_pol[ndx]));
77  }
78 
79  edm::FileInPath user_decay = pset.getParameter<edm::FileInPath>("user_decay_file");
80  decay_table_s = decay_table.fullPath();
81  pdt_s = pdt.fullPath();
82  user_decay_s = user_decay.fullPath();
83 
84  //-->pythia_params = pset.getParameter< std::vector<std::string> >("processParameters");
85 
86 
87  // any number of alias names for forced decays can be specified using dynamic std vector of strings
88  forced_names = pset.getParameter< std::vector<std::string> >("list_forced_decays");
89 
90 
91  // fill up default list of particles to be declared stable in the "master generator"
92  // these are assumed to be PDG ID's
93  // in case of combo with Pythia6, translation is done in Pythia6Hadronizer
94  //
95  // Note: Pythia6's kc=43, 44, and 84 commented out because they're obsolete (per S.Mrenna)
96  //
97  m_PDGs.push_back( 300553 ) ;
98  m_PDGs.push_back( 511 ) ;
99  m_PDGs.push_back( 521 ) ;
100  m_PDGs.push_back( 523 ) ;
101  m_PDGs.push_back( 513 ) ;
102  m_PDGs.push_back( 533 ) ;
103  m_PDGs.push_back( 531 ) ;
104 
105  m_PDGs.push_back( 15 ) ;
106 
107  m_PDGs.push_back( 413 ) ;
108  m_PDGs.push_back( 423 ) ;
109  m_PDGs.push_back( 433 ) ;
110  m_PDGs.push_back( 411 ) ;
111  m_PDGs.push_back( 421 ) ;
112  m_PDGs.push_back( 431 ) ;
113  m_PDGs.push_back( 10411 );
114  m_PDGs.push_back( 10421 );
115  m_PDGs.push_back( 10413 );
116  m_PDGs.push_back( 10423 );
117  m_PDGs.push_back( 20413 );
118  m_PDGs.push_back( 20423 );
119 
120  m_PDGs.push_back( 415 );
121  m_PDGs.push_back( 425 );
122  m_PDGs.push_back( 10431 );
123  m_PDGs.push_back( 20433 );
124  m_PDGs.push_back( 10433 );
125  m_PDGs.push_back( 435 );
126 
127  m_PDGs.push_back( 310 );
128  m_PDGs.push_back( 311 );
129  m_PDGs.push_back( 313 );
130  m_PDGs.push_back( 323 );
131  m_PDGs.push_back( 10321 );
132  m_PDGs.push_back( 10311 );
133  m_PDGs.push_back( 10313 );
134  m_PDGs.push_back( 10323 );
135  m_PDGs.push_back( 20323 );
136  m_PDGs.push_back( 20313 );
137  m_PDGs.push_back( 325 );
138  m_PDGs.push_back( 315 );
139 
140  m_PDGs.push_back( 100313 );
141  m_PDGs.push_back( 100323 );
142  m_PDGs.push_back( 30313 );
143  m_PDGs.push_back( 30323 );
144  m_PDGs.push_back( 30343 );
145  m_PDGs.push_back( 30353 );
146  m_PDGs.push_back( 30363 );
147 
148  m_PDGs.push_back( 111 );
149  m_PDGs.push_back( 221 );
150  m_PDGs.push_back( 113 );
151  m_PDGs.push_back( 213 );
152  m_PDGs.push_back( 223 );
153  m_PDGs.push_back( 331 );
154  m_PDGs.push_back( 333 );
155  m_PDGs.push_back( 20213 );
156  m_PDGs.push_back( 20113 );
157  m_PDGs.push_back( 215 );
158  m_PDGs.push_back( 115 );
159  m_PDGs.push_back( 10213 );
160  m_PDGs.push_back( 10113 );
161  m_PDGs.push_back( 9000111 ); // PDG ID = 9000111, Pythia6 ID = 10111
162  m_PDGs.push_back( 9000211 ); // PDG ID = 9000211, Pythia6 ID = 10211
163  m_PDGs.push_back( 9010221 ); // PDG ID = 9010211, Pythia6 ID = ???
164  m_PDGs.push_back( 10221 );
165  m_PDGs.push_back( 20223 );
166  m_PDGs.push_back( 20333 );
167  m_PDGs.push_back( 225 );
168  m_PDGs.push_back( 9020221 ); // PDG ID = 9020211, Pythia6 ID = ???
169  m_PDGs.push_back( 335 );
170  m_PDGs.push_back( 10223 );
171  m_PDGs.push_back( 10333 );
172  m_PDGs.push_back( 100213 );
173  m_PDGs.push_back( 100113 );
174 
175  m_PDGs.push_back( 441 );
176  m_PDGs.push_back( 100441 );
177  m_PDGs.push_back( 443 );
178  m_PDGs.push_back( 100443 );
179  m_PDGs.push_back( 9000443 );
180  m_PDGs.push_back( 9010443 );
181  m_PDGs.push_back( 9020443 );
182  m_PDGs.push_back( 10441 );
183  m_PDGs.push_back( 20443 );
184  m_PDGs.push_back( 445 );
185 
186  m_PDGs.push_back( 30443 );
187  m_PDGs.push_back( 551 );
188  m_PDGs.push_back( 553 );
189  m_PDGs.push_back( 100553 );
190  m_PDGs.push_back( 200553 );
191  m_PDGs.push_back( 10551 );
192  m_PDGs.push_back( 20553 );
193  m_PDGs.push_back( 555 );
194  m_PDGs.push_back( 10553 );
195 
196  m_PDGs.push_back( 110551 );
197  m_PDGs.push_back( 120553 );
198  m_PDGs.push_back( 100555 );
199  m_PDGs.push_back( 210551 );
200  m_PDGs.push_back( 220553 );
201  m_PDGs.push_back( 200555 );
202  m_PDGs.push_back( 30553 );
203  m_PDGs.push_back( 20555 );
204 
205  m_PDGs.push_back( 557 );
206  m_PDGs.push_back( 130553 );
207  m_PDGs.push_back( 120555 );
208  m_PDGs.push_back( 100557 );
209  m_PDGs.push_back( 110553 );
210  m_PDGs.push_back( 210553 );
211  m_PDGs.push_back( 10555 );
212  m_PDGs.push_back( 110555 );
213 
214  m_PDGs.push_back( 4122 );
215  m_PDGs.push_back( 4132 );
216  // m_PDGs.push_back( 84 ); // obsolete
217  m_PDGs.push_back( 4112 );
218  m_PDGs.push_back( 4212 );
219  m_PDGs.push_back( 4232 );
220  m_PDGs.push_back( 4222 );
221  m_PDGs.push_back( 4322 );
222  m_PDGs.push_back( 4312 );
223 
224  m_PDGs.push_back( 13122 );
225  m_PDGs.push_back( 13124 );
226  m_PDGs.push_back( 23122 );
227  m_PDGs.push_back( 33122 );
228  m_PDGs.push_back( 43122 );
229  m_PDGs.push_back( 53122 );
230  m_PDGs.push_back( 13126 );
231  m_PDGs.push_back( 13212 );
232  m_PDGs.push_back( 13241 );
233 
234  m_PDGs.push_back( 3126 );
235  m_PDGs.push_back( 3124 );
236  m_PDGs.push_back( 3122 );
237  m_PDGs.push_back( 3222 );
238  m_PDGs.push_back( 2214 );
239  m_PDGs.push_back( 2224 );
240  m_PDGs.push_back( 3324 );
241  m_PDGs.push_back( 2114 );
242  m_PDGs.push_back( 1114 );
243  m_PDGs.push_back( 3112 );
244  m_PDGs.push_back( 3212 );
245  m_PDGs.push_back( 3114 );
246  m_PDGs.push_back( 3224 );
247  m_PDGs.push_back( 3214 );
248  m_PDGs.push_back( 3216 );
249  m_PDGs.push_back( 3322 );
250  m_PDGs.push_back( 3312 );
251  m_PDGs.push_back( 3314 );
252  m_PDGs.push_back( 3334 );
253 
254  m_PDGs.push_back( 4114 );
255  m_PDGs.push_back( 4214 );
256  m_PDGs.push_back( 4224 );
257  m_PDGs.push_back( 4314 );
258  m_PDGs.push_back( 4324 );
259  m_PDGs.push_back( 4332 );
260  m_PDGs.push_back( 4334 );
261  //m_PDGs.push_back( 43 ); // obsolete (?)
262  //m_PDGs.push_back( 44 ); // obsolete (?)
263  m_PDGs.push_back( 10443 );
264 
265  m_PDGs.push_back( 5122 );
266  m_PDGs.push_back( 5132 );
267  m_PDGs.push_back( 5232 );
268  m_PDGs.push_back( 5332 );
269  m_PDGs.push_back( 5222 );
270  m_PDGs.push_back( 5112 );
271  m_PDGs.push_back( 5212 );
272  m_PDGs.push_back( 541 );
273  m_PDGs.push_back( 14122 );
274  m_PDGs.push_back( 14124 );
275  m_PDGs.push_back( 5312 );
276  m_PDGs.push_back( 5322 );
277  m_PDGs.push_back( 10521 );
278  m_PDGs.push_back( 20523 );
279  m_PDGs.push_back( 10523 );
280 
281  m_PDGs.push_back( 525 );
282  m_PDGs.push_back( 10511 );
283  m_PDGs.push_back( 20513 );
284  m_PDGs.push_back( 10513 );
285  m_PDGs.push_back( 515 );
286  m_PDGs.push_back( 10531 );
287  m_PDGs.push_back( 20533 );
288  m_PDGs.push_back( 10533 );
289  m_PDGs.push_back( 535 );
290  m_PDGs.push_back( 543 );
291  m_PDGs.push_back( 545 );
292  m_PDGs.push_back( 5114 );
293  m_PDGs.push_back( 5224 );
294  m_PDGs.push_back( 5214 );
295  m_PDGs.push_back( 5314 );
296  m_PDGs.push_back( 5324 );
297  m_PDGs.push_back( 5334 );
298  m_PDGs.push_back( 10541 );
299  m_PDGs.push_back( 10543 );
300  m_PDGs.push_back( 20543 );
301 
302  m_PDGs.push_back( 4424 );
303  m_PDGs.push_back( 4422 );
304  m_PDGs.push_back( 4414 );
305  m_PDGs.push_back( 4412 );
306  m_PDGs.push_back( 4432 );
307  m_PDGs.push_back( 4434 );
308 
309  m_PDGs.push_back( 130 );
310 
311  // now check if we need to override default list of particles/IDs
312  if ( pset.exists("operates_on_particles") )
313  {
314  std::vector<int> tmpPIDs = pset.getParameter< std::vector<int> >("operates_on_particles");
315  if ( tmpPIDs.size() > 0 )
316  {
317  if ( tmpPIDs[0] > 0 ) // 0 means default !!!
318  {
319  m_PDGs.clear();
320  m_PDGs = tmpPIDs;
321  }
322  }
323  }
324  m_Py6Service = new Pythia6Service;
325 }
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
gen::EvtGenInterface::forced_names
std::vector< std::string > forced_names
Definition: EvtGenInterface.h:94
edm::ParameterSet::exists
bool exists(std::string const &parameterName) const
checks if a parameter exists
Definition: ParameterSet.cc:759
gen::EvtGenInterface::polarize_ids
std::vector< int > polarize_ids
Definition: EvtGenInterface.h:80
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:67
gen::EvtGenInterface::decay_table_s
std::string decay_table_s
Definition: EvtGenInterface.h:91
gen::EvtGenInterface::polarizations
std::map< int, float > polarizations
Definition: EvtGenInterface.h:82
gen::EvtGenInterface::polarize_pol
std::vector< double > polarize_pol
Definition: EvtGenInterface.h:81
gen::EvtGenInterface::the_engine
myEvtRandomEngine * the_engine
Definition: EvtGenInterface.h:88
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
edm::FileInPath::fullPath
std::string fullPath() const
Definition: FileInPath.cc:165
gen::EvtGenInterface::user_decay_s
std::string user_decay_s
Definition: EvtGenInterface.h:93
EvtGenInterface::~EvtGenInterface ( )

Definition at line 327 of file EvtGenInterface.cc.

References gather_cfg::cout.

328 {
329  std::cout << " EvtGenProducer terminating ... " << std::endl;
330  delete m_Py6Service;
331 }
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:67
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121

Member Function Documentation

void EvtGenInterface::addToHepMC ( HepMC::GenParticle *  partHep,
EvtId  idEvt,
HepMC::GenEvent *  theEvent,
bool  del_daug 
)

Definition at line 451 of file EvtGenInterface.cc.

References gather_cfg::cout, Vector3DBase< T, FrameTag >::cross(), configurableAnalysis::GenParticle, groupFilesInBlocks::ntotal, edm::second(), sistrip::STRING, Vector3DBase< T, FrameTag >::unit(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

452 {
453  // Set spin type
454  EvtSpinType::spintype stype = EvtPDL::getSpinType(idEvt);
455  EvtParticle* partEvt;
456  switch (stype){
457  case EvtSpinType::SCALAR:
458  partEvt = new EvtScalarParticle();
459  break;
460  case EvtSpinType::STRING:
461  partEvt = new EvtStringParticle();
462  break;
463  case EvtSpinType::DIRAC:
464  partEvt = new EvtDiracParticle();
465  break;
466  case EvtSpinType::VECTOR:
467  partEvt = new EvtVectorParticle();
468  break;
469  case EvtSpinType::RARITASCHWINGER:
470  partEvt = new EvtRaritaSchwingerParticle();
471  break;
472  case EvtSpinType::TENSOR:
473  partEvt = new EvtTensorParticle();
474  break;
475  case EvtSpinType::SPIN5HALF: case EvtSpinType::SPIN3: case EvtSpinType::SPIN7HALF: case EvtSpinType::SPIN4:
476  partEvt = new EvtHighSpinParticle();
477  break;
478  default:
479  std::cout << "Unknown spintype in EvtSpinType!" << std::endl;
480  return;
481  }
482 
483  // Generate decay
484  EvtVector4R momEvt;
485  HepMC::FourVector momHep = partHep->momentum();
486  momEvt.set(momHep.t(),momHep.x(),momHep.y(),momHep.z());
487  EvtVector4R posEvt;
488  HepMC::GenVertex* initVert = partHep->production_vertex();
489  HepMC::FourVector posHep = initVert->position();
490  posEvt.set(posHep.t(),posHep.x(),posHep.y(),posHep.z());
491  partEvt->init(idEvt,momEvt);
492  if (stype == EvtSpinType::DIRAC
493  && polarizations.find(partHep->pdg_id()) != polarizations.end()) {
494  // std::cout << "Polarize particle" << std::endl;
495  //Particle is spin 1/2, so we can polarize it.
496  //Check polarizations map for particle, grab its polarization if it exists
497  // and make the spin density matrix
498  float pol = polarizations.find(partHep->pdg_id())->second;
499  GlobalVector pPart(momHep.x(), momHep.y(), momHep.z());
500  //std::cout << "Polarizing particle with PDG ID "
501  // << partHep->pdg_id()
502  // << " at " << pol*100 << "%" << std::endl;
503  GlobalVector zHat(0., 0., 1.);
504  GlobalVector zCrossP = zHat.cross(pPart);
505  GlobalVector polVec = pol * zCrossP.unit();
506 
507  EvtSpinDensity theSpinDensity;
508  theSpinDensity.SetDim(2);
509  theSpinDensity.Set(0, 0, EvtComplex(1./2. + polVec.z()/2., 0.));
510  theSpinDensity.Set(0, 1, EvtComplex(polVec.x()/2., -polVec.y()/2.));
511  theSpinDensity.Set(1, 0, EvtComplex(polVec.x()/2., polVec.y()/2.));
512  theSpinDensity.Set(1, 1, EvtComplex(1./2. - polVec.z()/2., 0.));
513 
514  partEvt->setSpinDensityForwardHelicityBasis(theSpinDensity);
515 
516  } else {
517  partEvt->setDiagonalSpinDensity();
518  }
519  partEvt->decay();
520 
521  // extend the search of candidates to be forced to EvtGen decay products and delete their daughters **
522  // otherwise they wouldn't get their chance to take part in the forced decay lottery **
523  if (del_daug) go_through_daughters(partEvt); // recursive function go_through_daughters will do **
524 
525  // Change particle in stdHEP format
526  static EvtStdHep evtstdhep;
527 
528  evtstdhep.init();
529  partEvt->makeStdHep(evtstdhep);
530 
531  ntotal++;
532  partEvt->deleteTree();
533 
534  // ********* Now add to the HepMC Event **********
535 
536  // Then loop on evtstdhep to add vertexes...
537  HepMC::GenVertex* theVerts[200];
538  for (int ivert = 0; ivert < 200; ivert++) {
539  theVerts[ivert] = 0;
540  }
541 
542  for (int ipart = 0; ipart < evtstdhep.getNPart(); ipart++) {
543  int theMum = evtstdhep.getFirstMother(ipart);
544  if (theMum != -1 && !theVerts[theMum]) {
545  EvtVector4R theVpos = evtstdhep.getX4(ipart) + posEvt;
546  theVerts[theMum] =
547  new HepMC::GenVertex(HepMC::FourVector(theVpos.get(1),
548  theVpos.get(2),
549  theVpos.get(3),
550  theVpos.get(0)),0);
551  }
552  }
553 
554  // ...then particles
555  partHep->set_status(2);
556  if (theVerts[0]) theVerts[0]->add_particle_in( partHep );
557 
558  for (int ipart2 = 1; ipart2 < evtstdhep.getNPart(); ipart2++) {
559  int idHep = evtstdhep.getStdHepID(ipart2);
560  HepMC::GenParticle* thePart =
561  new HepMC::GenParticle( HepMC::FourVector(evtstdhep.getP4(ipart2).get(1),
562  evtstdhep.getP4(ipart2).get(2),
563  evtstdhep.getP4(ipart2).get(3),
564  evtstdhep.getP4(ipart2).get(0)),
565  idHep,
566  evtstdhep.getIStat(ipart2));
567  npartial++;
568  thePart->suggest_barcode(npartial + nPythia);
569  int theMum2 = evtstdhep.getFirstMother(ipart2);
570  if (theMum2 != -1 && theVerts[theMum2]) theVerts[theMum2]->add_particle_out( thePart );
571  if (theVerts[ipart2]) theVerts[ipart2]->add_particle_in( thePart );
572 
573  }
574 
575  for (int ipart3 = 0; ipart3 < evtstdhep.getNPart(); ipart3++) {
576  if (theVerts[ipart3]) theEvent->add_vertex( theVerts[ipart3] );
577  }
578 
579 }
PV3DBase::y
T y() const
Definition: PV3DBase.h:63
gen::EvtGenInterface::go_through_daughters
void go_through_daughters(EvtParticle *part)
Definition: EvtGenInterface.cc:592
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:250
Vector3DBase::cross
Vector3DBase< typename PreciseFloatType< T, U >::Type, FrameTag > cross(const Vector3DBase< U, FrameTag > &v) const
Definition: Vector3DBase.h:119
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
gen::EvtGenInterface::polarizations
std::map< int, float > polarizations
Definition: EvtGenInterface.h:82
Vector3DBase::unit
Vector3DBase unit() const
Definition: Vector3DBase.h:57
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
PV3DBase::x
T x() const
Definition: PV3DBase.h:62
HepMC::GenEvent * EvtGenInterface::decay ( HepMC::GenEvent *  evt)
virtual

Reimplemented from gen::EvtGenInterfaceBase.

Definition at line 370 of file EvtGenInterface.cc.

References edm::hlt::Exception, i, cmsHarvester::index, gen::k, edm::errors::LogicError, nevent, gen::p, and ntuplemaker::status.

371 {
372  if(the_engine->engine() == nullptr) {
373  throw edm::Exception(edm::errors::LogicError)
374  << "The EvtGen code attempted to use a random number engine while\n"
375  << "the engine pointer was null in EvtGenInterface::decay. This might\n"
376  << "mean that the code was modified to generate a random number outside\n"
377  << "the event and beginLuminosityBlock methods, which is not allowed.\n";
378  }
379  CLHEP::RandFlat m_flat(*the_engine->engine(), 0., 1.);
380 
381  Pythia6Service::InstanceWrapper guard(m_Py6Service); // grab Py6 instance
382 
383  nevent++;
384  npartial = 0;
385  // std::cout << "nevent = " << nevent << std::endl ;
386 
387  int idHep,ipart,status;
388  EvtId idEvt;
389 
390  nPythia = evt->particles_size();
391  // FIX A MEMORY LEAK (RC)
392  // HepMC::GenEvent* newEvt = new HepMC::GenEvent( *evt );
393 
394  // First pass through undecayed Pythia particles to decay particles known to EvtGen left stable by Pythia
395  // except candidates to be forced which will be searched later to include EvtGen decay products
396  nlist = 0;
397 
398  // Notice dynamical use of evt
399  for (HepMC::GenEvent::particle_const_iterator p= evt->particles_begin(); p != evt->particles_end(); ++p)
400  {
401  status = (*p)->status();
402 
403  if(status==1) { // only not decayed (status = 1) particles
404 
405 
406  idHep = (*p)->pdg_id();
407  int do_force=0;
408  for(int i=0;i<nforced;i++) // First check if part with forced decay
409  { // In that case do not decay immediately
410  if(idHep == forced_Hep[i]) // (only 1 per event will be forced)
411  { // Fill list
412  update_candlist(i,*p);
413  do_force=1;
414  }
415  }
416  if(do_force==0) // particles with decays not forced are decayed immediately
417  {
418  idEvt = EvtPDL::evtIdFromStdHep(idHep);
419  ipart = idEvt.getId();
420  if (ipart==-1) continue; // particle not known to EvtGen
421  if (EvtDecayTable::getNMode(ipart)==0) continue; // particles stable for EvtGen
422  addToHepMC(*p,idEvt,evt,true); // generate decay
423  }
424  }
425  }
426 
427  if(nlist!=0)
428  {
429  // decide randomly which one to decay as alias
430  int which = (int)(nlist * m_flat.fire());
431  if (which == nlist) which = nlist-1;
432 
433  for(int k=0;k < nlist; k++)
434  {
435  if(k == which || !usePythia) {
436  addToHepMC(listp[k],forced_Evt[index[k]],evt,false); // decay as alias
437  }
438  else
439  {
440  int id_non_alias = forced_Evt[index[k]].getId();
441  EvtId non_alias(id_non_alias,id_non_alias); // create new EvtId with id = alias
442  addToHepMC(listp[k],non_alias,evt,false); // decay as standard (non alias)
443  }
444  }
445  }
446 
447  return evt;
448 
449 }
i
int i
Definition: DBlmapReader.cc:9
gen::EvtGenInterface::addToHepMC
void addToHepMC(HepMC::GenParticle *partHep, EvtId idEvt, HepMC::GenEvent *theEvent, bool del_daug)
Definition: EvtGenInterface.cc:451
gen::EvtGenInterface::update_candlist
void update_candlist(int theIndex, HepMC::GenParticle *thePart)
Definition: EvtGenInterface.cc:617
gen::EvtGenInterface::forced_Evt
std::vector< EvtId > forced_Evt
Definition: EvtGenInterface.h:70
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:67
gen::p
double p[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:80
myEvtRandomEngine::engine
CLHEP::HepRandomEngine * engine() const
Definition: myEvtRandomEngine.h:40
gen::EvtGenInterface::listp
HepMC::GenParticle * listp[10]
Definition: EvtGenInterface.h:85
gen::k
int k[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:79
gen::EvtGenInterface::the_engine
myEvtRandomEngine * the_engine
Definition: EvtGenInterface.h:88
gen::EvtGenInterface::forced_Hep
std::vector< int > forced_Hep
Definition: EvtGenInterface.h:71
ntuplemaker.status
tuple status
Definition: ntuplemaker.py:245
double EvtGenInterface::flat ( )
static

Definition at line 648 of file EvtGenInterface.cc.

References edm::hlt::Exception.

Referenced by phoran_().

648  {
649  if ( !fRandomEngine ) {
650  throw cms::Exception("LogicError")
651  << "EvtGenInterface::flat: Attempt to generate random number when engine pointer is null\n"
652  << "This might mean that the code was modified to generate a random number outside the\n"
653  << "event and beginLuminosityBlock methods, which is not allowed.\n";
654  }
655  return fRandomEngine->flat();
656 }
gen::EvtGenInterface::fRandomEngine
static CLHEP::HepRandomEngine * fRandomEngine
Definition: EvtGenInterface.h:96
void EvtGenInterface::go_through_daughters ( EvtParticle *  part)

Definition at line 592 of file EvtGenInterface.cc.

References newFWLiteAna::found, i, and gen::k.

592  {
593 
594  int NDaug=part->getNDaug();
595  if(NDaug)
596  {
597  EvtParticle* Daughter;
598  for (int i=0;i<NDaug;i++)
599  {
600  Daughter=part->getDaug(i);
601  int idHep = EvtPDL::getStdHep(Daughter->getId());
602  int found=0;
603  for(int k=0;k<nforced;k++)
604  {
605  if(idHep == forced_Hep[k])
606  {
607  found = 1;
608  Daughter->deleteDaughters();
609  }
610  }
611  if (!found) go_through_daughters(Daughter);
612  }
613  }
614 }
i
int i
Definition: DBlmapReader.cc:9
gen::EvtGenInterface::go_through_daughters
void go_through_daughters(EvtParticle *part)
Definition: EvtGenInterface.cc:592
gen::k
int k[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:79
part
part
Definition: HCALResponse.h:20
gen::EvtGenInterface::forced_Hep
std::vector< int > forced_Hep
Definition: EvtGenInterface.h:71
void EvtGenInterface::init ( void  )
virtual

Reimplemented from gen::EvtGenInterfaceBase.

Definition at line 333 of file EvtGenInterface.cc.

References edm::hlt::Exception, newFWLiteAna::found, and getId().

334 {
335  m_EvtGen = new EvtGen (decay_table_s.c_str(), pdt_s.c_str(), the_engine);
336  // 4th parameter should be rad cor - set to PHOTOS (default)
337 
338  if (!useDefault) m_EvtGen->readUDecay( user_decay_s.c_str() );
339 
340  nforced = 0;
341  for (auto const& forced_name : forced_names) {
342  ++nforced;
343  EvtId found = EvtPDL::getId(forced_name);
344  if (found.getId() == -1) {
345  throw cms::Exception("Configuration")
346  << "name in part list for forced decays not found: " << forced_name;
347  }
348  if (found.getId() == found.getAlias()) {
349  throw cms::Exception("Configuration")
350  << "name in part list for forced decays is not an alias: " << forced_name;
351  }
352  forced_Evt.push_back(found); // forced_Evt is the list of EvtId's
353  forced_Hep.push_back(EvtPDL::getStdHep(found)); // forced_Hep is the list of stdhep codes
354  }
355 
356  Pythia6Service::InstanceWrapper guard(m_Py6Service); // grab Py6 instance
357 
358  // Do here initialization of EvtPythia then restore original settings
359  if (usePythia) EvtPythia::pythiaInit(0);
360 
361 // no need - will be done via Pythia6Hadronizer::declareStableParticles
362 //
363 // for( std::vector<std::string>::const_iterator itPar = pythia_params.begin(); itPar != pythia_params.end(); ++itPar ) {
364 // call_pygive(*itPar);
365 // }
366 
367  return;
368 }
gen::EvtGenInterface::forced_names
std::vector< std::string > forced_names
Definition: EvtGenInterface.h:94
getId
static unsigned int getId(void)
Definition: DQMStoreStats.h:110
gen::EvtGenInterface::forced_Evt
std::vector< EvtId > forced_Evt
Definition: EvtGenInterface.h:70
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:67
gen::EvtGenInterface::decay_table_s
std::string decay_table_s
Definition: EvtGenInterface.h:91
gen::EvtGenInterface::the_engine
myEvtRandomEngine * the_engine
Definition: EvtGenInterface.h:88
gen::EvtGenInterface::forced_Hep
std::vector< int > forced_Hep
Definition: EvtGenInterface.h:71
gen::EvtGenInterface::user_decay_s
std::string user_decay_s
Definition: EvtGenInterface.h:93
const std::vector<int>& gen::EvtGenInterface::operatesOnParticles ( )
inlinevirtual

Reimplemented from gen::EvtGenInterfaceBase.

Definition at line 54 of file EvtGenInterface.h.

References gen::EvtGenInterfaceBase::m_PDGs.

54 { return m_PDGs; }
void EvtGenInterface::setRandomEngine ( CLHEP::HepRandomEngine *  v)
virtual

Implements gen::EvtGenInterfaceBase.

Definition at line 642 of file EvtGenInterface.cc.

References gen::v.

642  {
643  the_engine->setRandomEngine(v);
644  fRandomEngine=v;
645  m_Py6Service->setRandomEngine(v);
646 }
gen::v
double v[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:80
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:67
gen::EvtGenInterface::fRandomEngine
static CLHEP::HepRandomEngine * fRandomEngine
Definition: EvtGenInterface.h:96
gen::EvtGenInterface::the_engine
myEvtRandomEngine * the_engine
Definition: EvtGenInterface.h:88
gen::Pythia6Service::setRandomEngine
void setRandomEngine(CLHEP::HepRandomEngine *v)
Definition: Pythia6Service.h:48
myEvtRandomEngine::setRandomEngine
void setRandomEngine(CLHEP::HepRandomEngine *v)
Definition: myEvtRandomEngine.h:38
void EvtGenInterface::update_candlist ( int  theIndex,
HepMC::GenParticle *  thePart 
)

Definition at line 617 of file EvtGenInterface.cc.

References CastorDataFrameFilter_impl::check(), edm::hlt::Exception, cmsHarvester::index, and theIndex().

618 {
619  if(nlist<10) // not nice ... but is 10 a reasonable maximum?
620  {
621  bool isThere = false;
622  if (nlist) {
623  for (int check=0; check < nlist; check++) {
624  if (listp[check]->barcode() == thePart->barcode()) isThere = true;
625  }
626  }
627  if (!isThere) {
628  listp[nlist] = thePart;
629  index[nlist++] = theIndex;
630  }
631  }
632  else
633  {
634  throw cms::Exception("runtime")
635  << "more than 10 candidates to be forced ";
636  return;
637  }
638  return;
639 }
theIndex
theIndex(0)
CastorDataFrameFilter_impl::check
bool check(const DataFrame &df, bool capcheck, bool dvercheck)
Definition: CastorDataFrameFilter.cc:6
gen::EvtGenInterface::listp
HepMC::GenParticle * listp[10]
Definition: EvtGenInterface.h:85

Member Data Documentation

std::string gen::EvtGenInterface::decay_table_s
private

Definition at line 91 of file EvtGenInterface.h.

std::vector<EvtId> gen::EvtGenInterface::forced_Evt
private

Definition at line 70 of file EvtGenInterface.h.

std::vector<int> gen::EvtGenInterface::forced_Hep
private

Definition at line 71 of file EvtGenInterface.h.

std::vector<std::string> gen::EvtGenInterface::forced_names
private

Definition at line 94 of file EvtGenInterface.h.

CLHEP::HepRandomEngine * EvtGenInterface::fRandomEngine
staticprivate

Definition at line 96 of file EvtGenInterface.h.

int gen::EvtGenInterface::index[10]
private

Definition at line 86 of file EvtGenInterface.h.

Referenced by BeautifulSoup.PageElement::_invert().

HepMC::GenParticle* gen::EvtGenInterface::listp[10]
private

Definition at line 85 of file EvtGenInterface.h.

EvtGen* gen::EvtGenInterface::m_EvtGen
private

Definition at line 69 of file EvtGenInterface.h.

Pythia6Service* gen::EvtGenInterface::m_Py6Service
private

Definition at line 67 of file EvtGenInterface.h.

int gen::EvtGenInterface::nevent
private

Definition at line 73 of file EvtGenInterface.h.

int gen::EvtGenInterface::nforced
private

Definition at line 72 of file EvtGenInterface.h.

int gen::EvtGenInterface::nlist
private

Definition at line 84 of file EvtGenInterface.h.

int gen::EvtGenInterface::npartial
private

Definition at line 73 of file EvtGenInterface.h.

int gen::EvtGenInterface::nPythia
private

Definition at line 75 of file EvtGenInterface.h.

int gen::EvtGenInterface::ntotal
private

Definition at line 73 of file EvtGenInterface.h.

std::string gen::EvtGenInterface::pdt_s
private

Definition at line 92 of file EvtGenInterface.h.

std::map<int, float> gen::EvtGenInterface::polarizations
private

Definition at line 82 of file EvtGenInterface.h.

std::vector<int> gen::EvtGenInterface::polarize_ids
private

Definition at line 80 of file EvtGenInterface.h.

std::vector<double> gen::EvtGenInterface::polarize_pol
private

Definition at line 81 of file EvtGenInterface.h.

myEvtRandomEngine* gen::EvtGenInterface::the_engine
private

Definition at line 88 of file EvtGenInterface.h.

bool gen::EvtGenInterface::useDefault
private

Definition at line 90 of file EvtGenInterface.h.

bool gen::EvtGenInterface::usePythia
private

Definition at line 76 of file EvtGenInterface.h.

std::string gen::EvtGenInterface::user_decay_s
private

Definition at line 93 of file EvtGenInterface.h.

Generated for CMSSW Reference Manual by doxygen 1.8.5