ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4_TTL_EIC.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4_TTL_EIC.C
1 #ifndef MACRO_G4TTLEIC_C
2 #define MACRO_G4TTLEIC_C
3 
4 #include "GlobalVariables.C"
5 
8 
9 #include <g4main/PHG4Reco.h>
10 
11 #include <string>
12 
13 R__LOAD_LIBRARY(libg4detectors.so)
14 
15 int make_forward_station(string name, PHG4Reco *g4Reco, double zpos, double Rmin,
16  double Rmax,double tSilicon, double xoffset=0);
17 int make_forward_station_basic(string name, PHG4Reco *g4Reco, double zpos, double Rmin,
18  double Rmax,double tSilicon);
19 int make_barrel_layer_basic(string name, PHG4Reco *g4Reco,
20  double radius, double halflength, double tSilicon, double zOffset);
21 int make_barrel_layer_LYSO_basic(string name, PHG4Reco *g4Reco,
22  double radius, double halflength, double tSilicon, double zOffset);
23 int make_barrel_layer(string name, PHG4Reco *g4Reco,
24  double radius, double halflength, double tSilicon, double zOffset);
25 
26 //-----------------------------------------------------------------------------------//
27 namespace Enable
28 {
29  bool FTTL = false;
30  bool ETTL = false;
31  bool CTTL = false;
32 }
33 
34 namespace G4TTL
35 {
36  int layer[3] = { 2, 1, 2};
37  double positionToVtx[3][3] = { {-169., -172., -309.5}, {80., 114.7, 0. }, { 287., 289., 340.} };
38  double minExtension[3][3] = { {8, 8, 15.3}, {218, 180, 0 }, {11.62, 11.7, 13.8 } };
39  double maxExtension[3][3] = { {61., 61. , 200}, {-40, 0, 0 }, {170., 170., 250 } };
40  double xoffsetFTTLIP6[3] = { -6., -6., -6.};
41  double xoffsetFTTLIP8[3] = { 8.4, 8.4, 8.4};
42  namespace SETTING
43  {
44  bool optionCEMC = false;
45  bool optionEEMCH = true;
46  bool optionBasicGeo = false;
47  bool optionLYSO = false;
48  int optionDR = 0;
49  int optionGeo = 7;
50  int optionGran = 1;
51  } // namespace SETTING
52 
53 
54 
55  // 2, LGAD based ToF by Wei Li:
56  // Present the recent simulation studies and detector layout for the LGAD based ToF,
57  // which can be used as an outer tracker. The current detector can reach around 30 micron spatial
58  // resolution with 500 micron AC-LGAD sensor.
59  // This detector can provide better coverage and provide further constraints in track
60  // reconstruction within high pseudorapidity regions.
61 
62  const double PositionResolution(30e-4);
63 
64 } // namespace G4TTL
65 
66 
67 //-----------------------------------------------------------------------------------//
68 void TTL_Init()
69 {
72 
74  G4TTL::positionToVtx[0][0] = -155.5;
75  G4TTL::positionToVtx[0][1] = -158.5;
76  }
77 
79  G4TTL::maxExtension[0][0] = 80.;
80  G4TTL::maxExtension[0][1] = 80.;
81  G4TTL::positionToVtx[1][0] = 92.;
83  }
84 
86  G4TTL::minExtension[1][0] = 180;
87  G4TTL::maxExtension[1][0] = -25;
88  }
89 
90  if (G4TTL::SETTING::optionGeo == 1){
91  cout << "TTL setup infront of ECals with 2 layers fwd/bwd & 1 layer barrel" << endl;
92  }
93  if (G4TTL::SETTING::optionGeo == 2){
94  cout << "TTL setup infront of ECals with 2 layers fwd/bwd & 1 layer barrel, lower barrel layer" << endl;
95  G4TTL::positionToVtx[1][0] = 50.;
97  G4TTL::minExtension[1][0] = 100.;
98  G4TTL::maxExtension[1][0] = 0.;
99  }
101  cout << "TTL setup infront of ECals with 1 layers fwd/bwd & 1 layer barrel, lower barrel layer" << endl;
102  G4TTL::positionToVtx[1][0] = 50.;
104  G4TTL::minExtension[1][0] = 100.;
105  G4TTL::maxExtension[1][0] = 0.;
106  G4TTL::layer[0] = 1;
107  G4TTL::layer[2] = 1;
114  }
115  if (G4TTL::SETTING::optionGeo == 4){
116  cout << "TTL setup infront of ECals with 2 layers fwd/bwd & 1 layer barrel, 1 layer before HCals everywhere" << endl;
117  G4TTL::layer[0] = 3;
118  G4TTL::layer[1] = 2;
120  G4TTL::layer[2] = 3;
121  }
122  if (G4TTL::SETTING::optionGeo == 5){
123  // Option 5 is 1 layer fwd/bwd, with LYSO in the central barrel. We use the geometry for option 3 since it's the same.
124  // However, we create a separate option to simplify setting up the configuration.
125  cout << "TTL setup using LYSO in central barrel" << endl;
127  }
128  if(G4TTL::SETTING::optionGeo == 7){
129  cout << "TTL one forward disk in front of dRICH and one backward disk in front of EEMC, barrel CTTL center at radius 64cm" << endl;
130  // single disk in front of dRICH (full eta)
131  G4TTL::layer[2] = 1;
132  G4TTL::minExtension[2][0] = 7.0;
133  G4TTL::maxExtension[2][0] = 87;
134  G4TTL::positionToVtx[2][0] = 182.;
135  G4TTL::xoffsetFTTLIP6[0] = -2.7;
136  G4TTL::xoffsetFTTLIP8[0] = 3.0;
137 
138  // single disk in front of EEMC
139  G4TTL::layer[0] = 1;
140  G4TTL::maxExtension[0][0] = 64;
141 
142  // barrel layer at 64cm
143  G4TTL::positionToVtx[1][0] = 64.;
144  G4TTL::minExtension[1][0] = 140;
145  G4TTL::maxExtension[1][0] = 0;
146  }
147  if(G4TTL::SETTING::optionGeo == 8){
148  G4TTL::layer[2] = 3;
149  cout << "TTL forward disk 1 in front of dRICH" << endl;
150  G4TTL::minExtension[2][0] = 7.0;
151  G4TTL::maxExtension[2][0] = 87;
152  G4TTL::positionToVtx[2][0] = 182.;
153  cout << "additional two small TTL disks in front of FEMC" << endl;
154  G4TTL::minExtension[2][1] = 11.62;
155  G4TTL::minExtension[2][2] = 11.7;
156  G4TTL::maxExtension[2][1] = 60.;
157  G4TTL::maxExtension[2][2] = 60.;
158  G4TTL::positionToVtx[2][1] = 287.;
159  G4TTL::positionToVtx[2][2] = 289.;
160 
161  G4TTL::layer[0] = 1;
162  G4TTL::maxExtension[0][0] = 64;
163 
164  G4TTL::positionToVtx[1][0] = 64.;
165  G4TTL::minExtension[1][0] = 140;
166  G4TTL::maxExtension[1][0] = 0;
167  }
168 
170  cout << "TTL setup infront of ECals with 2 layers fwd/bwd & 1 layer barrel, 1 layer before HCals everywhere choosen!" << endl;
171  cout << "conflicting in forward region with DR calo, reducing by 1 layer!" << endl;
172  G4TTL::layer[2] = 2;
173  } else if (G4TTL::SETTING::optionDR == 1 && G4TTL::SETTING::optionGeo == 4 ){
174  cout << "TTL setup infront of ECals with 2 layers fwd/bwd & 1 layer barrel, 1 layer before HCals everywhere choosen!" << endl;
175  cout << "conflicting in forward region with DR calo, reducing adding larger cutout!" << endl;
176  G4TTL::minExtension[2][2] = 2.5;
177  }
178 
179 }
180 //-----------------------------------------------------------------------------------//
181 void FTTLSetup(PHG4Reco *g4Reco, TString fttloption = "")
182 {
183  const double cm = PHG4Sector::Sector_Geometry::Unit_cm();
184  const double mm = .1 * cm;
185  const double um = 1e-3 * mm;
186 
187  for (Int_t i = 0; i < G4TTL::layer[2]; i++){
188  cout << G4TTL::positionToVtx[2][i] << "\t" << G4TTL::minExtension[2][i] << "\t" << G4TTL::maxExtension[2][i] << endl;
191  } else {
192  make_forward_station_basic(Form("FTTL_%d", i), g4Reco, G4TTL::positionToVtx[2][i], G4TTL::minExtension[2][i], G4TTL::maxExtension[2][i], 85*um);
193  }
194  }
195 }
196 
197 
198 //-----------------------------------------------------------------------------------//
199 void ETTLSetup(PHG4Reco *g4Reco, TString ettloption = "")
200 {
201  const double cm = PHG4Sector::Sector_Geometry::Unit_cm();
202  const double mm = .1 * cm;
203  const double um = 1e-3 * mm;
204  for (Int_t i = 0; i < G4TTL::layer[0]; i++){
205  cout << G4TTL::positionToVtx[0][i] << "\t" << G4TTL::minExtension[0][i] << "\t" << G4TTL::maxExtension[0][i] << endl;
207  make_forward_station(Form("ETTL_%d", i), g4Reco, G4TTL::positionToVtx[0][i], G4TTL::minExtension[0][i], G4TTL::maxExtension[0][i], 85*um, 0);
208  } else {
209  make_forward_station_basic(Form("ETTL_%d", i), g4Reco, G4TTL::positionToVtx[0][i], G4TTL::minExtension[0][i], G4TTL::maxExtension[0][i], 85*um);
210  }
211  }
212 }
213 
214 //-----------------------------------------------------------------------------------//
215 void CTTLSetup(PHG4Reco *g4Reco, TString cttloption = "")
216 {
217  const double cm = PHG4Sector::Sector_Geometry::Unit_cm();
218  const double mm = .1 * cm;
219  const double um = 1e-3 * mm;
220 
221  for (Int_t i = 0; i < G4TTL::layer[1]; i++){
222  cout << "Radius: " << G4TTL::positionToVtx[1][i] << "\tLength: " << G4TTL::minExtension[1][i] << "\tz-Offset: " << G4TTL::maxExtension[1][i] << endl;
224  make_barrel_layer_LYSO_basic(Form("CTTL_%d",i), g4Reco, G4TTL::positionToVtx[1][i], G4TTL::minExtension[1][i], 85*um, G4TTL::maxExtension[1][i]);
226  make_barrel_layer_basic(Form("CTTL_%d",i), g4Reco, G4TTL::positionToVtx[1][i], G4TTL::minExtension[1][i], 85*um, G4TTL::maxExtension[1][i]);
227  } else {
228  make_barrel_layer(Form("CTTL_%d",i), g4Reco, G4TTL::positionToVtx[1][i], G4TTL::minExtension[1][i], 85*um, G4TTL::maxExtension[1][i]);
229  }
230  }
231 }
232 
233 
234 //-----------------------------------------------------------------------------------//
235 int make_forward_station(string name, PHG4Reco *g4Reco,
236  double zpos, double rMin, double rMax,
237  double tSilicon, //silicon thickness
238  double xoffset = 0 )
239 {
240  cout << "r min: " << rMin << "\t r max: " << rMax << "\t z: " << zpos << endl;
241 
242  // always facing the interaction point
243  double polar_angle = 0;
244  double place_z(zpos);
245  if (place_z < 0){
246  place_z = -place_z;
247  polar_angle = M_PI;
248  }
250  ttl = new PHG4TTLSubsystem(name);
251  ttl->SetDetailed(false);
252  ttl->SuperDetector(name);
253  ttl->set_double_param("polar_angle", polar_angle); //
254  ttl->set_int_param("isForward", 1); //
255  ttl->set_double_param("place_z", place_z * cm); //
256  ttl->set_double_param("rMin", rMin * cm); //
257  ttl->set_double_param("rMax", rMax * cm); //
258  ttl->set_double_param("offset_x", xoffset * cm); //
259  ttl->set_double_param("tSilicon", tSilicon); //
260  // ttl->OverlapCheck(true);
262 
263  g4Reco->registerSubsystem(ttl);
264 
265 
267  {
268  TRACKING::FastKalmanFilter->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
269  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
270  G4TTL::PositionResolution, // const float radres,
271  G4TTL::PositionResolution, // const float phires,
272  tSilicon / sqrt(12.), // const float lonres, *ignored in plane detector*
273  1, // const float eff,
274  0); // const float noise
276  TRACKING::ProjectionNames.insert(name);
277  }
279  {
280  TRACKING::FastKalmanFilterInnerTrack->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
281  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
282  G4TTL::PositionResolution, // const float radres,
283  G4TTL::PositionResolution, // const float phires,
284  tSilicon / sqrt(12.), // const float lonres, *ignored in plane detector*
285  1, // const float eff,
286  0); // const float noise
287  }
288 
289  return 0;
290 }
291 
292 
293 //-----------------------------------------------------------------------------------//
295  double zpos, double rmin, double rmax,
296  double tSilicon) //silicon thickness
297 {
298  // cout
299  // << "make_GEM_station - GEM construction with PHG4SectorSubsystem - make_GEM_station_EdgeReadout of "
300  // << name << endl;
301 
302  // always facing the interaction point
303  double etamin = -TMath::Log(TMath::Tan(rmin/TMath::Abs(zpos)/2));
304  double etamax = -TMath::Log(TMath::Tan(rmax/TMath::Abs(zpos)/2));
305 
306  double polar_angle = 0;
307  if (zpos < 0){
308  zpos = -zpos;
309  polar_angle = M_PI;
310  etamin = -etamin;
311  etamax = -etamax;
312  }
313  if (etamax < etamin){
314  double t = etamax;
315  etamax = etamin;
316  etamin = t;
317  }
318 
320  ttl = new PHG4SectorSubsystem(name);
321 
322  ttl->SuperDetector(name);
323 
324  ttl->get_geometry().set_normal_polar_angle(polar_angle);
330  ttl->get_geometry().set_N_Sector(1);
331  ttl->get_geometry().set_material("G4_AIR");
333 
334  const double cm = PHG4Sector::Sector_Geometry::Unit_cm();
335  const double mm = .1 * cm;
336  const double um = 1e-3 * mm;
337  // build up layers
338 
339  ttl->get_geometry().AddLayer("SiliconSensor", "G4_Si", tSilicon, true, 100);
340  ttl->get_geometry().AddLayer("Metalconnection", "G4_Al", 100 * um, false, 100);
341  ttl->get_geometry().AddLayer("HDI", "G4_KAPTON", 20 * um, false, 100);
342  ttl->get_geometry().AddLayer("Cooling", "G4_WATER", 100 * um, false, 100);
343  ttl->get_geometry().AddLayer("Support", "G4_GRAPHITE", 50 * um, false, 100);
344  ttl->get_geometry().AddLayer("Support_Gap", "G4_AIR", 1 * cm, false, 100);
345  ttl->get_geometry().AddLayer("Support2", "G4_GRAPHITE", 50 * um, false, 100);
346 
347  g4Reco->registerSubsystem(ttl);
348 
349 
351  {
352  TRACKING::FastKalmanFilter->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
353  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
354  G4TTL::PositionResolution, // const float radres,
355  G4TTL::PositionResolution, // const float phires,
356  tSilicon / sqrt(12.), // const float lonres, *ignored in plane detector*
357  1, // const float eff,
358  0); // const float noise
360  TRACKING::ProjectionNames.insert(name);
361  }
363  {
364  TRACKING::FastKalmanFilterInnerTrack->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
365  PHG4TrackFastSim::Vertical_Plane, // const DETECTOR_TYPE phg4dettype,
366  G4TTL::PositionResolution, // const float radres,
367  G4TTL::PositionResolution, // const float phires,
368  tSilicon / sqrt(12.), // const float lonres, *ignored in plane detector*
369  1, // const float eff,
370  0); // const float noise
371  }
372 
373  return 0;
374 }
375 
376 //-----------------------------------------------------------------------------------//
377 int make_barrel_layer(string name, PHG4Reco *g4Reco,
378  double radius, double halflength, double tSilicon, double zOffset )
379 {
380  // cout << "r min: " << rMin << "\t r max: " << rMax << "\t z: " << zpos << endl;
382  cout << "The improved barrel TTL layer is placed at a fixed radius of rMin = 80 * cm and is meant to only be used with the BECAL!" << endl;
383  }
385  ttl = new PHG4TTLSubsystem(name);
386  ttl->SetDetailed(false);
387  ttl->SuperDetector(name);
388  ttl->set_int_param("isForward", 0); //
389  ttl->set_double_param("place_z", zOffset * cm); //
390  ttl->set_double_param("rMin", radius * cm); //
391  ttl->set_double_param("length", 2.0 * halflength * cm);
392  ttl->set_double_param("tSilicon", tSilicon); //
394 
395  g4Reco->registerSubsystem(ttl);
396 
397 
399  {
400  TRACKING::FastKalmanFilter->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
401  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
402  tSilicon / sqrt(12.), // const float radres,
403  G4TTL::PositionResolution, // const float phires,
404  G4TTL::PositionResolution, // const float lonres,
405  1, // const float eff,
406  0); // const float noise
408 
409  TRACKING::ProjectionNames.insert(name);
410  }
412  {
413  TRACKING::FastKalmanFilterInnerTrack->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
414  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
415  tSilicon / sqrt(12.), // const float radres,
416  G4TTL::PositionResolution, // const float phires,
417  G4TTL::PositionResolution, // const float lonres,
418  1, // const float eff,
419  0); // const float noise
420 
421  }
422 
423  return 0;
424 }
425 
426 
427 
428 
429 //-----------------------------------------------------------------------------------//
430 int make_barrel_layer_basic(string name, PHG4Reco *g4Reco,
431  double radius, double halflength, double tSilicon, double zOffset){
432 
433  //---------------------------------
434  //build barrel layer
435  //---------------------------------
436  const int nSubLayer = 7;
437 
438  string layerName[nSubLayer] = {"SiliconSensor", "Metalconnection", "HDI", "Cooling",
439  "Support1", "Support_Gap", "Support2"};
440  string material[nSubLayer] = {"G4_Si", "G4_Al", "G4_KAPTON", "G4_WATER",
441  "G4_GRAPHITE", "G4_AIR", "G4_GRAPHITE"};
442  double thickness[nSubLayer] = {tSilicon , 15 * um, 20 * um, 100 * um,
443  50 * um, 1, 50 * um};
444 
445  double max_bh_radius = 0.;
447 // cout << "started to create cylinder layer: " << name << endl;
448 
449  double currRadius = radius;
450 // cout << currRadius << endl;
451  for (int l = 0; l < nSubLayer; l++) {
452 // cout << name <<"_"<< layerName[l] << endl;
453  cyl = new PHG4CylinderSubsystem(name + "_" + layerName[l],l);
454  cyl->SuperDetector(name);
455  cyl->set_double_param("radius", currRadius);
456  cyl->set_double_param("length", 2.0 * halflength);
457  cyl->set_string_param("material", material[l]);
458  cyl->set_double_param("thickness", thickness[l]);
459  cyl->set_double_param("place_x", 0.);
460  cyl->set_double_param("place_y", 0.);
461  cyl->set_double_param("place_z", zOffset);
462 
463  if (l == 0) cyl->SetActive(); //only the Silicon Sensor is active
464  cyl->OverlapCheck(true);
465  g4Reco->registerSubsystem(cyl);
466  currRadius = currRadius+thickness[l];
467 // cout << currRadius << endl;
468  }
469 
470 
472  {
473  TRACKING::FastKalmanFilter->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
474  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
475  tSilicon / sqrt(12.), // const float radres,
476  G4TTL::PositionResolution, // const float phires,
477  G4TTL::PositionResolution, // const float lonres,
478  1, // const float eff,
479  0); // const float noise
481 
482  TRACKING::ProjectionNames.insert(name);
483  }
485  {
486  TRACKING::FastKalmanFilterInnerTrack->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
487  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
488  tSilicon / sqrt(12.), // const float radres,
489  G4TTL::PositionResolution, // const float phires,
490  G4TTL::PositionResolution, // const float lonres,
491  1, // const float eff,
492  0); // const float noise
493 
494  }
495 
496  return 0;
497 }
498 
499 //-----------------------------------------------------------------------------------//
501  double radius, double halflength, double tSilicon, double zOffset){
502 
503  //---------------------------------
504  //build barrel layer
505  //---------------------------------
506  const int nSubLayer = 4;
507 
508  string layerName[nSubLayer] = {"Cooling", "Crystal","SiliconSensor", "Motherboard"};
509  string material[nSubLayer] = {"G4_Al", "LSO", "G4_Si", "FR4"};
510  double thickness[nSubLayer] = {0.031 * cm, 0.029 * cm, tSilicon, 0.033 * cm};
511 
512  double max_bh_radius = 0.;
514 // cout << "started to create cylinder layer: " << name << endl;
515 
516  double currRadius = radius;
517 // cout << currRadius << endl;
518  for (int l = 0; l < nSubLayer; l++) {
519 // cout << name <<"_"<< layerName[l] << endl;
520  cyl = new PHG4CylinderSubsystem(name + "_" + layerName[l],l);
521  cyl->SuperDetector(name);
522  cyl->set_double_param("radius", currRadius);
523  cyl->set_double_param("length", 2.0 * halflength);
524  cyl->set_string_param("material", material[l]);
525  cyl->set_double_param("thickness", thickness[l]);
526  cyl->set_double_param("place_x", 0.);
527  cyl->set_double_param("place_y", 0.);
528  cyl->set_double_param("place_z", zOffset);
529 
530  if (l == 2) cyl->SetActive(); //only the Silicon Sensor is active
532  g4Reco->registerSubsystem(cyl);
533  currRadius = currRadius+thickness[l];
534 // cout << currRadius << endl;
535  }
536 
538  {
539  float resLGAD_barrel = G4TTL::PositionResolution;
541  resLGAD_barrel = 35e-1; // https://cds.cern.ch/record/2667167/files/CMS-TDR-020.pdf page 33 bottom
542  }
543  TRACKING::FastKalmanFilter->add_phg4hits(string("G4HIT_") + name, // const std::string& phg4hitsNames,
544  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
545  999, // const float radres,
546  resLGAD_barrel, // const float phires,
547  resLGAD_barrel, // const float lonres,
548  1, // const float eff, NOTE: Different from 0.95 in Modular
549  0); // const float noise
551 
552  TRACKING::ProjectionNames.insert(name);
553  }
554 
555  return 0;
556 }
557 
558 #endif
559 
560 //-----------------------------------------------------------------------------------//