ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4_mRwell_EIC.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4_mRwell_EIC.C
1 #ifndef MACRO_G4mRWELL_C
2 #define MACRO_G4mRWELL_C
3 
4 #include <GlobalVariables.C>
5 
6 #include <fun4all/Fun4AllServer.h>
7 #include <g4detectors/PHG4CylinderSubsystem.h>
8 #include <g4detectors/PHG4SectorSubsystem.h>
9 #include <g4main/PHG4Reco.h>
10 #include <g4trackfastsim/PHG4TrackFastSim.h>
11 
12 R__LOAD_LIBRARY(libfun4all.so)
13 R__LOAD_LIBRARY(libg4testbench.so)
14 R__LOAD_LIBRARY(libg4detectors.so)
15 
16 namespace Enable
17 {
18  bool RWELL = false;
19  bool RWELL_OVERLAPCHECK = false;
20 } // namespace Enable
21 
22 namespace RWELL
23 {
24  //All units specified in cm unless stated otherwise
25  // const int n_layer = 2; //tracker layers
26  // //const double nom_radius[RWELL::n_layer] = {79.5,90.0}; //77 to not hit DIRC
27  // // const double nom_radius[RWELL::n_layer] = {78.67, 90.0}; //77 to not hit DIRC
28  // const double nom_radius[RWELL::n_layer] = {69 - 1.6 - 2.6, 78.67}; //77 to not hit DIRC
29  // const double nom_driftgap[RWELL::n_layer] = {0.4, 0.4};
30  // const double nom_length[RWELL::n_layer] = {300.0, 300.0};
31 
32  const int n_layer = 2; //tracker layers
33  const double nom_radius[RWELL::n_layer] = {51., 77.0175};
34  double e_length_uRwell[RWELL::n_layer] = {106., 197};
35  double h_length_uRwell[RWELL::n_layer] = { e_length_uRwell[0], 145};
36  const double nom_driftgap[RWELL::n_layer] = {0.4, 0.4};
37  const double nom_length[RWELL::n_layer] = {2*e_length_uRwell[0], 342.0};
38 
39 /*
40  //ECCE Proposal Values
41  const int n_layer = 3; //tracker layers
42  const double nom_radius[RWELL::n_layer] = {33.14, 51., 77.0175}; //ECCE Proposal values
43  double e_length_uRwell[RWELL::n_layer] = {40.08, 106., 197};
44  double h_length_uRwell[RWELL::n_layer] = {e_length_uRwell[0], e_length_uRwell[1], 145};
45  const double nom_driftgap[RWELL::n_layer] = {0.4, 0.4, 0.4};
46  const double nom_length[RWELL::n_layer] = {2*e_length_uRwell[0], 2*e_length_uRwell[1], 342.0};
47 */
48  int subsysID = 0;
49 } //namespace RWELL
50 
51 void RWellInit(int verbosity = 0)
52 {
53  BlackHoleGeometry::max_radius = std::max(BlackHoleGeometry::max_radius, RWELL::nom_radius[RWELL::n_layer - 1] / 10. + 0.7);
54  BlackHoleGeometry::max_z = std::max(BlackHoleGeometry::max_z, RWELL::nom_length[RWELL::n_layer - 1] / 2.0);
55  BlackHoleGeometry::min_z = std::min(BlackHoleGeometry::min_z, -RWELL::nom_length[RWELL::n_layer - 1] / 2.0);
56 }
57 
58 double Build_G4_RWell_Bare(PHG4Reco* g4Reco,
59  double rwellrad = 80.0,
60  double driftgap = 1.83,
61  double length = 200.0,
62  int index = 0)
63 {
64  bool OverlapCheck = Enable::OVERLAPCHECK || Enable::RWELL_OVERLAPCHECK;
65 
66  gSystem->Load("libfun4all");
67  gSystem->Load("libg4detectors.so");
68  gSystem->Load("libg4testbench.so");
69  gSystem->Load("libg4trackfastsim.so");
70 
71  //double driftgap = 3.0;
72  //double rwellrad = radius;
73  //double length = length;
74  double rsum = 0.0;
75 
76  //MPGD parameters (units in cm)
77  double kapton_thickness = 0.0175; //cm
78  double cu_thickness = 0.002; //cm
79  double pcb_thickness = 0.010; //cm
80  double prepreg_thickness = 0.005; //cm
81 
82  //Suppoort parameters (units in cm)
83  //string supMat = "PEEK"; //support material
84  //string supMat = "G4_Galactic"; //support material
85  // making carbon fiber epoxy
86  // G4Material *cfrp_intt = new G4Material("CFRP_INTT", density = 1.69 * g / cm3, ncomponents = 3);
87  string supMat = "CFRP_INTT";
88  //inner tube
89  double support_01_thickness = 0.50;
90  double support_01_length = 7.2;
91  //inner ring
92  double support_02_thickness = 1.6;
93  double support_02_length = 1.2;
94  //outer ring
95  double support_03_thickness = 0.50;
96  double support_03_length = 1.2;
97 
98  PHG4CylinderSubsystem* rwell_cyl;
99 
100  double Shift = (RWELL::h_length_uRwell[index] - RWELL::e_length_uRwell[index])/2.0;
101  // here is our uRwell:
102  //Gass layer
103  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_%d", index), RWELL::subsysID);
104  rwell_cyl->set_double_param("radius", rwellrad);
105  rwell_cyl->set_string_param("material", "G4_Ar");
106  rwell_cyl->set_double_param("thickness", driftgap);
107  rwell_cyl->set_int_param("lengthviarapidity", 0);
108  rwell_cyl->set_double_param("place_z", Shift);
109  rwell_cyl->set_double_param("length", length);
110  rwell_cyl->SuperDetector("RWELL");
111  rwell_cyl->SetActive(1);
112  rwell_cyl->OverlapCheck(OverlapCheck);
113  g4Reco->registerSubsystem(rwell_cyl);
114 
115  //Kapton
116  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Kapton_%d", index), RWELL::subsysID);
117  rwell_cyl->set_double_param("radius", rwellrad - kapton_thickness);
118  rwell_cyl->set_string_param("material", "G4_KAPTON");
119  rwell_cyl->set_double_param("thickness", kapton_thickness);
120  rwell_cyl->set_int_param("lengthviarapidity", 0);
121  rwell_cyl->set_double_param("place_z", Shift);
122  rwell_cyl->set_double_param("length", length);
123  rwell_cyl->SuperDetector("RWELL");
124  rwell_cyl->SetActive(0);
125  rwell_cyl->OverlapCheck(OverlapCheck);
126  g4Reco->registerSubsystem(rwell_cyl);
127  //Cu
128  rsum = rwellrad + driftgap;
129  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Cu_%d", index), RWELL::subsysID);
130  rwell_cyl->set_double_param("radius", rsum);
131  rwell_cyl->set_string_param("material", "G4_Cu");
132  rwell_cyl->set_double_param("thickness", cu_thickness);
133  rwell_cyl->set_int_param("lengthviarapidity", 0);
134  rwell_cyl->set_double_param("place_z", Shift);
135  rwell_cyl->set_double_param("length", length);
136  rwell_cyl->SuperDetector("RWELL");
137  rwell_cyl->SetActive(0);
138  rwell_cyl->OverlapCheck(OverlapCheck);
139  g4Reco->registerSubsystem(rwell_cyl);
140 
141  //Prepreg
142  rsum += cu_thickness;
143  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_PrePreg_%d", index), RWELL::subsysID);
144  rwell_cyl->set_double_param("radius", rsum);
145  rwell_cyl->set_string_param("material", "NOMEX");
146  rwell_cyl->set_double_param("thickness", prepreg_thickness);
147  rwell_cyl->set_int_param("lengthviarapidity", 0);
148  rwell_cyl->set_double_param("place_z", Shift);
149  rwell_cyl->set_double_param("length", length);
150  rwell_cyl->SuperDetector("RWELL");
151  rwell_cyl->SetActive(0);
152  g4Reco->registerSubsystem(rwell_cyl);
153  rwell_cyl->OverlapCheck(OverlapCheck);
154 
155  //PCB
156  rsum += prepreg_thickness;
157  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_PCB_%d", index), RWELL::subsysID);
158  rwell_cyl->set_double_param("radius", rsum);
159  rwell_cyl->set_string_param("material", "FR4");
160  rwell_cyl->set_double_param("thickness", pcb_thickness);
161  rwell_cyl->set_int_param("lengthviarapidity", 0);
162  rwell_cyl->set_double_param("place_z", Shift);
163  rwell_cyl->set_double_param("length", length);
164  rwell_cyl->SuperDetector("RWELL");
165  rwell_cyl->SetActive(0);
166  rwell_cyl->OverlapCheck(OverlapCheck);
167  g4Reco->registerSubsystem(rwell_cyl);
168 
169  return rwellrad;
170 }
171 
172 double Build_G4_RWell_Sup01(PHG4Reco* g4Reco,
173  double rwellrad = 80.0,
174  double driftgap = 1.5,
175  double length = 200.0,
176  int index = 0)
177 {
178  gSystem->Load("libfun4all");
179  gSystem->Load("libg4detectors.so");
180 
181  bool OverlapCheck = Enable::OVERLAPCHECK || Enable::RWELL_OVERLAPCHECK;
182 
183  //double driftgap = gap;
184  //double rwellrad = radius;
185  //double length = length;
186  double rsum = 0.0;
187 
188  //MPGD parameters (units in cm)
189  double kapton_thickness = 0.0175; //cm
190  double cu_thickness = 0.002; //cm
191  double pcb_thickness = 0.010; //cm
192  double prepreg_thickness = 0.005; //cm
193 
194  //Suppoort parameters (units in cm)
195  //string supMat = "PEEK"; //support material
196  //string supMat = "G4_Galactic"; //support material
197  // making carbon fiber epoxy
198  // G4Material *cfrp_intt = new G4Material("CFRP_INTT", density = 1.69 * g / cm3, ncomponents = 3);
199  string supMat = "CFRP_INTT";
200  //inner tube
201  double support_01_thickness = 0.50;
202  double support_01_length = 7.2;
203  //inner ring
204  double support_02_thickness = 1.6;
205  double support_02_length = 1.2;
206  //outer ring
207  double support_03_thickness = 0.50;
208  double support_03_length = 1.2;
209 
210  PHG4CylinderSubsystem* rwell_cyl;
211 
212  double Shift = (RWELL::h_length_uRwell[index] - RWELL::e_length_uRwell[index])/2.0;
213  // here is our uRwell:
214  //Gass layer
215  rwell_cyl = new PHG4CylinderSubsystem("RWELL", index);
216  rwell_cyl->set_double_param("radius", rwellrad);
217  rwell_cyl->set_string_param("material", "G4_METHANE");
218  rwell_cyl->set_double_param("thickness", driftgap);
219  rwell_cyl->set_int_param("lengthviarapidity", 0);
220  rwell_cyl->set_double_param("place_z", Shift);
221  rwell_cyl->set_double_param("length", length);
222  // rwell_cyl->SuperDetector("RWELL"); // breakout RWELL into individual layers
223  rwell_cyl->SetActive(1);
224  rwell_cyl->OverlapCheck(OverlapCheck);
225  g4Reco->registerSubsystem(rwell_cyl);
226  ++RWELL::subsysID;
227 
228  //Kapton
229  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Kapton_%d", index), RWELL::subsysID);
230  rwell_cyl->set_double_param("radius", rwellrad - kapton_thickness);
231  rwell_cyl->set_string_param("material", "G4_KAPTON");
232  rwell_cyl->set_double_param("thickness", kapton_thickness);
233  rwell_cyl->set_int_param("lengthviarapidity", 0);
234  rwell_cyl->set_double_param("place_z", Shift);
235  rwell_cyl->set_double_param("length", length);
236  rwell_cyl->SuperDetector("RWELL");
237  rwell_cyl->SetActive(0);
238  rwell_cyl->OverlapCheck(OverlapCheck);
239  g4Reco->registerSubsystem(rwell_cyl);
240  ++RWELL::subsysID;
241  //Cu
242  rsum = rwellrad + driftgap;
243  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Cu_%d", index), RWELL::subsysID);
244  rwell_cyl->set_double_param("radius", rsum);
245  rwell_cyl->set_string_param("material", "G4_Cu");
246  rwell_cyl->set_double_param("thickness", cu_thickness);
247  rwell_cyl->set_int_param("lengthviarapidity", 0);
248  rwell_cyl->set_double_param("place_z", Shift);
249  rwell_cyl->set_double_param("length", length);
250  rwell_cyl->SuperDetector("RWELL");
251  rwell_cyl->SetActive(0);
252  g4Reco->registerSubsystem(rwell_cyl);
253  rwell_cyl->OverlapCheck(OverlapCheck);
254  ++RWELL::subsysID;
255 
256  //Prepreg
257  rsum += cu_thickness;
258  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_PrePreg_%d", index), RWELL::subsysID);
259  rwell_cyl->set_double_param("radius", rsum);
260  rwell_cyl->set_string_param("material", "NOMEX");
261  rwell_cyl->set_double_param("thickness", prepreg_thickness);
262  rwell_cyl->set_int_param("lengthviarapidity", 0);
263  rwell_cyl->set_double_param("length", length);
264  rwell_cyl->SuperDetector("RWELL");
265  rwell_cyl->SetActive(0);
266  rwell_cyl->OverlapCheck(OverlapCheck);
267  g4Reco->registerSubsystem(rwell_cyl);
268  ++RWELL::subsysID;
269 
270  //PCB
271  rsum += prepreg_thickness;
272  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_PCB_%d", index), RWELL::subsysID);
273  rwell_cyl->set_double_param("radius", rsum);
274  rwell_cyl->set_string_param("material", "FR4");
275  rwell_cyl->set_double_param("thickness", pcb_thickness);
276  rwell_cyl->set_int_param("lengthviarapidity", 0);
277  rwell_cyl->set_double_param("place_z", Shift);
278  rwell_cyl->set_double_param("length", length);
279  rwell_cyl->SuperDetector("RWELL");
280  rwell_cyl->SetActive(0);
281  rwell_cyl->OverlapCheck(OverlapCheck);
282  g4Reco->registerSubsystem(rwell_cyl);
283  ++RWELL::subsysID;
284 
285  //---Support structure--
286  //tube
287  rsum += pcb_thickness;
288  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Support01_0_%d", index), RWELL::subsysID); //RWELL_<support type>_<location>_<index>
289  //Support 01 = tube
290  //Support 02 = inner ring
291  //Support 03 = outer ring
292  //location 0 = z < 0
293  //location 1 = z > 0
294  rwell_cyl->set_double_param("radius", rsum);
295  rwell_cyl->set_string_param("material", supMat);
296  rwell_cyl->set_double_param("thickness", support_01_thickness);
297  rwell_cyl->set_int_param("lengthviarapidity", 0);
298  rwell_cyl->set_double_param("place_z", -length / 2 + Shift + 0.5);
299  rwell_cyl->set_double_param("length", support_01_length);
300  rwell_cyl->SuperDetector("RWELL");
301  rwell_cyl->SetActive(0);
302  rwell_cyl->OverlapCheck(OverlapCheck);
303  g4Reco->registerSubsystem(rwell_cyl);
304  ++RWELL::subsysID;
305  //tube 2
306  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Support01_1_%d", index), RWELL::subsysID);
307  rwell_cyl->set_double_param("radius", rsum);
308  rwell_cyl->set_string_param("material", supMat);
309  rwell_cyl->set_double_param("thickness", support_01_thickness);
310  rwell_cyl->set_int_param("lengthviarapidity", 0);
311  rwell_cyl->set_double_param("place_z", length / 2 + Shift);
312  rwell_cyl->set_double_param("length", support_01_length);
313  rwell_cyl->SuperDetector("RWELL");
314  rwell_cyl->SetActive(0);
315  rwell_cyl->OverlapCheck(OverlapCheck);
316  g4Reco->registerSubsystem(rwell_cyl);
317  ++RWELL::subsysID;
318  //inner ring
319  rsum += support_01_thickness;
320  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Support02_0_%d", index), RWELL::subsysID);
321  rwell_cyl->set_double_param("radius", rsum);
322  rwell_cyl->set_string_param("material", supMat);
323  rwell_cyl->set_double_param("thickness", support_02_thickness);
324  rwell_cyl->set_int_param("lengthviarapidity", 0);
325  rwell_cyl->set_double_param("place_z", -length / 2 + support_01_length / 2 + Shift);
326  rwell_cyl->set_double_param("length", support_02_length);
327  rwell_cyl->SuperDetector("RWELL");
328  rwell_cyl->SetActive(0);
329  rwell_cyl->OverlapCheck(OverlapCheck);
330  g4Reco->registerSubsystem(rwell_cyl);
331  ++RWELL::subsysID;
332  //inner ring 2
333  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Support02_1_%d", index), RWELL::subsysID);
334  rwell_cyl->set_double_param("radius", rsum);
335  rwell_cyl->set_string_param("material", supMat);
336  rwell_cyl->set_double_param("thickness", support_02_thickness);
337  rwell_cyl->set_int_param("lengthviarapidity", 0);
338  rwell_cyl->set_double_param("place_z", length / 2 + Shift - support_01_length / 2);
339  rwell_cyl->set_double_param("length", support_02_length);
340  rwell_cyl->SuperDetector("RWELL");
341  rwell_cyl->SetActive(0);
342  rwell_cyl->OverlapCheck(OverlapCheck);
343  g4Reco->registerSubsystem(rwell_cyl);
344  ++RWELL::subsysID;
345  //outer ring
346  rsum += support_02_thickness;
347  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Support03_0_%d", index), RWELL::subsysID);
348  rwell_cyl->set_double_param("radius", rsum);
349  rwell_cyl->set_string_param("material", supMat);
350  rwell_cyl->set_double_param("thickness", support_03_thickness);
351  rwell_cyl->set_int_param("lengthviarapidity", 0);
352  rwell_cyl->set_double_param("length", support_03_length);
353  rwell_cyl->set_double_param("place_z", -length / 2 + Shift + support_01_length / 2);
354  rwell_cyl->SuperDetector("RWELL");
355  rwell_cyl->SetActive(0);
356  rwell_cyl->OverlapCheck(OverlapCheck);
357  g4Reco->registerSubsystem(rwell_cyl);
358  ++RWELL::subsysID;
359  //outer ring 2
360  rwell_cyl = new PHG4CylinderSubsystem(Form("RWELL_Support03_1_%d", index), RWELL::subsysID);
361  rwell_cyl->set_double_param("radius", rsum);
362  rwell_cyl->set_string_param("material", supMat);
363  rwell_cyl->set_double_param("thickness", support_03_thickness);
364  rwell_cyl->set_int_param("lengthviarapidity", 0);
365  rwell_cyl->set_double_param("length", support_03_length);
366  rwell_cyl->set_double_param("place_z", length / 2 + Shift - support_01_length / 2);
367  rwell_cyl->SuperDetector("RWELL");
368  rwell_cyl->SetActive(0);
369  rwell_cyl->OverlapCheck(OverlapCheck);
370  g4Reco->registerSubsystem(rwell_cyl);
371  ++RWELL::subsysID;
372 
373  return rwellrad;
374 }
375 
377 // 0: bare RWell (no support structure)
378 // 1: Implimentation of FIT support rings
379 double RWellSetup(PHG4Reco* g4Reco,
380  int type = 1)
381 {
382  double radius = 0;
383 
384  for (int ilyr = 0; ilyr < RWELL::n_layer; ilyr++) //RWELL trackers are registered in Build_RWELL macro
385  {
386  if (type == 0)
387  {
388  radius = Build_G4_RWell_Bare(g4Reco, //returns RWELL radiaus
389  RWELL::nom_radius[ilyr], //radius
390  RWELL::nom_driftgap[ilyr], //driftgap,
391  RWELL::nom_length[ilyr], //length
392  ilyr); //index
393  }
394  if (type == 1)
395  {
396  radius = Build_G4_RWell_Sup01(g4Reco, //returns RWELL radiaus
397  RWELL::nom_radius[ilyr], //radius
398  RWELL::nom_driftgap[ilyr], //driftgap,
399  RWELL::nom_length[ilyr], //length
400  ilyr); //index
401  }
402 
403 // // sourav: For spatial resolution the mRwell I will use about 55 microns
404 // // (usually it is in between 40-60 microns depending on the angle of incidence of
405 // // primary tracks when mRwell are used in microTPC mode i.e drift gap of 3-4 mm) .
406  if (TRACKING::FastKalmanFilter)
407  {
408  TRACKING::FastKalmanFilter->add_phg4hits(string("G4HIT_") + string(Form("RWELL_%d", ilyr)), // const std::string& phg4hitsNames,
409  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
410  1. / sqrt(12.), // const float radres,
411  55e-4, // const float phires,
412  55e-4, // const float lonres,
413  1, // const float eff,
414  0); // const float noise
415  TRACKING::FastKalmanFilter->add_cylinder_state(Form("RWELL_%d", ilyr), RWELL::nom_radius[ilyr]);
416  TRACKING::ProjectionNames.insert(Form("RWELL_%d", ilyr));
417  }
418 
419  // FastKalmanFilterInnerTrack has everything within DIRC
420  if (TRACKING::FastKalmanFilterInnerTrack and RWELL::nom_radius[ilyr] < 70)
421  TRACKING::FastKalmanFilterInnerTrack->add_phg4hits(string("G4HIT_") + string(Form("RWELL_%d", ilyr)), // const std::string& phg4hitsNames,
422  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
423  1. / sqrt(12.), // const float radres,
424  55e-4, // const float phires,
425  55e-4, // const float lonres,
426  1, // const float eff,
427  0); // const float noise
428  // only for layers that is close to the silicon tracker system, use in FastKalmanFilterSiliconTrack
429  if (TRACKING::FastKalmanFilterSiliconTrack and RWELL::nom_radius[ilyr] < 50)
430  TRACKING::FastKalmanFilterSiliconTrack->add_phg4hits(string("G4HIT_") + string(Form("RWELL_%d", ilyr)), // const std::string& phg4hitsNames,
431  PHG4TrackFastSim::Cylinder, // const DETECTOR_TYPE phg4dettype,
432  1. / sqrt(12.), // const float radres,
433  55e-4, // const float phires,
434  55e-4, // const float lonres,
435  1, // const float eff,
436  0); // const float noise
437  }
438  return radius; //cm
439 }
440 
441 // Central detector cell reco is disabled as EIC setup use the fast tracking sim for now
442 void RWell_Cells(int verbosity = 0)
443 {
444  // runs the cellularization of the energy deposits (g4hits)
445  // into detector hits (g4cells)
446 
447  //---------------
448  // Load libraries
449  //---------------
450 
451  gSystem->Load("libfun4all.so");
452  gSystem->Load("libg4detectors.so");
453 
454  //---------------
455  // Fun4All server
456  //---------------
457 
458  Fun4AllServer* se = Fun4AllServer::instance();
459 
460  //-----------
461  // SVTX cells
462  //-----------
463 
464  return;
465 }
466 // Central detector reco is disabled as EIC setup use the fast tracking sim for now
467 void RWell_Reco(int verbosity = 0)
468 {
469  //---------------
470  // Load libraries
471  //---------------
472  gSystem->Load("libfun4all.so");
473 
474  //---------------
475  // Fun4All server
476  //---------------
477 
478  Fun4AllServer* se = Fun4AllServer::instance();
479 
480  return;
481 }
482 
483 #endif