CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
FlatRandomEGunProducer.cc
Go to the documentation of this file.
1 /*
2  * \author Julia Yarba
3  */
4 
5 #include <ostream>
6 
8 
11 
14 
15 using namespace edm;
16 using namespace std;
17 
20 {
21 
22  ParameterSet defpset ;
23  // ParameterSet pgun_params = pset.getParameter<ParameterSet>("PGunParameters") ;
24  ParameterSet pgun_params =
25  pset.getParameter<ParameterSet>("PGunParameters") ;
26 
27  // doesn't seem necessary to check if pset is empty - if this
28  // is the case, default values will be taken for params
29  fMinE = pgun_params.getParameter<double>("MinE");
30  fMaxE = pgun_params.getParameter<double>("MaxE");
31 
32  produces<HepMCProduct>();
33  produces<GenEventInfoProduct>();
34 
35  cout << "Internal FlatRandomEGun is initialzed" << endl ;
36 // cout << "It is going to generate " << remainingEvents() << "events" << endl ;
37 
38 }
39 
41 {
42  // no need to cleanup fEvt since it's done in HepMCProduct
43 }
44 
46 {
47 
48  if ( fVerbosity > 0 )
49  {
50  cout << " FlatRandomEGunProducer : Begin New Event Generation" << endl ;
51  }
52 
53  // event loop (well, another step in it...)
54 
55  // no need to clean up GenEvent memory - done in HepMCProduct
56 
57  // here re-create fEvt (memory)
58  //
59  fEvt = new HepMC::GenEvent() ;
60 
61  // now actualy, cook up the event from PDGTable and gun parameters
62  //
63 
64  // 1st, primary vertex
65  //
66  HepMC::GenVertex* Vtx = new HepMC::GenVertex( HepMC::FourVector(0.,0.,0.));
67 
68  // loop over particles
69  //
70  int barcode = 1;
71  for (unsigned int ip=0; ip<fPartIDs.size(); ip++)
72  {
73  double energy = fRandomGenerator->fire(fMinE, fMaxE) ;
74  double eta = fRandomGenerator->fire(fMinEta, fMaxEta) ;
75  double phi = fRandomGenerator->fire(fMinPhi, fMaxPhi) ;
76  int PartID = fPartIDs[ip] ;
77  const HepPDT::ParticleData*
78  PData = fPDGTable->particle(HepPDT::ParticleID(abs(PartID))) ;
79  double mass = PData->mass().value() ;
80  double mom2 = energy*energy - mass*mass ;
81  double mom = 0. ;
82  if (mom2 > 0.)
83  {
84  mom = sqrt(mom2) ;
85  }
86  else
87  {
88  mom = 0. ;
89  }
90  double theta = 2.*atan(exp(-eta)) ;
91  double px = mom*sin(theta)*cos(phi) ;
92  double py = mom*sin(theta)*sin(phi) ;
93  double pz = mom*cos(theta) ;
94 
95  HepMC::FourVector p(px,py,pz,energy) ;
96  HepMC::GenParticle* Part =
97  new HepMC::GenParticle(p,PartID,1);
98  Part->suggest_barcode( barcode ) ;
99  barcode++ ;
100  Vtx->add_particle_out(Part);
101 
102  if ( fAddAntiParticle )
103  {
104  HepMC::FourVector ap(-px,-py,-pz,energy) ;
105  int APartID = -PartID ;
106  if ( PartID == 22 || PartID == 23 )
107  {
108  APartID = PartID ;
109  }
110  HepMC::GenParticle* APart =
111  new HepMC::GenParticle(ap,APartID,1);
112  APart->suggest_barcode( barcode ) ;
113  barcode++ ;
114  Vtx->add_particle_out(APart) ;
115  }
116 
117  }
118  fEvt->add_vertex(Vtx) ;
119  fEvt->set_event_number(e.id().event()) ;
120  fEvt->set_signal_process_id(20) ;
121 
122 
123  if ( fVerbosity > 0 )
124  {
125  fEvt->print() ;
126  }
127 
128  auto_ptr<HepMCProduct> BProduct(new HepMCProduct()) ;
129  BProduct->addHepMCData( fEvt );
130  e.put(BProduct);
131 
132  auto_ptr<GenEventInfoProduct> genEventInfo(new GenEventInfoProduct(fEvt));
133  e.put(genEventInfo);
134 
135  if ( fVerbosity > 0 )
136  {
137  // for testing purpose only
138  //fEvt->print() ; // for some strange reason, it prints NO event info
139  // after it's been put into edm::Event...
140  cout << " FlatRandomEGunProducer : Event Generation Done " << endl;
141  }
142 }
FlatRandomEGunProducer(const ParameterSet &pset)
T getParameter(std::string const &) const
EventNumber_t event() const
Definition: EventID.h:44
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Geom::Theta< T > theta() const
T eta() const
ESHandle< HepPDT::ParticleDataTable > fPDGTable
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:116
T sqrt(T t)
Definition: SSEVec.h:48
CLHEP::RandFlat * fRandomGenerator
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::vector< int > fPartIDs
HepPDT::ParticleData ParticleData
edm::EventID id() const
Definition: EventBase.h:56
tuple cout
Definition: gather_cfg.py:121
virtual void produce(Event &e, const EventSetup &es) override
Definition: DDAxes.h:10