ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Fun4All_G4_Momentum_Projection.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file Fun4All_G4_Momentum_Projection.C
1 #ifndef FUN4ALL_G4_MOMENTUM_PROJECTION_C
2 #define FUN4ALL_G4_MOMENTUM_PROJECTION_C
3 
5 
8 
10 #include <g4main/PHG4ParticleGun.h>
11 #include <g4main/PHG4Reco.h>
13 
18 #include <fun4all/Fun4AllServer.h>
19 #include <fun4all/SubsysReco.h>
20 
21 #include <phool/recoConsts.h>
22 
23 R__LOAD_LIBRARY(libfun4all.so)
24 R__LOAD_LIBRARY(libg4testbench.so)
25 R__LOAD_LIBRARY(libg4detectors.so)
26 R__LOAD_LIBRARY(libg4trackfastsim.so)
27 
28 int Fun4All_G4_Momentum_Projection(const int nEvents = 1000, const string &evalfile = "FastTrackingEval.root", const string &outfile = "")
29 {
31  // Make the Server
34  se->Verbosity(0);
35 
37  // if you want to use a fixed seed for reproducible results
38  //rc->set_IntFlag("RANDOMSEED", 12345); // if you want to use a fixed seed
39  // PHG4ParticleGenerator generates particle
40  // distributions in eta/phi/mom range
41  PHG4ParticleGenerator *gen = new PHG4ParticleGenerator("PGENERATOR");
42  gen->set_name("pi-");
43  gen->set_vtx(0, 0, 0);
44  gen->set_eta_range(0.5, 0.5);
45  gen->set_mom_range(2, 2); // GeV/c
46  gen->set_phi_range(0., 90. / 180. * M_PI); // 0-90 deg
47 
48  se->registerSubsystem(gen);
49 
50  PHG4Reco *g4Reco = new PHG4Reco();
51  g4Reco->set_field(1.5); // 1.5 T solenoidal field
52 
53  double si_thickness[6] = {0.02, 0.02, 0.0625, 0.032, 0.032, 0.032};
54  double svxrad[6] = {2.71, 4.63, 11.765, 25.46, 41.38, 63.66};
55  double length[6] = {20., 20., 36., -1., -1., -1.}; // -1 use eta coverage to determine length
57  // here is our silicon:
58  for (int ilayer = 0; ilayer < 6; ilayer++)
59  {
60  cyl = new PHG4CylinderSubsystem("SVTX", ilayer);
61  cyl->set_double_param("radius", svxrad[ilayer]);
62  cyl->set_string_param("material", "G4_Si");
63  cyl->set_double_param("thickness", si_thickness[ilayer]);
64  cyl->SetActive();
65  cyl->SuperDetector("SVTX");
66  if (length[ilayer] > 0)
67  {
68  cyl->set_double_param("length", length[ilayer]);
69  }
70  g4Reco->registerSubsystem(cyl);
71  }
72 
73  // Black hole swallows everything - prevent loopers from returning
74  // to inner detectors
75  cyl = new PHG4CylinderSubsystem("BlackHole", 0);
76  cyl->set_double_param("radius", 80); // 80 cm
77  cyl->set_double_param("thickness", 0.1); // does not matter (but > 0)
78  cyl->SetActive();
79  cyl->BlackHole(); // eats everything
80  g4Reco->registerSubsystem(cyl);
81 
83  g4Reco->registerSubsystem(truth);
84 
85  se->registerSubsystem(g4Reco);
86 
87  //---------------------------
88  // fast pattern recognition and full Kalman filter
89  // output evaluation file for truth track and reco tracks are PHG4TruthInfoContainer
90  //---------------------------
91  PHG4TrackFastSim *kalman = new PHG4TrackFastSim("PHG4TrackFastSim");
92  kalman->set_use_vertex_in_fitting(false);
93  kalman->set_sub_top_node_name("SVTX");
94  kalman->set_trackmap_out_name("SvtxTrackMap");
95 
96  // add Si Trtacker
97  kalman->add_phg4hits(
98  "G4HIT_SVTX", // const std::string& phg4hitsNames,
99  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
100  300e-4, // radial-resolution [cm]
101  30e-4, // azimuthal-resolution [cm]
102  1, // z-resolution [cm]
103  1, // efficiency,
104  0 // noise hits
105  );
106 
107  kalman->add_cylinder_state("MyCylinder1", 2.); // projection onto cylinder with radius = 2cm
108  kalman->add_cylinder_state("MyCylinder2", 70.); // projection onto cylinder with radius = 70cm
109  kalman->add_zplane_state("MyPlane1", 100.); // projection onto z-plane at 100cm
110 
111  se->registerSubsystem(kalman);
112 
113  PHG4TrackFastSimEval *fast_sim_eval = new PHG4TrackFastSimEval("FastTrackingEval");
114  fast_sim_eval->set_filename(evalfile);
115  // Add the above projections to evaluation ntuple
116  fast_sim_eval->AddProjection("MyCylinder1");
117  fast_sim_eval->AddProjection("MyCylinder2");
118  fast_sim_eval->AddProjection("MyPlane1");
119 
120  se->registerSubsystem(fast_sim_eval);
121  //---------------------------
122 
123  //---------------------------
124  // output DST file for further offlien analysis
125  //---------------------------
126  if (!outfile.empty())
127  {
128  Fun4AllOutputManager *out = new Fun4AllDstOutputManager("DSTOUT", outfile);
129  se->registerOutputManager(out);
130  }
132  se->registerInputManager(in);
133 
134  if (nEvents > 0)
135  {
136  se->run(nEvents);
137  // finish job - close and save output files
138  se->End();
139  std::cout << "All done" << std::endl;
140 
141  // cleanup - delete the server and exit
142  delete se;
143  gSystem->Exit(0);
144  }
145  return 0;
146 }
147 
148 PHG4ParticleGenerator *get_gen(const char *name = "PGENERATOR")
149 {
152  return pgun;
153 }
154 
155 #endif