Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2026-04-02 08:15:26

0001 #ifndef MACRO_G4SETUPSPHENIX_C
0002 #define MACRO_G4SETUPSPHENIX_C
0003 
0004 #include <GlobalVariables.C>
0005 
0006 #include <G4_BeamLine.C>
0007 #include <G4_BlackHole.C>
0008 #include <G4_CEmc_Albedo.C>
0009 #include <G4_CEmc_Spacal.C>
0010 #include <G4_EPD.C>
0011 #include <G4_HcalIn_ref.C>
0012 #include <G4_HcalOut_ref.C>
0013 #include <G4_Magnet.C>
0014 #include <G4_Mbd.C>
0015 #include <G4_PSTOF.C>
0016 #include <G4_Pipe.C>
0017 #include <G4_PlugDoor.C>
0018 #include <G4_TrkrSimulation.C>
0019 #include <G4_User.C>
0020 #include <G4_World.C>
0021 #include <G4_ZDC.C>
0022 
0023 #include <g4detectors/PHG4CylinderSubsystem.h>
0024 
0025 #include <g4eval/PHG4DstCompressReco.h>
0026 
0027 #include <g4main/PHG4Reco.h>
0028 #include <g4main/PHG4TruthSubsystem.h>
0029 
0030 #include <phfield/PHFieldConfig.h>
0031 
0032 #include <g4decayer/EDecayType.hh>
0033 
0034 #include <fun4all/Fun4AllDstOutputManager.h>
0035 #include <fun4all/Fun4AllServer.h>
0036 
0037 #include <sstream>
0038 
0039 R__LOAD_LIBRARY(libg4decayer.so)
0040 R__LOAD_LIBRARY(libg4detectors.so)
0041 
0042 void G4Init()
0043 {
0044   // Check on invalid combinations
0045   if (Enable::CEMC && Enable::CEMCALBEDO)
0046   {
0047     std::cout << "Enable::CEMCALBEDO and Enable::CEMC cannot be set simultanously" << std::endl;
0048     gSystem->Exit(1);
0049   }
0050   // load detector/material macros and execute Init() function
0051 
0052   if (Enable::PIPE) PipeInit();
0053   if (Enable::MVTX) MvtxInit();
0054   if (Enable::INTT) InttInit();
0055   if (Enable::TPC) TPCInit();
0056   if (Enable::MICROMEGAS) MicromegasInit();
0057   if (Enable::MBD) MbdInit();
0058   if (Enable::CEMCALBEDO) CEmcAlbedoInit();
0059   if (Enable::CEMC) CEmcInit();
0060   if (Enable::HCALIN) HCalInnerInit();
0061   if (Enable::MAGNET) MagnetInit();
0062   MagnetFieldInit();  // We want the field - even if the magnet volume is disabled
0063   if (Enable::HCALOUT) HCalOuterInit();
0064   if (Enable::PLUGDOOR) PlugDoorInit();
0065   if (Enable::EPD) EPDInit();
0066   if (Enable::BEAMLINE)
0067   {
0068     BeamLineInit();
0069     if (Enable::ZDC) ZDCInit();
0070   }
0071   if (Enable::USER) UserInit();
0072   if (Enable::BLACKHOLE) BlackHoleInit();
0073 }
0074 
0075 int G4Setup()
0076 {
0077   //---------------
0078   // Fun4All server
0079   //---------------
0080 
0081   Fun4AllServer *se = Fun4AllServer::instance();
0082 
0083   PHG4Reco *g4Reco = new PHG4Reco();
0084   g4Reco->set_rapidity_coverage(1.1);  // according to drawings
0085   WorldInit(g4Reco);
0086   // PYTHIA 6
0087   if (G4P6DECAYER::decayType != EDecayType::kAll)
0088   {
0089     g4Reco->set_force_decay(G4P6DECAYER::decayType);
0090   }
0091   // EvtGen
0092   g4Reco->CustomizeEvtGenDecay(EVTGENDECAYER::DecayFile);
0093 
0094   double fieldstrength;
0095   std::istringstream stringline(G4MAGNET::magfield);
0096   stringline >> fieldstrength;
0097   if (stringline.fail())
0098   {  // conversion to double fails -> we have a string
0099 
0100     if (G4MAGNET::magfield.find("sphenix3dbigmapxyz") != std::string::npos ||
0101         G4MAGNET::magfield.find(".root") == std::string::npos)
0102     {
0103       g4Reco->set_field_map(G4MAGNET::magfield, PHFieldConfig::Field3DCartesian);
0104     }
0105     else
0106     {
0107       g4Reco->set_field_map(G4MAGNET::magfield, PHFieldConfig::kField2D);
0108     }
0109   }
0110   else
0111   {
0112     g4Reco->set_field(fieldstrength);                  // use const soleniodal field
0113     G4MAGNET::magfield_tracking = G4MAGNET::magfield;  // set tracking fieldmap to value
0114   }
0115   g4Reco->set_field_rescale(G4MAGNET::magfield_rescale);
0116 
0117   // the radius is an older protection against overlaps, it is not
0118   // clear how well this works nowadays but it doesn't hurt either
0119   double radius = 0.;
0120 
0121   if (Enable::PIPE) radius = Pipe(g4Reco, radius);
0122   if (Enable::MVTX) radius = Mvtx(g4Reco, radius);
0123   if (Enable::INTT) radius = Intt(g4Reco, radius);
0124   if (Enable::TPC) radius = TPC(g4Reco, radius);
0125   if (Enable::MICROMEGAS) Micromegas(g4Reco);
0126   if (Enable::MBD) Mbd(g4Reco);
0127   if (Enable::CEMCALBEDO) CEmcAlbedo(g4Reco);
0128   if (Enable::CEMC) radius = CEmc(g4Reco, radius, 8);
0129   if (Enable::HCALIN) radius = HCalInner(g4Reco, radius, 4);
0130   if (Enable::MAGNET) radius = Magnet(g4Reco, radius);
0131   if (Enable::HCALOUT) radius = HCalOuter(g4Reco, radius, 4);
0132   if (Enable::PLUGDOOR) PlugDoor(g4Reco);
0133   if (Enable::EPD) EPD(g4Reco);
0134   if (Enable::BEAMLINE)
0135   {
0136     BeamLineDefineMagnets(g4Reco);
0137     BeamLineDefineBeamPipe(g4Reco);
0138     if (Enable::ZDC) ZDCSetup(g4Reco);
0139   }
0140   if (Enable::USER) UserDetector(g4Reco);
0141 
0142   //----------------------------------------
0143   // BLACKHOLE
0144 
0145   if (Enable::BLACKHOLE) BlackHole(g4Reco, radius);
0146 
0147   PHG4TruthSubsystem *truth = new PHG4TruthSubsystem();
0148   g4Reco->registerSubsystem(truth);
0149 
0150   // finally adjust the world size in case the default is too small
0151   WorldSize(g4Reco, radius);
0152 
0153   se->registerSubsystem(g4Reco);
0154   return 0;
0155 }
0156 
0157 void ShowerCompress(int /*verbosity*/ = 0)
0158 {
0159   Fun4AllServer *se = Fun4AllServer::instance();
0160 
0161   PHG4DstCompressReco *compress = new PHG4DstCompressReco("PHG4DstCompressReco");
0162   compress->AddHitContainer("G4HIT_PIPE");
0163   compress->AddHitContainer("G4HIT_SVTXSUPPORT");
0164   compress->AddHitContainer("G4HIT_CEMC_ELECTRONICS");
0165   compress->AddHitContainer("G4HIT_CEMC");
0166   compress->AddHitContainer("G4HIT_ABSORBER_CEMC");
0167   compress->AddHitContainer("G4HIT_CEMC_SPT");
0168   compress->AddHitContainer("G4HIT_ABSORBER_HCALIN");
0169   compress->AddHitContainer("G4HIT_HCALIN");
0170   compress->AddHitContainer("G4HIT_HCALIN_SPT");
0171   compress->AddHitContainer("G4HIT_MAGNET");
0172   compress->AddHitContainer("G4HIT_ABSORBER_HCALOUT");
0173   compress->AddHitContainer("G4HIT_HCALOUT");
0174   compress->AddHitContainer("G4HIT_BH_1");
0175   compress->AddHitContainer("G4HIT_BH_FORWARD_PLUS");
0176   compress->AddHitContainer("G4HIT_BH_FORWARD_NEG");
0177   compress->AddHitContainer("G4HIT_BBC");
0178   compress->AddCellContainer("G4CELL_CEMC");
0179   compress->AddCellContainer("G4CELL_HCALIN");
0180   compress->AddCellContainer("G4CELL_HCALOUT");
0181   compress->AddTowerContainer("TOWER_SIM_CEMC");
0182   compress->AddTowerContainer("TOWER_RAW_CEMC");
0183   compress->AddTowerContainer("TOWER_CALIB_CEMC");
0184   compress->AddTowerContainer("TOWER_SIM_HCALIN");
0185   compress->AddTowerContainer("TOWER_RAW_HCALIN");
0186   compress->AddTowerContainer("TOWER_CALIB_HCALIN");
0187   compress->AddTowerContainer("TOWER_SIM_HCALOUT");
0188   compress->AddTowerContainer("TOWER_RAW_HCALOUT");
0189   compress->AddTowerContainer("TOWER_CALIB_HCALOUT");
0190   se->registerSubsystem(compress);
0191 
0192   return;
0193 }
0194 
0195 void DstCompress(Fun4AllDstOutputManager *out)
0196 {
0197   if (out)
0198   {
0199     out->StripNode("G4HIT_PIPE");
0200     out->StripNode("G4HIT_SVTXSUPPORT");
0201     out->StripNode("G4HIT_CEMC_ELECTRONICS");
0202     out->StripNode("G4HIT_CEMC");
0203     out->StripNode("G4HIT_ABSORBER_CEMC");
0204     out->StripNode("G4HIT_CEMC_SPT");
0205     out->StripNode("G4HIT_ABSORBER_HCALIN");
0206     out->StripNode("G4HIT_HCALIN");
0207     out->StripNode("G4HIT_HCALIN_SPT");
0208     out->StripNode("G4HIT_MAGNET");
0209     out->StripNode("G4HIT_ABSORBER_HCALOUT");
0210     out->StripNode("G4HIT_HCALOUT");
0211     out->StripNode("G4HIT_BH_1");
0212     out->StripNode("G4HIT_BH_FORWARD_PLUS");
0213     out->StripNode("G4HIT_BH_FORWARD_NEG");
0214     out->StripNode("G4HIT_BBC");
0215     out->StripNode("G4CELL_CEMC");
0216     out->StripNode("G4CELL_HCALIN");
0217     out->StripNode("G4CELL_HCALOUT");
0218   }
0219 }
0220 #endif