src_JetMedianTree/macro/Fun4All_JetRhoMedian.C

0001 #ifndef MACRO_FUN4ALLJETANA_C
0002 #define MACRO_FUN4ALLJETANA_C
0003
0004 #include <fun4all/Fun4AllDstInputManager.h>
0005 #include <fun4all/Fun4AllInputManager.h>
0006 #include <fun4all/Fun4AllServer.h>
0007 #include <fun4all/SubsysReco.h>
0008 #include <g4centrality/PHG4CentralityReco.h>
0009 #include <jetbase/FastJetAlgo.h>
0010 #include <jetbase/JetReco.h>
0011 #include <jetbase/TowerJetInput.h>
0012 #include <g4jets/TruthJetInput.h>
0013 /* #include <g4vertex/GlobalVertexReco.h> */
0014 #include <jetbackground/FastJetAlgoSub.h>
0015 #include <fun4all/PHTFileServer.h>
0016 #include <jetbackground/RetowerCEMC.h>
0017 #include <HIJetReco.C>
0018
0019 // here you need your package name (set in configure.ac)
0020 #include <fastjetmedianbkg/JetRhoMedian.h>
0021 #include <fastjetmedianbkg/PrintTowers.h>
0022 R__LOAD_LIBRARY(libfun4all.so)
0023 R__LOAD_LIBRARY(libg4centrality.so)
0024 R__LOAD_LIBRARY(libg4jets.so)
0025 R__LOAD_LIBRARY(libjetbase.so)
0026 R__LOAD_LIBRARY(libjetbackground.so)
0027 R__LOAD_LIBRARY(libjetrhomedian.so)
0028 R__LOAD_LIBRARY(libjetbackground.so)
0029 /* R__LOAD_LIBRARY(libg4vertex.so) */
0030
0031 void Fun4All_JetRhoMedian(
0032     const int     nevnt = 10
0033   , const float   jet_R = 0.4
0034   , const string& out_name = "test.root"
0035   , const float   min_lead_truth_pt = 30
0036   , const string& list_calo_cluster=""
0037   , const string& list_truth_jet=""
0038   , const string& list_bbc=""
0039   , const string& list_global=""
0040   , const string& list_truth_g4hit=""
0041   , const int     verbosity=0
0042   )
0043 {
0044   gSystem->Load("libjetrhomedian");
0045   gSystem->Load("libg4dst");
0046
0047   Fun4AllServer *se = Fun4AllServer::instance();
0048
0049   if (list_bbc!="") {
0050     PHG4CentralityReco *cent = new PHG4CentralityReco();
0051     cent->Verbosity(0);
0052     cent->GetCalibrationParameters().ReadFromFile("centrality", "xml", 0, 0, string(getenv("CALIBRATIONROOT")) + string("/Centrality/"));
0053     se->registerSubsystem( cent );
0054
0055     Fun4AllInputManager *inp_bbc = new Fun4AllDstInputManager("BBC_info");
0056     inp_bbc->AddListFile(list_bbc);
0057     se->registerInputManager(inp_bbc);
0058   } else {
0059     cout << " Fatal: missing bbc input file. Exiting." << endl;
0060     gSystem->Exit(0);
0061   }
0062
0063   // add the CEMC retowering
0064   RetowerCEMC *rcemc = new RetowerCEMC();
0065   /* rcemc->Verbosity(verbosity); */
0066   /* rcemc->set_towerinfo(true); */
0067   se->registerSubsystem(rcemc);
0068   //
0069   // do the background subtraction
0070   Enable::HIJETS_TRUTH=false;
0071   HIJetReco();
0072
0073   std::string sub1_jet_name = Form("AntiKt_Tower_r0%i_Sub1", ((int)(jet_R*10)));
0074   std::string truth_jet_name = Form("AntiKt_Truth_r0%i", ((int)(jet_R*10)));
0075
0076   JetRhoMedian *jetRhoMedian = new JetRhoMedian(out_name, jet_R, truth_jet_name, sub1_jet_name, min_lead_truth_pt);
0077   jetRhoMedian->add_input(new TowerJetInput(Jet::CEMC_TOWER_RETOWER));
0078   jetRhoMedian->add_input(new TowerJetInput(Jet::HCALIN_TOWER));
0079   jetRhoMedian->add_input(new TowerJetInput(Jet::HCALOUT_TOWER));
0080   jetRhoMedian->Verbosity(0);
0081   se->registerSubsystem(jetRhoMedian);
0082
0083
0084   if (list_truth_jet!="") {
0085     Fun4AllInputManager *inp_truth_jet = new Fun4AllDstInputManager("DSTtruth");
0086     inp_truth_jet->AddListFile(list_truth_jet,1); // adding the option "1" confirms to use, even if file is large
0087     se->registerInputManager(inp_truth_jet);
0088   } else {
0089     cout << " Fatal: missing truth jet list file. Exiting." << endl;
0090     gSystem->Exit(0);
0091   }
0092
0093   if (list_calo_cluster!="") {
0094     Fun4AllInputManager *inp_calo_cluster = new Fun4AllDstInputManager("DSTcalocluster");
0095     inp_calo_cluster->AddListFile(list_calo_cluster,1);
0096     se->registerInputManager(inp_calo_cluster);
0097   } else {
0098     cout << " Fatal: missing calo cluster list file. Exiting." << endl;
0099     gSystem->Exit(0);
0100   }
0101
0102
0103   if (list_global!="") {
0104     Fun4AllInputManager *inp_global = new Fun4AllDstInputManager("DSTglobal");
0105     inp_global->AddListFile(list_global,1); // adding the option "1" confirms to use, even if file is large
0106     se->registerInputManager(inp_global);
0107   } else {
0108     cout << " Fatal: missing global dst list file. Exiting." << endl;
0109     gSystem->Exit(0);
0110   }
0111
0112   if (list_truth_g4hit!="") {
0113     Fun4AllInputManager *inp_truth_g4hit = new Fun4AllDstInputManager("DST_truth_g4hit");
0114     inp_truth_g4hit->AddListFile(list_truth_g4hit,1); // adding the option "1" confirms to use, even if file is large
0115     se->registerInputManager(inp_truth_g4hit);
0116   } else {
0117     cout << " Fatal: missing truth g4hit list file. Exiting." << endl;
0118     gSystem->Exit(0);
0119   }
0120
0121
0122   se->run(nevnt);
0123   se->End();
0124   delete se;
0125   cout << " Done in Fun4All_JetRhoMedian.C " << endl;
0126   gSystem->Exit(0);
0127
0128 }
0129
0130 #endif