VandySkimmer/macros/Fun4All_VandySkimmer.C

0001 #include <string>
0002 #include <iostream>
0003 #include <format>
0004
0005 #include <fun4all/Fun4AllBase.h>
0006 #include <fun4all/Fun4AllUtils.h>
0007 #include <fun4all/Fun4AllServer.h>
0008 #include <fun4all/Fun4AllInputManager.h>
0009 #include <fun4all/Fun4AllDstInputManager.h>
0010 #include <fun4all/SubsysReco.h>
0011 #include <ffamodules/CDBInterface.h>
0012 #include <Calo_Calib.C>
0013
0014 #include <caloreco/RawClusterBuilderTopo.h>
0015
0016 #include <jetbase/JetReco.h>
0017 #include <jetbase/TowerJetInput.h>
0018 #include <jetbackground/FastJetAlgoSub.h>
0019
0020 #include <g4jets/TruthJetInput.h>
0021
0022 //#include <globalvertex/GlobalVertexReco.h>
0023 //#include <GlobalVertex.h>
0024 //#include <MbdDigitization.h>
0025 //#include <MbdReco.h>
0026
0027 #include <jetbackground/TimingCut.h>
0028
0029 #include <vandyskimmer/VandyJetDSTSkimmer.h>
0030 //#include <jetbackground/RetowerCEMC.h>
0031
0032
0033 R__LOAD_LIBRARY(libfun4all.so)
0034 R__LOAD_LIBRARY(libfun4allraw.so)
0035 R__LOAD_LIBRARY(libcalo_io.so)
0036 R__LOAD_LIBRARY(libffamodules.so)
0037 R__LOAD_LIBRARY(libVandySkimmer.so)
0038 R__LOAD_LIBRARY(libjetbase.so)
0039 R__LOAD_LIBRARY(libjetbackground.so)
0040
0041 void Fun4All_VandySkimmer(const std::string caloDSTlist, const std::string jetDSTlist, const std::string outDir = "/sphenix/tg/tg01/jets/bkimelman/wEEC/")
0042 {
0043
0044   bool doSim = false;
0045
0046   Fun4AllServer* se=Fun4AllServer::instance();
0047   se->Verbosity(0);
0048
0049   std::ifstream ifs(caloDSTlist);
0050   std::string filepath;
0051
0052   int runnumber = 0;
0053
0054   std::getline(ifs,filepath);
0055   std::pair<int,int> runseg = Fun4AllUtils::GetRunSegment(filepath);
0056   runnumber = runseg.first;
0057
0058   Fun4AllInputManager *inCalo = new Fun4AllDstInputManager("InputManagerCalo");
0059   inCalo->AddListFile(caloDSTlist);
0060   se->registerInputManager(inCalo);
0061
0062   Fun4AllInputManager *inJet = new Fun4AllDstInputManager("InputManagerJet");
0063   inJet->AddListFile(jetDSTlist);
0064   se->registerInputManager(inJet);
0065
0066   auto rc = recoConsts::instance();
0067   rc->set_StringFlag("CDB_GLOBALTAG", "ProdA_2024");
0068   rc->set_uint64Flag("TIMESTAMP", runnumber);
0069   CDBInterface::instance()->Verbosity(0);
0070   Process_Calo_Calib();
0071
0072   RawClusterBuilderTopo* ClusterBuilder = new RawClusterBuilderTopo("HcalRawClusterBuilderTopo");
0073   ClusterBuilder->Verbosity(0);
0074   ClusterBuilder->set_nodename("TOPOCLUSTER_ALLCALO");
0075   ClusterBuilder->set_enable_HCal(true);
0076   ClusterBuilder->set_enable_EMCal(true);
0077   ClusterBuilder->set_noise(0.0053, 0.0351, 0.0684); // 3sigma of pedestal noise
0078   ClusterBuilder->set_significance(4.0, 2.0, 1.0);
0079   ClusterBuilder->allow_corner_neighbor(true);
0080   ClusterBuilder->set_do_split(true);
0081   ClusterBuilder->set_minE_local_max(1.0, 2.0, 0.5);
0082   ClusterBuilder->set_R_shower(0.025);
0083   ClusterBuilder->set_use_only_good_towers(true);
0084   ClusterBuilder->set_absE(true);
0085   se->registerSubsystem(ClusterBuilder);
0086
0087   JetReco *towerjetreco = new JetReco("towerJetReco");
0088   TowerJetInput* emtji = new TowerJetInput(Jet::CEMC_TOWERINFO_RETOWER,"TOWERINFO_CALIB_RETOWER");
0089   TowerJetInput* ohtji = new TowerJetInput(Jet::HCALIN_TOWERINFO,"TOWERINFO_CALIB");
0090   TowerJetInput* ihtji = new TowerJetInput(Jet::HCALOUT_TOWERINFO,"TOWERINFO_CALIB");
0091   emtji->set_GlobalVertexType(GlobalVertex::VTXTYPE::MBD);
0092   zzohtji->set_GlobalVertexType(GlobalVertex::VTXTYPE::MBD);
0093   ihtji->set_GlobalVertexType(GlobalVertex::VTXTYPE::MBD);
0094   towerjetreco->add_input(emtji);
0095   towerjetreco->add_input(ohtji);
0096   towerjetreco->add_input(ihtji);
0097   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.2), "AntiKt_r02");
0098   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.3), "AntiKt_r03");
0099   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.4), "AntiKt_r04");
0100   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.5), "AntiKt_r05");
0101   towerjetreco->set_algo_node("ANTIKT");
0102   towerjetreco->set_input_node("TOWER");
0103   towerjetreco->Verbosity(0);
0104   se->registerSubsystem(towerjetreco);
0105
0106   if(doSim)
0107   {
0108     JetReco *truthjetreco = new JetReco("truthJetReco");
0109     truthjetreco->add_input(new TruthJetInput(Jet::PARTICLE));
0110     truthjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.2), "AntiKt_Truth_r02");
0111     truthjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.3), "AntiKt_Truth_r03");
0112     truthjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.4), "AntiKt_Truth_r04");
0113     truthjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.5), "AntiKt_Truth_r05");
0114     truthjetreco->set_algo_node("ANTIKT");
0115     truthjetreco->set_input_node("TRUTH");
0116     truthjetreco->Verbosity(0);
0117     se->registerSubsystem(truthjetreco);
0118   }
0119
0120   TimingCut* tc = new TimingCut("AntiKt_r04","TimingCutModule",false); //arguments: 1. jet node name to use, 2. name for fun4all module, 3. doAbort. The last of these tells the module whether or not to abort events that fail the cuts.
0121   tc->Verbosity(0);
0122   se->registerSubsystem(tc);
0123
0124   VandyJetDSTSkimmer *vs = new VandyJetDSTSkimmer("VandyJetDSTSkimmer");
0125   vs->SetRunnumber(runnumber);
0126   vs->SetOutfileName(std::format("{}/VandyDSTs_run2pp_ana521_2025p007_v001-{}.root",outDir, runnumber).c_str());
0127   vs->SetDoSim(doSim);
0128   se->registerSubsystem(vs);
0129
0130   se->run(-1);
0131   se->End();
0132   se->PrintTimer();
0133
0134   delete se;
0135   gSystem->Exit(0);
0136 }