UE_in_pp/macro/HIJetReco.C

0001 #ifndef MACRO_HIJETRECO_C
0002 #define MACRO_HIJETRECO_C
0003
0004 #include <GlobalVariables.C>
0005
0006 #include <jetbase/FastJetAlgo.h>
0007 #include <jetbase/JetReco.h>
0008 #include <jetbase/TowerJetInput.h>
0009 #include <g4jets/TruthJetInput.h>
0010
0011 #include <jetbackground/CopyAndSubtractJets.h>
0012 #include <jetbackground/DetermineTowerBackground.h>
0013 #include <jetbackground/DetermineTowerRho.h>
0014 #include <jetbackground/FastJetAlgoSub.h>
0015 #include <jetbackground/RetowerCEMC.h>
0016 #include <jetbackground/SubtractTowers.h>
0017 #include <jetbackground/SubtractTowersCS.h>
0018 #include <jetbackground/TowerRho.h>
0019
0020 #include <globalvertex/GlobalVertex.h>
0021
0022 #include <fun4all/Fun4AllServer.h>
0023
0024 R__LOAD_LIBRARY(libjetbase.so)
0025 R__LOAD_LIBRARY(libg4jets.so)
0026 R__LOAD_LIBRARY(libjetbackground.so)
0027 R__LOAD_LIBRARY(libglobalvertex.so)
0028
0029 namespace Enable
0030 {
0031   bool HIJETS = false;
0032   int HIJETS_VERBOSITY = 0;
0033   bool HIJETS_MC = false;
0034   bool HIJETS_TRUTH = false;
0035 }  // namespace Enable
0036
0037 namespace HIJETS
0038 {
0039   bool do_flow = false; // should be set to true once the EPD event plane correction is implemented
0040   bool do_CS = false;
0041   bool is_pp = true;  // turn off functionality only relevant for nucleon collisions
0042   std::string tower_prefix = "TOWERINFO_CALIB";
0043   bool do_vertex_type = true;
0044   GlobalVertex::VTXTYPE vertex_type = GlobalVertex::MBD;
0045 }  // namespace HIJETS
0046
0047
0048 void HIJetReco()
0049 {
0050   int verbosity = std::max(Enable::VERBOSITY, Enable::HIJETS_VERBOSITY);
0051
0052
0053   //---------------
0054   // Fun4All server
0055   //---------------
0056   Fun4AllServer *se = Fun4AllServer::instance();
0057
0058   if (Enable::HIJETS_MC && Enable::HIJETS_TRUTH)
0059     {
0060       JetReco *truthjetreco = new JetReco();
0061       TruthJetInput *tji = new TruthJetInput(Jet::PARTICLE);
0062       tji->add_embedding_flag(0);  // changes depending on signal vs. embedded
0063       truthjetreco->add_input(tji);
0064       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.2), "AntiKt_Truth_r02");
0065       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.3), "AntiKt_Truth_r03");
0066       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.4), "AntiKt_Truth_r04");
0067       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.5), "AntiKt_Truth_r05");
0068       truthjetreco->set_algo_node("ANTIKT");
0069       truthjetreco->set_input_node("TRUTH");
0070       truthjetreco->Verbosity(verbosity);
0071       se->registerSubsystem(truthjetreco);
0072     }
0073
0074   RetowerCEMC *rcemc = new RetowerCEMC();
0075   rcemc->Verbosity(verbosity);
0076   rcemc->set_towerinfo(true);
0077   rcemc->set_frac_cut(0.5); //fraction of retower that must be masked to mask the full retower
0078   rcemc->set_towerNodePrefix(HIJETS::tower_prefix);
0079   se->registerSubsystem(rcemc);
0080
0081
0082   JetReco *towerjetreco = new JetReco();
0083   TowerJetInput *incemc = new TowerJetInput(Jet::CEMC_TOWERINFO_RETOWER,HIJETS::tower_prefix);
0084   TowerJetInput *inihcal = new TowerJetInput(Jet::HCALIN_TOWERINFO,HIJETS::tower_prefix);
0085   TowerJetInput *inohcal = new TowerJetInput(Jet::HCALOUT_TOWERINFO,HIJETS::tower_prefix);
0086   if (HIJETS::do_vertex_type)
0087   {
0088     incemc->set_GlobalVertexType(HIJETS::vertex_type);
0089     inihcal->set_GlobalVertexType(HIJETS::vertex_type);
0090     inohcal->set_GlobalVertexType(HIJETS::vertex_type);
0091   }
0092   towerjetreco->add_input(incemc);
0093   towerjetreco->add_input(inihcal);
0094   towerjetreco->add_input(inohcal);
0095   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.2), "AntiKt_TowerInfo_HIRecoSeedsRaw_r02");
0096   towerjetreco->set_algo_node("ANTIKT");
0097   towerjetreco->set_input_node("TOWER");
0098   towerjetreco->Verbosity(verbosity);
0099   se->registerSubsystem(towerjetreco);
0100
0101   DetermineTowerBackground *dtb = new DetermineTowerBackground();
0102   dtb->SetBackgroundOutputName("TowerInfoBackground_Sub1");
0103   dtb->SetFlow(HIJETS::do_flow);
0104   dtb->SetSeedType(0);
0105   dtb->SetSeedJetD(3);
0106   dtb->set_towerinfo(true);
0107   dtb->Verbosity(verbosity);
0108   dtb->set_towerNodePrefix(HIJETS::tower_prefix);
0109   se->registerSubsystem(dtb);
0110
0111   CopyAndSubtractJets *casj = new CopyAndSubtractJets();
0112   casj->SetFlowModulation(HIJETS::do_flow);
0113   casj->Verbosity(verbosity);
0114   casj->set_towerinfo(true);
0115   casj->set_towerNodePrefix(HIJETS::tower_prefix);
0116   se->registerSubsystem(casj);
0117
0118
0119   DetermineTowerBackground *dtb2 = new DetermineTowerBackground();
0120   dtb2->SetBackgroundOutputName("TowerInfoBackground_Sub2");
0121   dtb2->SetFlow(HIJETS::do_flow);
0122   dtb2->SetSeedType(1);
0123   dtb2->SetSeedJetPt(7);
0124   dtb2->Verbosity(verbosity);
0125   dtb2->set_towerinfo(true);
0126   dtb2->set_towerNodePrefix(HIJETS::tower_prefix);
0127   se->registerSubsystem(dtb2);
0128
0129
0130   SubtractTowers *st = new SubtractTowers();
0131   st->SetFlowModulation(HIJETS::do_flow);
0132   st->Verbosity(verbosity);
0133   st->set_towerinfo(true);
0134   st->set_towerNodePrefix(HIJETS::tower_prefix);
0135   se->registerSubsystem(st);
0136
0137   towerjetreco = new JetReco();
0138   incemc = new TowerJetInput(Jet::CEMC_TOWERINFO_SUB1,HIJETS::tower_prefix);
0139   inihcal = new TowerJetInput(Jet::HCALIN_TOWERINFO_SUB1,HIJETS::tower_prefix);
0140   inohcal = new TowerJetInput(Jet::HCALOUT_TOWERINFO_SUB1,HIJETS::tower_prefix);
0141   if (HIJETS::do_vertex_type)
0142   {
0143     incemc->set_GlobalVertexType(HIJETS::vertex_type);
0144     inihcal->set_GlobalVertexType(HIJETS::vertex_type);
0145     inohcal->set_GlobalVertexType(HIJETS::vertex_type);
0146   }
0147   towerjetreco->add_input(incemc);
0148   towerjetreco->add_input(inihcal);
0149   towerjetreco->add_input(inohcal);
0150   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.2, verbosity), "AntiKt_Tower_r02_Sub1");
0151   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.3, verbosity), "AntiKt_Tower_r03_Sub1");
0152   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.4, verbosity), "AntiKt_Tower_r04_Sub1");
0153   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.5, verbosity), "AntiKt_Tower_r05_Sub1");
0154   towerjetreco->set_algo_node("ANTIKT");
0155   towerjetreco->set_input_node("TOWER");
0156   towerjetreco->Verbosity(verbosity);
0157   se->registerSubsystem(towerjetreco);
0158
0159   return;
0160
0161 }
0162
0163
0164 // ----------------------------------------------------------------------------
0165 //! Determine rho from tower input to jet reco (necessary for jet QA)
0166 // ----------------------------------------------------------------------------
0167 void DoRhoCalculation()
0168 {
0169
0170   // set verbosity
0171   int verbosity = std::max(Enable::VERBOSITY, Enable::HIJETS_VERBOSITY);
0172
0173   //---------------
0174   // Fun4All server
0175   //---------------
0176   Fun4AllServer* se = Fun4AllServer::instance();
0177
0178   // run rho calculations w/ default parameters
0179   DetermineTowerRho* towRhoCalc = new DetermineTowerRho();
0180   towRhoCalc -> add_method(TowerRho::Method::AREA);
0181   towRhoCalc -> add_method(TowerRho::Method::MULT);
0182   towRhoCalc -> add_tower_input( new TowerJetInput(Jet::CEMC_TOWERINFO_RETOWER) );
0183   towRhoCalc -> add_tower_input( new TowerJetInput(Jet::HCALIN_TOWERINFO) );
0184   towRhoCalc -> add_tower_input( new TowerJetInput(Jet::HCALOUT_TOWERINFO) );
0185   se -> registerSubsystem( towRhoCalc );
0186
0187   // exit back to main macro
0188   return;
0189
0190 }
0191
0192 #endif