Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4_HIJetReco.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4_HIJetReco.C
1 #ifndef MACRO_G4HIJETRECO_C
2 #define MACRO_G4HIJETRECO_C
3 
4 #include <GlobalVariables.C>
5 
6 #include <g4jets/FastJetAlgo.h>
7 #include <g4jets/JetReco.h>
8 #include <g4jets/TowerJetInput.h>
9 #include <g4jets/TruthJetInput.h>
10 #include <g4centrality/PHG4CentralityReco.h>
11 
12 #include <jetbackground/CopyAndSubtractJets.h>
13 #include <jetbackground/DetermineTowerBackground.h>
14 #include <jetbackground/FastJetAlgoSub.h>
15 #include <jetbackground/RetowerCEMC.h>
16 #include <jetbackground/SubtractTowers.h>
17 #include <jetbackground/SubtractTowersCS.h>
18 
19 #include <fun4all/Fun4AllServer.h>
20 #include <fun4all/SubsysReco.h>
21 #include <fun4all/Fun4AllInputManager.h>
22 #include <fun4all/Fun4AllDstInputManager.h>
23 
24 #include <fun4all/Fun4AllDstOutputManager.h>
25 #include <fun4all/Fun4AllOutputManager.h>
26 
27 #include <phool/PHRandomSeed.h>
28 #include <phool/recoConsts.h>
29 
30 #include <myjetanalysis/MyJetAnalysis.h>
31 #include <negativetowerremover/NegativeTowerRemover.h>
32 
33 #include <fastjet/PseudoJet.hh>
34 
35 R__LOAD_LIBRARY(libfun4all.so)
36 R__LOAD_LIBRARY(libg4jets.so)
37 R__LOAD_LIBRARY(libjetbackground.so)
38 R__LOAD_LIBRARY(libmyjetanalysis.so)
39 R__LOAD_LIBRARY(libg4centrality.so)
40 
41 namespace Enable
42 {
43  bool HIJETS = false;
44  int HIJETS_VERBOSITY = 0;
45 } // namespace Enable
46 
47 namespace G4HIJETS
48 {
49  bool do_flow = false;
50  bool do_CS = false;
51 } // namespace G4HIJETS
52 
53 void G4_HIJetReco(const char *filelisttruth = "dst_truth.list",
54  const char *filelisttruthjet = "dst_truth_jet.list",
55  const char *filelisthits = "dst_bbc_g4hit.list",
56  const char *filelistcalo = "dst_calo_cluster.list",
57  const string &outname = "output_000002",
58  int n_skip = 2,
59  int n_event = 5
60  )
61  {
62  gSystem->Load("libg4dst");
63  gSystem->Load("libmyjetanalysis");
64 
65  int verbosity = 0;//std::max(Enable::VERBOSITY, Enable::HIJETS_VERBOSITY);
66 
67  //---------------
68  // Fun4All server
69  //--------------1-
70 
71  Fun4AllServer *se = Fun4AllServer::instance();
72  se->Verbosity(verbosity);
73 
74  PHG4CentralityReco *cent = new PHG4CentralityReco();
75  cent->Verbosity(0);
76  cent->GetCalibrationParameters().ReadFromFile("centrality", "xml", 0, 0, string(getenv("CALIBRATIONROOT")) + string("/Centrality/"));
77  se->registerSubsystem( cent );
78 
79 //JetRecolysis::InitRun" << end`l;
80 // //m_jetEvalStack = shared_ptr<JetEvalStack>(new JetEvalStack(topNode, m_recoJetName, m_truthJetName));
81 // //m_jetEvalStack->get_stvx_eval_stack()->set_use_initial_vertex(initial_vertex);
82 // return Fun4AllReturnCodes::EVENT_OK;
83 // *truthjetreco = new JetReco();
84  //TruthJetInput *tji = new TruthJetInput(Jet::PARTICLE);
85  //tji->add_embedding_flag(1); // changes depending on signal vs. embedded
86  //truthjetreco->add_input(tji);
87  //truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.2), "AntiKt_Truth_r02_sub");
88  //truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.3), "AntiKt_Truth_r03_sub");
89  //truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.4), "AntiKt_Truth_r04_sub");
90  //truthjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.5), "AntiKt_Truth_r05_sub");
92  //truthjetreco->set_input_node("TRUTH");
93  //truthjetreco->Verbosity(verbosity);
94  //se->registerSubsystem(truthjetreco);
95 
96  RetowerCEMC *rcemc = new RetowerCEMC();
97  rcemc->Verbosity(verbosity);
98  se->registerSubsystem(rcemc);
99 
100  JetReco *towerjetreco = new JetReco();
101  //towerjetreco->add_input(new TowerJetInput(Jet::CEMC_TOWER));
102  towerjetreco->add_input(new TowerJetInput(Jet::CEMC_TOWER_RETOWER));
103  towerjetreco->add_input(new TowerJetInput(Jet::HCALIN_TOWER));
104  towerjetreco->add_input(new TowerJetInput(Jet::HCALOUT_TOWER));
105  towerjetreco->add_algo(new FastJetAlgo(Jet::ANTIKT, 0.2), "AntiKt_Tower_HIRecoSeedsRaw_r02");
106  towerjetreco->set_algo_node("ANTIKT");
107  towerjetreco->set_input_node("TOWER");
108  towerjetreco->Verbosity(verbosity);
109  se->registerSubsystem(towerjetreco);
110 
111  DetermineTowerBackground *dtb = new DetermineTowerBackground();
112  dtb->SetBackgroundOutputName("TowerBackground_Sub1");
113  dtb->SetFlow(G4HIJETS::do_flow);
114  dtb->SetSeedType(0);
115  dtb->SetSeedJetD(3);
116  dtb->Verbosity(verbosity);
117  se->registerSubsystem(dtb);
118 
119  CopyAndSubtractJets *casj = new CopyAndSubtractJets();
120  casj->SetFlowModulation(G4HIJETS::do_flow);
121  //casj->SetBackgroundInputName("TowerBackground_Sub1");
122  //casj->SetJetInputName("AntiKt_Tower_HIRecoSeedsRaw_r02");
123  //casj->SetJetOutputName("AntiKt_Tower_HIRecoSeedsSub_r02");
124  casj->Verbosity(verbosity);
125  se->registerSubsystem(casj);
126 
127  DetermineTowerBackground *dtb2 = new DetermineTowerBackground();
128  dtb2->SetBackgroundOutputName("TowerBackground_Sub2");
129  dtb2->SetFlow(G4HIJETS::do_flow);
130  dtb2->SetSeedType(1);
131  dtb2->SetSeedJetPt(7);
132  dtb2->Verbosity(verbosity);
133  se->registerSubsystem(dtb2);
134 
135  //CopyAndSubtractJets *casj_final = new CopyAndSubtractJets();
136  //casj_final->SetFlowModulation(G4HIJETS::do_flow);
137  //casj_final->SetBackgroundInputName("TowerBackground_Sub2");
138  //casj_final->SetJetInputName("AntiKt_Tower_HIRecoSeedsRaw_r02");
139  //casj_final->SetJetOutputName("AntiKt_Tower_r02_Sub1_V2");
140  //casj_final->Verbosity(verbosity);
141  //se->registerSubsystem(casj_final);
142 
143  //cout << "Started Subtracting Towers" << endl;
144 
145  SubtractTowers *st = new SubtractTowers();
146  st->SetFlowModulation(G4HIJETS::do_flow);
147  st->Verbosity(verbosity);
148  se->registerSubsystem(st);
149  /*
150  NegativeTowerRemover *ntr = new NegativeTowerRemover();
151  ntr->SetFlowModulation(G4HIJETS::do_flow);
152  ntr->Verbosity(verbosity);
153  se->registerSubsystem(ntr);
154  */
155 
156  //JetReco *towerjetreco = new JetReco();
157  towerjetreco = new JetReco();
158  towerjetreco->add_input(new TowerJetInput(Jet::CEMC_TOWER_SUB1));
159  towerjetreco->add_input(new TowerJetInput(Jet::HCALIN_TOWER_SUB1));
160  towerjetreco->add_input(new TowerJetInput(Jet::HCALOUT_TOWER_SUB1));
161  towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.2, verbosity), "AntiKt_Tower_r02_Sub1");
162  //towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.3, verbosity), "AntiKt_Tower_r03_Sub1");
163  //towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.4, verbosity), "AntiKt_Tower_r04_Sub1");
164  //towerjetreco->add_algo(new FastJetAlgoSub(Jet::ANTIKT, 0.5, verbosity), "AntiKt_Tower_r05_Sub1");
165 
166  towerjetreco->set_algo_node("ANTIKT");
167  towerjetreco->set_input_node("TOWER");
168  towerjetreco->Verbosity(verbosity);
169  se->registerSubsystem(towerjetreco);
170 
171  MyJetAnalysis *myJetAnalysis = new MyJetAnalysis("AntiKt_Tower_r02_Sub1", "AntiKt_Truth_r02", "myjet_ppPileup_R02_" + outname + ".root");
172  myJetAnalysis->setPtRange(15, 200);
173  myJetAnalysis->setEtaRange(-1.1, 1.1);
174  myJetAnalysis->setMindR(0.2);
175  //myJetAnalysis->doECCprint(false);
176  se->registerSubsystem(myJetAnalysis);
177 
178  //cout << "Finished Constructing the Jets" << endl;
179 
180  //MyJetAnalysis *myJetAnalysis_R04 = new MyJetAnalysis("AntiKt_Tower_r04", "AntiKt_Truth_r04", "output/myjet_ppPileup_R04_" + outname + ".root");
181  //myJetAnalysis_R04->setPtRange(15, 200);
182  //myJetAnalysis_R04->setEtaRange(-1.1, 1.1);
183  //myJetAnalysis_R04->setMindR(0.2);
184  //myJetAnalysis_R04->doECCprint(false);
185  //se->registerSubsystem(myJetAnalysis_R04);
186 
187  string inputFile = "DST_TRUTH_JET_pythia8_Jet30_sHijing_0_20fm_50kHz_bkg_0_20fm-0000000062-00000.root";
188  string inputFile0 = "DST_CALO_CLUSTER_pythia8_Jet30_sHijing_0_20fm_50kHz_bkg_0_20fm-0000000062-00000.root";
189  string inputFile1 = "DST_TRUTH_pythia8_Jet30_sHijing_0_20fm_50kHz_bkg_0_20fm-0000000062-00000.root";
190  string inputFile2 = "DST_BBC_G4HIT_pythia8_Jet30_sHijing_0_20fm_50kHz_bkg_0_20fm-0000000062-00000.root";
191 
192  //cout << "Got Strings" << endl;
193 
194  Fun4AllInputManager *intrue2 = new Fun4AllDstInputManager("DSTtruth2");
195  intrue2->AddListFile(filelisttruthjet,1);
196  se->registerInputManager(intrue2);
197 
198  Fun4AllInputManager *intrue = new Fun4AllDstInputManager("DSTtruth");
199  intrue->AddListFile(filelisttruth,1);
200  se->registerInputManager(intrue);
201 
202  Fun4AllInputManager *in = new Fun4AllDstInputManager("DSTcalo");
203  in->AddListFile(filelistcalo,1);
204  se->registerInputManager(in);
205 
206  Fun4AllInputManager *in2 = new Fun4AllDstInputManager("DSThit");
207  in2->AddListFile(filelisthits,1);
208  se->registerInputManager(in2);
209 
210  //se->skip(134);
211 
212  //se->run(-1);
213  //cout << "Ran" << endl;
214  se->run(n_event);
215  se->skip(n_skip);
216 
217  se->End();
218 
219  gSystem->Exit(0);
220  return 0;
221 }
222 #endif