StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StarGenerator.h
1 #ifndef __StarGenerator_h__
2 #define __StarGenerator_h__
3 
4 #include "StMaker.h"
5 #include "TString.h"
6 #include "TDatabasePDG.h"
7 #include "StarCallf77.h"
8 #include "TLorentzVector.h"
9 
10 #include "StarGenerator/EVENT/StarGenEvent.h"
11 #include "StarGenerator/EVENT/StarGenStats.h"
12 #include "StarGenerator/UTIL/StarParticleData.h"
13 
14 
34 class StarGenerator : public StMaker
35 {
36  public:
37 
38  StarGenerator( const Char_t *name="" );
39  ~StarGenerator(){ /* nada */ };
40 
41  protected:
43  const Double_t mm;
45  const Double_t cm;
46 
47  public:
51  virtual Int_t Init() = 0;
52 
54  Int_t InitRun( Int_t runnumber ) { return kStOK; };
55 
57  virtual Int_t PreGenerate(){ return kStOK; }
59  virtual Int_t PostGenerate(){ return kStOK; };
60 
66  virtual Int_t Generate() = 0;
67 
69  virtual void FillPP( StarGenEvent *event ){ };
71  virtual void FillAA( StarGenEvent *event ){ };
73  virtual void FillEP( StarGenEvent *event ){ };
75  virtual void FillEA( StarGenEvent *event ){ };
77  virtual void FillUser( StarGenEvent *event ){ };
78 
80  Int_t PreGenerateHook();
82  Int_t PostGenerateHook();
84 
86  Int_t Make() { return Generate(); }
87 
89  void Clear( const Option_t *opts="" );
90 
93  Int_t Finalize();
94 
95 
96  /********************************************************************************
97  *
98  * Configuration Options
99  *
100  *******************************************************************************/
101 
102 
104  void SetBlue( const Char_t *b );// { mBlue = b; }
106  void SetYell( const Char_t *y ); //{ mYell = y; }
110  void SetFrame( const Char_t *frame, const Double_t val );
115  void SetFrame( const Char_t *frame, const Double_t *pblue, const Double_t *pyell );
116 
118  void SetImpact( Double_t bmin, Double_t bmax )
119  {
120  mImpactMin = bmin;
121  mImpactMax = bmax;
122  }
123 
125  Int_t GetNumberOfParticles(){ return mNumberOfParticles; }
126 
130  void SetPileup( Bool_t p, Double_t prob=0.10 ){ mIsPileup = p; if ( p ) mPileup=prob; else mPileup=1.0; }
131 
132 
134  Bool_t IsPileup(){ return mIsPileup; }
135 
137  Double_t GetPileup(){ return mPileup; }
138 
139 
142  void SetOutputTree ( TTree *tree );
143 
147  void SetInputTree ( TTree *tree, const Char_t *bname=0 );
148 
149  void SetInputFile ( const Char_t *name, const Char_t *treename="genevents", const Char_t *bname=0 );
150 
152 
153 
154  Int_t IOmode(){ return mIOMode; }
155 
157  StarGenEvent *Event(){ return mEvent; }
158 
160  virtual StarGenStats Stats() { return StarGenStats(); }
161 
162 
163  virtual const char *GetCVS() const
164  {static const char cvs[]="Tag $Name: $ $Id: StarGenerator.h,v 1.6 2015/04/24 18:21:11 jwebb Exp $ built " __DATE__ " " __TIME__ ; return cvs;}
165 
166 
167  private:
168 
169  // Generator ID
170  Int_t mId;
171  // Set generator ID
172  void SetId( Int_t id ) { mId = id; }
173 
174  protected:
175  // Number of particles in this event
176  Int_t mNumberOfParticles;
177 
180 
182  TString mBlue;
184  TString mYell;
186  TString mFrame;
188  Double_t mRootS;
190  Double_t mDirect;
191 
193  Double_t mImpactMin;
195  Double_t mImpactMax;
196 
198  TLorentzVector mBlueMomentum;
200  TLorentzVector mYellMomentum;
202  Double_t mBlueMass;
204  Double_t mYellMass;
205 
207  Double_t mVertex[3];
209  Double_t mSpread[3];
210 
212  Int_t mIOMode;
213  TTree *mInputTree;
215  TTree *mOutputTree;
217 
219  Bool_t mIsPileup;
221  Double_t mPileup;
222 
225 
227  TString mLibrary;
228 
229  friend class StarPrimaryMaker;
230 
231  ClassDef(StarGenerator,1);
232 
233 };
234 #endif
Double_t mRootS
CMS energy or incident beam momentum for fixed target collisions.
TString mLibrary
Concrete event generator library.
virtual Int_t Generate()=0
void SetOutputTree(TTree *tree)
TString mYell
Name of the yellow beam particle (-z)
Bool_t IsPileup()
Returns true if this is a pileup generator.
Int_t IOmode()
Returns the ID of this event generator.
StarGenEvent * Event()
Retrieves the event record.
void SetFrame(const Char_t *frame, const Double_t val)
Double_t mVertex[3]
Position of the interaction vertex.
Int_t Make()
If called, make simply aliases to the Generate() function.
Definition: StarGenerator.h:86
void SetImpact(Double_t bmin, Double_t bmax)
Set the minimum and maximum impact parameters for the collision (HI generators)
Bool_t mIsPileup
Flag which indicates whether this is a pileup simulation.
TLorentzVector mYellMomentum
4-Momentum of the yellow beam
Int_t mIOMode
IO flag. 0=no IO, 1=write, 2=read.
virtual void FillEP(StarGenEvent *event)
(Optional) Method to fill a DIS event
Definition: StarGenerator.h:73
ABC for defining event generator interfaces.
Definition: StarGenerator.h:34
void Clear(const Option_t *opts="")
Clear the event.
void SetInputTree(TTree *tree, const Char_t *bname=0)
StarGenEvent * mEvent
Generated event.
Double_t mImpactMin
Minimum impact parameter in a HI collision.
virtual void FillEA(StarGenEvent *event)
(Optional) Method to fill a eA event
Definition: StarGenerator.h:75
TTree * mOutputTree
Pointer to the tree for writing out events.
Interface to PDG information.
Int_t InitRun(Int_t runnumber)
Developer may not provide run initialization.
Definition: StarGenerator.h:54
Double_t mImpactMax
Maximum impact parameter in a HI collision.
void SetBlue(const Char_t *b)
Sets the particle species for the blue beam.
Double_t mDirect
Direction (+1 = W, -1 = E) of the beam in fixted target mode.
Int_t GetNumberOfParticles()
Returns the number of particles in the event.
TLorentzVector mBlueMomentum
4-Momentum of the blue beam
virtual Int_t PostGenerate()
Developers may provide a post-generate method which will execute after Generate().
Definition: StarGenerator.h:59
Double_t mSpread[3]
Width of the interaction vertex.
Double_t mPileup
Indicates the probability that this generator should be run in pileup mode.
const Double_t cm
If an event generator measures distances in cm, then multiply by cm.
Definition: StarGenerator.h:45
End of run statistics for event generators.
Definition: StarGenStats.h:21
TTree * mInputTree
Pointer to the tree for reading in events.
const Double_t mm
If an event generator measures distances in mm, then multiply by mm.
Definition: StarGenerator.h:39
Base class for event records.
Definition: StarGenEvent.h:81
Definition: Stypes.h:40
virtual void FillUser(StarGenEvent *event)
(Optional) Method to fill a user-defined event
Definition: StarGenerator.h:77
void SetPileup(Bool_t p, Double_t prob=0.10)
virtual StarGenStats Stats()
Create and retrieve end-of-run statistics.
Main steering class for event generation.
Double_t mBlueMass
Mass of the blue beam.
void SetYell(const Char_t *y)
Sets the particle species for the yellow beam.
Double_t mYellMass
Mass of the yellow beam.
StarParticleData * mParticleDb
Database for particle information (based on TDatabasePDG)
virtual Int_t Init()=0
TString mFrame
Frame of the collision, i.e. CMS, FIXT, 3MOM, 4MOM, 5MOM.
TString mBlue
Name of the blue beam particle (+z)
virtual void FillPP(StarGenEvent *event)
(Optional) Method to fill a PP event
Definition: StarGenerator.h:69
virtual Int_t PreGenerate()
Developers may provide a pre-generate method which will execute before Generate().
Definition: StarGenerator.h:57
virtual void FillAA(StarGenEvent *event)
(Optional) Method to fill a AA event
Definition: StarGenerator.h:71
Double_t GetPileup()
Returns the pileup probability.