ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4_EEMC_hybrid.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4_EEMC_hybrid.C
1 #ifndef MACRO_G4EEMCHYBRID_C
2 #define MACRO_G4EEMCHYBRID_C
3 
4 #include <GlobalVariables.C>
5 
6 #include <g4calo/RawTowerBuilderByHitIndex.h>
7 #include <g4calo/RawTowerDigitizer.h>
8 #include <g4eiccalos/PHG4HybridHomogeneousCalorimeterSubsystem.h>
9 #include <g4eiccalos/PHG4ForwardCalCellReco.h>
10 
11 #include <g4eval/CaloEvaluator.h>
12 
13 #include <g4main/PHG4Reco.h>
14 
15 #include <eiccaloreco/RawClusterBuilderkMA.h>
16 #include <eiccaloreco/RawClusterBuilderHelper.h>
17 
18 #include <caloreco/RawClusterBuilderFwd.h>
19 #include <caloreco/RawClusterBuilderTemplate.h>
20 #include <caloreco/RawTowerCalibration.h>
21 
22 #include <fun4all/Fun4AllServer.h>
23 
24 R__LOAD_LIBRARY(libcalo_reco.so)
25 R__LOAD_LIBRARY(libg4eiccalos.so)
26 R__LOAD_LIBRARY(libg4detectors.so)
27 R__LOAD_LIBRARY(libg4eval.so)
28 
29 namespace Enable
30 {
31  bool EEMCH = false;
32  bool EEMCH_ABSORBER = false;
33  bool EEMCH_CELL = false;
34  bool EEMCH_TOWER = false;
35  bool EEMCH_CLUSTER = false;
36  bool EEMCH_EVAL = false;
37  bool EEMCH_OVERLAPCHECK = false;
38  int EEMCH_VERBOSITY = 0;
39 } // namespace Enable
40 
41 
42 namespace G4EEMCH
43 {
44  int use_projective_geometry = 0;
45 
46  // double Gdz = 18. + 0.0001; // These 2 paras are only served as the dimension of the black hole
47  // double Gz0 = -170.;
48  double Gdz = 20. + 0.1;
49  double Gz0 = -180.;
50 
51  namespace SETTING
52  {
53  bool USEHYBRID = false;
54  bool USECEMCGeo = false;
55  bool USECUSTOMMAPNOCARBON = false;
56  bool USECUSTOMMAP30CM = false;
57  bool USECUSTOMMAPCARBON = false;
58  bool USECUSTOMMAPUPDATED = false;
59  } // namespace SETTING
60 
61  // Digitization (default photon digi):
62  RawTowerDigitizer::enu_digi_algorithm TowerDigi = RawTowerDigitizer::kSimple_photon_digitization;
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 
70  enum enu_Eemc_clusterizer
71  {
72  kEemcGraphClusterizer,
73  kEemcTemplateClusterizer
74  };
75  //default template clusterizer, as developed by Sasha Bazilevsky
76  enu_Eemc_clusterizer Eemc_clusterizer = kEemcTemplateClusterizer;
77  // graph clusterizer
78  //enu_Eemc_clusterizer Eemc_clusterizer = kEemcGraphClusterizer;
79 
80 } // namespace G4EEMC
81 
82 
83 void EEMCHInit()
84 {
85  BlackHoleGeometry::max_radius = std::max(BlackHoleGeometry::max_radius, 90.);
86  // from towerMap_EEMC_v006.txt
87  BlackHoleGeometry::min_z = std::min(BlackHoleGeometry::min_z, G4EEMCH::Gz0 - G4EEMCH::Gdz / 2.);
88  BlackHoleGeometry::max_z = std::max(BlackHoleGeometry::max_z, 30.);
89 
90 }
91 
92 
93 void EEMCHSetup(PHG4Reco *g4Reco)
94 {
95  bool AbsorberActive = Enable::ABSORBER || Enable::EEMCH_ABSORBER;
96  bool OverlapCheck = Enable::OVERLAPCHECK || Enable::EEMCH_OVERLAPCHECK;
97  int verbosity = std::max(Enable::VERBOSITY, Enable::EEMCH_VERBOSITY);
98 
100  ostringstream mapping_eemc_1, mapping_eemc_2;
101 
102  cout << "hybrid: " << G4EEMCH::SETTING::USEHYBRID << "\t CEMC:" << G4EEMCH::SETTING::USECEMCGeo << endl;
103 
104  PHG4HybridHomogeneousCalorimeterSubsystem *eemc_crystal = new PHG4HybridHomogeneousCalorimeterSubsystem("EEMC");
105  eemc_crystal->SuperDetector("EEMC");
106  eemc_crystal->SetActive();
107  if (AbsorberActive)
108  eemc_crystal->SetAbsorberActive();
109 
110  if(G4EEMCH::SETTING::USECUSTOMMAPNOCARBON){
111  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185_noCarbon.txt";
112  } else if (G4EEMCH::SETTING::USECUSTOMMAP30CM){
113  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185_30cm.txt";
114  } else if (G4EEMCH::SETTING::USECUSTOMMAPCARBON){
115  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185cm_EEEMCcarbon.txt";
116  } else if (G4EEMCH::SETTING::USECUSTOMMAPUPDATED){
117  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185cm_updatedGeo.txt";
118  } else if (G4EEMCH::SETTING::USEHYBRID && !G4EEMCH::SETTING::USECEMCGeo)
119  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_crystal_200cm_SciGlassBarrel.txt";
120  else if (G4EEMCH::SETTING::USEHYBRID && G4EEMCH::SETTING::USECEMCGeo)
121  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_crystal_200cm_CEMCBarrel.txt";
122  else if (!G4EEMCH::SETTING::USEHYBRID && !G4EEMCH::SETTING::USECEMCGeo)
123 // mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_200cm_SciGlassBarrel.txt";
124  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185cm.txt";
125  else if (!G4EEMCH::SETTING::USEHYBRID && G4EEMCH::SETTING::USECEMCGeo)
126  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_200cm_CEMCBarrel.txt";
127  else {
128  cout << "*******************************************************************************" << endl;
129  cout << "****** ATTENTION no EEMC set as your settings aren't correct ******" << endl;
130  cout << "*******************************************************************************" << endl;
131  return;
132  }
133  cout << "setting EEMC crystal mapping: " << mapping_eemc_1.str() << endl;
134  eemc_crystal->set_string_param("mappingtower", mapping_eemc_1.str());
135  eemc_crystal->OverlapCheck(OverlapCheck);
136 
137  g4Reco->registerSubsystem(eemc_crystal);
138 
139  if (G4EEMCH::SETTING::USEHYBRID){
140  if (!G4EEMCH::SETTING::USECEMCGeo)
141  mapping_eemc_2 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_glass_200cm_SciGlassBarrel.txt";
142  else if ( G4EEMCH::SETTING::USECEMCGeo)
143  mapping_eemc_2 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_glass_200cm_CEMCBarrel.txt";
144  else {
145  cout << "*******************************************************************************" << endl;
146  cout << "****** requested hybrid option but no glass mapping set ******" << endl;
147  cout << "*******************************************************************************" << endl;
148  return;
149  }
150 
151  cout << "setting EEMC glass mapping: " << mapping_eemc_2.str() << endl;
152  PHG4HybridHomogeneousCalorimeterSubsystem *eemc_glass = new PHG4HybridHomogeneousCalorimeterSubsystem("EEMC_glass");
153  eemc_glass->SuperDetector("EEMC_glass");
154  eemc_glass->SetActive();
155  if (AbsorberActive)
156  eemc_glass->SetAbsorberActive();
157 
158  eemc_glass->set_string_param("mappingtower", mapping_eemc_2.str());
159  eemc_glass->OverlapCheck(OverlapCheck);
160  g4Reco->registerSubsystem(eemc_glass);
161 
162  }
163 
164 }
165 
166 void EEMCH_Cells()
167 {}
168 
169 void EEMCH_Towers()
170 {
171  int verbosity = std::max(Enable::VERBOSITY, Enable::EEMCH_VERBOSITY);
172 
173  Fun4AllServer *se = Fun4AllServer::instance();
174 
175  ostringstream mapping_eemc_1, mapping_eemc_2;
176 
177  if(G4EEMCH::SETTING::USECUSTOMMAPNOCARBON){
178  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185_noCarbon.txt";
179  } else if (G4EEMCH::SETTING::USECUSTOMMAP30CM){
180  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185_30cm.txt";
181  } else if (G4EEMCH::SETTING::USECUSTOMMAPCARBON){
182  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185cm_EEEMCcarbon.txt";
183  } else if (G4EEMCH::SETTING::USECUSTOMMAPUPDATED){
184  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185cm_updatedGeo_TemporaryTowerFix.txt";
185  } else if (G4EEMCH::SETTING::USEHYBRID && !G4EEMCH::SETTING::USECEMCGeo)
186  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_crystal_200cm_SciGlassBarrel.txt";
187  else if (G4EEMCH::SETTING::USEHYBRID && G4EEMCH::SETTING::USECEMCGeo)
188  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_crystal_200cm_CEMCBarrel.txt";
189  else if (!G4EEMCH::SETTING::USEHYBRID && !G4EEMCH::SETTING::USECEMCGeo)
190 // mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_200cm_SciGlassBarrel.txt";
191  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_185cm.txt";
192  else if (!G4EEMCH::SETTING::USEHYBRID && G4EEMCH::SETTING::USECEMCGeo)
193  mapping_eemc_1 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_purecrystal_200cm_CEMCBarrel.txt";
194 
195  // CMS lead tungstate barrel ECAL at 18 degree centrigrade: 4.5 photoelectrons per MeV
196  // lead tungsten test in Orsay is 15~20 p.e. per MeV, sci-glass is 5 p.e. per MeV
197  const double EEMC_photoelectron_per_GeV_crystal = 15000;
198  const double EEMC_photoelectron_per_GeV_glass = 5000;
199 
200  //the original values are [8, 16], no noise case[0, 0], really high case[80, 160]
201  const double crystal_pedestal_ADC = 0, crystal_zero_suppression_ADC = 0;
202  const double glass_pedestal_ADC = 0, glass_zero_suppression_ADC = 0;
203 
204  RawTowerBuilderByHitIndex *tower_EEMC_crystal = new RawTowerBuilderByHitIndex("TowerBuilder_EEMC_crystal");
205  tower_EEMC_crystal->Detector("EEMC");
206  tower_EEMC_crystal->set_sim_tower_node_prefix("SIM");
207  tower_EEMC_crystal->GeometryTableFile(mapping_eemc_1.str());
208  se->registerSubsystem(tower_EEMC_crystal);
209 
210  // Calorimeter digitization
211  RawTowerDigitizer *TowerDigitizer_EEMC_crystal = new RawTowerDigitizer("EEMCRawTowerDigitizer_crystal");
212  TowerDigitizer_EEMC_crystal->Detector("EEMC");
213  TowerDigitizer_EEMC_crystal->Verbosity(verbosity);
214  TowerDigitizer_EEMC_crystal->set_raw_tower_node_prefix("RAW");
215  TowerDigitizer_EEMC_crystal->set_digi_algorithm(G4EEMCH::TowerDigi);
216  TowerDigitizer_EEMC_crystal->set_pedstal_central_ADC(0);
217  TowerDigitizer_EEMC_crystal->set_pedstal_width_ADC(crystal_pedestal_ADC); // eRD1 test beam setting
218  TowerDigitizer_EEMC_crystal->set_photonelec_ADC(1); //not simulating ADC discretization error
219  TowerDigitizer_EEMC_crystal->set_photonelec_yield_visible_GeV(EEMC_photoelectron_per_GeV_crystal);
220  TowerDigitizer_EEMC_crystal->set_zero_suppression_ADC(crystal_zero_suppression_ADC); // eRD1 test beam setting
221  se->registerSubsystem(TowerDigitizer_EEMC_crystal);
222 
223  // Calorimeter calibration
224  RawTowerCalibration *TowerCalibration_EEMC_crystal = new RawTowerCalibration("EEMCRawTowerCalibration_crystal");
225  TowerCalibration_EEMC_crystal->Detector("EEMC");
226  TowerCalibration_EEMC_crystal->Verbosity(verbosity);
227  TowerCalibration_EEMC_crystal->set_calib_algorithm(RawTowerCalibration::kSimple_linear_calibration);
228  if (G4EEMCH::TowerDigi == RawTowerDigitizer::kNo_digitization)
229  TowerCalibration_EEMC_crystal->set_calib_const_GeV_ADC(1.);
230  else
231  TowerCalibration_EEMC_crystal->set_calib_const_GeV_ADC(1. / EEMC_photoelectron_per_GeV_crystal);
232  TowerCalibration_EEMC_crystal->set_pedstal_ADC(0);
233  se->registerSubsystem(TowerCalibration_EEMC_crystal);
234 
235  if (G4EEMCH::SETTING::USEHYBRID){
236  if (!G4EEMCH::SETTING::USECEMCGeo)
237  mapping_eemc_2 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_glass_200cm_SciGlassBarrel.txt";
238  else if ( G4EEMCH::SETTING::USECEMCGeo)
239  mapping_eemc_2 << getenv("CALIBRATIONROOT") << "/CrystalCalorimeter/mapping/crystal_mapping/tower_map_glass_200cm_CEMCBarrel.txt";
240 
241  RawTowerBuilderByHitIndex *tower_EEMC_glass = new RawTowerBuilderByHitIndex("TowerBuilder_EEMC_glass");
242  tower_EEMC_glass->Detector("EEMC_glass");
243  tower_EEMC_glass->set_sim_tower_node_prefix("SIM");
244  tower_EEMC_glass->GeometryTableFile(mapping_eemc_2.str());
245  se->registerSubsystem(tower_EEMC_glass);
246 
247  RawTowerDigitizer *TowerDigitizer_EEMC_glass = new RawTowerDigitizer("EEMCRawTowerDigitizer_glass");
248  TowerDigitizer_EEMC_glass->Detector("EEMC_glass");
249  TowerDigitizer_EEMC_glass->Verbosity(verbosity);
250  TowerDigitizer_EEMC_glass->set_raw_tower_node_prefix("RAW");
251  TowerDigitizer_EEMC_glass->set_digi_algorithm(G4EEMCH::TowerDigi);
252  TowerDigitizer_EEMC_glass->set_pedstal_central_ADC(0);
253  TowerDigitizer_EEMC_glass->set_pedstal_width_ADC(glass_pedestal_ADC); // eRD1 test beam setting
254  TowerDigitizer_EEMC_glass->set_photonelec_ADC(1); //not simulating ADC discretization error
255  TowerDigitizer_EEMC_glass->set_photonelec_yield_visible_GeV(EEMC_photoelectron_per_GeV_glass);
256  TowerDigitizer_EEMC_glass->set_zero_suppression_ADC(glass_zero_suppression_ADC); // eRD1 test beam setting
257  se->registerSubsystem(TowerDigitizer_EEMC_glass);
258 
259  RawTowerCalibration *TowerCalibration_EEMC_glass = new RawTowerCalibration("EEMCRawTowerCalibration_glass");
260  TowerCalibration_EEMC_glass->Detector("EEMC_glass");
261  TowerCalibration_EEMC_glass->Verbosity(verbosity);
262  TowerCalibration_EEMC_glass->set_calib_algorithm(RawTowerCalibration::kSimple_linear_calibration);
263  if (G4EEMCH::TowerDigi == RawTowerDigitizer::kNo_digitization)
264  TowerCalibration_EEMC_glass->set_calib_const_GeV_ADC(1.);
265  else
266  TowerCalibration_EEMC_glass->set_calib_const_GeV_ADC(1. / EEMC_photoelectron_per_GeV_glass);
267  TowerCalibration_EEMC_glass->set_pedstal_ADC(0);
268  se->registerSubsystem(TowerCalibration_EEMC_glass);
269  }
270 }
271 
272 
273 void EEMCH_Clusters()
274 {
275  int verbosity = std::max(Enable::VERBOSITY, Enable::EEMCH_VERBOSITY);
276  Fun4AllServer *se = Fun4AllServer::instance();
277 
278  if (G4EEMCH::Eemc_clusterizer == G4EEMCH::kEemcTemplateClusterizer)
279  {
280  RawClusterBuilderHelper *ClusterBuilder = new RawClusterBuilderkMA("EEMCRawClusterBuilderkMA");
281  ClusterBuilder->Detector("EEMC");
282  ClusterBuilder->set_seed_e(0.1);
283  ClusterBuilder->set_agg_e(0.001);
284  se->registerSubsystem(ClusterBuilder);
285  /*
286  RawClusterBuilderTemplate *ClusterBuilder_crystal = new RawClusterBuilderTemplate("EEMCRawClusterBuilderTemplate_crystal");
287  ClusterBuilder_crystal->Detector("EEMC");
288  ClusterBuilder_crystal->Verbosity(2);
289  se->registerSubsystem(ClusterBuilder_crystal);
290 
291  if (G4EEMCH::SETTING::USEHYBRID){
292  RawClusterBuilderTemplate *ClusterBuilder_glass = new RawClusterBuilderTemplate("EEMCRawClusterBuilderTemplate_glass");
293  ClusterBuilder_glass->Detector("EEMC_glass");
294  ClusterBuilder_glass->Verbosity(verbosity);
295  se->registerSubsystem(ClusterBuilder_glass);
296  }
297  */
298  }
299  else if (G4EEMCH::Eemc_clusterizer == G4EEMCH::kEemcGraphClusterizer)
300  {
301 
302  RawClusterBuilderFwd *ClusterBuilder_crystal = new RawClusterBuilderFwd("EEMCRawClusterBuilderFwd_crystal");
303  ClusterBuilder_crystal->Detector("EEMC");
304  ClusterBuilder_crystal->Verbosity(verbosity);
305  ClusterBuilder_crystal->Verbosity(2);
306  se->registerSubsystem(ClusterBuilder_crystal);
307 
308  if (G4EEMCH::SETTING::USEHYBRID){
309  RawClusterBuilderFwd *ClusterBuilder_glass = new RawClusterBuilderFwd("EEMCRawClusterBuilderFwd_glass");
310  ClusterBuilder_glass->Detector("EEMC_glass");
311  ClusterBuilder_glass->Verbosity(verbosity);
312  se->registerSubsystem(ClusterBuilder_glass);
313  }
314  }
315  else
316  {
317  cout << "EEMC_Clusters - unknown clusterizer setting " << G4EEMCH::Eemc_clusterizer << endl;
318  gSystem->Exit(1);
319  }
320  return;
321 }
322 
323 
324 void EEMCH_Eval(const std::string &outputfile)
325 {
326  int verbosity = std::max(Enable::VERBOSITY, Enable::EEMCH_VERBOSITY);
327 
328  Fun4AllServer *se = Fun4AllServer::instance();
329 
330  string outputroot = outputfile;
331  string remove_this = ".root";
332  size_t pos = outputroot.find(remove_this);
333  if (pos != string::npos){
334  outputroot.erase(pos, remove_this.length());
335  }
336  string outputrootc = outputroot+"_crystal.root";
337  string outputrootg = outputroot+"_glass.root";
338 
339  CaloEvaluator *eval_crystal = new CaloEvaluator("EEMCEVALUATOR", "EEMC", outputrootc.c_str());
340  eval_crystal->Verbosity(verbosity);
341  eval_crystal->set_do_cluster_eval(true);
342  se->registerSubsystem(eval_crystal);
343 
344  if (G4EEMCH::SETTING::USEHYBRID){
345  CaloEvaluator *eval_glass = new CaloEvaluator("EEMCGLASSEVALUATOR", "EEMC_glass", outputrootg.c_str());
346  eval_glass->Verbosity(verbosity);
347  eval_glass->set_do_cluster_eval(true);
348  se->registerSubsystem(eval_glass);
349  }
350  return;
351 }
352 #endif