macros/leptoquarks_fun4all/G4_HcalIn_ref.C

0001 //Inner HCal construction macro
0002
0003 //Inner HCal absorber material selector:
0004 //false - Default, absorber material is SS310
0005 //true - Choose if you want Aluminum
0006 const bool inner_hcal_material_Al = false;
0007
0008
0009 enum enu_HCalIn_clusterizer
0010 {
0011   kHCalInGraphClusterizer,
0012
0013   kHCalInTemplateClusterizer
0014 };
0015
0016 //! template clusterizer, RawClusterBuilderTemplate, as developed by Sasha Bazilevsky
0017 enu_HCalIn_clusterizer HCalIn_clusterizer = kHCalInTemplateClusterizer;
0018 //! graph clusterizer, RawClusterBuilderGraph
0019 //enu_HCalIn_clusterizer HCalIn_clusterizer = kHCalInGraphClusterizer;
0020
0021
0022 // Init is called by G4Setup.C
0023 void HCalInnerInit() {}
0024
0025 double HCalInner(PHG4Reco* g4Reco,
0026          double radius,
0027          const int crossings,
0028          const int absorberactive = 0,
0029          int verbosity = 0) {
0030   // all sizes are in cm!
0031
0032   gSystem->Load("libg4detectors.so");
0033   gSystem->Load("libg4testbench.so");
0034
0035
0036   PHG4InnerHcalSubsystem *hcal = new PHG4InnerHcalSubsystem("HCALIN");
0037   // these are the parameters you can change with their default settings
0038   // hcal->set_string_param("material","SS310");
0039   if(inner_hcal_material_Al)
0040     {
0041       cout <<"HCalInner - construct inner HCal absorber with G4_Al"<<endl;
0042       hcal->set_string_param("material","G4_Al");
0043     }
0044   // hcal->set_double_param("inner_radius", 117.27);
0045   //-----------------------------------------
0046   // the light correction can be set in a single call
0047   // hcal->set_double_param("light_balance_inner_corr", NAN);
0048   // hcal->set_double_param("light_balance_inner_radius", NAN);
0049   // hcal->set_double_param("light_balance_outer_corr", NAN);
0050   // hcal->set_double_param("light_balance_outer_radius", NAN);
0051   // hcal->SetLightCorrection(NAN,NAN,NAN,NAN);
0052   //-----------------------------------------
0053   // hcal->set_double_param("outer_radius", 134.42);
0054   // hcal->set_double_param("place_x", 0.);
0055   // hcal->set_double_param("place_y", 0.);
0056   // hcal->set_double_param("place_z", 0.);
0057   // hcal->set_double_param("rot_x", 0.);
0058   // hcal->set_double_param("rot_y", 0.);
0059   // hcal->set_double_param("rot_z", 0.);
0060   // hcal->set_double_param("scinti_eta_coverage", 1.1);
0061   // hcal->set_double_param("scinti_gap_neighbor", 0.1);
0062   // hcal->set_double_param("scinti_inner_gap", 0.85);
0063   // hcal->set_double_param("scinti_outer_gap", 1.22 * (5.0 / 4.0));
0064   // hcal->set_double_param("scinti_outer_radius", 133.3);
0065   // hcal->set_double_param("scinti_tile_thickness", 0.7);
0066   // hcal->set_double_param("size_z", 175.94 * 2);
0067   // hcal->set_double_param("steplimits", NAN);
0068   // hcal->set_double_param("tilt_angle", 36.15);
0069
0070   // hcal->set_int_param("light_scint_model", 1);
0071   // hcal->set_int_param("ncross", 0);
0072   // hcal->set_int_param("n_towers", 64);
0073   // hcal->set_int_param("n_scinti_plates_per_tower", 4);
0074   // hcal->set_int_param("n_scinti_tiles", 12);
0075
0076   // hcal->set_string_param("material", "SS310");
0077
0078   hcal->SetActive();
0079   hcal->SuperDetector("HCALIN");
0080   if (absorberactive)
0081     {
0082       hcal->SetAbsorberActive();
0083     }
0084   hcal->OverlapCheck(overlapcheck);
0085
0086   g4Reco->registerSubsystem( hcal );
0087
0088   radius = hcal->get_double_param("outer_radius");
0089
0090   HCalInner_SupportRing(g4Reco,absorberactive);
0091
0092   radius += no_overlapp;
0093   return radius;
0094 }
0095
0096 //! A rough version of the inner HCal support ring, from Richie's CAD drawing. - Jin
0097 void HCalInner_SupportRing(PHG4Reco* g4Reco,
0098                const int absorberactive = 0) {
0099
0100   gSystem->Load("libg4detectors.so");
0101   gSystem->Load("libg4testbench.so");
0102
0103   const double z_ring1 = (2025 + 2050) / 2. / 10.;
0104   const double z_ring2 = (2150 + 2175) / 2. / 10.;
0105   const double dz = 25. / 10.;
0106   const double innerradius = 116.;
0107   const double maxradius = 178.0 - 0.001; // avoid touching the outer HCal envelop volumne
0108
0109   const double z_rings[] =
0110     { -z_ring2, -z_ring1, z_ring1, z_ring2, 0, 0, 0, 0 };
0111
0112   PHG4CylinderSubsystem *cyl;
0113
0114   for (int i = 0; i < 4; i++)
0115     {
0116       cyl = new PHG4CylinderSubsystem("HCALIN_SPT_N1", i);
0117       cyl->set_double_param("place_z",z_rings[i]);
0118       cyl->SuperDetector("HCALIN_SPT");
0119       cyl->set_double_param("radius",innerradius);
0120       cyl->set_int_param("lengthviarapidity",0);
0121       cyl->set_double_param("length",dz);
0122       cyl->set_string_param("material","SS310");
0123       cyl->set_double_param("thickness",maxradius - innerradius);
0124       if (absorberactive)
0125     {
0126       cyl->SetActive();
0127     }
0128       g4Reco->registerSubsystem(cyl);
0129     }
0130
0131   return;
0132 }
0133
0134
0135 void HCALInner_Cells(int verbosity = 0) {
0136
0137   gSystem->Load("libfun4all.so");
0138   gSystem->Load("libg4detectors.so");
0139   Fun4AllServer *se = Fun4AllServer::instance();
0140
0141   PHG4HcalCellReco *hc = new PHG4HcalCellReco("HCALIN_CELLRECO");
0142   hc->Detector("HCALIN");
0143   //  hc->Verbosity(2);
0144   // check for energy conservation - needs modified "infinite" timing cuts
0145   // 0-999999999
0146   //  hc->checkenergy();
0147   // timing cuts with their default settings
0148   // hc->set_double_param("tmin",0.);
0149   // hc->set_double_param("tmax",60.0);
0150   // or all at once:
0151   // hc->set_timing_window(0.0,60.0);
0152   se->registerSubsystem(hc);
0153
0154   return;
0155 }
0156
0157 void HCALInner_Towers(int verbosity = 0) {
0158
0159   gSystem->Load("libg4calo.so");
0160   gSystem->Load("libcalo_reco.so");
0161   Fun4AllServer *se = Fun4AllServer::instance();
0162
0163   HcalRawTowerBuilder *TowerBuilder = new HcalRawTowerBuilder("HcalInRawTowerBuilder");
0164   TowerBuilder->Detector("HCALIN");
0165   TowerBuilder->set_sim_tower_node_prefix("SIM");
0166   TowerBuilder->Verbosity(verbosity);
0167   se->registerSubsystem( TowerBuilder );
0168
0169   // From 2016 Test beam sim
0170   RawTowerDigitizer *TowerDigitizer = new RawTowerDigitizer("HcalInRawTowerDigitizer");
0171   TowerDigitizer->Detector("HCALIN");
0172 //  TowerDigitizer->set_raw_tower_node_prefix("RAW_LG");
0173   TowerDigitizer->set_digi_algorithm(
0174        RawTowerDigitizer::kSimple_photon_digitalization);
0175   TowerDigitizer->set_pedstal_central_ADC(0);
0176   TowerDigitizer->set_pedstal_width_ADC(1); // From Jin's guess. No EMCal High Gain data yet! TODO: update
0177   TowerDigitizer->set_photonelec_ADC(32. / 5.);
0178   TowerDigitizer->set_photonelec_yield_visible_GeV(32. / 5 / (0.4e-3));
0179   TowerDigitizer->set_zero_suppression_ADC(-0); // no-zero suppression
0180   se->registerSubsystem(TowerDigitizer);
0181
0182   //Default sampling fraction for SS310
0183   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
0184
0185   if(inner_hcal_material_Al) visible_sample_fraction_HCALIN = 0.162166; //for "G4_Al", Abhisek Sen <sen.abhisek@gmail.com>
0186
0187   RawTowerCalibration *TowerCalibration = new RawTowerCalibration("HcalInRawTowerCalibration");
0188   TowerCalibration->Detector("HCALIN");
0189 //  TowerCalibration->set_raw_tower_node_prefix("RAW_LG");
0190 //  TowerCalibration->set_calib_tower_node_prefix("CALIB_LG");
0191   TowerCalibration->set_calib_algorithm(RawTowerCalibration::kSimple_linear_calibration);
0192   TowerCalibration->set_calib_const_GeV_ADC(0.4e-3 / visible_sample_fraction_HCALIN);
0193   TowerCalibration->set_pedstal_ADC(0);
0194   TowerCalibration->set_zero_suppression_GeV(-1); // no-zero suppression
0195   se->registerSubsystem(TowerCalibration);
0196
0197   return;
0198 }
0199
0200 void HCALInner_Clusters(int verbosity = 0) {
0201   gSystem->Load("libcalo_reco.so");
0202
0203   Fun4AllServer *se = Fun4AllServer::instance();
0204
0205   if (HCalIn_clusterizer == kHCalInTemplateClusterizer)
0206   {
0207     RawClusterBuilderTemplate* ClusterBuilder = new RawClusterBuilderTemplate("HcalInRawClusterBuilderTemplate");
0208     ClusterBuilder->Detector("HCALIN");
0209     ClusterBuilder->Verbosity(verbosity);
0210     se->registerSubsystem( ClusterBuilder );
0211
0212   }
0213   else if (HCalIn_clusterizer == kHCalInGraphClusterizer)
0214   {
0215     RawClusterBuilderGraph* ClusterBuilder = new RawClusterBuilderGraph("HcalInRawClusterBuilderGraph");
0216     ClusterBuilder->Detector("HCALIN");
0217     ClusterBuilder->Verbosity(verbosity);
0218     se->registerSubsystem( ClusterBuilder );
0219
0220   }
0221   else
0222   {
0223     cout <<"HCalIn_Clusters - unknown clusterizer setting!"<<endl;
0224     exit(1);
0225   }
0226   return;
0227 }
0228
0229 void HCALInner_Eval(std::string outputfile, int verbosity = 0) {
0230
0231   gSystem->Load("libfun4all.so");
0232   gSystem->Load("libg4eval.so");
0233   Fun4AllServer *se = Fun4AllServer::instance();
0234
0235   CaloEvaluator* eval = new CaloEvaluator("HCALINEVALUATOR", "HCALIN", outputfile.c_str());
0236   eval->Verbosity(verbosity);
0237   se->registerSubsystem( eval );
0238
0239   return;
0240 }