EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4_HcalIn_ref.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4_HcalIn_ref.C
1 //Inner HCal reconstruction macro
2 #ifndef MACRO_G4HCALINREF_C
3 #define MACRO_G4HCALINREF_C
4 
5 #include <GlobalVariables.C>
6 #include <QA.C>
7 
10 
14 
15 #include <g4eval/CaloEvaluator.h>
16 
17 #include <g4main/PHG4Reco.h>
18 
19 #include <caloreco/RawClusterBuilderGraph.h>
20 #include <caloreco/RawClusterBuilderTemplate.h>
21 #include <caloreco/RawTowerCalibration.h>
22 #include <qa_modules/QAG4SimulationCalorimeter.h>
23 
24 #include <fun4all/Fun4AllServer.h>
25 
26 R__LOAD_LIBRARY(libcalo_reco.so)
27 R__LOAD_LIBRARY(libg4calo.so)
28 R__LOAD_LIBRARY(libg4detectors.so)
29 R__LOAD_LIBRARY(libg4eval.so)
30 R__LOAD_LIBRARY(libqa_modules.so)
31 
32 void HCalInner_SupportRing(PHG4Reco *g4Reco);
33 
34 namespace Enable
35 {
36  bool HCALIN = false;
37  bool HCALIN_ABSORBER = false;
38  bool HCALIN_OVERLAPCHECK = false;
39  bool HCALIN_CELL = false;
40  bool HCALIN_TOWER = false;
41  bool HCALIN_CLUSTER = false;
42  bool HCALIN_EVAL = false;
43  bool HCALIN_QA = false;
44  bool HCALIN_SUPPORT = false;
45  int HCALIN_VERBOSITY = 0;
46 } // namespace Enable
47 
48 namespace G4HCALIN
49 {
50  double support_ring_outer_radius = 178.0 - 0.001;
51  double support_ring_z_ring2 = (2150 + 2175) / 2. / 10.;
52  double dz = 25. / 10.;
53 
54  //Inner HCal absorber material selector:
55  //false - old version, absorber material is SS310
56  //true - default Choose if you want Aluminum
57  bool inner_hcal_material_Al = true;
58 
59  int inner_hcal_eic = 0;
60 
61  // Digitization (default photon digi):
63  // directly pass the energy of sim tower to digitized tower
64  // kNo_digitization
65  // simple digitization with photon statistics, single amplitude ADC conversion and pedestal
66  // kSimple_photon_digitization
67  // digitization with photon statistics on SiPM with an effective pixel N, ADC conversion and pedestal
68  // kSiPM_photon_digitization
69 
71  {
73 
75  };
76 
80  //enu_HCalIn_clusterizer HCalIn_clusterizer = kHCalInGraphClusterizer;
81 } // namespace G4HCALIN
82 
83 // Init is called by G4Setup.C
84 void HCalInnerInit(const int iflag = 0)
85 {
89  if (iflag == 1)
90  {
92  }
93 }
94 
95 double HCalInner(PHG4Reco *g4Reco,
96  double radius,
97  const int crossings)
98 {
99  bool AbsorberActive = Enable::ABSORBER || Enable::HCALIN_ABSORBER;
100  bool OverlapCheck = Enable::OVERLAPCHECK || Enable::HCALIN_OVERLAPCHECK;
102 
103  // all sizes are in cm!
105  // these are the parameters you can change with their default settings
106  // hcal->set_string_param("material","SS310");
108  {
109  if (verbosity > 0)
110  {
111  cout << "HCalInner - construct inner HCal absorber with G4_Al" << endl;
112  }
113  hcal->set_string_param("material", "G4_Al");
114  }
115  else
116  {
117  if (verbosity > 0)
118  {
119  cout << "HCalInner - construct inner HCal absorber with SS310" << endl;
120  }
121  hcal->set_string_param("material", "SS310");
122  }
123  // hcal->set_double_param("inner_radius", 117.27);
124  //-----------------------------------------
125  // the light correction can be set in a single call
126  // hcal->set_double_param("light_balance_inner_corr", NAN);
127  // hcal->set_double_param("light_balance_inner_radius", NAN);
128  // hcal->set_double_param("light_balance_outer_corr", NAN);
129  // hcal->set_double_param("light_balance_outer_radius", NAN);
130  // hcal->SetLightCorrection(NAN,NAN,NAN,NAN);
131  //-----------------------------------------
132  // hcal->set_double_param("outer_radius", 134.42);
133  // hcal->set_double_param("place_x", 0.);
134  // hcal->set_double_param("place_y", 0.);
135  // hcal->set_double_param("place_z", 0.);
136  // hcal->set_double_param("rot_x", 0.);
137  // hcal->set_double_param("rot_y", 0.);
138  // hcal->set_double_param("rot_z", 0.);
139  // hcal->set_double_param("scinti_eta_coverage", 1.1);
140  // hcal->set_double_param("scinti_gap_neighbor", 0.1);
141  // hcal->set_double_param("scinti_inner_gap", 0.85);
142  // hcal->set_double_param("scinti_outer_gap", 1.22 * (5.0 / 4.0));
143  // hcal->set_double_param("scinti_outer_radius", 133.3);
144  // hcal->set_double_param("scinti_tile_thickness", 0.7);
145  // hcal->set_double_param("size_z", 175.94 * 2);
146  // hcal->set_double_param("steplimits", NAN);
147  // hcal->set_double_param("tilt_angle", 36.15);
148 
149  // hcal->set_int_param("light_scint_model", 1);
150  // hcal->set_int_param("ncross", 0);
151  // hcal->set_int_param("n_towers", 64);
152  // hcal->set_int_param("n_scinti_plates_per_tower", 4);
153  // hcal->set_int_param("n_scinti_tiles", 12);
154 
155  // hcal->set_string_param("material", "SS310");
156 
157  hcal->SetActive();
158  hcal->SuperDetector("HCALIN");
159  if (AbsorberActive)
160  {
161  hcal->SetAbsorberActive();
162  }
163  hcal->OverlapCheck(OverlapCheck);
164 
165  g4Reco->registerSubsystem(hcal);
166 
167  radius = hcal->get_double_param("outer_radius");
168 
169  HCalInner_SupportRing(g4Reco);
170 
171  radius += no_overlapp;
172  return radius;
173 }
174 
176 void HCalInner_SupportRing(PHG4Reco *g4Reco)
177 {
178  bool AbsorberActive = Enable::SUPPORT || Enable::HCALIN_SUPPORT;
179 
180  const double z_ring1 = (2025 + 2050) / 2. / 10.;
181  const double innerradius_sphenix = 116.;
182  const double innerradius_ephenix_hadronside = 138.;
183  const double z_rings[] =
185 
187 
188  for (int i = 0; i < 4; i++)
189  {
190  double innerradius = innerradius_sphenix;
191  if (z_rings[i] > 0 && G4HCALIN::inner_hcal_eic == 1)
192  {
193  innerradius = innerradius_ephenix_hadronside;
194  }
195  cyl = new PHG4CylinderSubsystem("HCALIN_SPT_N1", i);
196  cyl->set_double_param("place_z", z_rings[i]);
197  cyl->SuperDetector("HCALIN_SPT");
198  cyl->set_double_param("radius", innerradius);
199  cyl->set_int_param("lengthviarapidity", 0);
200  cyl->set_double_param("length", G4HCALIN::dz);
201  cyl->set_string_param("material", "SS310");
202  cyl->set_double_param("thickness", G4HCALIN::support_ring_outer_radius - innerradius);
203  if (AbsorberActive)
204  {
205  cyl->SetActive();
206  }
207  g4Reco->registerSubsystem(cyl);
208  }
209 
210  return;
211 }
212 
214 {
216 
218 
219  PHG4HcalCellReco *hc = new PHG4HcalCellReco("HCALIN_CELLRECO");
220  hc->Detector("HCALIN");
221  // hc->Verbosity(2);
222  // check for energy conservation - needs modified "infinite" timing cuts
223  // 0-999999999
224  // hc->checkenergy();
225  // timing cuts with their default settings
226  // hc->set_double_param("tmin",0.);
227  // hc->set_double_param("tmax",60.0);
228  // or all at once:
229  // hc->set_timing_window(0.0,60.0);
230  // this sets all cells to a fixed energy for debugging
231  // hc->set_fixed_energy(1.);
232  se->registerSubsystem(hc);
233 
234  return;
235 }
236 
238 {
241 
242  HcalRawTowerBuilder *TowerBuilder = new HcalRawTowerBuilder("HcalInRawTowerBuilder");
243  TowerBuilder->Detector("HCALIN");
244  TowerBuilder->set_sim_tower_node_prefix("SIM");
245  // this sets specific decalibration factors
246  // for a given cell
247  // TowerBuilder->set_cell_decal_factor(1,10,0.1);
248  // for a whole tower
249  // TowerBuilder->set_tower_decal_factor(0,10,0.2);
250  TowerBuilder->Verbosity(verbosity);
251  se->registerSubsystem(TowerBuilder);
252 
253  // From 2016 Test beam sim
254  RawTowerDigitizer *TowerDigitizer = new RawTowerDigitizer("HcalInRawTowerDigitizer");
255  TowerDigitizer->Detector("HCALIN");
256  // TowerDigitizer->set_raw_tower_node_prefix("RAW_LG");
257  TowerDigitizer->set_digi_algorithm(G4HCALIN::TowerDigi);
258  TowerDigitizer->set_pedstal_central_ADC(0);
259  TowerDigitizer->set_pedstal_width_ADC(1); // From Jin's guess. No EMCal High Gain data yet! TODO: update
260  TowerDigitizer->set_photonelec_ADC(32. / 5.);
261  TowerDigitizer->set_photonelec_yield_visible_GeV(32. / 5 / (0.4e-3));
262  TowerDigitizer->set_zero_suppression_ADC(-0); // no-zero suppression
263  se->registerSubsystem(TowerDigitizer);
264 
265  //Default sampling fraction for SS310
266  double visible_sample_fraction_HCALIN = 0.0631283; //, /gpfs/mnt/gpfs04/sphenix/user/jinhuang/prod_analysis/hadron_shower_res_nightly/./G4Hits_sPHENIX_pi-_eta0_16GeV-0000.root_qa.rootQA_Draw_HCALIN_G4Hit.pdf
267 
268  if (G4HCALIN::inner_hcal_material_Al) visible_sample_fraction_HCALIN = 0.162166; //for "G4_Al", Abhisek Sen <sen.abhisek@gmail.com>
269 
270  RawTowerCalibration *TowerCalibration = new RawTowerCalibration("HcalInRawTowerCalibration");
271  TowerCalibration->Detector("HCALIN");
272  // TowerCalibration->set_raw_tower_node_prefix("RAW_LG");
273  // TowerCalibration->set_calib_tower_node_prefix("CALIB_LG");
276  {
277  // 0.176 extracted from electron sims (edep(scintillator)/edep(total))
278  TowerCalibration->set_calib_const_GeV_ADC(1. / 0.176);
279  }
280  else
281  {
282  TowerCalibration->set_calib_const_GeV_ADC(0.4e-3 / visible_sample_fraction_HCALIN);
283  }
284  TowerCalibration->set_pedstal_ADC(0);
285  se->registerSubsystem(TowerCalibration);
286 
287  return;
288 }
289 
291 {
293 
295 
297  {
298  RawClusterBuilderTemplate *ClusterBuilder = new RawClusterBuilderTemplate("HcalInRawClusterBuilderTemplate");
299  ClusterBuilder->Detector("HCALIN");
300  ClusterBuilder->SetCylindricalGeometry(); // has to be called after Detector()
301  ClusterBuilder->Verbosity(verbosity);
302  se->registerSubsystem(ClusterBuilder);
303  }
305  {
306  RawClusterBuilderGraph *ClusterBuilder = new RawClusterBuilderGraph("HcalInRawClusterBuilderGraph");
307  ClusterBuilder->Detector("HCALIN");
308  ClusterBuilder->Verbosity(verbosity);
309  se->registerSubsystem(ClusterBuilder);
310  }
311  else
312  {
313  cout << "HCalIn_Clusters - unknown clusterizer setting!" << endl;
314  exit(1);
315  }
316  return;
317 }
318 
319 void HCALInner_Eval(const std::string &outputfile)
320 {
323 
324  CaloEvaluator *eval = new CaloEvaluator("HCALINEVALUATOR", "HCALIN", outputfile);
325  eval->Verbosity(verbosity);
326  se->registerSubsystem(eval);
327 
328  return;
329 }
330 
332 {
334 
337  qa->Verbosity(verbosity);
338  se->registerSubsystem(qa);
339 
340  return;
341 }
342 
343 #endif