Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:19:42

0001 #ifndef MACRO_G4CEMCSPACAL_C
0002 #define MACRO_G4CEMCSPACAL_C
0003 
0004 #include <GlobalVariables.C>
0005 #include <QA.C>
0006 
0007 #include <g4detectors/PHG4CylinderCellReco.h>
0008 #include <g4detectors/PHG4CylinderGeom_Spacalv1.h>
0009 #include <g4detectors/PHG4CylinderSubsystem.h>
0010 #include <g4detectors/PHG4FullProjSpacalCellReco.h>
0011 #include <g4detectors/PHG4SpacalSubsystem.h>
0012 
0013 #include <g4calo/RawTowerBuilder.h>
0014 #include <g4calo/RawTowerDigitizer.h>
0015 
0016 #include <g4eval/CaloEvaluator.h>
0017 
0018 #include <g4main/PHG4Reco.h>
0019 #include <g4main/PHG4Utils.h>
0020 
0021 #include <caloreco/RawClusterBuilderGraph.h>
0022 #include <caloreco/RawClusterBuilderTemplate.h>
0023 #include <caloreco/RawClusterPositionCorrection.h>
0024 #include <caloreco/RawTowerCalibration.h>
0025 #include <qa_modules/QAG4SimulationCalorimeter.h>
0026 
0027 #include <fun4all/Fun4AllServer.h>
0028 
0029 double
0030 CEmc_1DProjectiveSpacal(PHG4Reco *g4Reco, double radius, const int crossings);
0031 
0032 double
0033 CEmc_2DProjectiveSpacal(PHG4Reco *g4Reco, double radius, const int crossings);
0034 
0035 R__LOAD_LIBRARY(libcalo_reco.so)
0036 R__LOAD_LIBRARY(libg4calo.so)
0037 R__LOAD_LIBRARY(libg4detectors.so)
0038 R__LOAD_LIBRARY(libg4eval.so)
0039 R__LOAD_LIBRARY(libqa_modules.so)
0040 
0041 namespace Enable
0042 {
0043   bool CEMC = false;
0044   bool CEMC_ABSORBER = false;
0045   bool CEMC_OVERLAPCHECK = false;
0046   bool CEMC_CELL = false;
0047   bool CEMC_TOWER = false;
0048   bool CEMC_CLUSTER = false;
0049   bool CEMC_EVAL = false;
0050   bool CEMC_QA = false;
0051   int CEMC_VERBOSITY = 0;
0052 }  // namespace Enable
0053 
0054 namespace G4CEMC
0055 {
0056   int Min_cemc_layer = 1;
0057   int Max_cemc_layer = 1;
0058 
0059   // Digitization (default photon digi):
0060   RawTowerDigitizer::enu_digi_algorithm TowerDigi = RawTowerDigitizer::kSimple_photon_digitization;
0061   //  RawTowerDigitizer::enu_digi_algorithm TowerDigi = RawTowerDigitizer::kNo_digitization;
0062   // directly pass the energy of sim tower to digitized tower
0063   // kNo_digitization
0064   // simple digitization with photon statistics, single amplitude ADC conversion and pedestal
0065   // kSimple_photon_digitization
0066   // digitization with photon statistics on SiPM with an effective pixel N, ADC conversion and pedestal
0067   // kSiPM_photon_digitization
0068 
0069   // set a default value for SPACAL configuration
0070   //  // 1D azimuthal projective SPACAL (fast)
0071   //int Cemc_spacal_configuration = PHG4CylinderGeom_Spacalv1::k1DProjectiveSpacal;
0072   //   2D azimuthal projective SPACAL (slow)
0073   int Cemc_spacal_configuration = PHG4CylinderGeom_Spacalv1::k2DProjectiveSpacal;
0074 
0075   enum enu_Cemc_clusterizer
0076   {
0077     kCemcGraphClusterizer,
0078 
0079     kCemcTemplateClusterizer
0080   };
0081 
0082   //! template clusterizer, RawClusterBuilderTemplate, as developed by Sasha Bazilevsky
0083   enu_Cemc_clusterizer Cemc_clusterizer = kCemcTemplateClusterizer;
0084   //! graph clusterizer, RawClusterBuilderGraph
0085   //enu_Cemc_clusterizer Cemc_clusterizer = kCemcGraphClusterizer;
0086 
0087 }  // namespace G4CEMC
0088 
0089 // black hole parameters are set in CEmc function
0090 // needs a dummy argument to play with current G4Setup_sPHENIX.C
0091 void CEmcInit(const int i = 0)
0092 {
0093 }
0094 
0095 //! EMCal main setup macro
0096 double
0097 CEmc(PHG4Reco *g4Reco, double radius, const int crossings)
0098 {
0099   if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k1DProjectiveSpacal)
0100   {
0101     return CEmc_1DProjectiveSpacal(/*PHG4Reco**/ g4Reco, /*double*/ radius, /*const int */ crossings);
0102   }
0103   else if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k2DProjectiveSpacal)
0104   {
0105     return CEmc_2DProjectiveSpacal(/*PHG4Reco**/ g4Reco, /*double*/ radius, /*const int */ crossings);
0106   }
0107   else
0108   {
0109     std::cout
0110         << "G4_CEmc_Spacal.C::CEmc - Fatal Error - unrecognized SPACAL configuration #"
0111         << G4CEMC::Cemc_spacal_configuration << ". Force exiting..." << std::endl;
0112     exit(-1);
0113     return 0;
0114   }
0115 }
0116 
0117 //! EMCal setup macro - 1D azimuthal projective SPACAL
0118 double
0119 CEmc_1DProjectiveSpacal(PHG4Reco *g4Reco, double radius, const int crossings)
0120 {
0121   bool AbsorberActive = Enable::ABSORBER || Enable::CEMC_ABSORBER;
0122   bool OverlapCheck = Enable::OVERLAPCHECK || Enable::CEMC_OVERLAPCHECK;
0123 
0124   double emc_inner_radius = 95.;  // emc inner radius from engineering drawing
0125   double cemcthickness = 12.7;
0126   double emc_outer_radius = emc_inner_radius + cemcthickness;  // outer radius
0127 
0128   if (radius > emc_inner_radius)
0129   {
0130     cout << "inconsistency: pstof outer radius: " << radius
0131          << " larger than emc inner radius: " << emc_inner_radius
0132          << endl;
0133     gSystem->Exit(-1);
0134   }
0135 
0136   //  boundary check
0137   if (radius > emc_inner_radius - 1.5 - no_overlapp)
0138   {
0139     cout << "G4_CEmc_Spacal.C::CEmc() - expect radius < " << emc_inner_radius - 1.5 - no_overlapp << " to install SPACAL" << endl;
0140     exit(1);
0141   }
0142   radius = emc_inner_radius - 1.5 - no_overlapp;
0143 
0144   // 1.5cm thick teflon as an approximation for EMCAl light collection + electronics (10% X0 total estimated)
0145   PHG4CylinderSubsystem *cyl = new PHG4CylinderSubsystem("CEMC_ELECTRONICS", 0);
0146   cyl->SuperDetector("CEMC_ELECTRONICS");
0147   cyl->set_double_param("radius", radius);
0148   cyl->set_string_param("material", "G4_TEFLON");
0149   cyl->set_double_param("thickness", 1.5);
0150   if (AbsorberActive) cyl->SetActive();
0151   g4Reco->registerSubsystem(cyl);
0152 
0153   radius += 1.5;
0154   radius += no_overlapp;
0155 
0156   int ilayer = G4CEMC::Min_cemc_layer;
0157   PHG4SpacalSubsystem *cemc = new PHG4SpacalSubsystem("CEMC", ilayer);
0158   cemc->set_double_param("radius", emc_inner_radius);
0159   cemc->set_double_param("thickness", cemcthickness);
0160 
0161   cemc->SetActive();
0162   cemc->SuperDetector("CEMC");
0163   if (AbsorberActive) cemc->SetAbsorberActive();
0164   cemc->OverlapCheck(OverlapCheck);
0165 
0166   g4Reco->registerSubsystem(cemc);
0167 
0168   if (ilayer > G4CEMC::Max_cemc_layer)
0169   {
0170     cout << "layer discrepancy, current layer " << ilayer
0171          << " max cemc layer: " << G4CEMC::Max_cemc_layer << endl;
0172   }
0173 
0174   radius += cemcthickness;
0175   radius += no_overlapp;
0176 
0177   // 0.5cm thick Stainless Steel as an approximation for EMCAl support system
0178   cyl = new PHG4CylinderSubsystem("CEMC_SPT", 0);
0179   cyl->SuperDetector("CEMC_SPT");
0180   cyl->set_double_param("radius", radius);
0181   cyl->set_string_param("material", "SS310");  // SS310 Stainless Steel
0182   cyl->set_double_param("thickness", 0.5);
0183   if (AbsorberActive) cyl->SetActive();
0184   g4Reco->registerSubsystem(cyl);
0185   radius += 0.5;
0186   // this is the z extend and outer radius of the support structure and therefore the z extend
0187   // and radius of the surrounding black holes
0188   BlackHoleGeometry::max_z = std::max(BlackHoleGeometry::max_z, 149.47);
0189   BlackHoleGeometry::min_z = std::min(BlackHoleGeometry::min_z, -149.47);
0190   BlackHoleGeometry::max_radius = std::max(BlackHoleGeometry::max_radius, radius);
0191   radius += no_overlapp;
0192 
0193   return radius;
0194 }
0195 
0196 //! 2D full projective SPACAL
0197 double
0198 CEmc_2DProjectiveSpacal(PHG4Reco *g4Reco, double radius, const int crossings)
0199 {
0200   bool AbsorberActive = Enable::ABSORBER || Enable::CEMC_ABSORBER;
0201   bool OverlapCheck = Enable::OVERLAPCHECK || Enable::CEMC_OVERLAPCHECK;
0202 
0203   double emc_inner_radius = 92;  // emc inner radius from engineering drawing
0204   double cemcthickness = 24.00000 - no_overlapp;
0205 
0206   //max radius is 116 cm;
0207   double emc_outer_radius = emc_inner_radius + cemcthickness;  // outer radius
0208   assert(emc_outer_radius < 116);
0209 
0210   if (radius > emc_inner_radius)
0211   {
0212     cout << "inconsistency: preshower radius+thickness: " << radius
0213          << " larger than emc inner radius: " << emc_inner_radius << endl;
0214     gSystem->Exit(-1);
0215   }
0216 
0217   // the radii are only to determined the thickness of the cemc
0218   radius = emc_inner_radius;
0219 
0220   // 1.5cm thick teflon as an approximation for EMCAl light collection + electronics (10% X0 total estimated)
0221   PHG4CylinderSubsystem *cyl = new PHG4CylinderSubsystem("CEMC_ELECTRONICS", 0);
0222   cyl->set_double_param("radius", radius);
0223   cyl->set_string_param("material", "G4_TEFLON");
0224   cyl->set_double_param("thickness", 1.5 - no_overlapp);
0225   cyl->SuperDetector("CEMC_ELECTRONICS");
0226   cyl->OverlapCheck(OverlapCheck);
0227   if (AbsorberActive) cyl->SetActive();
0228   g4Reco->registerSubsystem(cyl);
0229 
0230   radius += 1.5;
0231   cemcthickness -= 1.5 + no_overlapp;
0232 
0233   // 0.5cm thick Stainless Steel as an approximation for EMCAl support system
0234   cyl = new PHG4CylinderSubsystem("CEMC_SPT", 0);
0235   cyl->SuperDetector("CEMC_SPT");
0236   cyl->set_double_param("radius", radius + cemcthickness - 0.5);
0237   cyl->set_string_param("material", "SS310");  // SS310 Stainless Steel
0238   cyl->set_double_param("thickness", 0.5 - no_overlapp);
0239   cyl->OverlapCheck(OverlapCheck);
0240   if (AbsorberActive) cyl->SetActive();
0241   g4Reco->registerSubsystem(cyl);
0242 
0243   // this is the z extend and outer radius of the support structure and therefore the z extend
0244   // and radius of the surrounding black holes
0245   double sptlen = PHG4Utils::GetLengthForRapidityCoverage(radius + cemcthickness);
0246   BlackHoleGeometry::max_z = std::max(BlackHoleGeometry::max_z, sptlen);
0247   BlackHoleGeometry::min_z = std::min(BlackHoleGeometry::min_z, -sptlen);
0248   BlackHoleGeometry::max_radius = std::max(BlackHoleGeometry::max_radius, radius + cemcthickness);
0249 
0250   cemcthickness -= 0.5 + no_overlapp;
0251 
0252   int ilayer = 0;
0253   PHG4SpacalSubsystem *cemc;
0254 
0255   cemc = new PHG4SpacalSubsystem("CEMC", ilayer);
0256 
0257   cemc->set_int_param("virualize_fiber", 0);
0258   cemc->set_int_param("azimuthal_seg_visible", 1);
0259   cemc->set_int_param("construction_verbose", 0);
0260   cemc->Verbosity(0);
0261 
0262   cemc->UseCalibFiles(PHG4DetectorSubsystem::xml);
0263   cemc->SetCalibrationFileDir(string(getenv("CALIBRATIONROOT")) + string("/CEMC/Geometry_2018ProjTilted/"));
0264   cemc->set_double_param("radius", radius);            // overwrite minimal radius
0265   cemc->set_double_param("thickness", cemcthickness);  // overwrite thickness
0266 
0267   cemc->SetActive();
0268   cemc->SuperDetector("CEMC");
0269   if (AbsorberActive) cemc->SetAbsorberActive();
0270   cemc->OverlapCheck(OverlapCheck);
0271 
0272   g4Reco->registerSubsystem(cemc);
0273 
0274   if (ilayer > G4CEMC::Max_cemc_layer)
0275   {
0276     cout << "layer discrepancy, current layer " << ilayer
0277          << " max cemc layer: " << G4CEMC::Max_cemc_layer << endl;
0278   }
0279 
0280   radius += cemcthickness;
0281   radius += no_overlapp;
0282 
0283   return radius;
0284 }
0285 
0286 void CEMC_Cells()
0287 {
0288   int verbosity = std::max(Enable::VERBOSITY, Enable::CEMC_VERBOSITY);
0289 
0290   Fun4AllServer *se = Fun4AllServer::instance();
0291 
0292   if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k1DProjectiveSpacal)
0293   {
0294     PHG4CylinderCellReco *cemc_cells = new PHG4CylinderCellReco("CEMCCYLCELLRECO");
0295     cemc_cells->Detector("CEMC");
0296     cemc_cells->Verbosity(verbosity);
0297     for (int i = G4CEMC::Min_cemc_layer; i <= G4CEMC::Max_cemc_layer; i++)
0298     {
0299       //          cemc_cells->etaphisize(i, 0.024, 0.024);
0300       const double radius = 95;
0301       cemc_cells->cellsize(i, 2 * M_PI / 256. * radius, 2 * M_PI / 256. * radius);
0302     }
0303     se->registerSubsystem(cemc_cells);
0304   }
0305   else if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k2DProjectiveSpacal)
0306   {
0307     PHG4FullProjSpacalCellReco *cemc_cells = new PHG4FullProjSpacalCellReco("CEMCCYLCELLRECO");
0308     cemc_cells->Detector("CEMC");
0309     cemc_cells->Verbosity(verbosity);
0310     cemc_cells->get_light_collection_model().load_data_file(
0311         string(getenv("CALIBRATIONROOT")) + string("/CEMC/LightCollection/Prototype3Module.xml"),
0312         "data_grid_light_guide_efficiency", "data_grid_fiber_trans");
0313     se->registerSubsystem(cemc_cells);
0314   }
0315   else
0316   {
0317     cout << "G4_CEmc_Spacal.C::CEmc - Fatal Error - unrecognized SPACAL configuration #"
0318          << G4CEMC::Cemc_spacal_configuration << ". Force exiting..." << endl;
0319     gSystem->Exit(-1);
0320     return;
0321   }
0322 
0323   return;
0324 }
0325 
0326 void CEMC_Towers()
0327 {
0328   int verbosity = std::max(Enable::VERBOSITY, Enable::CEMC_VERBOSITY);
0329 
0330   Fun4AllServer *se = Fun4AllServer::instance();
0331 
0332   RawTowerBuilder *TowerBuilder = new RawTowerBuilder("EmcRawTowerBuilder");
0333   TowerBuilder->Detector("CEMC");
0334   TowerBuilder->set_sim_tower_node_prefix("SIM");
0335   TowerBuilder->Verbosity(verbosity);
0336   se->registerSubsystem(TowerBuilder);
0337 
0338   double sampling_fraction = 1;
0339   if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k1DProjectiveSpacal)
0340   {
0341     sampling_fraction = 0.0234335;  //from production:/gpfs02/phenix/prod/sPHENIX/preCDR/pro.1-beta.3/single_particle/spacal1d/zerofield/G4Hits_sPHENIX_e-_eta0_8GeV.root
0342   }
0343   else if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k2DProjectiveSpacal)
0344   {
0345     //      sampling_fraction = 0.02244; //from production: /gpfs02/phenix/prod/sPHENIX/preCDR/pro.1-beta.3/single_particle/spacal2d/zerofield/G4Hits_sPHENIX_e-_eta0_8GeV.root
0346     //    sampling_fraction = 2.36081e-02;  //from production: /gpfs02/phenix/prod/sPHENIX/preCDR/pro.1-beta.5/single_particle/spacal2d/zerofield/G4Hits_sPHENIX_e-_eta0_8GeV.root
0347     //    sampling_fraction = 1.90951e-02; // 2017 Tilt porjective SPACAL, 8 GeV photon, eta = 0.3 - 0.4
0348     sampling_fraction = 2e-02;  // 2017 Tilt porjective SPACAL, tower-by-tower calibration
0349   }
0350   else
0351   {
0352     std::cout
0353         << "G4_CEmc_Spacal.C::CEMC_Towers - Fatal Error - unrecognized SPACAL configuration #"
0354         << G4CEMC::Cemc_spacal_configuration << ". Force exiting..." << std::endl;
0355     exit(-1);
0356     return;
0357   }
0358 
0359   const double photoelectron_per_GeV = 500;  //500 photon per total GeV deposition
0360 
0361   bool doSimple = true;
0362 
0363   RawTowerDigitizer *TowerDigitizer = new RawTowerDigitizer("EmcRawTowerDigitizer");
0364   TowerDigitizer->Detector("CEMC");
0365   TowerDigitizer->Verbosity(verbosity);
0366   TowerDigitizer->set_digi_algorithm(G4CEMC::TowerDigi); 
0367   TowerDigitizer->set_variable_pedestal(true);  //read ped central and width from calibrations file comment next 2 lines if true
0368                                                 //  TowerDigitizer->set_pedstal_central_ADC(0);
0369                                                 //  TowerDigitizer->set_pedstal_width_ADC(8);  // eRD1 test beam setting
0370   TowerDigitizer->set_photonelec_ADC(1);        //not simulating ADC discretization error
0371   TowerDigitizer->set_photonelec_yield_visible_GeV(photoelectron_per_GeV / sampling_fraction);
0372   TowerDigitizer->set_variable_zero_suppression(true);  //read zs values from calibrations file comment next line if true
0373                                                         //  TowerDigitizer->set_zero_suppression_ADC(16);  // eRD1 test beam setting
0374   TowerDigitizer->GetParameters().ReadFromFile("CEMC", "xml", 0, 0,
0375                                                string(getenv("CALIBRATIONROOT")) + string("/CEMC/TowerCalibCombinedParams_2020/"));  // calibration database
0376 
0377 
0378   // tower digitizer settings for doing decalibration -JEF May '22
0379 
0380   TowerDigitizer->set_UseConditionsDB(false);
0381 
0382   //---------------
0383   // if conditions db is enabled in the line above, how to handle filename
0384   // needs decided see below examples  (TBD by Chris P)
0385   //------------------
0386   // Some standard decal files specification where full decal is happening
0387   // uses the same db file accessor api as for TowerCalib and thus format 
0388   // (format can change along with accessor internals, but user
0389   // needs to know/give a readable format file).
0390   // Decal tower by tow. factors in file are apply as a multiplicative 
0391   // factor to raw energy/adc see below about adc level
0392   //---
0393   //  TowerDigitizer->set_DoTowerDecal(true,"emcal_corr_sec12bands.root",false);
0394   //  TowerDigitizer->set_DoTowerDecal(true,"emcal_corr1_29.root",false);
0395   TowerDigitizer->set_DoTowerDecal(true,"emcal_newPatternCinco.root",false);
0396 
0397   // third parameter (doInverse) specifies if you want 
0398   //to instead apply the reciprocal  of the TowbyTow factors 
0399   // in the db file as a multiplicative factor 
0400   // here is an example of that
0401   //TowerDigitizer->set_DoTowerDecal(true,"emcal_newPatternCinco.root",true);
0402 
0403   // in actuality we do not actually suffer from digitization
0404   // effects on the entire pulse amplitude but rather sample
0405   //  ~12-14 times (pulse/pedestal itself about ~8 times)
0406   // and both the pedestal and pulse are extracted as continuous 
0407   // (or near continous) quantities.  In the near future the 
0408   // simple ("old fashioned") digitization scheme needs to implement
0409   // a full sim of the pulse extraction procedure. 
0410   // therefore for now when running in decal mode, as a temporary
0411   // more realistic approximation of this procedure, we change the 
0412   // energy response to continuous adc/energy values rather than digitized.
0413   // this behavior currently only occurs if running in the decal mode
0414   // ....
0415   // If you want to apply this non-digitizing part of the code
0416   // WITHOUT mod'ing the energy (ie w/o decal), call the DoDecal function without 
0417   // specifiying a filename as in the following example
0418   // 
0419   //  TowerDigitizer->set_DoTowerDecal(true, "",false);
0420 
0421   se->registerSubsystem(TowerDigitizer);
0422   
0423 
0424   RawTowerCalibration *TowerCalibration = new RawTowerCalibration("EmcRawTowerCalibration");
0425   TowerCalibration->Detector("CEMC");
0426   TowerCalibration->Verbosity(verbosity);
0427   
0428 
0429 
0430 
0431   if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k1DProjectiveSpacal)
0432     {
0433     if (G4CEMC::TowerDigi == RawTowerDigitizer::kNo_digitization)
0434     {
0435       // just use sampling fraction set previously
0436       TowerCalibration->set_calib_const_GeV_ADC(1.0 / sampling_fraction);
0437     }
0438     else
0439     {
0440       TowerCalibration->set_calib_algorithm(RawTowerCalibration::kSimple_linear_calibration);
0441       TowerCalibration->set_calib_const_GeV_ADC(1. / photoelectron_per_GeV);
0442       TowerCalibration->set_pedstal_ADC(0);
0443     }
0444   }
0445   else if (G4CEMC::Cemc_spacal_configuration == PHG4CylinderGeom_Spacalv1::k2DProjectiveSpacal)
0446   {
0447     if (G4CEMC::TowerDigi == RawTowerDigitizer::kNo_digitization)
0448     {
0449       // just use sampling fraction set previously
0450       TowerCalibration->set_calib_const_GeV_ADC(1.0 / sampling_fraction);
0451     }
0452     else
0453     {
0454       
0455       if (!doSimple) 
0456     {
0457 
0458 
0459       //for tower by tower cal
0460       TowerCalibration->set_calib_algorithm(RawTowerCalibration::kTower_by_tower_calibration);
0461       TowerCalibration->GetCalibrationParameters().ReadFromFile("CEMC", "xml", 0, 0,
0462                                 string(getenv("CALIBRATIONROOT")) + string("/CEMC/TowerCalibCombinedParams_2020/"));  // calibration database
0463       TowerCalibration->set_variable_GeV_ADC(true);                                                                                                   //read GeV per ADC from calibrations file comment next line if true
0464       //    TowerCalibration->set_calib_const_GeV_ADC(1. / photoelectron_per_GeV / 0.9715);                                                             // overall energy scale based on 4-GeV photon simulations
0465       TowerCalibration->set_variable_pedestal(true);                                                                                                  //read pedestals from calibrations file comment next line if true
0466       //    TowerCalibration->set_pedstal_ADC(0);
0467       ///////////////////////////
0468 
0469     }
0470       else // dosimple
0471     {
0472 
0473       TowerCalibration->set_calib_algorithm(RawTowerCalibration::kDbfile_tbt_gain_corr);
0474       TowerCalibration->set_UseConditionsDB(false);
0475 
0476       // Some standard db calo calibration files specification 
0477       //  
0478       // uses the db file accessor api and thus format 
0479       // (format can change along with accessor internals, but user
0480       // needs to know/give a readable format file).
0481       // Decal tower by tow. factors in file are apply as a multiplicative 
0482 
0483       //TowerCalibration->set_CalibrationFileName("emcal_corr1_29.root");
0484       //TowerCalibration->set_CalibrationFileName("inv_emcal_corr_sec12bands.root");
0485       TowerCalibration->set_CalibrationFileName("emcal_corr1_00.root");
0486       //TowerCalibration->set_CalibrationFileName("emcal_newPatternCinco.root");
0487 
0488       // since for this loop we avert the tower by tower improvements 
0489       // in the kTower_by_tower xml-file based cemc calibration 
0490       // which can be reimplemented in the new tower by tower dbfile format
0491       // but for now we make a single overall reduction factor of 0.87
0492       // that matches the same average calibration correction
0493       // changing the following line to the one after it.  -JEF
0494       //          TowerCalibration->set_calib_const_GeV_ADC(1. / photoelectron_per_GeV / 0.9715);                                                             // overall energy scale based on 4-GeV photon simulations
0495           TowerCalibration->set_calib_const_GeV_ADC(0.87 * 1./ photoelectron_per_GeV / 0.9715);                                                             // overall energy scale based on 4-GeV photon simulations
0496       TowerCalibration->set_pedstal_ADC(0);
0497       // note that in the TowerCalibration object, the pedestal subtraction is no
0498       // longer applied, the above line simply follows suit with all the other 
0499       // "calib_algorithms" on that object
0500 
0501       ///////////////////////////
0502     }
0503       
0504     }
0505   }
0506   else
0507   {
0508     cout << "G4_CEmc_Spacal.C::CEMC_Towers - Fatal Error - unrecognized SPACAL configuration #"
0509          << G4CEMC::Cemc_spacal_configuration << ". Force exiting..." << endl;
0510     gSystem->Exit(-1);
0511     return;
0512   }
0513   se->registerSubsystem(TowerCalibration);
0514 
0515   return;
0516 }
0517 
0518 void CEMC_Clusters()
0519 {
0520   int verbosity = std::max(Enable::VERBOSITY, Enable::CEMC_VERBOSITY);
0521 
0522   Fun4AllServer *se = Fun4AllServer::instance();
0523 
0524   if (G4CEMC::Cemc_clusterizer == G4CEMC::kCemcTemplateClusterizer)
0525   {
0526     RawClusterBuilderTemplate *ClusterBuilder = new RawClusterBuilderTemplate("EmcRawClusterBuilderTemplate");
0527     ClusterBuilder->Detector("CEMC");
0528     ClusterBuilder->Verbosity(verbosity);
0529     ClusterBuilder->set_threshold_energy(0.030);  // This threshold should be the same as in CEMCprof_Thresh**.root file below
0530     std::string emc_prof = getenv("CALIBRATIONROOT");
0531     emc_prof += "/EmcProfile/CEMCprof_Thresh30MeV.root";
0532     ClusterBuilder->LoadProfile(emc_prof);
0533     se->registerSubsystem(ClusterBuilder);
0534   }
0535   else if (G4CEMC::Cemc_clusterizer == G4CEMC::kCemcGraphClusterizer)
0536   {
0537     RawClusterBuilderGraph *ClusterBuilder = new RawClusterBuilderGraph("EmcRawClusterBuilderGraph");
0538     ClusterBuilder->Detector("CEMC");
0539     ClusterBuilder->Verbosity(verbosity);
0540     se->registerSubsystem(ClusterBuilder);
0541   }
0542   else
0543   {
0544     cout << "CEMC_Clusters - unknown clusterizer setting!" << endl;
0545     exit(1);
0546   }
0547 
0548   RawClusterPositionCorrection *clusterCorrection = new RawClusterPositionCorrection("CEMC");
0549 
0550   clusterCorrection->Get_eclus_CalibrationParameters().ReadFromFile("CEMC_RECALIB", "xml", 0, 0,
0551                                                                     //raw location
0552                                                                     string(getenv("CALIBRATIONROOT")) + string("/CEMC/PositionRecalibration_EMCal_9deg_tilt/"));
0553 
0554   clusterCorrection->Get_ecore_CalibrationParameters().ReadFromFile("CEMC_ECORE_RECALIB", "xml", 0, 0,
0555                                                                     //raw location
0556                                                                     string(getenv("CALIBRATIONROOT")) + string("/CEMC/PositionRecalibration_EMCal_9deg_tilt/"));
0557 
0558   clusterCorrection->Verbosity(verbosity);
0559   se->registerSubsystem(clusterCorrection);
0560 
0561   return;
0562 }
0563 void CEMC_Eval(const std::string &outputfile)
0564 {
0565   int verbosity = std::max(Enable::VERBOSITY, Enable::CEMC_VERBOSITY);
0566 
0567   Fun4AllServer *se = Fun4AllServer::instance();
0568 
0569   CaloEvaluator *eval = new CaloEvaluator("CEMCEVALUATOR", "CEMC", outputfile);
0570   eval->Verbosity(verbosity);
0571   se->registerSubsystem(eval);
0572 
0573   return;
0574 }
0575 
0576 void CEMC_QA()
0577 {
0578   int verbosity = std::max(Enable::QA_VERBOSITY, Enable::CEMC_VERBOSITY);
0579 
0580   Fun4AllServer *se = Fun4AllServer::instance();
0581   QAG4SimulationCalorimeter *qa = new QAG4SimulationCalorimeter("CEMC");
0582   qa->Verbosity(verbosity);
0583   se->registerSubsystem(qa);
0584 
0585   return;
0586 }
0587 
0588 #endif