EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4_Tracking_Modular.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4_Tracking_Modular.C
1 #ifndef MACRO_G4TRACKINGMODULAR_C
2 #define MACRO_G4TRACKINGMODULAR_C
3 
4 #include <GlobalVariables.C>
5 
6 #include <G4_CEmc_EIC.C>
7 #include <G4_FEMC_EIC.C>
8 #include <G4_FHCAL.C>
9 #include <G4_GEM_EIC.C>
10 #include <G4_Mvtx_EIC.C>
11 #include <G4_TPC_EIC.C>
12 #include <G4_AllSilicon.C>
13 #include <G4_TTL_EIC.C>
14 
15 #include <g4trackfastsim/PHG4TrackFastSim.h>
16 #include <trackreco/PHRaveVertexing.h>
17 #include <g4trackfastsim/PHG4TrackFastSimEval.h>
18 #include <fun4all/Fun4AllServer.h>
19 
20 #include <vector>
21 
22 R__LOAD_LIBRARY(libtrack_reco.so)
23 R__LOAD_LIBRARY(libg4trackfastsim.so)
24 
25 namespace Enable
26 {
27  bool TRACKING = false;
28  bool TRACKING_EVAL = false;
29  bool TRACKING_EVAL_DETAILED = false;
30  bool TRACKING_INNER = false;
31  int TRACKING_VERBOSITY = 0;
32 } // namespace Enable
33 
34 namespace G4TRACKING
35 {
36  bool DISPLACED_VERTEX = false;
37  bool PROJECTION_EEMC = false;
38  bool PROJECTION_CEMC = false;
39  bool PROJECTION_FEMC = false;
40  bool PROJECTION_FHCAL = false;
41  bool PROJECTION_EHCAL = false;
42  bool PROJECTION_DRCALO = false;
43 } // namespace G4TRACKING
44 
45 namespace TRACKING
46 {
47  std::string TrackNodeNameInner = "InnerTrackMap";
48 } // namespace TRACKING
49 
50 //-----------------------------------------------------------------------------//
51 void TrackingInit()
52 {
53  TRACKING::TrackNodeName = "TrackMap";
54  TRACKING::TrackNodeNameInner = "InnerTrackMap";
55 }
56 //-----------------------------------------------------------------------------//
57 void Tracking_Reco(TString specialSetting = "")
58 {
59  int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
60  //---------------
61  // Fun4All server
62  //---------------
63 
64  Fun4AllServer *se = Fun4AllServer::instance();
65 
66  PHG4TrackFastSim *kalman = new PHG4TrackFastSim("PHG4TrackFastSim");
67  kalman->Verbosity(verbosity);
68  // kalman->Smearing(false);
69 
70  // Store a list of FastSim objects. This way, we can uniformly access multiple objects,
71  // setting the same objets on each when appropriate.
72  std::vector<PHG4TrackFastSim *> kalmanTrackers = {kalman};
73 
74  // Option to perform an additional, separate tracking using only the inner detectors.
75  PHG4TrackFastSim *kalmanInnerTracking = nullptr;
76  if (Enable::TRACKING_INNER) {
77  kalmanInnerTracking = new PHG4TrackFastSim("PHG4TrackFastSimInnerTracking");
78  kalmanInnerTracking->Verbosity(verbosity);
79  // kalmanInnerTracking->Smearing(false);
80  // Add to the list of trackers so we can uniformly access it as needed.
81  kalmanTrackers.push_back(kalmanInnerTracking);
82  }
83  for (auto k : kalmanTrackers) {
84  if (G4TRACKING::DISPLACED_VERTEX){
85  // do not use truth vertex in the track fitting,
86  // which would lead to worse momentum resolution for prompt tracks
87  // but this allows displaced track analysis including DCA and vertex finding
88  k->set_use_vertex_in_fitting(false);
89  k->set_vertex_xy_resolution(0); // do not smear the vertex used in the built-in DCA calculation
90  k->set_vertex_z_resolution(0); // do not smear the vertex used in the built-in DCA calculation
91  k->enable_vertexing(true); // enable vertex finding and fitting
92  } else {
93  // constraint to a primary vertex and use it as part of the fitting level arm
94  k->set_use_vertex_in_fitting(true);
95  k->set_vertex_xy_resolution(50e-4);
96  k->set_vertex_z_resolution(50e-4);
97  }
98  }
99 
100  kalman->set_sub_top_node_name("TRACKS");
101  kalman->set_trackmap_out_name(TRACKING::TrackNodeName);
102  if (Enable::TRACKING_INNER) {
103  kalmanInnerTracking->set_sub_top_node_name("TRACKS");
104  kalmanInnerTracking->set_trackmap_out_name(TRACKING::TrackNodeNameInner);
105  }
106 
107  //-------------------------
108  // Barrel upgrade (LANL)
109  //-------------------------
110  // Different Barrel versions documented in arXiv:2009.0288
111  if (Enable::BARREL){
112  int nLayer = 5;
113  if (G4BARREL::SETTING::BARRELV6)
114  nLayer = 4;
115  float_t pitch = 20e-4; // default pitch size
116  double r[6] = { 3.64, 4.81, 5.98, 16.0, 22.0, -1}; //cm
117  if (specialSetting.Contains("BARRELV4")){
118  nLayer = 6;
119  r[3] = 9.2;
120  r[4] = 17.;
121  r[5] = 27.;
122  }
123 
124  for (Int_t i = 0; i < nLayer; i++){
125  for (auto k : kalmanTrackers) {
126  k->add_phg4hits(Form("G4HIT_BARREL_%d", i), // const std::string& phg4hitsNames,
127  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
128  999., // const float radres, (not used in cylindrical geom)
129  pitch/sqrt(12), // const float phires,
130  pitch/sqrt(12), // const float lonres,
131  1, // const float eff,
132  0); // const float noise
133  if (Enable::TRACKING_EVAL_DETAILED){
134  k->add_cylinder_state(Form("BARREL_%d", i), r[i]);
135  }
136  }
137  }
138  }
139 
140  //-------------------------
141  // FST (LANL)
142  //-------------------------
143  // Different FST versions documented in arXiv:2009.0288
144  if (Enable::FST) {
145  // FORWARD DISKS
146  Int_t nDisks = 5;
147  float zFWDdisks[6] = {35, 53, 77, 101, 125, 270};
148  if (specialSetting.Contains("FSTV4"))
149  nDisks = 6;
150  if (specialSetting.Contains("FSTV2"))
151  zFWDdisks[4] = 270;
152  Float_t pitch = 20e-4; // default pitch size
153  Int_t llargerPitch = -1; // layer number above which pitch size is higher, if -1 all initialized with small pitch
154  if (specialSetting.Contains("FSTV3") || specialSetting.Contains("FSTV41"))
155  llargerPitch = 3;
156  if (specialSetting.Contains("FSTV42"))
157  llargerPitch = 4;
158 
159  // intializing hits for different layers
160  for (int i = 0; i < nDisks; i++) {
161  if (llargerPitch != -1 && i >= llargerPitch)
162  pitch = 36.4e-4;
163  for (auto k : kalmanTrackers) {
164  k->add_phg4hits(Form("G4HIT_FST_%d", i), // const std::string& phg4hitsNames,
165  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
166  pitch/sqrt(12), // const float radres,
167  pitch/sqrt(12), // const float phires,
168  999., // const float lonres, number not used in vertical plane geometry
169  1, // const float eff,
170  0); // const float noise
171  if (Enable::TRACKING_EVAL_DETAILED){
172  k->add_zplane_state(Form("FST_%d", i), zFWDdisks[i]);
173  }
174  }
175  }
176  }
177 
178  //-------------------------
179  // Barrel ALLSILICON version (LBL)
180  //-------------------------
181  char nodename[100];
182  if (Enable::ALLSILICON){
183  // CENTRAL BARREL
184  float_t pitch = 10e-4;
185  float rBarrel[6] = {3.3, 5.7, 21.0, 22.68, 39.30, 43.23};
186  for (int i = 10; i < 16; i++) {
187  sprintf(nodename, "G4HIT_LBLVTX_CENTRAL_%d", i);
188 
189  for (auto k : kalmanTrackers) {
190  k->add_phg4hits(
191  nodename, // const std::string& phg4hitsNames
192  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype
193  999., // radial-resolution [cm] (this number is not used in cylindrical geometry)
194  pitch/sqrt(12), // azimuthal (arc-length) resolution [cm]
195  pitch/sqrt(12), // longitudinal (z) resolution [cm]
196  1, // efficiency (fraction)
197  0 // hit noise
198  );
199  if (Enable::TRACKING_EVAL_DETAILED){
200  k->add_cylinder_state(Form("LBLVTX_CENTRAL_%d", i), rBarrel[i-10]);
201  }
202  }
203  }
204  // FORWARD DISKS
205  float zFWDdisks[5] = {25, 49, 73, 97, 121};
206  for (int i = 20; i < 25; i++) {
207  sprintf(nodename, "G4HIT_LBLVTX_FORWARD_%d", i);
208  for (auto k : kalmanTrackers) {
209  k->add_phg4hits(
210  nodename, // const std::string& phg4hitsNames
211  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype
212  pitch/sqrt(12), // radial-resolution [cm]
213  pitch/sqrt(12), // azimuthal (arc-length) resolution [cm]
214  999., // longitudinal (z) resolution [cm] (this number is not used in vertical plane geometry)
215  1, // efficiency (fraction)
216  0 // hit noise
217  );
218  if (Enable::TRACKING_EVAL_DETAILED){
219  k->add_zplane_state(Form("LBLVTX_FORWARD_%d", i), zFWDdisks[i-20]);
220  }
221  }
222  }
223 
224  // BACKWARD DISKS
225  float zBWDdisks[5] = {-25, -49, -73, -97, -121};
226  for (int i = 30; i < 35; i++) {
227  sprintf(nodename, "G4HIT_LBLVTX_BACKWARD_%d", i);
228  for (auto k : kalmanTrackers) {
229  k->add_phg4hits(
230  nodename, // const std::string& phg4hitsNames
231  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype
232  pitch/sqrt(12), // radial-resolution [cm]
233  pitch/sqrt(12), // azimuthal (arc-length) resolution [cm]
234  999., // longitudinal (z) resolution [cm] (this number is not used in vertical plane geometry)
235  1, // efficiency (fraction)
236  0 // hit noise
237  );
238  if (Enable::TRACKING_EVAL_DETAILED){
239  k->add_zplane_state(Form("LBLVTX_BACKWARD_%d", i), zBWDdisks[i-30]);
240  }
241  }
242  }
243  }
244 
245  //-------------------------
246  // Timing tracking Layer (TTL - ORNL/Rice)
247  //-------------------------
248  // position resol improvement
249  float posResImp = sqrt(12);
250  float pitchLGAD = 500e-4;
251  if (G4TTL::SETTING::optionGran == 2)
252  posResImp = sqrt(256);
253  else if (G4TTL::SETTING::optionGran == 3)
254  pitchLGAD = 1300e-4;
255  float resLGAD = pitchLGAD/posResImp;
256 
257  // central barrel
258  if (Enable::CTTL){
259  float resLGAD_barrel = resLGAD;
260  if(G4TTL::SETTING::optionLYSO){
261  resLGAD_barrel = 35e-1; // https://cds.cern.ch/record/2667167/files/CMS-TDR-020.pdf page 33 bottom
262  }
263  for (int i = 0; i < G4TTL::layer[1]; i++){
264  kalman->add_phg4hits(Form("G4HIT_CTTL_%d", i), // const std::string& phg4hitsNames,
265  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
266  999, // const float radres,
267  resLGAD_barrel, // const float phires,
268  resLGAD_barrel, // const float lonres, *ignored in plane detector*
269  0.95, // const float eff,
270  0); // const float noise
271  kalman -> add_cylinder_state(Form("CTTL_%d",i), G4TTL::positionToVtx[1][i]);
272  }
273  }
274 
275  // electron going direction
276  if (Enable::ETTL){
277  for (int i = 0; i < G4TTL::layer[0]; i++){
278  kalman->add_phg4hits(Form("G4HIT_ETTL_%d", i), // const std::string& phg4hitsNames,
279  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
280  resLGAD, // const float radres,
281  resLGAD, // const float phires,
282  999., // const float lonres, *ignored in plane detector*
283  0.95, // const float eff,
284  0); // const float noise
285  kalman -> add_zplane_state(Form("ETTL_%d",i), G4TTL::positionToVtx[0][i]);
286  }
287  }
288 
289  // forward hadron going direction
290  if (Enable::FTTL){
291  for (int i = 0; i < G4TTL::layer[2]; i++){
292  kalman->add_phg4hits(Form("G4HIT_FTTL_%d", i), // const std::string& phg4hitsNames,
293  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
294  resLGAD, // const float radres,
295  resLGAD, // const float phires,
296  999., // const float lonres, *ignored in plane detector*
297  0.95, // const float eff,
298  0); // const float noise
299  kalman -> add_zplane_state(Form("FTTL_%d",i), G4TTL::positionToVtx[2][i]);
300  }
301  }
302 
303  //-------------------------
304  // MVTX (sPHENIX) with better resol
305  //-------------------------
306  if (Enable::MVTX){
307  // MAPS
308  kalman->add_phg4hits(
309  "G4HIT_MVTX", // const std::string& phg4hitsNames,
310  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
311  5e-4, // const float radres,
312  5e-4, // const float phires,
313  5e-4, // const float lonres,
314  1, // const float eff,
315  0 // const float noise
316  );
317  }
318  //-------------------------
319  // TPC
320  //-------------------------
321  if (Enable::TPC){
322  kalman->add_phg4hits(
323  "G4HIT_TPC", // const std::string& phg4hitsNames,
324  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
325  1, // const float radres,
326  200e-4, // const float phires,
327  500e-4, // const float lonres,
328  1, // const float eff,
329  0 // const float noise
330  );
331  }
332 
333  //-------------------------
334  // EGEM
335  //-------------------------
336  if (Enable::EGEM){
337  // GEM, 70um azimuthal resolution, 1cm radial strips
338  int minEGEM = 0;
339  if (!Enable::EGEM_FULL) // all to disable inner most discs
340  minEGEM = 2;
341  for (int i = minEGEM; i < 4; i++){
342  kalman->add_phg4hits(
343  Form("G4HIT_EGEM_%d", i), // const std::string& phg4hitsNames,
344  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
345  1. / sqrt(12.), // const float radres,
346  70e-4, // const float phires,
347  100e-4, // const float lonres,
348  1, // const float eff,
349  0 // const float noise
350  );
351  }
352  }
353  //-------------------------
354  // FGEM
355  //-------------------------
356  if (Enable::FGEM) {
357  int minFGEM = 0;
358  if (Enable::FGEM_ORIG){
359  minFGEM = 0;
360  } else {
361  minFGEM = 2;
362  }
363  // GEM2, 70um azimuthal resolution, 1cm radial strips
364  for (int i = minFGEM; i < 5; i++) {
365  kalman->add_phg4hits(Form("G4HIT_FGEM_%d", i), // const std::string& phg4hitsNames,
366  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
367  1. / sqrt(12.), // const float radres,
368  70e-4, // const float phires,
369  100e-4, // const float lonres,
370  1, // const float eff,
371  0); // const float noise
372  }
373  }
374  //-------------------------
375  // FEMC
376  //-------------------------
377  // Saved track states (projections)
378  if (Enable::FEMC){
379  kalman->add_state_name("FEMC");
380  // kalman->add_zplane_state("FEMC", 310);
381  }
382 
383  //-------------------------
384  // DRCALO
385  //-------------------------
386  if (Enable::DRCALO) {
387  kalman -> add_zplane_state("DRCALO_0", 300);
388  }
389 
390  //-------------------------
391  // FHCAL
392  //-------------------------
393  if (Enable::FHCAL) {
394  kalman->add_state_name("FHCAL");
395  }
396 
397  //-------------------------
398  // LFHCAL
399  //-------------------------
400  if (Enable::LFHCAL) {
401  kalman->add_state_name("LFHCAL");
402  }
403  //-------------------------
404  // CEMC
405  //-------------------------
406 
407  if (Enable::CEMC){
408  kalman->add_state_name("CEMC");
409  }
410 
411  //-------------------------
412  // BECAL
413  //-------------------------
414 
415  if (Enable::BECAL){
416  kalman->add_state_name("BECAL");
417  }
418 
419  //-------------------------
420  // HCALIN
421  //-------------------------
422 
423  if (Enable::HCALIN){
424  kalman->add_state_name("HCALIN");
425  }
426 
427 
428  //-------------------------
429  // HCALOUT
430  //-------------------------
431 
432  if (Enable::HCALOUT){
433  kalman->add_state_name("HCALOUT");
434  }
435 
436  //-------------------------
437  // EEMC
438  //-------------------------
439  if ((Enable::EEMC || Enable::EEMCH ) )
440  {
441  kalman->add_state_name("EEMC");
442  }
443 
444 
445  //-------------------------
446  // EEMC
447  //-------------------------
448  if ((Enable::EEMCH ) )
449  {
450 // kalman->add_state_name("EEMCG");
451  kalman->add_zplane_state("EEMCG", -210);
452  }
453 
454  //-------------------------
455  // EEMC
456  //-------------------------
457  if ((Enable::EHCAL ) )
458  {
459  kalman->add_state_name("EHCAL");
460  }
461 
462  //-------------------------
463  // DIRC
464  //-------------------------
465  if ((Enable::DIRC ) )
466  {
467  // kalman->add_state_name("DIRC");
468  kalman->add_zplane_state("DIRC", 185);
469  }
470 
471 
472  for (auto k : kalmanTrackers) {
473  se->registerSubsystem(k);
474  }
475  return;
476 }
477 
478 //-----------------------------------------------------------------------------//
479 
480 void Tracking_Eval(const std::string &outputfile, TString specialSetting = "")
481 {
482  int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
483  //---------------
484  // Fun4All server
485  //---------------
486  Fun4AllServer *se = Fun4AllServer::instance();
487 
488  //----------------
489  // Fast Tracking evaluation
490  //----------------
491  PHG4TrackFastSimEval *fast_sim_eval = new PHG4TrackFastSimEval("FastTrackingEval");
492  fast_sim_eval->set_trackmapname(TRACKING::TrackNodeName);
493 
494  // save tracking planes in forward direction for single layer studies
495  if (Enable::TRACKING_EVAL_DETAILED){
496  if (Enable::FST){
497  Int_t nDisks = 5;
498  if (specialSetting.Contains("FSTV4"))
499  nDisks = 6;
500  for (Int_t i = 0; i < nDisks; i++){
501  fast_sim_eval->AddProjection(Form("FST_%d",i));
502  }
503  }
504  if (Enable::BARREL){
505  Int_t nLayer = 5;
506  if (specialSetting.Contains("BARRELV4"))
507  nLayer = 6;
508 
509  for (Int_t i = 0; i < nLayer; i++){
510  fast_sim_eval->AddProjection(Form("BARREL_%d",i));
511  }
512  }
513 
514  if (Enable::ALLSILICON){
515  for (int i = 10; i < 16; i++) {
516  fast_sim_eval->AddProjection(Form("LBLVTX_CENTRAL_%d", i));
517  }
518  for (int i = 20; i < 25; i++) {
519  fast_sim_eval->AddProjection(Form("LBLVTX_FORWARD_%d", i));
520  }
521  for (int i = 30; i < 35; i++) {
522  fast_sim_eval->AddProjection(Form("LBLVTX_BACKWARD_%d", i));
523  }
524  }
525  }
526 
527  // create projections on timing layers to read t and spacial coordinates
528  if (Enable::FTTL){
529  for (int l = 0; l < G4TTL::layer[2]; l++)
530  fast_sim_eval->AddProjection(Form("FTTL_%d",l));
531  }
532  if (Enable::ETTL){
533  for (int l = 0; l < G4TTL::layer[0]; l++)
534  fast_sim_eval->AddProjection(Form("ETTL_%d",l));
535  }
536  if (Enable::CTTL){
537  for (int l = 0; l < G4TTL::layer[1]; l++)
538  fast_sim_eval->AddProjection(Form("CTTL_%d",l));
539  }
540 
541  if(Enable::FHCAL && G4TRACKING::PROJECTION_FHCAL) fast_sim_eval->AddProjection("FHCAL");
542  if(Enable::FEMC && G4TRACKING::PROJECTION_FEMC) fast_sim_eval->AddProjection("FEMC");
543  // if(Enable::DRCALO && G4TRACKING::PROJECTION_DRCALO) fast_sim_eval->AddProjection("DRCALO_0");
544 
545  if(Enable::DIRC) fast_sim_eval->AddProjection("DIRC");
546 
547  // write to output file
548  fast_sim_eval->set_filename(outputfile);
549  se->registerSubsystem(fast_sim_eval);
550 }
551 #endif