Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
plot_ntp_track_out.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file plot_ntp_track_out.C
1 #include "TH2D.h"
2 #include "TH1D.h"
3 #include "TFile.h"
4 #include "TCanvas.h"
5 #include "TROOT.h"
6 #include "TStyle.h"
7 #include "TLine.h"
8 #include "TLegend.h"
9 #include "TF1.h"
10 
11 #include <iostream>
12 
13 void plot_ntp_track_out()
14 {
15  gROOT->SetStyle("Plain");
16  gStyle->SetOptStat(0);
17  gStyle->SetOptFit(1);
18  gStyle->SetOptTitle(0);
19 
20  double ptmax = 40.0;
21  double slice_low = 30.0;
22  double slice_high = 31.0;
23 
24  std::vector<std::string> finvec;
25  std::vector<std::string> legvec;
26  std::vector<int> col;
27 
28  //============
29 
30 
31  finvec.push_back("root_files/ntp_track_out_acts_MB_50khz_nomapscut.root");
32  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (MB+50, nomapscut)");
33  col.push_back(kRed);
34 
35  finvec.push_back("root_files/ntp_track_out_acts_MB_50khz_withmapscut.root");
36  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (MB+50, nmaps>1)");
37  col.push_back(kBlue);
38 
39 
40 
41  /*
42  finvec.push_back("root_files/ntp_track_out_lo_hough_genfitprop_genfit.root");
43  legvec.push_back( "Hough+GenfitTrkProp+GenfitTrkFitter (LO)");
44  col.push_back(kBlue);
45  */
46 
47 
48 
49  /*
50  finvec.push_back("root_files/ntp_track_out_stub_matcher_lo_4kevts.root");
51  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (lo)");
52  col.push_back(kRed);
53 
54  finvec.push_back("root_files/ntp_track_out_truth_assoc_lo_4kevts.root");
55  legvec.push_back( "TpcTracker+truthassoc+ActsTrkFit (lo)");
56  col.push_back(kBlack);
57  */
58 
59  /*
60  // repeat after merge with repo of finvec.push_back("root_files/ntp_track_out_4k_tune8_fixed_mag24_both_stub_matcher.root");
61  finvec.push_back("root_files/ntp_track_out.root");
62  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (0.02,0.015) (4K)");
63  col.push_back(kBlack);
64  */
65 
66  /*
67  finvec.push_back("root_files/ntp_track_out_4k_tune8_fixed_mag24_both_stub_matcher.root");
68  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (0.02,0.015) (4K)");
69  col.push_back(kBlack);
70  */
71 
72  /*
73  finvec.push_back("root_files/ntp_track_out_4k_tune7_fixed_mag24_both_stub_matcher.root");
74  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (0.02,0.02) (4K)");
75  col.push_back(kMagenta);
76  */
77 
78  /*
79  finvec.push_back("root_files/ntp_track_out_4k_tune1_fixed_mag24_both_stub_matcher.root");
80  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (0.03,0.02) (4K)");
81  col.push_back(kBlack);
82  */
83 
84  /*
85  finvec.push_back("root_files/ntp_track_out_4k_truth_assoc.root");
86  legvec.push_back( "TPCTracker+TruthSiliconAssoc+ActsTrkFit (4K)");
87  col.push_back(kRed);
88  */
89 
90  /*
91  finvec.push_back("root_files/ntp_track_out_4k_tune5_fixed_mag24_both_stub_matcher.root");
92  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (0.06,0.03) (4K)");
93  col.push_back(kBlack);
94  */
95 
96  /*
97  finvec.push_back("root_files/ntp_track_out_4k_hough_genfitprop_genfit.root");
98  legvec.push_back( "Hough+GenfitTrkProp+Genfit (4K)");
99  col.push_back(kBlue);
100  */
101 
102  /*
103  finvec.push_back("root_files/ntp_track_out_4k_tune2_fixed_mag24_both_stub_matcher.root");
104  legvec.push_back( "TpcTracker+stub matcher+ActsTrkFit (0.02,0.01) (4K)");
105  col.push_back(kMagenta);
106  */
107 
108  //==========
109 
110 
111  unsigned int nfiles = finvec.size();
112 
113  TCanvas *ctemp0 = new TCanvas("ctemp0","ctemp0",5,5,1600,800);
114  ctemp0->Divide(nfiles,1);
115  TCanvas *cslice = new TCanvas("cslice","cslice",5,5,1600,800);
116  cslice->Divide(nfiles,1);
117  TCanvas *cpt = new TCanvas("cpt","cpt",5,5,1500,800);
118  cpt->Divide(2,1);
119  cpt->cd(1); gPad->SetLeftMargin(0.2);
120  cpt->cd(2); gPad->SetLeftMargin(0.2);
121  TCanvas *ceff = new TCanvas("ceff","ceff",5,5,1200,800);
122  ceff->SetLeftMargin(0.18);
123  TCanvas *ctemp1 = new TCanvas("ctemp1","ctemp1",5,5,1200,800);
124  ctemp1->Divide(nfiles,1);
125  TCanvas *cslicexy = new TCanvas("cslicexy","cslicexy",5,5,1600,800);
126  cslicexy->Divide(nfiles,1);
127  TCanvas *cdcaxy = new TCanvas("cdcaxy","cdcaxy",5,5,1200,800);
128  cdcaxy->SetLeftMargin(0.2);
129  gPad->SetGrid();
130  TCanvas *ctemp2 = new TCanvas("ctemp2","ctemp2",5,5,1200,800);
131  ctemp2->Divide(nfiles,1);
132  TCanvas *csliceZ = new TCanvas("csliceZ","csliceZ",5,5,1600,800);
133  csliceZ->Divide(nfiles,1);
134  TCanvas *cdcaZ = new TCanvas("cdcaZ","cdcaZ",5,5,1200,800);
135  cdcaZ->SetLeftMargin(0.2);
136  gPad->SetGrid();
137  TCanvas *cmvtx= new TCanvas("cmvtx","cmvtx",5,5,1200,800);
138  cmvtx->SetLeftMargin(0.2);
139  TCanvas *ccomb= new TCanvas("ccomb","ccomb",5,5,1200,800);
140  ccomb->SetLeftMargin(0.2);
141 
142  TLegend *lpd = new TLegend(0.25, 0.6, 0.85, 0.75, "", "NDC");
143  lpd->SetBorderSize(1);
144  lpd->SetFillColor(kWhite);
145  lpd->SetFillStyle(1001);
146 
147  TF1 *fz = new TF1("fz","gaus");
148  fz->SetRange(-0.02, 0.02);
149 
150  for(unsigned int i=0; i<finvec.size();++i)
151  {
152  TFile *fin = new TFile(finvec[i].c_str());
153 
154  // pT resolution
155 
156  TH2D *hpt2d;
157  fin->GetObject("h1",hpt2d);
158  std::cout << hpt2d->GetEntries() << std::endl;
159  hpt2d->FitSlicesY();
160  TH1D *hptres = (TH1D*)gDirectory->Get("h1_2");
161  TH1D *hptcent = (TH1D*)gDirectory->Get("h1_1");
162  ctemp0->cd(i+1);
163  hpt2d->Draw("colz");
164 
165  int binlow = hpt2d->GetXaxis()->FindBin(slice_low);
166  int binhigh = hpt2d->GetXaxis()->FindBin(slice_high);
167  char hname[500];
168  sprintf(hname,"hptslice%i",i);
169  TH1D *hptslice = hpt2d->ProjectionY(hname,binlow, binhigh);
170  cslice->cd(i+1);
171  hptslice->DrawCopy();
172  hptslice->Fit("gaus");
173 
174  cpt->cd(1);
175  hptres->GetYaxis()->SetTitleOffset(2.1);
176  hptres->GetXaxis()->SetTitleOffset(1.2);
177  hptres->SetMarkerStyle(20);
178  hptres->SetMarkerSize(1.2);
179  hptres->SetMarkerColor(col[i]);
180  hptres->GetXaxis()->SetRangeUser(0, ptmax);
181  hptres->GetYaxis()->SetRangeUser(0.0, 0.08);
182  hptres->SetTitle(";p_{T} [GeV/c];#frac{#Delta p_{T}}{p_{T}} (resolution)");
183  if(i == 0)
184  hptres->DrawCopy("p");
185  else
186  hptres->DrawCopy("p same");
187 
188  lpd->AddEntry(hptres, legvec[i].c_str());
189 
190  cpt->cd(2);
191  hptcent->GetYaxis()->SetTitleOffset(2.1);
192  hptcent->GetXaxis()->SetTitleOffset(1.2);
193  hptcent->SetMarkerStyle(20);
194  hptcent->SetMarkerSize(1.2);
195  hptcent->SetMarkerColor(col[i]);
196  hptcent->GetXaxis()->SetRangeUser(0, ptmax);
197  hptcent->GetYaxis()->SetRangeUser(-0.1, 0.1);
198  hptcent->SetTitle(";p_{T} [GeV/c];#frac{#Delta p_{T}}{p_{T}} (offset)");
199  if(i==0)
200  hptcent->DrawCopy("p");
201  else
202  hptcent->DrawCopy("p same");
203 
204  // Efficiency
205  ceff->cd();
206  TH1D *hnum = 0;
207  TH1D *hden = 0;
208  fin->GetObject("h3_num",hnum);
209  fin->GetObject("h3_den",hden);
210 
211  TH1D* heff = (TH1D*)hden->Clone("heff");;
212 
213  for(int i=1;i<=hden->GetNbinsX();++i)
214  {
215  double pass = hnum->GetBinContent(i);
216  double all = hden->GetBinContent(i);
217  double eff = 0;
218  if(all > pass)
219  eff = pass/all;
220  else if(all > 0)
221  eff = 1.;
222  heff->SetBinContent(i, eff);
223  }
224 
225  heff->GetYaxis()->SetTitle("Efficiency");
226  heff->GetYaxis()->SetTitleOffset(1.5);
227  heff->GetXaxis()->SetTitle("p_{T} (GeV/c)");
228  heff->GetXaxis()->SetTitleOffset(1.2);
229  heff->SetMarkerStyle(20);
230  heff->SetMarkerSize(1.2);
231  heff->SetMarkerColor(col[i]);
232  heff->GetXaxis()->SetRangeUser(0.0, ptmax);
233  heff->GetYaxis()->SetRangeUser(0.0, 1.05);
234  if(i==0)
235  heff->DrawCopy("p");
236  else
237  heff->DrawCopy("p same");
238 
239  TLine *unit = new TLine(0,1.0,40,1.0);
240  unit->Draw();
241 
242  // dca xy resolution
243  cdcaxy->cd();
244  TH2D *hdca2d;
245  fin->GetObject("h2",hdca2d);
246  ctemp1->cd(i+1);
247  hdca2d->Draw("colz");
248 
249  binlow = hdca2d->GetXaxis()->FindBin(slice_low);
250  binhigh = hdca2d->GetXaxis()->FindBin(slice_high);
251  TH1D* hdcaxyslice = hdca2d->ProjectionY("hdcaxyslice",binlow, binhigh);
252  hdcaxyslice->SetLineColor(col[i]);
253  cslicexy->cd(i+1);
254  hdcaxyslice->Draw();
255  hdcaxyslice->Fit("gaus");
256 
257  cdcaxy->cd();
258  hdca2d->FitSlicesY();
259  TH1D*hdcares = (TH1D*)gDirectory->Get("h2_2");
260  hdcares->GetYaxis()->SetTitleOffset(2.1);
261  hdcares->GetYaxis()->SetTitle("DCA(r#phi) (cm)");
262  hdcares->GetXaxis()->SetTitleOffset(1.2);
263  hdcares->GetXaxis()->SetTitle("p_{T} (GeV/c)");
264  hdcares->SetMarkerStyle(20);
265  hdcares->SetMarkerSize(1.2);
266  hdcares->SetMarkerColor(col[i]);
267  hdcares->GetYaxis()->SetRangeUser(0, 0.008);
268  if(i==0)
269  hdcares->DrawCopy("p");
270  else
271  hdcares->DrawCopy("p same");
272 
273  // dca z resolution
274  TH2D *hdcaZ2d;
275  fin->GetObject("h3",hdcaZ2d);
276  ctemp2->cd(i+1);
277  hdcaZ2d->Draw("colz");
278 
279  binlow = hdcaZ2d->GetXaxis()->FindBin(slice_low);
280  binhigh = hdcaZ2d->GetXaxis()->FindBin(slice_high);
281  TH1D* hdcazslice = hdcaZ2d->ProjectionY("hdcazslice",binlow, binhigh);
282  hdcazslice->SetLineColor(col[i]);
283  csliceZ->cd(i+1);
284  hdcazslice->Draw();
285  hdcazslice->Fit("gaus");
286 
287  cdcaZ->cd();
288  hdcaZ2d->FitSlicesY(fz);
289  TH1D* hdcaZres = (TH1D*)gDirectory->Get("h3_2");
290 
291  hdcaZres->GetYaxis()->SetTitleOffset(2.1);
292  hdcaZres->GetYaxis()->SetTitle("DCA(Z) (cm)");
293  hdcaZres->GetXaxis()->SetTitleOffset(1.2);
294  hdcaZres->GetXaxis()->SetTitle("p_{T} (GeV/c");
295  hdcaZres->SetMarkerStyle(20);
296  hdcaZres->SetMarkerSize(1.2);
297  hdcaZres->SetMarkerColor(col[i]);
298  hdcaZres->GetYaxis()->SetRangeUser(0, 0.015);
299  if(i==0)
300  hdcaZres->DrawCopy("p");
301  else
302  hdcaZres->DrawCopy("p same");
303 
304  // Plot the fraction of tracks matched to the MVTX
305  cmvtx->cd();
306  TH2D *hnmvtx = 0;
307  fin->GetObject("h6",hnmvtx);
308  if(!hnmvtx) cout << " Did not get hnmvtx" << endl;
309 
310  //hnmvtx->DrawCopy("p");
311 
312  binlow = hnmvtx->GetYaxis()->FindBin(0.5);
313  binhigh = hnmvtx->GetYaxis()->FindBin(3.5);
314  cout << "binlow " << binlow << " binhigh " << binhigh << endl;
315  TH1D *hnmaps_hit = hnmvtx->ProjectionX("hnmaps_hit",binlow, binhigh);;
316 
317  binlow = hnmvtx->GetYaxis()->FindBin(-0.5);
318  binhigh = hnmvtx->GetYaxis()->FindBin(3.5);
319  cout << "binlow " << binlow << " binhigh " << binhigh << endl;
320  TH1D *hnmaps_all = hnmvtx->ProjectionX("hnmaps_all",binlow, binhigh);;
321 
322  TH1D *hnmaps = (TH1D*) hnmaps_all->Clone();
323  for(int i=1;i<=hnmaps_hit->GetNbinsX();++i)
324  {
325  double hit = hnmaps_hit->GetBinContent(i);
326  double all = hnmaps_all->GetBinContent(i);
327  double eff = 0;
328  if(all > hit)
329  eff = hit/all;
330  else if(all > 0)
331  eff = 1.;
332  hnmaps->SetBinContent(i, eff);
333  }
334 
335  hnmaps->GetXaxis()->SetTitle("p_{T} (GeV/c)");
336  hnmaps->GetYaxis()->SetTitle("MVTX association efficiency");
337  hnmaps->GetYaxis()->SetTitleOffset(1.5);
338  hnmaps->GetXaxis()->SetTitleOffset(1.2);
339  hnmaps->SetMinimum(0);
340  hnmaps->SetMaximum(1.1);
341  hnmaps->SetMarkerStyle(20);
342  hnmaps->SetMarkerSize(1.5);
343  hnmaps->SetMarkerColor(col[i]);
344 
345  if(i==0)
346  {
347  hnmaps->DrawCopy("p");
348  unit->Draw();
349  }
350  else
351  hnmaps->DrawCopy("p same");
352 
353  // Combine track efficiency and MVTX association efficiency
354 
355  ccomb->cd();
356 
357  cout << " hnmaps bins " << hnmaps->GetNbinsX() << " heff bins " << heff->GetNbinsX() << endl;
358  heff->Rebin(3);
359  cout << heff->GetNbinsX() << endl;
360 
361  TH1D *hcomb = (TH1D*) hnmaps->Clone();
362  for(int icomb=1; icomb< hnmaps->GetNbinsX(); ++icomb)
363  {
364  double vmaps = hnmaps->GetBinContent(icomb);
365  double vtrack = heff->GetBinContent(icomb) / 3;
366  double eff = 0;
367  eff = vmaps * vtrack;
368  //cout << " icomb " << icomb << " vmaps " << vmaps << " vtrack " << vtrack << " comb " << eff << endl;
369  hcomb->SetBinContent(icomb, eff);
370  }
371 
372  hcomb->GetXaxis()->SetRangeUser(0,39);
373  hcomb->GetYaxis()->SetTitle("Efficincy of MVTX associated tracks");
374 
375  TLine *unit2 = new TLine(0,1.0,39,1.0);
376 
377  if(i==0)
378  {
379  hcomb->DrawCopy("p");
380  unit2->Draw();
381  }
382  else
383  hcomb->DrawCopy("p same");
384 
385  }
386 
387  cpt->cd(2); lpd->Draw();
388 
389  cdcaxy->cd(); lpd->Draw();
390 
391  cdcaZ->cd(); lpd->Draw();
392 
393  ceff->cd(); lpd->Draw();
394 
395  cmvtx->cd(); lpd->Draw();
396 
397  ccomb->cd(); lpd->Draw();
398 
399 
400 
401 }
402