StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
starsim.decayer.Jpsi.C
1 // macro to instantiate the Geant3 from within
2 // STAR C++ framework and get the starsim prompt
3 // To use it do
4 // root4star starsim.C
5 
6 class St_geant_Maker;
7 St_geant_Maker *geant_maker = 0;
8 
9 class StarGenEvent;
10 StarGenEvent *event = 0;
11 
12 class StarPrimaryMaker;
13 StarPrimaryMaker *_primary = 0;
14 
15 class StarKinematics;
17 
18 TF1 *ptDist = 0;
19 TF1 *etaDist = 0;
20 
21 // ----------------------------------------------------------------------------
22 void geometry( TString tag, Bool_t agml=true )
23 {
24  TString cmd = "DETP GEOM "; cmd += tag;
25  if ( !geant_maker ) geant_maker = (St_geant_Maker *)chain->GetMaker("geant");
26  geant_maker -> LoadGeometry(cmd);
27  // if ( agml ) command("gexec $STAR_LIB/libxgeometry.so");
28 }
29 // ----------------------------------------------------------------------------
30 void command( TString cmd )
31 {
32  if ( !geant_maker ) geant_maker = (St_geant_Maker *)chain->GetMaker("geant");
33  geant_maker -> Do( cmd );
34 }
35 // ----------------------------------------------------------------------------
36 void trig( Int_t n=1 )
37 {
38  for ( Int_t i=0; i<n; i++ ) {
39 
40  // Clear the chain from the previous event
41  chain->Clear();
42  kinematics->Kine( 1, "J/psi", 5.0, 10.0, -2.0, 2.0 );
43  chain->Make();
44  // Print the event
45  _primary->event()->Print();
46  command("gprint kine");
47 
48  }
49 }
50 // ----------------------------------------------------------------------------
51 // ----------------------------------------------------------------------------
52 // ----------------------------------------------------------------------------
53 void Kinematics()
54 {
55 
56  // gSystem->Load( "libStarGeneratorPoolPythia6_4_23.so" );
57  gSystem->Load( "libKinematics.so");
58  kinematics = new StarKinematics();
59 
60  _primary->AddGenerator(kinematics);
61 }
62 // ----------------------------------------------------------------------------
63 // ----------------------------------------------------------------------------
64 // ----------------------------------------------------------------------------
65 void starsim( Int_t nevents=1, Int_t rngSeed=1234 )
66 {
67 
68  gROOT->ProcessLine(".L bfc.C");
69  {
70  TString simple = "y2012 geant gstar usexgeom agml debug"; bfc(0, simple );
71  }
72 
73  gSystem->Load( "libVMC.so");
74 
75  gSystem->Load( "StarGeneratorUtil.so" );
76  gSystem->Load( "StarGeneratorEvent.so" );
77  gSystem->Load( "StarGeneratorBase.so" );
78  gSystem->Load( "StarGeneratorDecay.so" );
79  gSystem->Load( "libPythia8_1_86.so");
80 
81  gSystem->Load( "libMathMore.so" );
82  gSystem->Load( "xgeometry.so" );
83 
84  // Setup RNG seed and map all ROOT TRandom here
85  StarRandom::seed( rngSeed );
87 
88 
89 
90  //
91  // Create the primary event generator and insert it
92  // before the geant maker
93  //
94  // StarPrimaryMaker *
95  _primary = new StarPrimaryMaker();
96  {
97  _primary -> SetFileName( "kinematics.starsim.root");
98  chain -> AddBefore( "geant", _primary );
99  }
100 
101  Kinematics();
102 
103 
104 
105  //
106  // Setup decay manager
107  //
108  StarDecayManager *decayMgr = AgUDecay::Manager();
109  StarPythia8Decayer *decayPy8 = new StarPythia8Decayer();
110  decayMgr->AddDecayer( 0, decayPy8 ); // Handle any decay requested
111  decayPy8->SetDebug(1);
112  decayPy8->Set("WeakSingleBoson:all = on");
113 
114 
115 
116  TString name;
117  double mass, lifetime, charge;
118  int tracktype, pdgcode, g3code;
119 
120  // Particle data
122  // One can add a particle to G3 using...
123  TParticlePDG* B_c_pl = data.GetParticle("B_c+");
124  TParticlePDG* B_c_mn = data.GetParticle("B_c-");
125  data.AddParticleToG3( B_c_pl, 50001 ); // pdg = 541
126  data.AddParticleToG3( B_c_mn, 50002 );
127  TParticlePDG* J_Psi = data.GetParticle("J/psi");
128  data.AddParticleToG3( J_Psi, 50003 ); // pdg = 443
129  decayPy8->Set("541:onMode = 0");
130  decayPy8->Set("541:onIfAll = 443 11 12"); // only the J/psi + e + nu modes
131  //decayPy8->Set("443:onMode = 0");
132 // decayPy8->Set("443:onIfAll = 11");
133 
134 
135 
136  //
137  // Initialize primary event generator and all sub makers
138  //
139  _primary -> Init();
140 
141  //
142  // Setup geometry and set starsim to use agusread for input
143  //
144  geometry("y2012");
145  command("gkine -4 0");
146  command("gfile o kinematics.starsim.fzd");
147 
148  //
149  // Trigger on nevents
150  //
151  trig( nevents );
152 
153  command("call agexit"); // Make sure that STARSIM exits properly
154 
155 }
156 // ----------------------------------------------------------------------------
157 
void Print(const Option_t *opts="head") const
Star Simple Kinematics Generator.
virtual void Clear(Option_t *option="")
User defined functions.
Definition: StChain.cxx:77
static StarParticleData & instance()
Returns a reference to the single instance of this class.
TParticlePDG * GetParticle(const Char_t *name)
Get a particle by name.
Interface to PDG information.
void AddGenerator(StarGenerator *gener)
Connects VMC to class(es) which handle particle decays.
Int_t Init()
Initialize generator.
virtual Int_t Make()
Definition: StChain.cxx:110
static void seed(UInt_t s)
Definition: StarRandom.cxx:119
Base class for event records.
Definition: StarGenEvent.h:81
StarGenEvent * event()
Return a pointer to the event.
Main steering class for event generation.
void AddDecayer(Int_t pdgid, TVirtualMCDecayer *decayer)
void SetDebug(int dbg=1)
Set the debug level.
static void capture()
Capture gRandom random number generator.
Definition: StarRandom.cxx:57
Sparse class to hold track kinematics.
void Kine(Int_t ntrack, const Char_t *type="pi+,pi-,K+,K-,proton,antiproton", Double_t ptlow=0.0, Double_t pthigh=500.0, Double_t ylow=-10.0, Double_t yhigh=+10.0, Double_t philow=0.0, Double_t phihigh=TMath::TwoPi())
TParticlePDG * AddParticleToG3(const char *name, const double mass, const double lifetime, const double charge, const int tracktype, const int pdgcode, const int g3code, const double *bratio=0, const int *mode=0)
void Set(const char *cmd)
Modify pythia8 behavior.