Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2026-04-04 08:15:42 

0001 #ifndef MACRO_FUN4ALLG4EMBED_C
0002 #define MACRO_FUN4ALLG4EMBED_C
0003 
0004 #include <GlobalVariables.C>
0005 
0006 #include <G4Setup_sPHENIX.C>
0007 #include <G4_Input.C>
0008 #include <G4_Production.C>
0009 
0010 #include <globalvertex/GlobalVertexReco.h>
0011 
0012 #include <phpythia8/PHPy8ParticleTrigger.h>
0013 #include <phpythia8/PHPy8JetTrigger.h>
0014 
0015 #include <caloembedding/HepMCCollisionVertex.h>
0016 #include <caloembedding/caloTowerEmbed.h>
0017 #include <caloembedding/CopyIODataNodes.h>
0018 #include <caloembedding/CombineTowerInfo.h>
0019 
0020 /*
0021 #include <eventselection/EventSelector.h>
0022 #include <eventselection/MinBiasCut.h>
0023 #include <eventselection/TowerChi2Cut.h>
0024 #include <eventselection/ZVertexCut.h>
0025 */
0026 
0027 #include <jetbase/JetReco.h>
0028 #include <jetbase/TowerJetInput.h>
0029 #include <jetbase/FastJetAlgo.h>
0030 
0031 #include <jetbackground/FastJetAlgoSub.h>
0032 #include <jetbackground/RetowerCEMC.h>
0033 
0034 #include <g4jets/TruthJetInput.h>
0035 
0036 #include <centrality/CentralityInfov2.h>
0037 
0038 #include <calotrigger/MinimumBiasInfov1.h>
0039 
0040 #include <ffamodules/CDBInterface.h>
0041 #include <ffamodules/FlagHandler.h>
0042 #include <ffamodules/SyncReco.h>
0043 
0044 #include <fun4all/Fun4AllDstOutputManager.h>
0045 #include <fun4all/Fun4AllOutputManager.h>
0046 #include <fun4all/Fun4AllServer.h>
0047 #include <fun4all/Fun4AllSyncManager.h>
0048 #include <fun4all/Fun4AllRunNodeInputManager.h>
0049 #include <fun4all/Fun4AllUtils.h>
0050 
0051 #include <phool/PHRandomSeed.h>
0052 #include <phool/recoConsts.h>
0053 
0054 #include <TRandom3.h>
0055 
0056 #include <format>
0057 
0058 R__LOAD_LIBRARY(libfun4all.so)
0059 R__LOAD_LIBRARY(libffamodules.so)
0060 R__LOAD_LIBRARY(libCaloEmbedding.so)
0061 // R__LOAD_LIBRARY( libeventselection.so )
0062 R__LOAD_LIBRARY(libglobalvertex.so)
0063 R__LOAD_LIBRARY(libcentrality.so)
0064 R__LOAD_LIBRARY(libcalotrigger.so)
0065 R__LOAD_LIBRARY(libjetbase.so)
0066 R__LOAD_LIBRARY(libjetbackground.so)
0067 R__LOAD_LIBRARY(libg4jets.so)
0068 // R__LOAD_LIBRARY(libg4dst.so)
0069 //  For HepMC Hijing
0070 //  try inputFile = /sphenix/sim/sim01/sphnxpro/sHijing_HepMC/sHijing_0-12fm.dat
0071 
0072 // namespace OUTPUTMANAGER
0073 //{
0074 // set<string> outfiles;
0075 // }
0076 
0077 void AddCommonNodes(Fun4AllOutputManager *out);
0078 
0079 int Fun4All_G4_Embed_data(
0080     const int nEvents = 2,
0081     const int segment = 00000,
0082     const std::string &embed_input_file0 = "DST_CALOFITTING-00054404-00002.root",
0083     const std::string &outdir = "./",
0084     const std::string &/*outnameEnd*/ = "embed_test.root",
0085     const std::string &jettrigger = "Jet10",
0086     const std::string &cdbtag = "MDC2")
0087 {
0088 
0089   std::cout << "segment: " << segment << std::endl;
0090 
0091   gSystem->Load("libg4dst.so");
0092 
0093   Fun4AllServer *se = Fun4AllServer::instance();
0094   se->Verbosity(1);
0095 
0096   int verbosity = 0;
0097   CDBInterface::instance()->Verbosity(1);
0098 
0099   std::pair<int, int> runseg = Fun4AllUtils::GetRunSegment(embed_input_file0);
0100   int dataRunNumber = runseg.first; 
0101   int dataSegment   = runseg.second; 
0102 
0103   int runnumber = 21;
0104 
0105   // Opt to print all random seed used for debugging reproducibility. Comment out to reduce stdout prints.
0106   PHRandomSeed::Verbosity(1);
0107 
0108   // just if we set some flags somewhere in this macro
0109   recoConsts *rc = recoConsts::instance();
0110   // By default every random number generator uses
0111   // PHRandomSeed() which reads /dev/urandom to get its seed
0112   // if the RANDOMSEED flag is set its value is taken as seed
0113   // You can either set this to a random value using PHRandomSeed()
0114   // which will make all seeds identical (not sure what the point of
0115   // this would be:
0116   //  rc->set_IntFlag("RANDOMSEED",PHRandomSeed());
0117   // or set it to a fixed value so you can debug your code
0118   //  rc->set_IntFlag("RANDOMSEED", 12345);
0119 
0120   //===============
0121   // conditions DB flags
0122   //===============
0123   Enable::CDB = true;
0124   // tag
0125   rc->set_StringFlag("CDB_GLOBALTAG", cdbtag);
0126   // 64 bit timestamp
0127   //rc->set_uint64Flag("TIMESTAMP", CDB::timestamp);
0128   rc->set_uint64Flag("TIMESTAMP", runnumber);
0129 
0130   //  pair<int, int> runseg = Fun4AllUtils::GetRunSegment(embed_input_file0);
0131   // int runnumber = runseg.first;
0132   // int segment = runseg.second;
0133   if (runnumber != 0)
0134   {
0135     rc->set_IntFlag("RUNNUMBER", runnumber);
0136     Fun4AllSyncManager *syncman = se->getSyncManager();
0137     syncman->SegmentNumber(segment);
0138   }
0139 
0140   // Sync Headers and Flags
0141   SyncReco *sync = new SyncReco();
0142   se->registerSubsystem(sync);
0143 
0144   //  G4CEMC::useTowerInfoV2 = false;
0145   //  G4HCALOUT::useTowerInfoV2 = false;
0146   //  G4HCALIN::useTowerInfoV2 = false;
0147 
0148   //  G4CEMC::useTowerInfoV4 = true;
0149   // G4HCALOUT::useTowerInfoV4 = true;
0150   // G4HCALIN::useTowerInfoV4 = true;
0151 
0152   //===============
0153   // Input options
0154   //===============
0155   // verbosity setting (applies to all input managers)
0156   Input::VERBOSITY = 0;
0157   // First enable the input generators
0158   // Either:
0159   // read previously generated g4-hits files, in this case it opens a DST and skips
0160   // the simulations step completely. The G4Setup macro is only loaded to get information
0161   // about the number of layers used for the cell reco code
0162   //  Input::READHITS = true;
0163   // INPUTREADHITS::filename[0] = inputFile;
0164   // if you use a filelist
0165   // INPUTREADHITS::listfile[0] = inputFile;
0166   // Or:
0167   // Use particle generator
0168   // And
0169   // Further choose to embed newly simulated events to a previous simulation. Not compatible with `readhits = true`
0170   // In case embedding into a production output, please double check your G4Setup_sPHENIX.C and G4_*.C consistent with those in the production macro folder
0171   // E.g. /sphenix/sim//sim01/production/2016-07-21/single_particle/spacal2d/
0172 
0173   Input::EMBED = false;
0174   INPUTEMBED::filename[0] = embed_input_file0;
0175   // no repeating of embedding background, stop processing when end of file reached
0176   INPUTEMBED::REPEAT = false;
0177 
0178   // if you use a filelist
0179   // INPUTEMBED::listfile[0] = embed_input_file;
0180 
0181   // Input::SIMPLE = true;
0182   // Input::SIMPLE_NUMBER = 2; // if you need 2 of them
0183   // Input::SIMPLE_VERBOSITY = 1;
0184 
0185   //  Input::PYTHIA6 = true;
0186   // Enable this is emulating the nominal pp/pA/AA collision vertex distribution
0187   Input::BEAM_CONFIGURATION = Input::AA_COLLISION; // for 2023 sims we want the AA geometry for no pileup sims
0188 
0189   Input::PYTHIA8 = true;
0190   if (Input::PYTHIA8)
0191   {
0192     std::string pythia8_config_file = std::string(getenv("CALIBRATIONROOT")) + "/Generators/JetStructure_TG/";
0193     std::cout << "pythia config path: " << pythia8_config_file << std::endl;
0194     if (jettrigger == "Jet10")
0195     {
0196       pythia8_config_file += "phpythia8_10GeV_JS_MDC2.cfg";
0197     }
0198     else if (jettrigger == "Jet30")
0199     {
0200       pythia8_config_file += "phpythia8_30GeV_JS_MDC2.cfg";
0201     }
0202     else if (jettrigger == "Jet40")
0203     {
0204       pythia8_config_file += "phpythia8_40GeV_JS_MDC2.cfg";
0205     }
0206     else if (jettrigger == "PhotonJet20")
0207     {
0208       pythia8_config_file += "phpythia8_15GeV_JS_MDC2.cfg";
0209     }
0210     else
0211     {
0212       std::cout << "invalid jettrigger: " << jettrigger << std::endl;
0213       gSystem->Exit(1);
0214     }
0215     PYTHIA8::config_file[0] = pythia8_config_file;
0216   }
0217 
0218   //-----------------
0219   // Initialize the selected Input/Event generation
0220   //-----------------
0221   InputInit();
0222 
0223   //--------------
0224   // Set generator specific options
0225   //--------------
0226   // can only be set after InputInit() is called
0227 
0228   if (Input::PYTHIA8)
0229   {
0230     //! apply sPHENIX nominal beam parameter with 2mrad crossing as defined in sPH-TRG-2020-001
0231     PHPy8JetTrigger *p8_js_signal_trigger = new PHPy8JetTrigger();
0232     p8_js_signal_trigger->SetEtaHighLow(1.5, -1.5); // Set eta acceptance for particles into the jet between +/- 1.5
0233     p8_js_signal_trigger->SetJetR(0.4);             // Set the radius for the trigger jet
0234 
0235     PHPy8ParticleTrigger *p8_photon_jet_trigger = new PHPy8ParticleTrigger();
0236     p8_photon_jet_trigger->SetStableParticleOnly(false);  // process unstable particles that include quarks
0237     p8_photon_jet_trigger->AddParticles(22);
0238     p8_photon_jet_trigger->SetEtaHighLow(1.5, -1.5);  // sample a rapidity range higher than the sPHENIX tracking pseudorapidity
0239     std::vector<int> partentsId{1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, -22, -21, -20, -19, -18, -17, -16, -15, -14, -13, -12, -11, -10, -9, -8, -7, -6, -5, -4, -3, -2, -1};
0240     p8_photon_jet_trigger->AddParents(partentsId);
0241 
0242 
0243     if (jettrigger == "Jet10")
0244     {
0245       p8_js_signal_trigger->SetMinJetPt(10);
0246     }
0247     else if (jettrigger == "Jet30")
0248     {
0249       p8_js_signal_trigger->SetMinJetPt(30);
0250     }
0251     else if (jettrigger == "Jet40")
0252     {
0253       p8_js_signal_trigger->SetMinJetPt(40);
0254     }
0255     
0256     if (jettrigger == "Jet10" ||
0257         jettrigger == "Jet30" ||
0258         jettrigger == "Jet40")
0259     {
0260       delete p8_photon_jet_trigger;
0261       p8_photon_jet_trigger = nullptr;
0262     }
0263     else if (jettrigger == "PhotonJet20")
0264     {
0265       delete p8_js_signal_trigger;
0266       p8_js_signal_trigger = nullptr;
0267       p8_photon_jet_trigger->SetPtLow(20);
0268     }
0269     else
0270     {
0271       std::cout << "invalid jettrigger: " << jettrigger << std::endl;
0272       gSystem->Exit(1);
0273     }
0274 
0275     if (p8_js_signal_trigger)
0276     {
0277       INPUTGENERATOR::Pythia8[0]->register_trigger(p8_js_signal_trigger);
0278       INPUTGENERATOR::Pythia8[0]->set_trigger_AND();
0279     }
0280     else if (p8_photon_jet_trigger)
0281     {
0282       INPUTGENERATOR::Pythia8[0]->register_trigger(p8_photon_jet_trigger);
0283       INPUTGENERATOR::Pythia8[0]->set_trigger_AND();
0284     }
0285 
0286     Input::ApplysPHENIXBeamParameter(INPUTGENERATOR::Pythia8[0]);
0287   }
0288 
0289   //--------------
0290   // Set Input Manager specific options
0291   //--------------
0292   // can only be set after InputInit() is called
0293 
0294   // register all input generators with Fun4All
0295 
0296   std::string DataTopNode = "TOPData";
0297   /*
0298   EventSelector *EventSelectorHandler = new EventSelector();
0299   EventSelectorHandler->Verbosity( 1 );
0300   auto mbc = new MinBiasCut( );
0301   mbc->SetNodeName("MinimumBiasInfo");
0302   EventSelectorHandler->AddCut(mbc);
0303 
0304   auto tcc = new TowerChi2Cut( );
0305   std::vector<std::string> TowerChi2Nodes = {
0306     "TOWERINFO_CALIB_CEMC",
0307     "TOWERINFO_CALIB_HCALIN",
0308     "TOWERINFO_CALIB_HCALOUT"
0309   };
0310   tcc->SetNodeNames(TowerChi2Nodes);
0311   EventSelectorHandler->AddCut( tcc );
0312 
0313   auto zvc =  new ZVertexCut( 20, -20 );
0314   zvc->SetNodeName( "GlobalVertexMap" );
0315   EventSelectorHandler->AddCut( zvc );
0316 
0317   EventSelectorHandler->PrintCuts();
0318   se->registerSubsystem( EventSelectorHandler, DataTopNode);
0319 
0320   */
0321 
0322   HepMCCollisionVertex *collvtx = new HepMCCollisionVertex();
0323   collvtx->Verbosity(0);
0324   se->registerSubsystem(collvtx);
0325 
0326   InputRegister();
0327 
0328   FlagHandler *flag = new FlagHandler();
0329   se->registerSubsystem(flag);
0330 
0331   // set up production relatedstuff
0332   Enable::PRODUCTION = true;
0333 
0334   //======================
0335   // Write the DST
0336   //======================
0337 
0338   Enable::DSTOUT = true;
0339   Enable::DSTOUT_COMPRESS = false;
0340   DstOut::OutputDir = outdir;
0341 
0342   if (Enable::PRODUCTION)
0343   {
0344     PRODUCTION::SaveOutputDir = DstOut::OutputDir;
0345     //    Production_CreateOutputDir();
0346   }
0347 
0348 
0349   //-----------------
0350   // Copy over RunHeader, EventHeader, centrality, globalvertex, and MinimumBiasInfo nodes
0351   //-----------------
0352   CopyIODataNodes *cp = new CopyIODataNodes();
0353   cp->set_CopyTowerInfo("TOWERINFO_CALIB_SEPD","TOWERINFO_CALIB_SEPD_data");
0354   cp->Verbosity(1);
0355   se->registerSubsystem(cp, DataTopNode);
0356 
0357   //======================
0358   // What to run
0359   //======================
0360 
0361   // Global options (enabled for all enables subsystems - if implemented)
0362   //  Enable::ABSORBER = true;
0363   //  Enable::OVERLAPCHECK = true;
0364   //  Enable::VERBOSITY = 1;
0365 
0366   Enable::MBD = true;
0367 
0368   // Enable::PIPE = true;
0369 
0370   // central tracking
0371   // Enable::MVTX = true;
0372 
0373   // Enable::INTT = true;
0374 
0375   // Enable::TPC = true;
0376 
0377   // Enable::MICROMEGAS = true;
0378 
0379   Enable::CEMC = true;
0380   // Enable::CEMC_ABSORBER = true;
0381   // Enable::CEMC_CELL = Enable::CEMC && true;
0382   // Enable::CEMC_TOWER = Enable::CEMC_CELL && true;
0383 
0384   Enable::HCALIN = true;
0385   // Enable::HCALIN_ABSORBER = true;
0386   // Enable::HCALIN_CELL = Enable::HCALIN && true;
0387   // Enable::HCALIN_TOWER = Enable::HCALIN_CELL && true;
0388 
0389   Enable::MAGNET = true;
0390   // Enable::MAGNET_ABSORBER = true;
0391 
0392   Enable::HCALOUT = true;
0393   // Enable::HCALOUT_ABSORBER = true;
0394   // Enable::HCALOUT_CELL = Enable::HCALOUT && true;
0395   // Enable::HCALOUT_TOWER = Enable::HCALOUT_CELL && true;
0396 
0397   Enable::EPD = true;
0398   Enable::EPD_TILE = Enable::EPD && true;
0399 
0400   // new settings using Enable namespace in GlobalVariables.C
0401   Enable::BLACKHOLE = true;
0402   Enable::BLACKHOLE_FORWARD_SAVEHITS = false; // disable forward/backward hits
0403   // Enable::BLACKHOLE_SAVEHITS = false; // turn off saving of bh hits
0404   // BlackHoleGeometry::visible = true;
0405 
0406   G4MAGNET::magfield = "FIELDMAP_GAP";
0407   G4MAGNET::magfield_OHCAL_steel = "FIELDMAP_STEEL";
0408   G4MAGNET::magfield_tracking = "FIELDMAP_TRACKING";
0409 
0410   // G4MAGNET::magfield = "/cvmfs/sphenix.sdcc.bnl.gov/calibrations/sphnxpro/cdb/FIELDMAP_GAP/65/a9/65a930ed6de9c0e049cd0f3ef226e6b4_sphenix3dbigmapxyz_gap_rebuild_v2.root";
0411   // G4MAGNET::magfield_OHCAL_steel = "/cvmfs/sphenix.sdcc.bnl.gov/calibrations/sphnxpro/cdb/FIELDMAP_STEEL/ea/f0/eaf07c2a2def280428b9aeea99d14fd1_sphenix3dbigmapxyz_steel_rebuild_v2.root";
0412   // G4MAGNET::magfield_tracking = "/cvmfs/sphenix.sdcc.bnl.gov/calibrations/sphnxpro/cdb/FIELDMAP_TRACKING/8e/4d/8e4d6c3b1660540a658da3a275af2bde_sphenix3dtrackingmapxyz.root";
0413 
0414   // Initialize the selected subsystems
0415   G4Init();
0416 
0417   //---------------------
0418   // GEANT4 Detector description
0419   //---------------------
0420   if (!Input::READHITS)
0421   {
0422     G4Setup();
0423   }
0424 
0425   //--------------------
0426   // GEANT4 Calorimeters
0427   //--------------------
0428 
0429   if (Enable::CEMC_CELL)
0430     CEMC_Cells();
0431 
0432   if (Enable::HCALIN_CELL)
0433     HCALInner_Cells();
0434 
0435   if (Enable::HCALOUT_CELL)
0436     HCALOuter_Cells();
0437 
0438   if (Enable::CEMC_TOWER)
0439     CEMC_Towers();
0440 
0441   if (Enable::HCALIN_TOWER)
0442     HCALInner_Towers();
0443 
0444   if (Enable::HCALOUT_TOWER)
0445     HCALOuter_Towers();
0446 
0447   if (Enable::EPD_TILE) EPD_Tiles();
0448 
0449   //truth jets
0450   JetReco *truthjets2 = new JetReco("TRUTHJETRECO2");
0451   truthjets2->add_input(new TruthJetInput(Jet::PARTICLE));
0452   truthjets2->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.2), "AntiKt_Truth_r02");
0453   truthjets2->set_algo_node("ANTIKT");
0454   truthjets2->set_input_node("TRUTH");
0455   truthjets2->Verbosity(verbosity);
0456   se->registerSubsystem(truthjets2);
0457 
0458   JetReco *truthjets3 = new JetReco("TRUTHJETRECO3");
0459   truthjets3->add_input(new TruthJetInput(Jet::PARTICLE));
0460   truthjets3->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.3), "AntiKt_Truth_r03");
0461   truthjets3->set_algo_node("ANTIKT");
0462   truthjets3->set_input_node("TRUTH");
0463   truthjets3->Verbosity(verbosity);
0464   se->registerSubsystem(truthjets3);
0465 
0466   JetReco *truthjets4 = new JetReco("TRUTHJETRECO4");
0467   truthjets4->add_input(new TruthJetInput(Jet::PARTICLE));
0468   truthjets4->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.4), "AntiKt_Truth_r04");
0469   truthjets4->set_algo_node("ANTIKT");
0470   truthjets4->set_input_node("TRUTH");
0471   truthjets4->Verbosity(verbosity);
0472   se->registerSubsystem(truthjets4);
0473 
0474   JetReco *truthjets5 = new JetReco("TRUTHJETRECO5");
0475   truthjets5->add_input(new TruthJetInput(Jet::PARTICLE));
0476   truthjets5->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.5), "AntiKt_Truth_r05");
0477   truthjets5->set_algo_node("ANTIKT");
0478   truthjets5->set_input_node("TRUTH");
0479   truthjets5->Verbosity(verbosity);
0480   se->registerSubsystem(truthjets5);
0481 
0482   JetReco *truthjets6 = new JetReco("TRUTHJETRECO6");
0483   truthjets6->add_input(new TruthJetInput(Jet::PARTICLE));
0484   truthjets6->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.6), "AntiKt_Truth_r06");
0485   truthjets6->set_algo_node("ANTIKT");
0486   truthjets6->set_input_node("TRUTH");
0487   truthjets6->Verbosity(verbosity);
0488   se->registerSubsystem(truthjets6);
0489 
0490   JetReco *truthjets7 = new JetReco("TRUTHJETRECO7");
0491   truthjets7->add_input(new TruthJetInput(Jet::PARTICLE));
0492   truthjets7->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.7), "AntiKt_Truth_r07");
0493   truthjets7->set_algo_node("ANTIKT");
0494   truthjets7->set_input_node("TRUTH");
0495   truthjets7->Verbosity(verbosity);
0496   se->registerSubsystem(truthjets7);
0497 
0498   JetReco *truthjets8 = new JetReco("TRUTHJETRECO8");
0499   truthjets8->add_input(new TruthJetInput(Jet::PARTICLE));
0500   truthjets8->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.8), "AntiKt_Truth_r08");
0501   truthjets8->set_algo_node("ANTIKT");
0502   truthjets8->set_input_node("TRUTH");
0503   truthjets8->Verbosity(verbosity);
0504   se->registerSubsystem(truthjets8);
0505 
0506    std::string test = CDBInterface::instance()->getUrl("CEMC_meanTime");// calling this line somehow prevents CDB bug when switching global tag 
0507    std::cout << "line to avoid CDB bug " << test << std::endl;
0508 
0509   std::string save_globaltag = rc->get_StringFlag("CDB_GLOBALTAG");
0510   int save_timestamp    = rc->get_uint64Flag("TIMESTAMP");
0511 
0512   rc->set_StringFlag("CDB_GLOBALTAG", "ProdA_2024");
0513   rc->set_uint64Flag("TIMESTAMP", dataRunNumber);
0514 
0515   std::string cemc_datacalib = CDBInterface::instance()->getUrl("CEMC_calib_ADC_to_ETower");
0516   std::string ohcal_datacalib = CDBInterface::instance()->getUrl("HCALOUT_calib_ADC_to_ETower");
0517   std::string ihcal_datacalib = CDBInterface::instance()->getUrl("HCALIN_calib_ADC_to_ETower");
0518 
0519   std::cout << "using data calibration" << std::endl << cemc_datacalib  << std::endl << ohcal_datacalib << std::endl << ihcal_datacalib  << std::endl;
0520 
0521   rc->set_StringFlag("CDB_GLOBALTAG", save_globaltag);
0522   rc->set_uint64Flag("TIMESTAMP", save_timestamp);
0523  
0524   // waveform sim
0525   CaloWaveformSim *caloWaveformSim = new CaloWaveformSim();
0526   caloWaveformSim->set_detector_type(CaloTowerDefs::CEMC);
0527   caloWaveformSim->set_detector("CEMC");
0528   caloWaveformSim->set_nsamples(12);
0529   caloWaveformSim->set_pedestalsamples(12);
0530   caloWaveformSim->set_timewidth(0.2);
0531   caloWaveformSim->set_peakpos(6);
0532   caloWaveformSim->set_noise_type(CaloWaveformSim::NOISE_NONE);
0533   caloWaveformSim->set_directURL_calib(cemc_datacalib);
0534   se->registerSubsystem(caloWaveformSim);
0535 
0536   caloWaveformSim = new CaloWaveformSim();
0537   caloWaveformSim->set_detector_type(CaloTowerDefs::HCALIN);
0538   caloWaveformSim->set_detector("HCALIN");
0539   caloWaveformSim->set_nsamples(12);
0540   caloWaveformSim->set_pedestalsamples(12);
0541   caloWaveformSim->set_timewidth(0.2);
0542   caloWaveformSim->set_peakpos(6);
0543   caloWaveformSim->set_noise_type(CaloWaveformSim::NOISE_NONE);
0544   caloWaveformSim->set_directURL_calib(ihcal_datacalib);
0545   se->registerSubsystem(caloWaveformSim);
0546 
0547   caloWaveformSim = new CaloWaveformSim();
0548   caloWaveformSim->set_detector_type(CaloTowerDefs::HCALOUT);
0549   caloWaveformSim->set_detector("HCALOUT");
0550   caloWaveformSim->set_nsamples(12);
0551   caloWaveformSim->set_pedestalsamples(12);
0552   caloWaveformSim->set_timewidth(0.2);
0553   caloWaveformSim->set_peakpos(6);
0554   caloWaveformSim->set_noise_type(CaloWaveformSim::NOISE_NONE);
0555   caloWaveformSim->set_directURL_calib(ohcal_datacalib);
0556   se->registerSubsystem(caloWaveformSim);
0557 
0558   se->Print("NODETREE");
0559 
0560   //-----------------
0561   // Calorimeter embedding
0562   //-----------------
0563 
0564   caloTowerEmbed *embedder_CEMC = new caloTowerEmbed("embedder_CEMC");
0565   embedder_CEMC->set_detector_type(CaloTowerDefs::CEMC);
0566   embedder_CEMC->set_removeBadTowers(false);
0567   embedder_CEMC->set_embedwaveform(true);
0568   embedder_CEMC->set_inputNodePrefix("TOWERS_");
0569   embedder_CEMC->Verbosity(verbosity);
0570   se->registerSubsystem(embedder_CEMC);
0571 
0572   caloTowerEmbed *embedder_IHCAL = new caloTowerEmbed("embedder_IHCAL");
0573   embedder_IHCAL->set_detector_type(CaloTowerDefs::HCALIN);
0574   embedder_IHCAL->set_removeBadTowers(false);
0575   embedder_IHCAL->set_embedwaveform(true);
0576   embedder_IHCAL->set_inputNodePrefix("TOWERS_");
0577   embedder_IHCAL->Verbosity(verbosity);
0578   se->registerSubsystem(embedder_IHCAL);
0579 
0580   caloTowerEmbed *embedder_OHCAL = new caloTowerEmbed("embedder_OHCal");
0581   embedder_OHCAL->set_detector_type(CaloTowerDefs::HCALOUT);
0582   embedder_OHCAL->set_removeBadTowers(false);
0583   embedder_OHCAL->set_embedwaveform(true);
0584   embedder_OHCAL->set_inputNodePrefix("TOWERS_");
0585   embedder_OHCAL->Verbosity(verbosity);
0586   se->registerSubsystem(embedder_OHCAL);
0587 
0588   // waveform fitting
0589 
0590   CaloTowerBuilder *ca2 = new CaloTowerBuilder();
0591   ca2->set_detector_type(CaloTowerDefs::CEMC);
0592   ca2->set_nsamples(12);
0593   ca2->set_dataflag(false);
0594   ca2->set_processing_type(CaloWaveformProcessing::TEMPLATE);
0595   ca2->set_builder_type(CaloTowerDefs::kWaveformTowerSimv1);
0596   // match our current ZS threshold ~60ADC for emcal
0597   ca2->set_softwarezerosuppression(true, 60);
0598   se->registerSubsystem(ca2);
0599 
0600   ca2 = new CaloTowerBuilder();
0601   ca2->set_detector_type(CaloTowerDefs::HCALIN);
0602   ca2->set_nsamples(12);
0603   ca2->set_dataflag(false);
0604   ca2->set_processing_type(CaloWaveformProcessing::TEMPLATE);
0605   ca2->set_builder_type(CaloTowerDefs::kWaveformTowerSimv1);
0606   ca2->set_softwarezerosuppression(true, 30);
0607   se->registerSubsystem(ca2);
0608 
0609   ca2 = new CaloTowerBuilder();
0610   ca2->set_detector_type(CaloTowerDefs::HCALOUT);
0611   ca2->set_nsamples(12);
0612   ca2->set_dataflag(false);
0613   ca2->set_processing_type(CaloWaveformProcessing::TEMPLATE);
0614   ca2->set_builder_type(CaloTowerDefs::kWaveformTowerSimv1);
0615   ca2->set_softwarezerosuppression(true, 30);
0616   se->registerSubsystem(ca2);
0617 
0618 
0619   ///////////////////////
0620   // EPD stuff
0621   CombineTowerInfo* comb = new CombineTowerInfo("CombineHCAL");
0622   comb->set_inputNodeA("TOWERINFO_CALIB_SEPD_data");
0623   comb->set_inputNodeB("TOWERINFO_CALIB_EPD");
0624   comb->set_outputNode("TOWERINFO_COMBINED_SEPD");
0625   comb->set_detector("EPD");
0626   se->registerSubsystem(comb);
0627 
0628 
0629 
0630   //--------------
0631   // Set up Input Managers
0632   //--------------
0633 
0634   for (auto &iter : INPUTEMBED::filename)
0635   {
0636     std::string mgrname = "DSTin" + std::to_string(iter.first);
0637     Fun4AllInputManager *hitsin = new Fun4AllDstInputManager(mgrname, "DST", DataTopNode);
0638     hitsin->fileopen(iter.second);
0639     hitsin->Verbosity(Input::VERBOSITY);
0640     if (INPUTEMBED::REPEAT)
0641     {
0642       hitsin->Repeat();
0643     }
0644     se->registerInputManager(hitsin);
0645   }
0646   // input for pedestals
0647   {
0648     TRandom3 randGen;
0649     // get seed
0650     unsigned int seed = PHRandomSeed();
0651     randGen.SetSeed(seed);
0652     // a int from 0 to 3259
0653     int sequence = randGen.Integer(3260);
0654     // pad the name
0655     std::ostringstream opedname;
0656     opedname << "pedestal-54256-0" << std::setw(4) << std::setfill('0') << sequence << ".root";
0657 
0658     std::string pedestalname = opedname.str();
0659 
0660     Fun4AllInputManager *hitsin = new Fun4AllNoSyncDstInputManager("DST2");
0661     hitsin->AddFile(pedestalname);
0662     hitsin->Repeat();
0663     se->registerInputManager(hitsin);
0664   }
0665   // input geom files
0666   {
0667     Fun4AllInputManager *ingeo = new Fun4AllRunNodeInputManager("DST_GEO");
0668     std::string geoLocation = CDBInterface::instance()->getUrl("calo_geo");
0669     ingeo->AddFile(geoLocation);
0670     se->registerInputManager(ingeo);
0671   }
0672 
0673   //  InputManagers();
0674 
0675   std::string outnameEnd2 =
0676     std::format("{}-{:08}-", jettrigger, runnumber) +
0677     std::format("{}-{:08}-{:05}.root","data", dataRunNumber, dataSegment);
0678 
0679   std::string FullOutFile = outdir + "DST_TRUTH_G4HIT_" + outnameEnd2;
0680   Fun4AllOutputManager *out = new Fun4AllDstOutputManager("TRUTHOUT", FullOutFile);
0681   AddCommonNodes(out);
0682   out->AddNode("G4TruthInfo");
0683   out->AddNode("G4HIT_BH_1");
0684   out->AddNode("PHHepMCGenEventMap");
0685   out->Verbosity(0);
0686   se->registerOutputManager(out);
0687   // OUTPUTMANAGER::outfiles.insert(FullOutFile);
0688   se->registerOutputManager(out);
0689 
0690   FullOutFile = outdir + "DST_TRUTH_JET_" + outnameEnd2;
0691   out = new Fun4AllDstOutputManager("JETOUT", FullOutFile);
0692   AddCommonNodes(out);
0693   out->AddNode("AntiKt_Truth_r02");
0694   out->AddNode("AntiKt_Truth_r03");
0695   out->AddNode("AntiKt_Truth_r04");
0696   out->AddNode("AntiKt_Truth_r05");
0697   out->AddNode("AntiKt_Truth_r06");
0698   out->AddNode("AntiKt_Truth_r07");
0699   out->AddNode("AntiKt_Truth_r08");
0700   out->AddNode("AntiKt_Sim_r02");
0701   out->AddNode("AntiKt_Sim_r03");
0702   out->AddNode("AntiKt_Sim_r04");
0703   out->AddNode("AntiKt_Sim_r05");
0704   out->AddNode("AntiKt_Sim_r06");
0705   out->AddNode("AntiKt_Sim_r07");
0706   out->AddNode("AntiKt_Sim_r08");
0707   out->AddNode("AntiKt_Sim_Retower_r02");
0708   out->AddNode("AntiKt_Sim_Retower_r03");
0709   out->AddNode("AntiKt_Sim_Retower_r04");
0710   out->AddNode("AntiKt_Sim_Retower_r05");
0711   out->AddNode("AntiKt_Sim_Retower_r06");
0712   out->AddNode("AntiKt_Sim_Retower_r07");
0713   out->AddNode("AntiKt_Sim_Retower_r08");
0714   out->Verbosity(0);
0715   se->registerOutputManager(out);
0716   // OUTPUTMANAGER::outfiles.insert(FullOutFile);
0717 
0718   FullOutFile = outdir + "DST_CALO_" + outnameEnd2;
0719   out = new Fun4AllDstOutputManager("CALOPROCESSED", FullOutFile);
0720   AddCommonNodes(out);
0721   out->AddNode("TOWERS_CEMC");
0722   out->AddNode("TOWERS_HCALIN");
0723   out->AddNode("TOWERS_HCALOUT");
0724   out->Verbosity(0);
0725   se->registerOutputManager(out);
0726   // OUTPUTMANAGER::outfiles.insert(FullOutFile);
0727 
0728   FullOutFile = outdir + "DST_GLOBAL_" + outnameEnd2;
0729   out = new Fun4AllDstOutputManager("GLOBALOUT", FullOutFile);
0730   AddCommonNodes(out);
0731   out->AddNode("GlobalVertexMap");
0732   out->AddNode("MinimumBiasInfo");
0733   out->AddNode("CentralityInfo");
0734   out->AddNode("MbdOut");
0735   out->AddNode("MbdPmtContainer");
0736   out->AddNode("MinimumBiasInfo");
0737   out->AddNode("CentralityInfo");
0738   out->AddNode("TOWERINFO_CALIB_SEPD_data");
0739   out->AddNode("TOWERINFO_CALIB_EPD");
0740   out->AddNode("TOWERINFO_COMBINED_SEPD");
0741   out->Verbosity(0);
0742   se->registerOutputManager(out);
0743   // OUTPUTMANAGER::outfiles.insert(FullOutFile);
0744 
0745   //-----------------
0746   // Event processing
0747   //-----------------
0748   // if we use a negative number of events we go back to the command line here
0749   if (nEvents < 0)
0750   {
0751     return 0;
0752   }
0753   // if we run the particle generator and use 0 it'll run forever
0754   if (nEvents == 0 && !Input::HEPMC && !Input::READHITS && INPUTEMBED::REPEAT)
0755   {
0756     std::cout << "using 0 for number of events is a bad idea when using particle generators" << std::endl;
0757     std::cout << "it will run forever, so I just return without running anything" << std::endl;
0758     return 0;
0759   }
0760 
0761   se->run(nEvents);
0762   se->PrintTimer();
0763 
0764   //-----
0765   // Exit
0766   //-----
0767 
0768   CDBInterface::instance()->Print(); // print used DB files
0769   se->End();
0770   std::cout << "All done" << std::endl;
0771   delete se;
0772   gSystem->Exit(0);
0773   return 0;
0774 }
0775 
0776 void AddCommonNodes(Fun4AllOutputManager *out)
0777 {
0778   out->AddNode("Sync");
0779   out->AddNode("EventHeader");
0780   return;
0781 }
0782 
0783 #endif
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!