Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
plotTracklets.cxx
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file plotTracklets.cxx
1 #include <TCanvas.h>
2 #include <TCut.h>
3 #include <TFile.h>
4 #include <TH1F.h>
5 #include <TH2F.h>
6 #include <TMath.h>
7 #include <TObjString.h>
8 #include <TRandom3.h>
9 #include <TTree.h>
10 #include <TTreeIndex.h>
11 
12 #include <fstream>
13 #include <iostream>
14 #include <sstream>
15 #include <string>
16 #include <vector>
17 
18 using namespace std;
19 
20 // in python, do this: np.logspace(-4,np.log10(0.5),num=51)
21 const int NBins = 50;
22 float edges[NBins + 1] = {1.00000000e-04, 1.18571250e-04, 1.40591414e-04, 1.66700998e-04, 1.97659458e-04, 2.34367291e-04, 2.77892228e-04, 3.29500290e-04, 3.90692614e-04, 4.63249118e-04, 5.49280272e-04, 6.51288487e-04, 7.72240903e-04,
23  9.15655696e-04, 1.08570441e-03, 1.28733329e-03, 1.52640718e-03, 1.80988009e-03, 2.14599745e-03, 2.54453601e-03, 3.01708817e-03, 3.57739917e-03, 4.24176693e-03, 5.02951609e-03, 5.96356012e-03, 7.07106781e-03,
24  8.38425353e-03, 9.94131425e-03, 1.17875406e-02, 1.39766343e-02, 1.65722701e-02, 1.96499479e-02, 2.32991889e-02, 2.76261397e-02, 3.27566592e-02, 3.88399805e-02, 4.60530506e-02, 5.46056779e-02, 6.47466352e-02,
25  7.67708949e-02, 9.10282102e-02, 1.07933287e-01, 1.27977848e-01, 1.51744935e-01, 1.79925867e-01, 2.13340350e-01, 2.52960321e-01, 2.99938216e-01, 3.55640493e-01, 4.21687380e-01, 5.00000000e-01};
26 
27 void makehist(TString infname, TString outfname)
28 {
29  TFile *fout = new TFile(outfname, "RECREATE");
30 
31  TH1F *hM_clusphi = new TH1F("hM_clusphi", "hM_clusphi", 140, -3.5, 3.5);
32  TH1F *hM_cluseta = new TH1F("hM_cluseta", "hM_cluseta", 160, -4, 4);
33  TH1F *hM_clusphisize = new TH1F("hM_clusphisize", "hM_clusphisize", 80, 0, 80);
34  TH2F *hM_clusphi_clusphisize = new TH2F("hM_clusphi_clusphisize", "hM_clusphi_clusphisize", 140, -3.5, 3.5, 80, 0, 80);
35  TH2F *hM_cluseta_clusphisize = new TH2F("hM_cluseta_clusphisize", "hM_cluseta_clusphisize", 160, -4, 4, 80, 0, 80);
36 
37  TH1F *hM_dEta_proto = new TH1F("hM_dEta_proto", "hM_dEta_proto", 200, -5, 5);
38  TH1F *hM_dPhi_proto = new TH1F("hM_dPhi_proto", "hM_dPhi_proto", 100, -0.5, 0.5);
39  TH1F *hM_dPhi_proto_altrange = new TH1F("hM_dPhi_proto_altrange", "hM_dPhi_proto_altrange", 100, -0.05, 0.05);
40  TH1F *hM_dR_proto = new TH1F("hM_dR_proto", "hM_dR_proto", 100, 0, 0.5);
41  TH1F *hM_dR_proto_altrange = new TH1F("hM_dR_proto_altrange", "hM_dR_proto_altrange", 50, 0, 0.05);
42  TH1F *hM_dR_proto_LogX = new TH1F("hM_dR_proto_LogX", "hM_dR_proto_LogX", NBins, edges);
43  TH1F *hM_Eta_proto = new TH1F("hM_Eta_proto", "hM_Eta_proto", 160, -4, 4);
44  TH1F *hM_Phi_proto = new TH1F("hM_Phi_proto", "hM_Phi_proto", 140, -3.5, 3.5);
45  TH2F *hM_Eta_vtxZ_proto_incl = new TH2F("hM_Eta_vtxZ_proto_incl", "hM_Eta_vtxZ_proto_incl", 280, -3.5, 3.5, 300, -50, 10);
46 
47  TH1F *hM_dEta_reco = new TH1F("hM_dEta_reco", "hM_dEta_reco", 200, -3, 3);
48  TH1F *hM_dEta_reco_altrange = new TH1F("hM_dEta_reco_altrange", "hM_dEta_reco_altrange", 100, -0.5, 0.5);
49  TH1F *hM_dEta_reco_altrange_Nclusle4000 = new TH1F("hM_dEta_reco_altrange_Nclusle4000", "hM_dEta_reco_altrange_Nclusle4000", 100, -0.5, 0.5);
50  TH1F *hM_dEta_reco_altrange_Nclusgt4000 = new TH1F("hM_dEta_reco_altrange_Nclusgt4000", "hM_dEta_reco_altrange_Nclusgt4000", 100, -0.5, 0.5);
51  TH1F *hM_dPhi_reco = new TH1F("hM_dPhi_reco", "hM_dPhi_reco", 200, -0.5, 0.5);
52  TH1F *hM_dPhi_reco_altrange = new TH1F("hM_dPhi_reco_altrange", "hM_dPhi_reco_altrange", 200, -0.05, 0.05);
53  TH1F *hM_dPhi_reco_altrange_Nclusle4000 = new TH1F("hM_dPhi_reco_altrange_Nclusle4000", "hM_dPhi_reco_altrange_Nclusle4000", 200, -0.05, 0.05);
54  TH1F *hM_dPhi_reco_altrange_Nclusgt4000 = new TH1F("hM_dPhi_reco_altrange_Nclusgt4000", "hM_dPhi_reco_altrange_Nclusgt4000", 200, -0.05, 0.05);
55  TH1F *hM_dR_reco = new TH1F("hM_dR_reco", "hM_dR_reco", 100, 0, 0.5);
56  TH1F *hM_dR_reco_Nclusle4000 = new TH1F("hM_dR_reco_Nclusle4000", "hM_dR_reco_Nclusle4000", 100, 0, 0.5);
57  TH1F *hM_dR_reco_Nclusgt4000 = new TH1F("hM_dR_reco_Nclusgt4000", "hM_dR_reco_Nclusgt4000", 100, 0, 0.5);
58  TH1F *hM_dR_reco_altrange = new TH1F("hM_dR_reco_altrange", "hM_dR_reco_altrange", 50, 0, 0.05);
59  TH1F *hM_dR_reco_LogX = new TH1F("hM_dR_reco_LogX", "hM_dR_reco_LogX", NBins, edges);
60  TH1F *hM_Eta_reco = new TH1F("hM_Eta_reco", "hM_Eta_reco", 160, -4, 4);
61  TH1F *hM_Eta_reco_Nclusle4000 = new TH1F("hM_Eta_reco_Nclusle4000", "hM_Eta_reco_Nclusle4000", 160, -4, 4);
62  TH1F *hM_Eta_reco_Nclusgt4000 = new TH1F("hM_Eta_reco_Nclusgt4000", "hM_Eta_reco_Nclusgt4000", 160, -4, 4);
63  TH1F *hM_Phi_reco = new TH1F("hM_Phi_reco", "hM_Phi_reco", 140, -3.5, 3.5);
64  TH1F *hM_Phi_reco_Nclusle4000 = new TH1F("hM_Phi_reco_Nclusle4000", "hM_Phi_reco_Nclusle4000", 140, -3.5, 3.5);
65  TH1F *hM_Phi_reco_Nclusgt4000 = new TH1F("hM_Phi_reco_Nclusgt4000", "hM_Phi_reco_Nclusgt4000", 140, -3.5, 3.5);
66  TH2F *hM_Eta_vtxZ_reco_incl = new TH2F("hM_Eta_vtxZ_reco_incl", "hM_Eta_vtxZ_reco_incl", 280, -3.5, 3.5, 300, -50, 10);
67  TH2F *hM_Eta_vtxZ_reco_Nclusle4000 = new TH2F("hM_Eta_vtxZ_reco_Nclusle4000", "hM_Eta_vtxZ_reco_Nclusle4000", 280, -3.5, 3.5, 300, -50, 10);
68  TH2F *hM_Eta_vtxZ_reco_Nclusgt4000 = new TH2F("hM_Eta_vtxZ_reco_Nclusgt4000", "hM_Eta_vtxZ_reco_Nclusgt4000", 280, -3.5, 3.5, 300, -50, 10);
69 
70  TH1F *hM_NClusLayer1 = new TH1F("hM_NClusLayer1", "hM_NClusLayer1", 100, 0, 5000);
71  TH1F *hM_NPrototkl = new TH1F("hM_NPrototkl", "hM_NPrototkl", 100, 0, 10000);
72  TH1F *hM_NRecotkl_Raw = new TH1F("hM_NRecotkl_Raw", "hM_NRecotkl_Raw", 100, 0, 5000);
73 
74  TH1F *hM_RecoPVz = new TH1F("hM_RecoPVz", "hM_RecoPVz", 300, -50, 10);
75  TH1F *hM_RecoPVz_Nclusle4000 = new TH1F("hM_RecoPVz_Nclusle4000", "hM_RecoPVz_Nclusle4000", 300, -50, 10);
76  TH1F *hM_RecoPVz_Nclusgt4000 = new TH1F("hM_RecoPVz_Nclusgt4000", "hM_RecoPVz_Nclusgt4000", 300, -50, 10);
77 
78  TFile *f = new TFile(infname, "READ");
79  TTree *t = (TTree *)f->Get("minitree");
80  t->BuildIndex("event"); // Reference: https://root-forum.cern.ch/t/sort-ttree-entries/13138
81  TTreeIndex *index = (TTreeIndex *)t->GetTreeIndex();
82  int event, NClusLayer1, NPrototkl, NRecotkl_Raw, NRecotkl_GenMatched, NGenHadron;
83  float PV_z, TruthPV_z;
84  vector<float> *clusphi = 0, *cluseta = 0, *clusphisize = 0;
85  vector<float> *prototkl_eta = 0, *prototkl_phi = 0, *prototkl_deta = 0, *prototkl_dphi = 0, *prototkl_dR = 0;
86  vector<float> *recotklraw_eta = 0, *recotklraw_phi = 0, *recotklraw_deta = 0, *recotklraw_dphi = 0, *recotklraw_dR = 0;
87  t->SetBranchAddress("event", &event);
88  t->SetBranchAddress("NClusLayer1", &NClusLayer1);
89  t->SetBranchAddress("NPrototkl", &NPrototkl);
90  t->SetBranchAddress("NRecotkl_Raw", &NRecotkl_Raw);
91  t->SetBranchAddress("PV_z", &PV_z);
92  t->SetBranchAddress("TruthPV_z", &TruthPV_z);
93  t->SetBranchAddress("clusphi", &clusphi);
94  t->SetBranchAddress("cluseta", &cluseta);
95  t->SetBranchAddress("clusphisize", &clusphisize);
96  t->SetBranchAddress("prototkl_eta", &prototkl_eta);
97  t->SetBranchAddress("prototkl_phi", &prototkl_phi);
98  t->SetBranchAddress("prototkl_deta", &prototkl_deta);
99  t->SetBranchAddress("prototkl_dphi", &prototkl_dphi);
100  t->SetBranchAddress("prototkl_dR", &prototkl_dR);
101  t->SetBranchAddress("recotklraw_eta", &recotklraw_eta);
102  t->SetBranchAddress("recotklraw_phi", &recotklraw_phi);
103  t->SetBranchAddress("recotklraw_deta", &recotklraw_deta);
104  t->SetBranchAddress("recotklraw_dphi", &recotklraw_dphi);
105  t->SetBranchAddress("recotklraw_dR", &recotklraw_dR);
106 
107  for (int ev = 0; ev < t->GetEntriesFast(); ev++)
108  {
109  Long64_t local = t->LoadTree(index->GetIndex()[ev]);
110  t->GetEntry(local);
111 
112  if (PV_z < -99.)
113  continue;
114 
115  cout << "Event=" << event << "; NClusLayer1=" << NClusLayer1 << "; NPrototkl=" << NPrototkl << "; NRecotkl_Raw=" << NRecotkl_Raw << endl;
116 
117  hM_NClusLayer1->Fill(NClusLayer1);
118  hM_NPrototkl->Fill(NPrototkl);
119  hM_NRecotkl_Raw->Fill(NRecotkl_Raw);
120 
121  hM_RecoPVz->Fill(PV_z);
122  if (NClusLayer1 <= 4000)
123  hM_RecoPVz_Nclusle4000->Fill(PV_z);
124  else
125  hM_RecoPVz_Nclusgt4000->Fill(PV_z);
126 
127  for (size_t j = 0; j < clusphisize->size(); j++)
128  {
129  hM_clusphi->Fill(clusphi->at(j));
130  hM_cluseta->Fill(cluseta->at(j));
131  hM_clusphisize->Fill(clusphisize->at(j));
132  hM_clusphi_clusphisize->Fill(clusphi->at(j), clusphisize->at(j));
133  hM_cluseta_clusphisize->Fill(cluseta->at(j), clusphisize->at(j));
134  }
135 
136  for (size_t j = 0; j < prototkl_eta->size(); j++)
137  {
138  hM_dEta_proto->Fill(prototkl_deta->at(j));
139  hM_dPhi_proto->Fill(prototkl_dphi->at(j));
140  hM_dPhi_proto_altrange->Fill(prototkl_dphi->at(j));
141  hM_dR_proto->Fill(prototkl_dR->at(j));
142  hM_dR_proto_altrange->Fill(prototkl_dR->at(j));
143  hM_dR_proto_LogX->Fill(prototkl_dR->at(j));
144  hM_Eta_proto->Fill(prototkl_eta->at(j));
145  hM_Phi_proto->Fill(prototkl_phi->at(j));
146 
147  hM_Eta_vtxZ_proto_incl->Fill(prototkl_eta->at(j), PV_z);
148  }
149 
150  for (size_t j = 0; j < recotklraw_eta->size(); j++)
151  {
152  hM_dEta_reco->Fill(recotklraw_deta->at(j));
153  hM_dEta_reco_altrange->Fill(recotklraw_deta->at(j));
154  hM_dPhi_reco->Fill(recotklraw_dphi->at(j));
155  hM_dPhi_reco_altrange->Fill(recotklraw_dphi->at(j));
156  hM_dR_reco->Fill(recotklraw_dR->at(j));
157  hM_dR_reco_altrange->Fill(recotklraw_dR->at(j));
158  hM_dR_reco_LogX->Fill(recotklraw_dR->at(j));
159  hM_Eta_reco->Fill(recotklraw_eta->at(j));
160  hM_Phi_reco->Fill(recotklraw_phi->at(j));
161  hM_Eta_vtxZ_reco_incl->Fill(recotklraw_eta->at(j), PV_z);
162 
163  if (NClusLayer1 <= 4000)
164  {
165  hM_dEta_reco_altrange_Nclusle4000->Fill(recotklraw_deta->at(j));
166  hM_dPhi_reco_altrange_Nclusle4000->Fill(recotklraw_dphi->at(j));
167  hM_dR_reco_Nclusle4000->Fill(recotklraw_dR->at(j));
168  hM_Eta_reco_Nclusle4000->Fill(recotklraw_eta->at(j));
169  hM_Phi_reco_Nclusle4000->Fill(recotklraw_phi->at(j));
170  hM_Eta_vtxZ_reco_Nclusle4000->Fill(recotklraw_eta->at(j), PV_z);
171  }
172  else
173  {
174  hM_dEta_reco_altrange_Nclusgt4000->Fill(recotklraw_deta->at(j));
175  hM_dPhi_reco_altrange_Nclusgt4000->Fill(recotklraw_dphi->at(j));
176  hM_dR_reco_Nclusgt4000->Fill(recotklraw_dR->at(j));
177  hM_Eta_reco_Nclusgt4000->Fill(recotklraw_eta->at(j));
178  hM_Phi_reco_Nclusgt4000->Fill(recotklraw_phi->at(j));
179  hM_Eta_vtxZ_reco_Nclusgt4000->Fill(recotklraw_eta->at(j), PV_z);
180  }
181  }
182  }
183 
184  f->Close();
185 
186  fout->cd();
187  hM_NClusLayer1->Write();
188  hM_NPrototkl->Write();
189  hM_NRecotkl_Raw->Write();
190  hM_RecoPVz->Write();
191  hM_RecoPVz_Nclusle4000->Write();
192  hM_RecoPVz_Nclusgt4000->Write();
193  hM_clusphi->Write();
194  hM_cluseta->Write();
195  hM_clusphisize->Write();
196  hM_clusphi_clusphisize->Write();
197  hM_cluseta_clusphisize->Write();
198  hM_dEta_proto->Write();
199  hM_dEta_reco->Write();
200  hM_dEta_proto->Write();
201  hM_dPhi_proto->Write();
202  hM_dPhi_proto_altrange->Write();
203  hM_dR_proto->Write();
204  hM_dR_proto_altrange->Write();
205  hM_dR_proto_LogX->Write();
206  hM_Eta_proto->Write();
207  hM_Phi_proto->Write();
208  hM_dEta_reco->Write();
209  hM_dEta_reco_altrange->Write();
210  hM_dEta_reco_altrange_Nclusle4000->Write();
211  hM_dEta_reco_altrange_Nclusgt4000->Write();
212  hM_dPhi_reco->Write();
213  hM_dPhi_reco_altrange->Write();
214  hM_dPhi_reco_altrange_Nclusle4000->Write();
215  hM_dPhi_reco_altrange_Nclusgt4000->Write();
216  hM_dR_reco->Write();
217  hM_dR_reco_Nclusle4000->Write();
218  hM_dR_reco_Nclusgt4000->Write();
219  hM_dR_reco_altrange->Write();
220  hM_dR_reco_LogX->Write();
221  hM_Eta_reco->Write();
222  hM_Eta_reco_Nclusle4000->Write();
223  hM_Eta_reco_Nclusgt4000->Write();
224  hM_Phi_reco->Write();
225  hM_Phi_reco_Nclusle4000->Write();
226  hM_Phi_reco_Nclusgt4000->Write();
227  hM_Eta_vtxZ_proto_incl->Write();
228  hM_Eta_vtxZ_reco_incl->Write();
229  hM_Eta_vtxZ_reco_Nclusle4000->Write();
230  hM_Eta_vtxZ_reco_Nclusgt4000->Write();
231 
232  fout->Close();
233 }
234 
235 int main(int argc, char *argv[])
236 {
237  TString infname;
238  TString outfname;
239 
240  if (argc == 1)
241  {
242  infname = "/sphenix/user/hjheng/TrackletAna/minitree/INTT/TrackletMinitree_ana382_zvtx-20cm_dummyAlignParams/TrackletAna_minitree_Evt0to2000_dRcut0p5.root";
243  outfname = "/sphenix/user/hjheng/TrackletAna/analysis_INTT/plot/hists/ana382_zvtx-20cm_dummyAlignParams/Hists_RecoTracklets.root";
244  }
245  else if (argc == 3)
246  {
247  infname = argv[1];
248  outfname = argv[2];
249  }
250  else
251  {
252  std::cout << "Usage: ./plotTracklets [infile] [outfile]" << std::endl;
253  return 0;
254  }
255 
256  cout << "[Run Info] Input file = " << infname << endl << " Output file = " << outfname << endl << "-----------" << endl;
257 
258  makehist(infname, outfname);
259 
260  return 0;
261 }