Jet-Study/macros/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   // Fun4All server
0054   //---------------
0055   Fun4AllServer *se = Fun4AllServer::instance();
0056
0057   if (Enable::HIJETS_MC && Enable::HIJETS_TRUTH)
0058     {
0059       JetReco *truthjetreco = new JetReco();
0060       TruthJetInput *tji = new TruthJetInput(Jet::PARTICLE);
0061       tji->add_embedding_flag(0);  // changes depending on signal vs. embedded
0062       truthjetreco->add_input(tji);
0063       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.2), "AntiKt_Truth_r02");
0064       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.3), "AntiKt_Truth_r03");
0065       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.4), "AntiKt_Truth_r04");
0066       truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.5), "AntiKt_Truth_r05");
0067       truthjetreco->set_algo_node("ANTIKT");
0068       truthjetreco->set_input_node("TRUTH");
0069       truthjetreco->Verbosity(verbosity);
0070       se->registerSubsystem(truthjetreco);
0071     }
0072
0073   RetowerCEMC *rcemc = new RetowerCEMC();
0074   rcemc->Verbosity(verbosity);
0075   rcemc->set_towerinfo(true);
0076   rcemc->set_frac_cut(0.5); //fraction of retower that must be masked to mask the full retower
0077   rcemc->set_towerNodePrefix(HIJETS::tower_prefix);
0078   se->registerSubsystem(rcemc);
0079
0080   JetReco *towerjetreco = new JetReco();
0081   TowerJetInput *incemc = new TowerJetInput(Jet::CEMC_TOWERINFO_RETOWER,HIJETS::tower_prefix);
0082   TowerJetInput *inihcal = new TowerJetInput(Jet::HCALIN_TOWERINFO,HIJETS::tower_prefix);
0083   TowerJetInput *inohcal = new TowerJetInput(Jet::HCALOUT_TOWERINFO,HIJETS::tower_prefix);
0084   if (HIJETS::do_vertex_type)
0085   {
0086     incemc->set_GlobalVertexType(HIJETS::vertex_type);
0087     inihcal->set_GlobalVertexType(HIJETS::vertex_type);
0088     inohcal->set_GlobalVertexType(HIJETS::vertex_type);
0089   }
0090   towerjetreco->add_input(incemc);
0091   towerjetreco->add_input(inihcal);
0092   towerjetreco->add_input(inohcal);
0093   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.2), "AntiKt_TowerInfo_HIRecoSeedsRaw_r02");
0094   towerjetreco->set_algo_node("ANTIKT");
0095   towerjetreco->set_input_node("TOWER");
0096   towerjetreco->Verbosity(verbosity);
0097   se->registerSubsystem(towerjetreco);
0098
0099   DetermineTowerBackground *dtb = new DetermineTowerBackground();
0100   dtb->SetBackgroundOutputName("TowerInfoBackground_Sub1");
0101   dtb->SetFlow(HIJETS::do_flow);
0102   dtb->SetSeedType(0);
0103   dtb->SetSeedJetD(3);
0104   dtb->set_towerinfo(true);
0105   dtb->Verbosity(verbosity);
0106   dtb->set_towerNodePrefix(HIJETS::tower_prefix);
0107   se->registerSubsystem(dtb);
0108
0109   CopyAndSubtractJets *casj = new CopyAndSubtractJets();
0110   casj->SetFlowModulation(HIJETS::do_flow);
0111   casj->Verbosity(verbosity);
0112   casj->set_towerinfo(true);
0113   casj->set_towerNodePrefix(HIJETS::tower_prefix);
0114   se->registerSubsystem(casj);
0115
0116   DetermineTowerBackground *dtb2 = new DetermineTowerBackground();
0117   dtb2->SetBackgroundOutputName("TowerInfoBackground_Sub2");
0118   dtb2->SetFlow(HIJETS::do_flow);
0119   dtb2->SetSeedType(1);
0120   dtb2->SetSeedJetPt(7);
0121   dtb2->Verbosity(verbosity);
0122   dtb2->set_towerinfo(true);
0123   dtb2->set_towerNodePrefix(HIJETS::tower_prefix);
0124   se->registerSubsystem(dtb2);
0125
0126   SubtractTowers *st = new SubtractTowers();
0127   st->SetFlowModulation(HIJETS::do_flow);
0128   st->Verbosity(verbosity);
0129   st->set_towerinfo(true);
0130   st->set_towerNodePrefix(HIJETS::tower_prefix);
0131   se->registerSubsystem(st);
0132
0133   towerjetreco = new JetReco();
0134   incemc = new TowerJetInput(Jet::CEMC_TOWERINFO_SUB1,HIJETS::tower_prefix);
0135   inihcal = new TowerJetInput(Jet::HCALIN_TOWERINFO_SUB1,HIJETS::tower_prefix);
0136   inohcal = new TowerJetInput(Jet::HCALOUT_TOWERINFO_SUB1,HIJETS::tower_prefix);
0137   if (HIJETS::do_vertex_type)
0138   {
0139     incemc->set_GlobalVertexType(HIJETS::vertex_type);
0140     inihcal->set_GlobalVertexType(HIJETS::vertex_type);
0141     inohcal->set_GlobalVertexType(HIJETS::vertex_type);
0142   }
0143   towerjetreco->add_input(incemc);
0144   towerjetreco->add_input(inihcal);
0145   towerjetreco->add_input(inohcal);
0146   // towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.2, verbosity), "AntiKt_Tower_r02_Sub1");
0147   // towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.3, verbosity), "AntiKt_Tower_r03_Sub1");
0148   towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.4, verbosity), "AntiKt_Tower_r04_Sub1");
0149   // towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.5, verbosity), "AntiKt_Tower_r05_Sub1");
0150   towerjetreco->set_algo_node("ANTIKT");
0151   towerjetreco->set_input_node("TOWER");
0152   towerjetreco->Verbosity(verbosity);
0153   se->registerSubsystem(towerjetreco);
0154
0155   return;
0156
0157 }
0158
0159
0160 // ----------------------------------------------------------------------------
0161 //! Determine rho from tower input to jet reco (necessary for jet QA)
0162 // ----------------------------------------------------------------------------
0163 void DoRhoCalculation()
0164 {
0165
0166   // set verbosity
0167   int verbosity = std::max(Enable::VERBOSITY, Enable::HIJETS_VERBOSITY);
0168
0169   //---------------
0170   // Fun4All server
0171   //---------------
0172   Fun4AllServer* se = Fun4AllServer::instance();
0173
0174   // run rho calculations w/ default parameters
0175   DetermineTowerRho* towRhoCalc = new DetermineTowerRho();
0176   towRhoCalc -> add_method(TowerRho::Method::AREA);
0177   towRhoCalc -> add_method(TowerRho::Method::MULT);
0178   towRhoCalc -> add_tower_input( new TowerJetInput(Jet::CEMC_TOWERINFO_RETOWER) );
0179   towRhoCalc -> add_tower_input( new TowerJetInput(Jet::HCALIN_TOWERINFO) );
0180   towRhoCalc -> add_tower_input( new TowerJetInput(Jet::HCALOUT_TOWERINFO) );
0181   se -> registerSubsystem( towRhoCalc );
0182
0183   // exit back to main macro
0184   return;
0185
0186 }
0187
0188 #endif