CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Private Attributes
Generator Class Reference

#include <Generator.h>

Public Member Functions

virtual const double eventWeight () const
 
 Generator (const edm::ParameterSet &p)
 
virtual const HepMC::GenEventgenEvent () const
 
virtual const math::XYZTLorentzVectorgenVertex () const
 
void HepMC2G4 (const HepMC::GenEvent *g, G4Event *e)
 
void nonCentralEvent2G4 (const HepMC::GenEvent *g, G4Event *e)
 
void setGenEvent (const HepMC::GenEvent *inpevt)
 
virtual ~Generator ()
 

Private Member Functions

bool isExotic (int pdgcode) const
 
bool isExoticNonDetectable (int pdgcode) const
 
bool IsInTheFilterList (int pdgcode) const
 
void particleAssignDaughters (G4PrimaryParticle *p, HepMC::GenParticle *hp, double length)
 
bool particlePassesPrimaryCuts (const G4ThreeVector &p) const
 
void setGenId (G4PrimaryParticle *p, int id) const
 

Private Attributes

HepMC::GenEventevt_
 
bool fEtaCuts
 
bool fFiductialCuts
 
LumiMonitorFilterfLumiFilter
 
bool fPCuts
 
bool fPDGFilter
 
bool fPhiCuts
 
bool fPtransCut
 
double maxZCentralCMS
 
std::vector< int > pdgFilter
 
bool pdgFilterSel
 
double theDecLenCut
 
double theDecRCut2
 
double theEtaCutForHector
 
double theMaxEtaCut
 
double theMaxPCut
 
double theMaxPhiCut
 
double theMinEtaCut
 
double theMinPCut
 
double theMinPhiCut
 
double theMinPtCut2
 
int verbose
 
math::XYZTLorentzVectorvtx_
 
double weight_
 
double Z_hector
 
double Z_lmax
 
double Z_lmin
 

Detailed Description

Definition at line 19 of file Generator.h.

Constructor & Destructor Documentation

◆ Generator()

Generator::Generator ( const edm::ParameterSet p)

Definition at line 23 of file Generator.cc.

References LumiMonitorFilter::Describe(), JetChargeProducer_cfi::exp, fEtaCuts, fFiductialCuts, fLumiFilter, fPCuts, fPDGFilter, fPhiCuts, fPtransCut, cuy::ii, BXlumiParameters_cfi::lumi, visualization-live-secondInstance_cfg::m, maxZCentralCMS, AlCaHLTBitMon_ParallelJobs::p, pdgFilter, pdgFilterSel, contentValuesCheck::ss, theDecLenCut, theDecRCut2, theEtaCutForHector, theMaxEtaCut, theMaxPCut, theMaxPhiCut, theMinEtaCut, theMinPCut, theMinPhiCut, theMinPtCut2, Z_hector, Z_lmax, and Z_lmin.

24  : fPCuts(p.getParameter<bool>("ApplyPCuts")),
25  fPtransCut(p.getParameter<bool>("ApplyPtransCut")),
26  fEtaCuts(p.getParameter<bool>("ApplyEtaCuts")),
27  fPhiCuts(p.getParameter<bool>("ApplyPhiCuts")),
28  theMinPhiCut(p.getParameter<double>("MinPhiCut")), // in radians (CMS standard)
29  theMaxPhiCut(p.getParameter<double>("MaxPhiCut")),
30  theMinEtaCut(p.getParameter<double>("MinEtaCut")),
31  theMaxEtaCut(p.getParameter<double>("MaxEtaCut")),
32  theMinPCut(p.getParameter<double>("MinPCut")), // in GeV (CMS standard)
33  theMaxPCut(p.getParameter<double>("MaxPCut")),
34  theEtaCutForHector(p.getParameter<double>("EtaCutForHector")),
35  verbose(p.getUntrackedParameter<int>("Verbosity", 0)),
36  fLumiFilter(nullptr),
37  evt_(nullptr),
38  vtx_(nullptr),
39  weight_(0),
40  Z_lmin(0),
41  Z_lmax(0),
42  Z_hector(0),
43  pdgFilterSel(false),
44  fPDGFilter(false) {
45  bool lumi = p.getParameter<bool>("ApplyLumiMonitorCuts");
46  if (lumi) {
48  }
49 
50  double theRDecLenCut = p.getParameter<double>("RDecLenCut") * CLHEP::cm;
51  theDecRCut2 = theRDecLenCut * theRDecLenCut;
52 
54 
55  double theDecLenCut = p.getParameter<double>("LDecLenCut") * CLHEP::cm;
56 
57  maxZCentralCMS = p.getParameter<double>("MaxZCentralCMS") * CLHEP::m;
58 
60 
61  pdgFilter.resize(0);
62  if (p.exists("PDGselection")) {
63  pdgFilterSel = (p.getParameter<edm::ParameterSet>("PDGselection")).getParameter<bool>("PDGfilterSel");
64  pdgFilter = (p.getParameter<edm::ParameterSet>("PDGselection")).getParameter<std::vector<int>>("PDGfilter");
65  if (!pdgFilter.empty()) {
66  fPDGFilter = true;
67  std::stringstream ss;
68  ss << "SimG4Core/Generator: ";
69  if (pdgFilterSel) {
70  ss << " Selecting only PDG ID = ";
71  } else {
72  ss << " Filtering out PDG ID = ";
73  }
74  for (unsigned int ii = 0; ii < pdgFilter.size(); ++ii) {
75  ss << pdgFilter[ii] << " ";
76  }
77  edm::LogVerbatim("SimG4CoreGenerator") << ss.str();
78  }
79  }
80 
81  if (fEtaCuts) {
82  Z_lmax = theRDecLenCut * ((1 - exp(-2 * theMaxEtaCut)) / (2 * exp(-theMaxEtaCut)));
83  Z_lmin = theRDecLenCut * ((1 - exp(-2 * theMinEtaCut)) / (2 * exp(-theMinEtaCut)));
84  }
85 
86  Z_hector = theRDecLenCut * ((1 - exp(-2 * theEtaCutForHector)) / (2 * exp(-theEtaCutForHector)));
87 
88  edm::LogVerbatim("SimG4CoreGenerator") << "SimG4Core/Generator: Rdecaycut= " << theRDecLenCut / CLHEP::cm
89  << " cm; Zdecaycut= " << theDecLenCut / CLHEP::cm
90  << "Z_min= " << Z_lmin / CLHEP::cm << " cm; Z_max= " << Z_lmax / CLHEP::cm
91  << " cm;\n"
92  << " MaxZCentralCMS = " << maxZCentralCMS / CLHEP::m
93  << " m;"
94  << " Z_hector = " << Z_hector / CLHEP::cm << " cm\n"
95  << " ApplyCuts: " << fFiductialCuts
96  << " PCuts: " << fPCuts << " PtransCut: " << fPtransCut
97  << " EtaCut: " << fEtaCuts << " PhiCut: " << fPhiCuts
98  << " LumiMonitorCut: " << lumi;
99  if (fFiductialCuts) {
100  edm::LogVerbatim("SimG4CoreGenerator")
101  << "SimG4Core/Generator: Pmin(GeV)= " << theMinPCut << "; Pmax(GeV)= " << theMaxPCut
102  << "; EtaMin= " << theMinEtaCut << "; EtaMax= " << theMaxEtaCut << "; PhiMin(rad)= " << theMinPhiCut
103  << "; PhiMax(rad)= " << theMaxPhiCut;
104  }
105  if (lumi) {
107  }
108 }
double theDecRCut2
Definition: Generator.h:55
Log< level::Info, true > LogVerbatim
bool fPhiCuts
Definition: Generator.h:46
double theMinPtCut2
Definition: Generator.h:53
double maxZCentralCMS
Definition: Generator.h:58
double theEtaCutForHector
Definition: Generator.h:56
double Z_lmax
Definition: Generator.h:64
double theMaxPhiCut
Definition: Generator.h:49
bool fPtransCut
Definition: Generator.h:44
double theDecLenCut
Definition: Generator.h:57
double theMinEtaCut
Definition: Generator.h:50
bool pdgFilterSel
Definition: Generator.h:66
double theMinPCut
Definition: Generator.h:52
LumiMonitorFilter * fLumiFilter
Definition: Generator.h:60
bool fFiductialCuts
Definition: Generator.h:47
bool fPDGFilter
Definition: Generator.h:67
double theMinPhiCut
Definition: Generator.h:48
double weight_
Definition: Generator.h:63
double theMaxEtaCut
Definition: Generator.h:51
ii
Definition: cuy.py:589
double theMaxPCut
Definition: Generator.h:54
math::XYZTLorentzVector * vtx_
Definition: Generator.h:62
HepMC::GenEvent * evt_
Definition: Generator.h:61
double Z_lmin
Definition: Generator.h:64
int verbose
Definition: Generator.h:59
std::vector< int > pdgFilter
Definition: Generator.h:65
bool fEtaCuts
Definition: Generator.h:45
double Z_hector
Definition: Generator.h:64
bool fPCuts
Definition: Generator.h:43

◆ ~Generator()

Generator::~Generator ( )
virtual

Definition at line 110 of file Generator.cc.

References fLumiFilter.

110 { delete fLumiFilter; }
LumiMonitorFilter * fLumiFilter
Definition: Generator.h:60

Member Function Documentation

◆ eventWeight()

virtual const double Generator::eventWeight ( ) const
inlinevirtual

Definition at line 32 of file Generator.h.

References weight_.

Referenced by RunManagerMTWorker::produce().

32 { return weight_; }
double weight_
Definition: Generator.h:63

◆ genEvent()

virtual const HepMC::GenEvent* Generator::genEvent ( ) const
inlinevirtual

Definition at line 30 of file Generator.h.

References evt_.

Referenced by RunManagerMTWorker::produce().

30 { return evt_; }
HepMC::GenEvent * evt_
Definition: Generator.h:61

◆ genVertex()

virtual const math::XYZTLorentzVector* Generator::genVertex ( ) const
inlinevirtual

Definition at line 31 of file Generator.h.

References vtx_.

Referenced by RunManagerMTWorker::produce().

31 { return vtx_; }
math::XYZTLorentzVector * vtx_
Definition: Generator.h:62

◆ HepMC2G4()

void Generator::HepMC2G4 ( const HepMC::GenEvent g,
G4Event *  e 
)

Definition at line 112 of file Generator.cc.

References funct::abs(), ALCARECOTkAlJpsiMuMu_cff::charge, class-composition::children, fEtaCuts, fFiductialCuts, fLumiFilter, fPCuts, fPDGFilter, fPhiCuts, fPtransCut, GenParticle::GenParticle, isExotic(), isExoticNonDetectable(), LumiMonitorFilter::isGoodForLumiMonitor(), IsInTheFilterList(), LogDebug, visualization-live-secondInstance_cfg::m, maxZCentralCMS, AlCaHLTBitMon_ParallelJobs::p, particleAssignDaughters(), dileptonTrigSettings_cff::pdg, pdgFilterSel, phi, multPhiCorr_741_25nsDY_cfi::px, multPhiCorr_741_25nsDY_cfi::py, setGenId(), mathSSE::sqrt(), contentValuesCheck::ss, mps_update::status, RandomServiceHelper::t1, theDecRCut2, theMaxPCut, theMaxPhiCut, theMinPtCut2, verbose, vtx_, weight_, testProducerWithPsetDescEmpty_cfi::x1, testProducerWithPsetDescEmpty_cfi::x2, protons_cff::xi, geometryCSVtoXML::xx, testProducerWithPsetDescEmpty_cfi::y1, testProducerWithPsetDescEmpty_cfi::y2, geometryCSVtoXML::yy, ppsModifySingularModes_cfi::z1, ppsModifySingularModes_cfi::z2, Z_hector, Z_lmax, and Z_lmin.

Referenced by RunManagerMTWorker::generateEvent().

112  {
113  HepMC::GenEvent *evt = new HepMC::GenEvent(*evt_orig);
114 
115  if (*(evt->vertices_begin()) == nullptr) {
116  std::stringstream ss;
117  ss << "SimG4Core/Generator: in event " << g4evt->GetEventID() << " Corrupted Event - GenEvent with no vertex \n";
118  throw SimG4Exception(ss.str());
119  }
120 
121  if (!evt->weights().empty()) {
122  weight_ = evt->weights()[0];
123  for (unsigned int iw = 1; iw < evt->weights().size(); ++iw) {
124  // terminate if the vector of weights contains a zero-weight
125  if (evt->weights()[iw] <= 0)
126  break;
127  weight_ *= evt->weights()[iw];
128  }
129  }
130 
131  if (vtx_ != nullptr) {
132  delete vtx_;
133  }
134  vtx_ = new math::XYZTLorentzVector((*(evt->vertices_begin()))->position().x(),
135  (*(evt->vertices_begin()))->position().y(),
136  (*(evt->vertices_begin()))->position().z(),
137  (*(evt->vertices_begin()))->position().t());
138 
139  edm::LogVerbatim("SimG4CoreGenerator") << "Generator: primary Vertex = (" << vtx_->x() << ", " << vtx_->y() << ", "
140  << vtx_->z() << ")";
141 
142  unsigned int ng4vtx = 0;
143  unsigned int ng4par = 0;
144 
145  for (HepMC::GenEvent::vertex_const_iterator vitr = evt->vertices_begin(); vitr != evt->vertices_end(); ++vitr) {
146  // loop for vertex, is it a real vertex?
147  // Set qvtx to true for any particles that should be propagated by GEANT,
148  // i.e., status 1 particles or status 2 particles that decay outside the
149  // beampipe.
150  G4bool qvtx = false;
151  HepMC::GenVertex::particle_iterator pitr;
152  for (pitr = (*vitr)->particles_begin(HepMC::children); pitr != (*vitr)->particles_end(HepMC::children); ++pitr) {
153  // For purposes of this function, the status is defined as follows:
154  // 1: particles are not decayed by generator
155  // 2: particles are decayed by generator but need to be propagated by GEANT
156  // 3: particles are decayed by generator and do not need to be propagated by GEANT
157  int status = (*pitr)->status();
158  int pdg = (*pitr)->pdg_id();
159  if (status > 3 && isExotic(pdg) && (!(isExoticNonDetectable(pdg)))) {
160  // In Pythia 8, there are many status codes besides 1, 2, 3.
161  // By setting the status to 2 for exotic particles, they will be
162  // checked: if its decay vertex is outside the beampipe, it will be
163  // propagated by GEANT. Some Standard Model particles, e.g., K0, cannot
164  // be propagated by GEANT, so do not change their status code.
165  status = 2;
166  }
167  if (status == 2 && abs(pdg) == 9900015) {
168  status = 3;
169  }
170 
171  // Particles which are not decayed by generator
172  if (status == 1) {
173  // filter out unwanted particles and vertices
175  continue;
176  }
177 
178  qvtx = true;
179  if (verbose > 2)
180  LogDebug("SimG4CoreGenerator") << "GenVertex barcode = " << (*vitr)->barcode()
181  << " selected for GenParticle barcode = " << (*pitr)->barcode();
182  break;
183  }
184  // The selection is made considering if the partcile with status = 2
185  // have the end_vertex with a radius greater than the radius of beampipe
186  // cylinder (no requirement on the Z of the vertex is applyed).
187  else if (status == 2) {
188  if ((*pitr)->end_vertex() != nullptr) {
189  double xx = (*pitr)->end_vertex()->position().x();
190  double yy = (*pitr)->end_vertex()->position().y();
191  double r_dd = xx * xx + yy * yy;
192  if (r_dd > theDecRCut2) {
193  qvtx = true;
194  if (verbose > 2)
195  LogDebug("SimG4CoreGenerator")
196  << "GenVertex barcode = " << (*vitr)->barcode()
197  << " selected for GenParticle barcode = " << (*pitr)->barcode() << " radius = " << std::sqrt(r_dd);
198  break;
199  }
200  } else {
201  // particles with status 2 without end_vertex are
202  // equivalent to stable
203  qvtx = true;
204  break;
205  }
206  }
207  }
208 
209  // if this vertex is inside fiductial volume inside the beam pipe
210  // and has no long-lived secondary the vertex is not saved
211  if (!qvtx) {
212  continue;
213  }
214 
215  double x1 = (*vitr)->position().x() * CLHEP::mm;
216  double y1 = (*vitr)->position().y() * CLHEP::mm;
217  double z1 = (*vitr)->position().z() * CLHEP::mm;
218  double t1 = (*vitr)->position().t() * CLHEP::mm / CLHEP::c_light;
219 
220  G4PrimaryVertex *g4vtx = new G4PrimaryVertex(x1, y1, z1, t1);
221 
222  for (pitr = (*vitr)->particles_begin(HepMC::children); pitr != (*vitr)->particles_end(HepMC::children); ++pitr) {
223  int status = (*pitr)->status();
224  int pdg = (*pitr)->pdg_id();
225  bool hasDecayVertex = (nullptr != (*pitr)->end_vertex());
226 
227  // Filter on allowed particle species if required
228  if (fPDGFilter) {
229  bool isInTheList = IsInTheFilterList(pdg);
230  if ((!pdgFilterSel && isInTheList) || (pdgFilterSel && !isInTheList)) {
231  if (0 < verbose)
232  edm::LogVerbatim("SimG4CoreGenerator")
233  << " Skiped GenParticle barcode= " << (*pitr)->barcode() << " PDGid= " << pdg << " status= " << status
234  << " isExotic: " << isExotic(pdg) << " isExoticNotDet: " << isExoticNonDetectable(pdg)
235  << " isInTheList: " << isInTheList << " hasDecayVertex: " << hasDecayVertex;
236  continue;
237  }
238  }
239 
240  if (0 < verbose) {
241  edm::LogVerbatim("SimG4CoreGenerator")
242  << "Generator: pdg= " << pdg << " status= " << status << " hasPreDefinedDecay: " << hasDecayVertex
243  << "\n isExotic: " << isExotic(pdg) << " isExoticNotDet: " << isExoticNonDetectable(pdg)
244  << " isInTheList: " << IsInTheFilterList(pdg) << "\n"
245  << " MaxZCentralCMS = " << maxZCentralCMS / CLHEP::m << " m; (x,y,z,t): (" << x1 << "," << y1 << "," << z1
246  << "," << t1 << ")";
247  }
248  if (status > 3 && isExotic(pdg) && (!(isExoticNonDetectable(pdg)))) {
249  status = hasDecayVertex ? 2 : 1;
250  }
251  if (status == 2 && abs(pdg) == 9900015) {
252  status = 3;
253  }
254 
255  // this particle has predefined decay but has no vertex
256  if (2 == status && !hasDecayVertex) {
257  edm::LogWarning("SimG4CoreGenerator: in event ")
258  << g4evt->GetEventID() << " a particle "
259  << " pdgid= " << pdg << " has status=2 but has no decay vertex, so will be fully tracked by Geant4";
260  status = 1;
261  }
262 
263  double x2 = x1;
264  double y2 = y1;
265  double z2 = z1;
266  double decay_length = 0.0;
267  if (2 == status) {
268  x2 = (*pitr)->end_vertex()->position().x();
269  y2 = (*pitr)->end_vertex()->position().y();
270  z2 = (*pitr)->end_vertex()->position().z();
271  decay_length = std::sqrt((x2 - x1) * (x2 - x1) + (y2 - y1) * (y2 - y1) + (z2 - z1) * (z2 - z1));
272  }
273 
274  bool toBeAdded = !fFiductialCuts;
275 
276  double px = (*pitr)->momentum().px();
277  double py = (*pitr)->momentum().py();
278  double pz = (*pitr)->momentum().pz();
279  double ptot = std::sqrt(px * px + py * py + pz * pz);
280  math::XYZTLorentzVector p(px, py, pz, (*pitr)->momentum().e());
281 
282  double ximpact = x1;
283  double yimpact = y1;
284  double zimpact = z1;
285 
286  // protection against numerical problems for extremely low momenta
287  // compute impact point at transition to Hector
288  const double minTan = 1.e-20;
289  if (std::abs(z1) < Z_hector && std::abs(pz) >= minTan * ptot) {
290  zimpact = (pz > 0.0) ? Z_hector : -Z_hector;
291  double del = (zimpact - z1) / pz;
292  ximpact += del * px;
293  yimpact += del * py;
294  }
295  double rimpact2 = ximpact * ximpact + yimpact * yimpact;
296 
297  if (verbose > 2)
298  LogDebug("SimG4CoreGenerator") << "Processing GenParticle barcode= " << (*pitr)->barcode() << " pdg= " << pdg
299  << " status= " << (*pitr)->status() << " st= " << status
300  << " rimpact(cm)= " << std::sqrt(rimpact2) / CLHEP::cm
301  << " zimpact(cm)= " << zimpact / CLHEP::cm << " ptot(GeV)= " << ptot
302  << " pz(GeV)= " << pz;
303 
304  // Particles of status 1 trasnported along the beam pipe
305  // HECTOR transport of protons are done in corresponding PPS producer
306  if (1 == status && std::abs(zimpact) >= Z_hector && rimpact2 <= theDecRCut2) {
307  // very forward n, nbar, gamma are allowed
308  toBeAdded = (2112 == std::abs(pdg) || 22 == pdg);
309  if (verbose > 1) {
310  edm::LogVerbatim("SimG4CoreGenerator")
311  << "GenParticle barcode = " << (*pitr)->barcode() << " very forward; to be added: " << toBeAdded;
312  }
313  } else {
314  // Standard case: particles not decayed by the generator and not forward
315  if (1 == status && (std::abs(zimpact) < Z_hector || rimpact2 > theDecRCut2)) {
316  // Ptot cut for all particles
317  if (fPCuts && (ptot < theMinPCut || ptot > theMaxPCut)) {
318  continue;
319  }
320  // phi cut is applied mainly for particle gun
321  if (fPhiCuts) {
322  double phi = p.phi();
323  if (phi < theMinPhiCut || phi > theMaxPhiCut) {
324  continue;
325  }
326  }
327  // eta cut is applied if position of the decay
328  // is within vacuum chamber and limited in Z
329  if (fEtaCuts) {
330  // eta cut
331  double xi = x1;
332  double yi = y1;
333  double zi = z1;
334 
335  // can be propagated along Z
336  if (std::abs(pz) >= minTan * ptot) {
337  if ((zi >= Z_lmax) & (pz < 0.0)) {
338  zi = Z_lmax;
339  } else if ((zi <= Z_lmin) & (pz > 0.0)) {
340  zi = Z_lmin;
341  } else {
342  if (pz > 0) {
343  zi = Z_lmax;
344  } else {
345  zi = Z_lmin;
346  }
347  }
348  double del = (zi - z1) / pz;
349  xi += del * px;
350  yi += del * py;
351  }
352  // check eta cut
353  if ((zi >= Z_lmin) & (zi <= Z_lmax) & (xi * xi + yi * yi < theDecRCut2)) {
354  continue;
355  }
356  }
357  if (fLumiFilter && !fLumiFilter->isGoodForLumiMonitor(*pitr)) {
358  continue;
359  }
360  toBeAdded = true;
361  if (verbose > 1)
362  edm::LogVerbatim("SimG4CoreGenerator")
363  << "GenParticle barcode = " << (*pitr)->barcode() << " passed case 1";
364 
365  // Decay chain outside the fiducial cylinder defined by theRDecLenCut
366  // are used for Geant4 tracking with predefined decay channel
367  // In the case of decay in vacuum particle is not tracked by Geant4
368  } else if (2 == status && x2 * x2 + y2 * y2 >= theDecRCut2 && std::abs(z2) < Z_hector) {
369  toBeAdded = true;
370  if (verbose > 1)
371  edm::LogVerbatim("SimG4CoreGenerator") << "GenParticle barcode = " << (*pitr)->barcode() << " passed case 2"
372  << " decay_length(cm)= " << decay_length / CLHEP::cm;
373  }
374  }
375  if (toBeAdded) {
376  G4PrimaryParticle *g4prim = new G4PrimaryParticle(pdg, px * CLHEP::GeV, py * CLHEP::GeV, pz * CLHEP::GeV);
377 
378  if (g4prim->GetG4code() != nullptr) {
379  g4prim->SetMass(g4prim->GetG4code()->GetPDGMass());
380  double charge = g4prim->GetG4code()->GetPDGCharge();
381 
382  // apply Pt cut
383  if (fPtransCut && 1 == status && 0.0 != charge && px * px + py * py < theMinPtCut2) {
384  delete g4prim;
385  continue;
386  }
387  g4prim->SetCharge(charge);
388  }
389 
390  // V.I. do not use SetWeight but the same code
391  // value of the code compute inside TrackWithHistory
392  // g4prim->SetWeight( 10000*(*vpitr)->barcode() ) ;
393  setGenId(g4prim, (*pitr)->barcode());
394 
395  if (2 == status) {
396  particleAssignDaughters(g4prim, (HepMC::GenParticle *)*pitr, decay_length);
397  }
398  if (verbose > 1)
399  g4prim->Print();
400 
401  ++ng4par;
402  g4vtx->SetPrimary(g4prim);
403  edm::LogVerbatim("SimG4CoreGenerator") << " " << ng4par << ". new Geant4 particle pdg= " << pdg
404  << " Ptot(GeV/c)= " << ptot << " Pt= " << std::sqrt(px * px + py * py)
405  << " status= " << status << "; dir= " << g4prim->GetMomentumDirection();
406  }
407  }
408 
409  if (verbose > 1)
410  g4vtx->Print();
411  g4evt->AddPrimaryVertex(g4vtx);
412  ++ng4vtx;
413  }
414 
415  // Add a protection for completely empty events (produced by LHCTransport):
416  // add a dummy vertex with no particle attached to it
417  if (ng4vtx == 0) {
418  G4PrimaryVertex *g4vtx = new G4PrimaryVertex(0.0, 0.0, 0.0, 0.0);
419  if (verbose > 1)
420  g4vtx->Print();
421 
422  g4evt->AddPrimaryVertex(g4vtx);
423  }
424 
425  edm::LogVerbatim("SimG4CoreGenerator") << "The list of Geant4 primaries includes " << ng4par << " particles in "
426  << ng4vtx << " vertex";
427  delete evt;
428 }
double theDecRCut2
Definition: Generator.h:55
Log< level::Info, true > LogVerbatim
bool fPhiCuts
Definition: Generator.h:46
bool isExoticNonDetectable(int pdgcode) const
Definition: Generator.cc:536
double theMinPtCut2
Definition: Generator.h:53
double maxZCentralCMS
Definition: Generator.h:58
double Z_lmax
Definition: Generator.h:64
double theMaxPhiCut
Definition: Generator.h:49
bool isExotic(int pdgcode) const
Definition: Generator.cc:528
bool fPtransCut
Definition: Generator.h:44
bool isGoodForLumiMonitor(const HepMC::GenParticle *) const
bool pdgFilterSel
Definition: Generator.h:66
void particleAssignDaughters(G4PrimaryParticle *p, HepMC::GenParticle *hp, double length)
Definition: Generator.cc:430
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
LumiMonitorFilter * fLumiFilter
Definition: Generator.h:60
bool fFiductialCuts
Definition: Generator.h:47
T sqrt(T t)
Definition: SSEVec.h:23
bool fPDGFilter
Definition: Generator.h:67
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double weight_
Definition: Generator.h:63
double theMaxPCut
Definition: Generator.h:54
math::XYZTLorentzVector * vtx_
Definition: Generator.h:62
double Z_lmin
Definition: Generator.h:64
int verbose
Definition: Generator.h:59
void setGenId(G4PrimaryParticle *p, int id) const
Definition: Generator.h:40
Log< level::Warning, false > LogWarning
bool IsInTheFilterList(int pdgcode) const
Definition: Generator.cc:543
bool fEtaCuts
Definition: Generator.h:45
double Z_hector
Definition: Generator.h:64
#define LogDebug(id)
bool fPCuts
Definition: Generator.h:43

◆ isExotic()

bool Generator::isExotic ( int  pdgcode) const
private

Definition at line 528 of file Generator.cc.

References funct::abs(), and EgammaValidation_cff::pdgid.

Referenced by HepMC2G4().

528  {
529  int pdgid = std::abs(pdgcode);
530  return ((pdgid >= 1000000 && pdgid < 4000000 && pdgid != 3000022) || // SUSY, R-hadron, and technicolor particles
531  pdgid == 17 || // 4th generation lepton
532  pdgid == 34 || // W-prime
533  pdgid == 37); // charged Higgs
534 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22

◆ isExoticNonDetectable()

bool Generator::isExoticNonDetectable ( int  pdgcode) const
private

Definition at line 536 of file Generator.cc.

References funct::abs(), ALCARECOTkAlJpsiMuMu_cff::charge, LHEGenericFilter_cfi::ParticleID, and EgammaValidation_cff::pdgid.

Referenced by HepMC2G4().

536  {
537  int pdgid = std::abs(pdgcode);
538  HepPDT::ParticleID pid(pdgcode);
539  int charge = pid.threeCharge();
540  return ((charge == 0) && (pdgid >= 1000000 && pdgid < 1000040)); // SUSY
541 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22

◆ IsInTheFilterList()

bool Generator::IsInTheFilterList ( int  pdgcode) const
private

Definition at line 543 of file Generator.cc.

References funct::abs(), pdgFilter, and EgammaValidation_cff::pdgid.

Referenced by HepMC2G4().

543  {
544  int pdgid = std::abs(pdgcode);
545  for (auto &pdg : pdgFilter) {
546  if (std::abs(pdg) == pdgid) {
547  return true;
548  }
549  }
550  return false;
551 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::vector< int > pdgFilter
Definition: Generator.h:65

◆ nonCentralEvent2G4()

void Generator::nonCentralEvent2G4 ( const HepMC::GenEvent g,
G4Event *  e 
)

Definition at line 553 of file Generator.cc.

References GenParticle::GenParticle, runTauDisplay::gp, mps_fire::i, ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::it, setGenId(), findQualityFiles::v, and verbose.

Referenced by RunManagerMTWorker::generateEvent().

553  {
554  int i = g4evt->GetNumberOfPrimaryVertex();
555  for (HepMC::GenEvent::particle_const_iterator it = evt->particles_begin(); it != evt->particles_end(); ++it) {
556  ++i;
557  HepMC::GenParticle *gp = (*it);
558 
559  // storing only particle with status == 1
560  if (gp->status() != 1)
561  continue;
562 
563  int pdg = gp->pdg_id();
564  G4PrimaryParticle *g4p = new G4PrimaryParticle(
565  pdg, gp->momentum().px() * CLHEP::GeV, gp->momentum().py() * CLHEP::GeV, gp->momentum().pz() * CLHEP::GeV);
566  if (g4p->GetG4code() != nullptr) {
567  g4p->SetMass(g4p->GetG4code()->GetPDGMass());
568  g4p->SetCharge(g4p->GetG4code()->GetPDGCharge());
569  }
570  setGenId(g4p, i);
571  G4PrimaryVertex *v = new G4PrimaryVertex(gp->production_vertex()->position().x() * CLHEP::mm,
572  gp->production_vertex()->position().y() * CLHEP::mm,
573  gp->production_vertex()->position().z() * CLHEP::mm,
574  gp->production_vertex()->position().t() * CLHEP::mm / CLHEP::c_light);
575  v->SetPrimary(g4p);
576  g4evt->AddPrimaryVertex(v);
577  if (verbose > 0)
578  v->Print();
579  } // end loop on HepMC particles
580 }
int verbose
Definition: Generator.h:59
void setGenId(G4PrimaryParticle *p, int id) const
Definition: Generator.h:40

◆ particleAssignDaughters()

void Generator::particleAssignDaughters ( G4PrimaryParticle *  p,
HepMC::GenParticle *  hp,
double  length 
)
private

Definition at line 430 of file Generator.cc.

References funct::abs(), class-composition::children, createTree::dd, LogDebug, AlCaHLTBitMon_ParallelJobs::p, setGenId(), mathSSE::sqrt(), mps_update::status, verbose, testProducerWithPsetDescEmpty_cfi::x1, testProducerWithPsetDescEmpty_cfi::x2, testProducerWithPsetDescEmpty_cfi::y1, testProducerWithPsetDescEmpty_cfi::y2, ppsModifySingularModes_cfi::z1, and ppsModifySingularModes_cfi::z2.

Referenced by HepMC2G4().

430  {
431  if (verbose > 1) {
432  LogDebug("SimG4CoreGenerator") << "Special case of long decay length \n"
433  << "Assign daughters with to mother with decaylength=" << decaylength / CLHEP::cm
434  << " cm";
435  }
436  math::XYZTLorentzVector p(vp->momentum().px(), vp->momentum().py(), vp->momentum().pz(), vp->momentum().e());
437 
438  // defined preassigned decay time
439  double proper_time = decaylength / (p.Beta() * p.Gamma() * c_light);
440  g4p->SetProperTime(proper_time);
441 
442  if (verbose > 2) {
443  LogDebug("SimG4CoreGenerator") << " px= " << p.px() << " py= " << p.py() << " pz= " << p.pz() << " e= " << p.e()
444  << " beta= " << p.Beta() << " gamma= " << p.Gamma()
445  << " Proper time= " << proper_time / CLHEP::ns << " ns";
446  }
447 
448  // the particle will decay after the same length if it
449  // has not interacted before
450  double x1 = vp->end_vertex()->position().x();
451  double y1 = vp->end_vertex()->position().y();
452  double z1 = vp->end_vertex()->position().z();
453 
454  for (HepMC::GenVertex::particle_iterator vpdec = vp->end_vertex()->particles_begin(HepMC::children);
455  vpdec != vp->end_vertex()->particles_end(HepMC::children);
456  ++vpdec) {
457  // transform decay products such that in the rest frame of mother
459  (*vpdec)->momentum().px(), (*vpdec)->momentum().py(), (*vpdec)->momentum().pz(), (*vpdec)->momentum().e());
460 
461  // children should only be taken into account once
462  G4PrimaryParticle *g4daught =
463  new G4PrimaryParticle((*vpdec)->pdg_id(), pdec.x() * CLHEP::GeV, pdec.y() * CLHEP::GeV, pdec.z() * CLHEP::GeV);
464 
465  if (g4daught->GetG4code() != nullptr) {
466  g4daught->SetMass(g4daught->GetG4code()->GetPDGMass());
467  g4daught->SetCharge(g4daught->GetG4code()->GetPDGCharge());
468  }
469 
470  // V.I. do not use SetWeight but the same code
471  // value of the code compute inside TrackWithHistory
472  setGenId(g4daught, (*vpdec)->barcode());
473 
474  int status = (*vpdec)->status();
475  if (verbose > 1)
476  LogDebug("SimG4CoreGenerator::::particleAssignDaughters")
477  << "Assigning a " << (*vpdec)->pdg_id() << " as daughter of a " << vp->pdg_id() << " status=" << status;
478 
479  if ((status == 2 || (status == 23 && std::abs(vp->pdg_id()) == 1000015) || (status > 50 && status < 100)) &&
480  (*vpdec)->end_vertex() != nullptr) {
481  double x2 = (*vpdec)->end_vertex()->position().x();
482  double y2 = (*vpdec)->end_vertex()->position().y();
483  double z2 = (*vpdec)->end_vertex()->position().z();
484  double dd = std::sqrt((x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2) + (z1 - z2) * (z1 - z2));
485  particleAssignDaughters(g4daught, *vpdec, dd);
486  }
487  (*vpdec)->set_status(1000 + status);
488  g4p->SetDaughter(g4daught);
489 
490  if (verbose > 1)
491  g4daught->Print();
492  }
493 }
string dd
Definition: createTree.py:154
void particleAssignDaughters(G4PrimaryParticle *p, HepMC::GenParticle *hp, double length)
Definition: Generator.cc:430
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
T sqrt(T t)
Definition: SSEVec.h:23
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int verbose
Definition: Generator.h:59
void setGenId(G4PrimaryParticle *p, int id) const
Definition: Generator.h:40
#define LogDebug(id)

◆ particlePassesPrimaryCuts()

bool Generator::particlePassesPrimaryCuts ( const G4ThreeVector &  p) const
private

Definition at line 496 of file Generator.cc.

References MillePedeFileConverter_cfg::e, PVValHelper::eta, fEtaCuts, RemoveAddSevLevel::flag, fPCuts, fPhiCuts, LogDebug, AlCaHLTBitMon_ParallelJobs::p, phi, theMaxEtaCut, theMaxPCut, theMaxPhiCut, and verbose.

496  {
497  bool flag = true;
498  double ptot = p.mag();
499  if (fPCuts && (ptot < theMinPCut * CLHEP::GeV || ptot > theMaxPCut * CLHEP::GeV)) {
500  flag = false;
501  }
502  if (fEtaCuts && flag) {
503  double pz = p.z();
504  if (ptot < pz + 1.e-10) {
505  flag = false;
506 
507  } else {
508  double eta = 0.5 * G4Log((ptot + pz) / (ptot - pz));
509  if (eta < theMinEtaCut || eta > theMaxEtaCut) {
510  flag = false;
511  }
512  }
513  }
514  if (fPhiCuts && flag) {
515  double phi = p.phi();
516  if (phi < theMinPhiCut || phi > theMaxPhiCut) {
517  flag = false;
518  }
519  }
520 
521  if (verbose > 2)
522  LogDebug("SimG4CoreGenerator") << "Generator ptot(GeV)= " << ptot / CLHEP::GeV << " eta= " << p.eta()
523  << " phi= " << p.phi() << " Flag= " << flag;
524 
525  return flag;
526 }
bool fPhiCuts
Definition: Generator.h:46
double theMaxPhiCut
Definition: Generator.h:49
double theMaxEtaCut
Definition: Generator.h:51
double theMaxPCut
Definition: Generator.h:54
int verbose
Definition: Generator.h:59
bool fEtaCuts
Definition: Generator.h:45
#define LogDebug(id)
bool fPCuts
Definition: Generator.h:43

◆ setGenEvent()

void Generator::setGenEvent ( const HepMC::GenEvent inpevt)
inline

Definition at line 24 of file Generator.h.

References evt_.

Referenced by RunManagerMTWorker::generateEvent().

24  {
25  evt_ = (HepMC::GenEvent *)inpevt;
26  return;
27  }
HepMC::GenEvent * evt_
Definition: Generator.h:61

◆ setGenId()

void Generator::setGenId ( G4PrimaryParticle *  p,
int  id 
) const
inlineprivate

Member Data Documentation

◆ evt_

HepMC::GenEvent* Generator::evt_
private

Definition at line 61 of file Generator.h.

Referenced by genEvent(), and setGenEvent().

◆ fEtaCuts

bool Generator::fEtaCuts
private

Definition at line 45 of file Generator.h.

Referenced by Generator(), HepMC2G4(), and particlePassesPrimaryCuts().

◆ fFiductialCuts

bool Generator::fFiductialCuts
private

Definition at line 47 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ fLumiFilter

LumiMonitorFilter* Generator::fLumiFilter
private

Definition at line 60 of file Generator.h.

Referenced by Generator(), HepMC2G4(), and ~Generator().

◆ fPCuts

bool Generator::fPCuts
private

Definition at line 43 of file Generator.h.

Referenced by Generator(), HepMC2G4(), and particlePassesPrimaryCuts().

◆ fPDGFilter

bool Generator::fPDGFilter
private

Definition at line 67 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ fPhiCuts

bool Generator::fPhiCuts
private

Definition at line 46 of file Generator.h.

Referenced by Generator(), HepMC2G4(), and particlePassesPrimaryCuts().

◆ fPtransCut

bool Generator::fPtransCut
private

Definition at line 44 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ maxZCentralCMS

double Generator::maxZCentralCMS
private

Definition at line 58 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ pdgFilter

std::vector<int> Generator::pdgFilter
private

Definition at line 65 of file Generator.h.

Referenced by Generator(), and IsInTheFilterList().

◆ pdgFilterSel

bool Generator::pdgFilterSel
private

Definition at line 66 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ theDecLenCut

double Generator::theDecLenCut
private

Definition at line 57 of file Generator.h.

Referenced by Generator().

◆ theDecRCut2

double Generator::theDecRCut2
private

Definition at line 55 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ theEtaCutForHector

double Generator::theEtaCutForHector
private

Definition at line 56 of file Generator.h.

Referenced by Generator().

◆ theMaxEtaCut

double Generator::theMaxEtaCut
private

Definition at line 51 of file Generator.h.

Referenced by Generator(), and particlePassesPrimaryCuts().

◆ theMaxPCut

double Generator::theMaxPCut
private

Definition at line 54 of file Generator.h.

Referenced by Generator(), HepMC2G4(), and particlePassesPrimaryCuts().

◆ theMaxPhiCut

double Generator::theMaxPhiCut
private

Definition at line 49 of file Generator.h.

Referenced by Generator(), HepMC2G4(), and particlePassesPrimaryCuts().

◆ theMinEtaCut

double Generator::theMinEtaCut
private

Definition at line 50 of file Generator.h.

Referenced by Generator().

◆ theMinPCut

double Generator::theMinPCut
private

Definition at line 52 of file Generator.h.

Referenced by Generator().

◆ theMinPhiCut

double Generator::theMinPhiCut
private

Definition at line 48 of file Generator.h.

Referenced by Generator().

◆ theMinPtCut2

double Generator::theMinPtCut2
private

Definition at line 53 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ verbose

int Generator::verbose
private

◆ vtx_

math::XYZTLorentzVector* Generator::vtx_
private

Definition at line 62 of file Generator.h.

Referenced by genVertex(), and HepMC2G4().

◆ weight_

double Generator::weight_
private

Definition at line 63 of file Generator.h.

Referenced by eventWeight(), and HepMC2G4().

◆ Z_hector

double Generator::Z_hector
private

Definition at line 64 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ Z_lmax

double Generator::Z_lmax
private

Definition at line 64 of file Generator.h.

Referenced by Generator(), and HepMC2G4().

◆ Z_lmin

double Generator::Z_lmin
private

Definition at line 64 of file Generator.h.

Referenced by Generator(), and HepMC2G4().