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

#include <EvtGenInterface.h>

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 update_candlist (int theIndex, HepMC::GenParticle *thePart)
 
 ~EvtGenInterface ()
 

Private Attributes

std::vector< EvtId > forced_Evt
 
std::vector< int > forced_Hep
 
int index [10]
 
HepMC::GenParticle * listp [10]
 
EvtGen * m_EvtGen
 
CLHEP::RandFlat * m_flat
 
std::vector< int > m_PDGs
 
Pythia6Service * m_Py6Service
 
int nevent
 
int nforced
 
int nlist
 
int npartial
 
int nPythia
 
int ntotal
 
std::map< int, float > polarizations
 
std::vector< int > polarize_ids
 
std::vector< double > polarize_pol
 
bool usePythia
 

Detailed Description

Definition at line 48 of file EvtGenInterface.h.

Constructor & Destructor Documentation

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

Definition at line 35 of file EvtGenInterface.cc.

References gather_cfg::cout, edm::hlt::Exception, edm::ParameterSet::exists(), newFWLiteAna::found, edm::FileInPath::fullPath(), edm::RandomNumberGenerator::getEngine(), getId(), edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), i, edm::Service< T >::isAvailable(), nevent, groupFilesInBlocks::ntotal, and AlCaHLTBitMon_QueryRunRegistry::string.

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

Definition at line 357 of file EvtGenInterface.cc.

References gather_cfg::cout.

358 {
359  std::cout << " EvtGenProducer terminating ... " << std::endl;
360  delete m_Py6Service;
361 }
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:68
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 454 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().

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

Definition at line 382 of file EvtGenInterface.cc.

References i, getHLTprescales::index, gen::k, nevent, gen::p, and ntuplemaker::status.

Referenced by gen::ExternalDecayDriver::decay().

383 {
384  Pythia6Service::InstanceWrapper guard(m_Py6Service); // grab Py6 instance
385 
386  nevent++;
387  npartial = 0;
388  // std::cout << "nevent = " << nevent << std::endl ;
389 
390  int idHep,ipart,status;
391  EvtId idEvt;
392 
393  nPythia = evt->particles_size();
394  // FIX A MEMORY LEAK (RC)
395  // HepMC::GenEvent* newEvt = new HepMC::GenEvent( *evt );
396 
397  // First pass through undecayed Pythia particles to decay particles known to EvtGen left stable by Pythia
398  // except candidates to be forced which will be searched later to include EvtGen decay products
399  nlist = 0;
400 
401  // Notice dynamical use of evt
402  for (HepMC::GenEvent::particle_const_iterator p= evt->particles_begin(); p != evt->particles_end(); ++p)
403  {
404  status = (*p)->status();
405 
406  if(status==1) { // only not decayed (status = 1) particles
407 
408 
409  idHep = (*p)->pdg_id();
410  int do_force=0;
411  for(int i=0;i<nforced;i++) // First check if part with forced decay
412  { // In that case do not decay immediately
413  if(idHep == forced_Hep[i]) // (only 1 per event will be forced)
414  { // Fill list
415  update_candlist(i,*p);
416  do_force=1;
417  }
418  }
419  if(do_force==0) // particles with decays not forced are decayed immediately
420  {
421  idEvt = EvtPDL::evtIdFromStdHep(idHep);
422  ipart = idEvt.getId();
423  if (ipart==-1) continue; // particle not known to EvtGen
424  if (EvtDecayTable::getNMode(ipart)==0) continue; // particles stable for EvtGen
425  addToHepMC(*p,idEvt,evt,true); // generate decay
426  }
427  }
428  }
429 
430  if(nlist!=0)
431  {
432  // decide randomly which one to decay as alias
433  int which = (int)(nlist*m_flat->fire());
434  if (which == nlist) which = nlist-1;
435 
436  for(int k=0;k < nlist; k++)
437  {
438  if(k == which || !usePythia) {
439  addToHepMC(listp[k],forced_Evt[index[k]],evt,false); // decay as alias
440  }
441  else
442  {
443  int id_non_alias = forced_Evt[index[k]].getId();
444  EvtId non_alias(id_non_alias,id_non_alias); // create new EvtId with id = alias
445  addToHepMC(listp[k],non_alias,evt,false); // decay as standard (non alias)
446  }
447  }
448  }
449 
450  return evt;
451 
452 }
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:454
gen::EvtGenInterface::update_candlist
void update_candlist(int theIndex, HepMC::GenParticle *thePart)
Definition: EvtGenInterface.cc:620
gen::EvtGenInterface::forced_Evt
std::vector< EvtId > forced_Evt
Definition: EvtGenInterface.h:74
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:68
gen::p
double p[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:82
gen::EvtGenInterface::listp
HepMC::GenParticle * listp[10]
Definition: EvtGenInterface.h:89
gen::k
int k[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:81
gen::EvtGenInterface::m_flat
CLHEP::RandFlat * m_flat
Definition: EvtGenInterface.h:72
gen::EvtGenInterface::forced_Hep
std::vector< int > forced_Hep
Definition: EvtGenInterface.h:75
ntuplemaker.status
tuple status
Definition: ntuplemaker.py:245
void EvtGenInterface::go_through_daughters ( EvtParticle *  part)

Definition at line 595 of file EvtGenInterface.cc.

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

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

Definition at line 363 of file EvtGenInterface.cc.

References hitfit::return.

Referenced by gen::ExternalDecayDriver::init().

364 {
365 
366  Pythia6Service::InstanceWrapper guard(m_Py6Service); // grab Py6 instance
367 
368  // Do here initialization of EvtPythia then restore original settings
369  if (usePythia) EvtPythia::pythiaInit(0);
370 
371 // no need - will be done via Pythia6Hadronizer::declareStableParticles
372 //
373 // for( std::vector<std::string>::const_iterator itPar = pythia_params.begin(); itPar != pythia_params.end(); ++itPar ) {
374 // call_pygive(*itPar);
375 // }
376 
377  return ;
378 
379 }
gen::EvtGenInterface::m_Py6Service
Pythia6Service * m_Py6Service
Definition: EvtGenInterface.h:68
const std::vector<int>& gen::EvtGenInterface::operatesOnParticles ( )
inline

Definition at line 57 of file EvtGenInterface.h.

References m_PDGs.

Referenced by gen::ExternalDecayDriver::init().

57 { return m_PDGs; }
gen::EvtGenInterface::m_PDGs
std::vector< int > m_PDGs
Definition: EvtGenInterface.h:70
void EvtGenInterface::update_candlist ( int  theIndex,
HepMC::GenParticle *  thePart 
)

Definition at line 620 of file EvtGenInterface.cc.

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

621 {
622  if(nlist<10) // not nice ... but is 10 a reasonable maximum?
623  {
624  bool isThere = false;
625  if (nlist) {
626  for (int check=0; check < nlist; check++) {
627  if (listp[check]->barcode() == thePart->barcode()) isThere = true;
628  }
629  }
630  if (!isThere) {
631  listp[nlist] = thePart;
632  index[nlist++] = theIndex;
633  }
634  }
635  else
636  {
637  throw cms::Exception("runtime")
638  << "more than 10 candidates to be forced ";
639  return;
640  }
641  return;
642 }
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:89

Member Data Documentation

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

Definition at line 74 of file EvtGenInterface.h.

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

Definition at line 75 of file EvtGenInterface.h.

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

Definition at line 90 of file EvtGenInterface.h.

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

Definition at line 89 of file EvtGenInterface.h.

EvtGen* gen::EvtGenInterface::m_EvtGen
private

Definition at line 73 of file EvtGenInterface.h.

CLHEP::RandFlat* gen::EvtGenInterface::m_flat
private

Definition at line 72 of file EvtGenInterface.h.

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

Definition at line 70 of file EvtGenInterface.h.

Referenced by operatesOnParticles().

Pythia6Service* gen::EvtGenInterface::m_Py6Service
private

Definition at line 68 of file EvtGenInterface.h.

int gen::EvtGenInterface::nevent
private

Definition at line 77 of file EvtGenInterface.h.

int gen::EvtGenInterface::nforced
private

Definition at line 76 of file EvtGenInterface.h.

int gen::EvtGenInterface::nlist
private

Definition at line 88 of file EvtGenInterface.h.

int gen::EvtGenInterface::npartial
private

Definition at line 77 of file EvtGenInterface.h.

int gen::EvtGenInterface::nPythia
private

Definition at line 79 of file EvtGenInterface.h.

int gen::EvtGenInterface::ntotal
private

Definition at line 77 of file EvtGenInterface.h.

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

Definition at line 86 of file EvtGenInterface.h.

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

Definition at line 84 of file EvtGenInterface.h.

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

Definition at line 85 of file EvtGenInterface.h.

bool gen::EvtGenInterface::usePythia
private

Definition at line 80 of file EvtGenInterface.h.

Generated for CMSSW Reference Manual by doxygen 1.8.5