ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/PixelOnlyTracking/MyPixelAnalysisRatios.cc
Revision: 1.14
Committed: Thu Jun 17 19:13:01 2010 UTC (14 years, 10 months ago) by dinardo
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.13: +0 -0 lines
State: FILE REMOVED
Error occurred while calculating annotation data.
Log Message:
*** empty log message ***

File Contents

# Content
1 // #############################################################
2 // # Program to compute efficiencies, purities and resolutions #
3 // # for pixel stand alone tracking and vertexing #
4 // #############################################################
5 // # Author: Mauro Dinardo #
6 // #############################################################
7
8 #include <TROOT.h>
9 #include <TF1.h>
10 #include <TH1.h>
11 #include <TH2.h>
12 #include <TFile.h>
13
14 #include <math.h>
15 #include <map>
16 #include <stdlib.h>
17 #include <sstream>
18 #include <string>
19
20
21 using namespace std;
22
23
24 // ####################
25 // # Histogram TRACKS #
26 // ####################
27
28 // @@@@@@ Pixel tracks @@@@@@
29 TH1F* H_pix_trk_pt;
30 TH1F* H_pix_trk_eta;
31
32 TH1F* H_pix_trk_doubleCounter_pt;
33 TH1F* H_pix_trk_doubleCounter_eta;
34
35 TH2F* H_pix_trk_normchi2_pt_eta[2];
36 TH2F* H_pix_trk_eta_phi;
37 TH2F* H_pix_trk_eta_pt;
38
39 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
40 // @ Histograms containing ratios @
41 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
42 TH1F* H_pix_trk_doubleCounter_pt_rel;
43 TH1F* H_pix_trk_doubleCounter_eta_rel;
44
45 // @@@@@@ Tracker tracks @@@@@@
46 TH1F* H_trk_trk_pt;
47 TH1F* H_trk_trk_eta;
48
49 TH1F* H_trk_trk_doubleCounter_pt;
50 TH1F* H_trk_trk_doubleCounter_eta;
51
52 TH2F* H_trk_trk_normchi2_pt_eta[2];
53 TH2F* H_trk_trk_eta_phi;
54 TH2F* H_trk_trk_eta_pt;
55
56 TH2F* H_trk_trk_eta_phi_algo;
57 TH2F* H_trk_trk_eta_pt_algo;
58
59 TH1F* H_trk_toteff_pt_num;
60 TH1F* H_trk_toteff_pt_den;
61 TH1F* H_trk_toteff_eta_num;
62 TH1F* H_trk_toteff_eta_den;
63 TH1F* H_trk_toteff_phi_num[3];
64 TH1F* H_trk_toteff_phi_den[3];
65
66 TH1F* H_trk_eff_pt_num;
67 TH1F* H_trk_eff_pt_den;
68 TH1F* H_trk_eff_eta_num;
69 TH1F* H_trk_eff_eta_den;
70 TH1F* H_trk_eff_phi_num[3];
71 TH1F* H_trk_eff_phi_den[3];
72
73 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
74 // @ Histograms containing ratios @
75 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
76 TH1F* H_trk_trk_doubleCounter_pt_rel;
77 TH1F* H_trk_trk_doubleCounter_eta_rel;
78 TH1F* H_trk_toteff_pt;
79 TH1F* H_trk_toteff_eta;
80 TH1F* H_trk_toteff_phi[3];
81 TH1F* H_trk_eff_pt;
82 TH1F* H_trk_eff_eta;
83 TH1F* H_trk_eff_phi[3];
84
85 // @@@@@@ Common pixel&tracker tracks @@@@@@
86 typedef map<string, map<string, TH1F*> > myMap;
87 myMap MapFit_d0; // Map to store d0 difference [pt, eta]
88 myMap MapFit_dz; // Map to store dz difference [pt, eta]
89 myMap MapFit_pt; // Map to store pt difference [pt, eta]
90 myMap MapFit_d0_pull; // Map to store d0 pulls [pt, eta]
91 myMap MapFit_dz_pull; // Map to store dz pulls [pt, eta]
92 myMap MapFit_pt_pull; // Map to store pt pulls [pt, eta]
93 myMap MapFit_d0err; // Map to store generalTrack d0err [pt, eta]
94 myMap MapFit_dzerr; // Map to store generalTrack dzerr [pt, eta]
95 myMap MapFit_pterr; // Map to store generalTrack pterr [pt, eta]
96
97 TH2F* H_trked_eta_phi_pur;
98 TH2F* H_trked_eta_phi_eff;
99 TH2F* H_trked_eta_pt_pur;
100 TH2F* H_trked_eta_pt_eff;
101 TH2F* H_trked_normchi2_pt_pur_eta[2];
102 TH2F* H_trked_normchi2_pt_eff_eta[2];
103
104 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
105 // @ Histograms containing ratios @
106 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
107 TH1F* H_ptres_eta_pt[5];
108 TH1F* H_ptpull_eta_pt[5];
109
110 TH1F* H_dzres_eta_pt[5];
111 TH1F* H_dzpull_eta_pt[5];
112
113 TH1F* H_d0res_eta_pt[5];
114 TH1F* H_d0pull_eta_pt[5];
115
116 TH1F* H_ptres_pt_eta[2];
117 TH1F* H_ptpull_pt_eta[2];
118
119 TH1F* H_dzres_pt_eta[2];
120 TH1F* H_dzpull_pt_eta[2];
121
122 TH1F* H_d0res_pt_eta[2];
123 TH1F* H_d0pull_pt_eta[2];
124
125 TH1F* H_skew_pt_eta[2];
126
127 TH1F* H_purity_trked_eta;
128 TH1F* H_purity_trked_pt_eta[2];
129 TH1F* H_purity_trked_pt_whole_eta;
130 TH1F* H_purity_trked_phi[3];
131
132 TH1F* H_efficiency_trked_eta;
133 TH1F* H_efficiency_trked_pt_eta[2];
134 TH1F* H_efficiency_trked_pt_whole_eta;
135 TH1F* H_efficiency_trked_phi[3];
136
137 TH1F* H_efficiency_trked_eta_algo;
138 TH1F* H_efficiency_trked_pt_eta_algo[2];
139 TH1F* H_efficiency_trked_pt_whole_eta_algo;
140 TH1F* H_efficiency_trked_phi_algo[3];
141
142 TH1F* H_tk_pur_eff;
143
144 TH2F* H_purity_trked_normchi2_pt_eta[2];
145 TH2F* H_efficiency_trked_normchi2_pt_eta[2];
146
147
148 // ######################
149 // # Histogram VERTICES #
150 // ######################
151
152 // @@@@@@ Common pixel&tracker vertices @@@@@@
153 TH1F* H_vx_counters;
154
155 map<string, TH1F*> map_H_vx_x_diff;
156 map<string, TH1F*> map_H_vx_y_diff;
157 map<string, TH1F*> map_H_vx_z_diff;
158 map<string, TH1F*> map_H_vx_x_diff_pull;
159 map<string, TH1F*> map_H_vx_y_diff_pull;
160 map<string, TH1F*> map_H_vx_z_diff_pull;
161 map<string, TH1F*> map_H_vx_x_diff_orig;
162 map<string, TH1F*> map_H_vx_y_diff_orig;
163 map<string, TH1F*> map_H_vx_z_diff_orig;
164 map<string, TH1F*> map_H_vx_x_err2;
165 map<string, TH1F*> map_H_vx_y_err2;
166 map<string, TH1F*> map_H_vx_z_err2;
167
168 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
169 // @ Histograms containing ratios @
170 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
171 TH1F* H_vx_x_diff_vs_trk;
172 TH1F* H_vx_y_diff_vs_trk;
173 TH1F* H_vx_z_diff_vs_trk;
174 TH1F* H_vx_x_diff_pull_vs_trk;
175 TH1F* H_vx_y_diff_pull_vs_trk;
176 TH1F* H_vx_z_diff_pull_vs_trk;
177 TH1F* H_vx_x_diff_orig_vs_trk;
178 TH1F* H_vx_y_diff_orig_vs_trk;
179 TH1F* H_vx_z_diff_orig_vs_trk;
180 TH1F* H_vx_x_diff_noerr_vs_trk;
181 TH1F* H_vx_y_diff_noerr_vs_trk;
182 TH1F* H_vx_z_diff_noerr_vs_trk;
183 TH1F* H_vx_pur_eff;
184
185
186 // @@@@@@ Duplicated variables @@@@@@
187 double EtaThr; // --> Must have at most 2 decimal digits
188 double EtaThrBin[3]; // --> Must have at most 2 decimal digits
189 double zRange; // Unit: [um]
190 double d0Range; // Unit: [cm]
191 double dzRange; // Unit: [cm]
192 double ptRange; // Unit: [GeV/c]
193 double d0errRange; // Unit: [um]^2
194 double dzerrRange; // Unit: [um]^2
195 double pterrRange; // Unit: [MeV/c]^2
196 double pullsRange;
197 double MaxPtEquiBinning;
198
199 // @@@@@@ Non duplicated variables @@@@@@
200 double PtThr[6]; // --> Must have at most 2 decimal digits
201 unsigned int MinEntriesEffPur; // Histogram minimum number of entries for efficiency and purity
202
203 // @@@@@@ Parameters from cfg file @@@@@@
204 double MaxEtaTkTrk; // Max eta that a track should have in order to be taken into account
205 double RangePt; // pt bins. Unit: [GeV/c] --> Must have at most 2 decimal digits
206 double PtStep; // pt bins. Unit: [GeV/c] --> Must have at most 2 decimal digits
207 double RangeEta; // eta bins --> Must have at most 2 decimal digits
208 double EtaStep; // eta bins --> Must have at most 2 decimal digits
209 double MinPtTk; // Mini Pt that a track must have in order to be taken into account. Unit: [GeV/c]
210 // --> Must have at most 2 decimal digits
211
212
213 void InitVariables ()
214 {
215 MinEntriesEffPur = 12;
216
217 MaxEtaTkTrk = 2.3;
218 RangePt = 8.0;
219 PtStep = 0.25;
220 RangeEta = 2.3;
221 EtaStep = 0.23;
222 MinPtTk = 0.75;
223
224 cout << "\n@@@@@@ CONFIGURATION PARAMETERS @@@@@@" << endl;
225
226 cout << "@@@@@@ GENERAL PARAMETERS @@@@@@" << endl;
227 cout << "MinEntriesEffPur --> " << MinEntriesEffPur << endl;
228
229 cout << "@@@@@@ TRACK PARAMETERS @@@@@@" << endl;
230 cout << "MaxEtaTkTrk --> " << MaxEtaTkTrk << endl;
231 cout << "RangePt --> " << RangePt << endl;
232 cout << "PtStep --> " << PtStep << endl;
233 cout << "RangeEta --> " << RangeEta << endl;
234 cout << "EtaStep --> " << EtaStep << endl;
235 cout << "MinPtTk --> " << MinPtTk << endl;
236
237 cout << "@@@@@@ VERTEX PARAMETERS @@@@@@" << endl;
238
239 MaxPtEquiBinning = 3.;
240
241 EtaThr = 1.5;
242 EtaThrBin[0] = 0;
243 EtaThrBin[1] = EtaThr;
244 EtaThrBin[2] = MaxEtaTkTrk;
245
246 PtThr[0] = MinPtTk;
247 for (int i = 1; i < 5; i++) PtThr[i] = PtThr[i-1] + 0.5;
248 PtThr[5] = RangePt;
249
250 zRange = 6000.; // Unit: [um]
251
252 d0Range = 0.6; // Unit: [cm]
253 dzRange = 0.6; // Unit: [cm]
254 ptRange = 12.0; // Unit: [GeV/c]
255 d0errRange = 10000.; // Unit: [um]^2
256 dzerrRange = 10000.; // Unit: [um]^2
257 pterrRange = 50000.; // Unit: [MeV/c]^2
258
259 pullsRange = 10.;
260 }
261
262
263 bool LoadHistos (string datatype, TFile* filename)
264 {
265 if (datatype.compare("Tk") == 0)
266 {
267 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
268 // @ Histograms containing numerators or denominators @
269 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
270 H_pix_trk_pt = (TH1F*)filename->Get("MyProcess/PixelTracks/Pixeltracks pt distribution");
271 H_pix_trk_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/Pixeltracks eta distribution");
272 H_pix_trk_doubleCounter_pt = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter pixeltracks pt distribution");
273 H_pix_trk_doubleCounter_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter pixeltracks eta distribution");
274
275 H_pix_trk_eta_phi = (TH2F*)filename->Get("MyProcess/PixelTracks/Pixeltracks distribution in eta and phi");
276 H_pix_trk_eta_pt = (TH2F*)filename->Get("MyProcess/PixelTracks/Pixeltracks distribution in eta and pt");
277
278 H_trk_trk_pt = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks pt distribution");
279 H_trk_trk_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks eta distribution");
280 H_trk_trk_doubleCounter_pt = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter trackertracks pt distribution");
281 H_trk_trk_doubleCounter_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter trackertracks eta distribution");
282
283 H_trk_trk_eta_phi = (TH2F*)filename->Get("MyProcess/PixelTracks/Trackertracks distribution in eta and phi");
284 H_trk_trk_eta_pt = (TH2F*)filename->Get("MyProcess/PixelTracks/Trackertracks distribution in eta and pt");
285
286 H_trk_trk_eta_phi_algo = (TH2F*)filename->Get("MyProcess/PixelTracks/Trackertracks distribution in eta and phi algo");
287 H_trk_trk_eta_pt_algo = (TH2F*)filename->Get("MyProcess/PixelTracks/Trackertracks distribution in eta and pt algo");
288
289 H_trk_eff_pt_num = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs pt numerator");
290 H_trk_eff_pt_den = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs pt denominator");
291 H_trk_eff_eta_num = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs eta numerator");
292 H_trk_eff_eta_den = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs eta denominator");
293
294 H_trk_toteff_pt_num = (TH1F*)filename->Get("MyProcess/PixelTracks/Total trackertracks efficiency vs pt numerator");
295 H_trk_toteff_pt_den = (TH1F*)filename->Get("MyProcess/PixelTracks/Total trackertracks efficiency vs pt denominator");
296 H_trk_toteff_eta_num = (TH1F*)filename->Get("MyProcess/PixelTracks/Total trackertracks efficiency vs eta numerator");
297 H_trk_toteff_eta_den = (TH1F*)filename->Get("MyProcess/PixelTracks/Total trackertracks efficiency vs eta denominator");
298
299 for (int i = 0; i < 3; i++)
300 {
301 stringstream i_str;
302 i_str << i + 1;
303 TString hName;
304
305 hName = "MyProcess/PixelTracks/Trackertracks efficiency vs phi numerator thr" + i_str.str();
306 H_trk_eff_phi_num[i] = (TH1F*)filename->Get(hName);
307
308 hName = "MyProcess/PixelTracks/Trackertracks efficiency vs phi denominator thr" + i_str.str();
309 H_trk_eff_phi_den[i] = (TH1F*)filename->Get(hName);
310
311 hName = "MyProcess/PixelTracks/Total trackertracks efficiency vs phi numerator thr" + i_str.str();
312 H_trk_toteff_phi_num[i] = (TH1F*)filename->Get(hName);
313
314 hName = "MyProcess/PixelTracks/Total trackertracks efficiency vs phi denominator thr" + i_str.str();
315 H_trk_toteff_phi_den[i] = (TH1F*)filename->Get(hName);
316 }
317
318 for (int i = 0; i < 2; i++)
319 {
320 stringstream i_str;
321 i_str << i + 1;
322 TString hName;
323
324 hName = "MyProcess/PixelTracks/Pixeltracks vs chi2DoF and pt cuts thr" + i_str.str();
325 H_pix_trk_normchi2_pt_eta[i] = (TH2F*)filename->Get(hName);
326
327 hName = "MyProcess/PixelTracks/Trackertracks vs chi2DoF and pt cuts thr" + i_str.str();
328 H_trk_trk_normchi2_pt_eta[i] = (TH2F*)filename->Get(hName);
329
330 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Matched pixeltracks vs chi2DoF and pt cuts pur thr" + i_str.str();
331 H_trked_normchi2_pt_pur_eta[i] = (TH2F*)filename->Get(hName);
332
333 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Matched pixeltracks vs chi2DoF and pt cuts eff thr" + i_str.str();
334 H_trked_normchi2_pt_eff_eta[i] = (TH2F*)filename->Get(hName);
335 }
336
337 H_trked_eta_phi_pur = (TH2F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Matched pixeltracks distribution in eta and phi pur");
338 H_trked_eta_phi_eff = (TH2F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Matched pixeltracks distribution in eta and phi eff");
339 H_trked_eta_pt_pur = (TH2F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Matched pixeltracks distribution in eta and pt pur");
340 H_trked_eta_pt_eff = (TH2F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Matched pixeltracks distribution in eta and pt eff");
341
342 stringstream PtCutStr;
343 stringstream EtaCutStr;
344 stringstream HistoName;
345 PtCutStr.precision(3);
346 EtaCutStr.precision(3);
347 HistoName.precision(3);
348
349 double PtCut;
350 double EtaCut;
351 for (PtCut = PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
352 {
353 PtCutStr.str("");
354 PtCutStr << PtCut;
355
356 for (EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
357 {
358 EtaCutStr.str("");
359 EtaCutStr << EtaCut;
360
361 HistoName.str("");
362 HistoName << "MyProcess/TrackFit/d0_pulls_pt_" << PtCut << "_eta_" << EtaCut;
363 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
364 MapFit_d0_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
365
366 HistoName.str("");
367 HistoName << "MyProcess/TrackFit/d0err_pt_" << PtCut << "_eta_" << EtaCut;
368 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
369 MapFit_d0err[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
370
371 HistoName.str("");
372 HistoName << "MyProcess/TrackFit/d0_pt_" << PtCut << "_eta_" << EtaCut;
373 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
374 MapFit_d0[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
375
376 HistoName.str("");
377 HistoName << "MyProcess/TrackFit/dz_pulls_pt_" << PtCut << "_eta_" << EtaCut;
378 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
379 MapFit_dz_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
380
381 HistoName.str("");
382 HistoName << "MyProcess/TrackFit/dzerr_pt_" << PtCut << "_eta_" << EtaCut;
383 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
384 MapFit_dzerr[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
385
386 HistoName.str("");
387 HistoName << "MyProcess/TrackFit/dz_pt_" << PtCut << "_eta_" << EtaCut;
388 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
389 MapFit_dz[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
390
391 HistoName.str("");
392 HistoName << "MyProcess/TrackFit/pt_pulls_pt_" << PtCut << "_eta_" << EtaCut;
393 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
394 MapFit_pt_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
395
396 HistoName.str("");
397 HistoName << "MyProcess/TrackFit/pterr_pt_" << PtCut << "_eta_" << EtaCut;
398 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
399 MapFit_pterr[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
400
401 HistoName.str("");
402 HistoName << "MyProcess/TrackFit/pt_pt_" << PtCut << "_eta_" << EtaCut;
403 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
404 MapFit_pt[PtCutStr.str().c_str()][EtaCutStr.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
405 }
406 }
407
408 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
409 // @ Histograms containing ratios @
410 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
411 H_pix_trk_doubleCounter_pt_rel = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter pixeltracks pt relative distribution");
412 H_pix_trk_doubleCounter_eta_rel = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter pixeltracks eta relative distribution");
413 H_pix_trk_doubleCounter_pt_rel->Reset();
414 H_pix_trk_doubleCounter_eta_rel->Reset();
415
416 H_trk_trk_doubleCounter_pt_rel = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter trackertracks pt relative distribution");
417 H_trk_trk_doubleCounter_eta_rel = (TH1F*)filename->Get("MyProcess/PixelTracks/Doublecounter trackertracks eta relative distribution");
418 H_trk_trk_doubleCounter_pt_rel->Reset();
419 H_trk_trk_doubleCounter_eta_rel->Reset();
420
421 H_trk_eff_pt = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs pt");
422 H_trk_eff_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs eta");
423 H_trk_eff_pt->Reset();
424 H_trk_eff_eta->Reset();
425
426 H_trk_toteff_pt = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs pt");
427 H_trk_toteff_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/Trackertracks efficiency vs eta");
428 H_trk_toteff_pt->Reset();
429 H_trk_toteff_eta->Reset();
430
431 for (int i = 0; i < 3; i++)
432 {
433 stringstream i_str;
434 i_str << i + 1;
435 TString hName;
436
437 hName = "MyProcess/PixelTracks/Trackertracks efficiency vs phi thr" + i_str.str();
438 H_trk_eff_phi[i] = (TH1F*)filename->Get(hName);
439 H_trk_eff_phi[i]->Reset();
440
441 hName = "MyProcess/PixelTracks/Total trackertracks efficiency vs phi thr" + i_str.str();
442 H_trk_toteff_phi[i] = (TH1F*)filename->Get(hName);
443 H_trk_toteff_phi[i]->Reset();
444 }
445
446 // #################################
447 // # Use this for constant binning #
448 // #################################
449 // const int nPtBins = rint(RangePt/PtStep);
450 // double* xPtBins = NULL;
451 // xPtBins = new double[nPtBins + 1];
452 // xPtBins[0] = 0.;
453 // for (int i = 1; i < nPtBins; i++) xPtBins[i] = xPtBins[i-1] + PtStep;
454 // xPtBins[nPtBins] = RangePt;
455
456 // #####################################
457 // # Use this for NON constant binning #
458 // #####################################
459 const int nPtBins = rint(MaxPtEquiBinning/PtStep) + 7;
460 double* xPtBins = NULL;
461 xPtBins = new double[nPtBins + 1];
462 xPtBins[0] = 0.;
463 for (int i = 1; i <= rint(MaxPtEquiBinning/PtStep); i++) xPtBins[i] = xPtBins[i-1] + PtStep;
464 xPtBins[(int)rint(MaxPtEquiBinning/PtStep) + 1] = 3.5;
465 xPtBins[(int)rint(MaxPtEquiBinning/PtStep) + 2] = 4.0;
466 xPtBins[(int)rint(MaxPtEquiBinning/PtStep) + 3] = 4.5;
467 xPtBins[(int)rint(MaxPtEquiBinning/PtStep) + 4] = 5.25;
468 xPtBins[(int)rint(MaxPtEquiBinning/PtStep) + 5] = 6.0;
469 xPtBins[(int)rint(MaxPtEquiBinning/PtStep) + 6] = 7.0;
470 xPtBins[(int)rint(MaxPtEquiBinning/PtStep) + 7] = RangePt;
471
472 for (int i = 0; i < 2; i++)
473 {
474 stringstream i_str;
475 i_str << i + 1;
476 TString hName;
477
478 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Purity pixeltracks vs chi2DoF and pt cuts thr" + i_str.str();
479 H_purity_trked_normchi2_pt_eta[i] = (TH2F*)filename->Get(hName);
480 H_purity_trked_normchi2_pt_eta[i]->Reset();
481
482 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Purity pixeltracks vs pt thr" + i_str.str();
483 H_purity_trked_pt_eta[i] = (TH1F*)filename->Get(hName);
484 H_purity_trked_pt_eta[i]->Reset();
485
486 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Efficiency pixeltracks vs chi2DoF and pt cuts thr" + i_str.str();
487 H_efficiency_trked_normchi2_pt_eta[i] = (TH2F*)filename->Get(hName);
488 H_efficiency_trked_normchi2_pt_eta[i]->Reset();
489
490 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Efficiency pixeltracks vs pt thr" + i_str.str();
491 H_efficiency_trked_pt_eta[i] = (TH1F*)filename->Get(hName);
492 H_efficiency_trked_pt_eta[i]->Reset();
493
494 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Algoefficiency pixeltracks vs pt thr" + i_str.str();
495 H_efficiency_trked_pt_eta_algo[i] = (TH1F*)filename->Get(hName);
496 H_efficiency_trked_pt_eta_algo[i]->Reset();
497
498 hName = "MyProcess/PixelTracks/CommonPxTkPlots/pt resolution vs pt thr" + i_str.str();
499 H_ptres_pt_eta[i] = (TH1F*)filename->Get(hName);
500 H_ptres_pt_eta[i]->Reset();
501 H_ptres_pt_eta[i] = (TH1F*)H_ptres_pt_eta[i]->Rebin(nPtBins, H_ptres_pt_eta[i]->GetName(), xPtBins);
502
503 hName = "MyProcess/PixelTracks/CommonPxTkPlots/pt pulls vs pt thr" + i_str.str();
504 H_ptpull_pt_eta[i] = (TH1F*)filename->Get(hName);
505 H_ptpull_pt_eta[i]->Reset();
506 H_ptpull_pt_eta[i] = (TH1F*)H_ptpull_pt_eta[i]->Rebin(nPtBins, H_ptpull_pt_eta[i]->GetName(), xPtBins);
507
508 hName = "MyProcess/PixelTracks/CommonPxTkPlots/dz resolution vs pt thr" + i_str.str();
509 H_dzres_pt_eta[i] = (TH1F*)filename->Get(hName);
510 H_dzres_pt_eta[i]->Reset();
511
512 hName = "MyProcess/PixelTracks/CommonPxTkPlots/dz pulls vs pt thr" + i_str.str();
513 H_dzpull_pt_eta[i] = (TH1F*)filename->Get(hName);
514 H_dzpull_pt_eta[i]->Reset();
515
516 hName = "MyProcess/PixelTracks/CommonPxTkPlots/d0 resolution vs pt thr" + i_str.str();
517 H_d0res_pt_eta[i] = (TH1F*)filename->Get(hName);
518 H_d0res_pt_eta[i]->Reset();
519
520 hName = "MyProcess/PixelTracks/CommonPxTkPlots/d0 pulls vs pt thr" + i_str.str();
521 H_d0pull_pt_eta[i] = (TH1F*)filename->Get(hName);
522 H_d0pull_pt_eta[i]->Reset();
523
524 hName = "MyProcess/PixelTracks/CommonPxTkPlots/SkewPixTrkpt vs pt thr" + i_str.str();
525 H_skew_pt_eta[i] = (TH1F*)filename->Get(hName);
526 H_skew_pt_eta[i]->Reset();
527 }
528
529 for (int i = 0; i < 5; i++)
530 {
531 stringstream i_str;
532 i_str << i + 1;
533 TString hName;
534
535 hName = "MyProcess/PixelTracks/CommonPxTkPlots/pt resolution vs eta thr" + i_str.str();
536 H_ptres_eta_pt[i] = (TH1F*)filename->Get(hName);
537 H_ptres_eta_pt[i]->Reset();
538
539 hName = "MyProcess/PixelTracks/CommonPxTkPlots/pt pulls vs eta thr" + i_str.str();
540 H_ptpull_eta_pt[i] = (TH1F*)filename->Get(hName);
541 H_ptpull_eta_pt[i]->Reset();
542
543 hName = "MyProcess/PixelTracks/CommonPxTkPlots/dz resolution vs eta thr" + i_str.str();
544 H_dzres_eta_pt[i] = (TH1F*)filename->Get(hName);
545 H_dzres_eta_pt[i]->Reset();
546
547 hName = "MyProcess/PixelTracks/CommonPxTkPlots/dz pulls vs eta thr" + i_str.str();
548 H_dzpull_eta_pt[i] = (TH1F*)filename->Get(hName);
549 H_dzpull_eta_pt[i]->Reset();
550
551 hName = "MyProcess/PixelTracks/CommonPxTkPlots/d0 resolution vs eta thr" + i_str.str();
552 H_d0res_eta_pt[i] = (TH1F*)filename->Get(hName);
553 H_d0res_eta_pt[i]->Reset();
554
555 hName = "MyProcess/PixelTracks/CommonPxTkPlots/d0 pulls vs eta thr" + i_str.str();
556 H_d0pull_eta_pt[i] = (TH1F*)filename->Get(hName);
557 H_d0pull_eta_pt[i]->Reset();
558 }
559
560 for (int i = 0; i < 3; i++)
561 {
562 stringstream i_str;
563 i_str << i + 1;
564 TString hName;
565
566 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Purity pixeltracks vs phi thr" + i_str.str();
567 H_purity_trked_phi[i] = (TH1F*)filename->Get(hName);
568 H_purity_trked_phi[i]->Reset();
569
570 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Efficiency pixeltracks vs phi thr" + i_str.str();
571 H_efficiency_trked_phi[i] = (TH1F*)filename->Get(hName);
572 H_efficiency_trked_phi[i]->Reset();
573
574 hName = "MyProcess/PixelTracks/CommonPxTkPlots/Algoefficiency pixeltracks vs phi thr" + i_str.str();
575 H_efficiency_trked_phi_algo[i] = (TH1F*)filename->Get(hName);
576 H_efficiency_trked_phi_algo[i]->Reset();
577 }
578
579 H_purity_trked_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Purity pixeltracks vs eta");
580 H_purity_trked_pt_whole_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Purity pixeltracks vs pt");
581 H_purity_trked_eta->Reset();
582 H_purity_trked_pt_whole_eta->Reset();
583
584 H_efficiency_trked_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Efficiency pixeltracks vs eta");
585 H_efficiency_trked_pt_whole_eta = (TH1F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Efficiency pixeltracks vs pt");
586 H_efficiency_trked_eta->Reset();
587 H_efficiency_trked_pt_whole_eta->Reset();
588
589 H_efficiency_trked_eta_algo = (TH1F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Algoefficiency pixeltracks vs eta");
590 H_efficiency_trked_pt_whole_eta_algo = (TH1F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Algoefficiency pixeltracks vs pt");
591 H_efficiency_trked_eta_algo->Reset();
592 H_efficiency_trked_pt_whole_eta_algo->Reset();
593
594 H_tk_pur_eff = (TH1F*)filename->Get("MyProcess/PixelTracks/CommonPxTkPlots/Purity Efficiency Tracks");
595 H_tk_pur_eff->Reset();
596
597 return true;
598 }
599 else if (datatype.compare("Vx") == 0)
600 {
601 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
602 // @ Histograms containing numerators or denominators @
603 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
604 H_vx_counters = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Counters Match Purity Efficiency for Vertices");
605
606 stringstream String;
607 stringstream HistoName;
608 String.precision(3);
609 HistoName.precision(3);
610
611 for (int i = 1; i < 100; i++)
612 {
613 String << i;
614
615 HistoName.str("");
616 HistoName << "MyProcess/VxNtrks/x_diff_Ntrk_" << String.str();
617 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
618 map_H_vx_x_diff[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
619
620 HistoName.str("");
621 HistoName << "MyProcess/VxNtrks/y_diff_Ntrk_" << String.str();
622 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
623 map_H_vx_y_diff[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
624
625 HistoName.str("");
626 HistoName << "MyProcess/VxNtrks/z_diff_Ntrk_" << String.str();
627 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
628 map_H_vx_z_diff[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
629
630 HistoName.str("");
631 HistoName << "MyProcess/VxNtrks/x_diff_pull_Ntrk_" << String.str();
632 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
633 map_H_vx_x_diff_pull[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
634
635 HistoName.str("");
636 HistoName << "MyProcess/VxNtrks/y_diff_pull_Ntrk_" << String.str();
637 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
638 map_H_vx_y_diff_pull[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
639
640 HistoName.str("");
641 HistoName << "MyProcess/VxNtrks/z_diff_pull_Ntrk_" << String.str();
642 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
643 map_H_vx_z_diff_pull[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
644
645 HistoName.str("");
646 HistoName << "MyProcess/VxNtrks/x_diff_orig_Ntrk_" << String.str();
647 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
648 map_H_vx_x_diff_orig[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
649
650 HistoName.str("");
651 HistoName << "MyProcess/VxNtrks/y_diff_orig_Ntrk_" << String.str();
652 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
653 map_H_vx_y_diff_orig[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
654
655 HistoName.str("");
656 HistoName << "MyProcess/VxNtrks/z_diff_orig_Ntrk_" << String.str();
657 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
658 map_H_vx_z_diff_orig[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
659
660 HistoName.str("");
661 HistoName << "MyProcess/VxNtrks/x_err2_Ntrk_" << String.str();
662 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
663 map_H_vx_x_err2[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
664
665 HistoName.str("");
666 HistoName << "MyProcess/VxNtrks/y_err2_Ntrk_" << String.str();
667 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
668 map_H_vx_y_err2[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
669
670 HistoName.str("");
671 HistoName << "MyProcess/VxNtrks/z_err2_Ntrk_" << String.str();
672 if ((TH1F*)filename->Get(HistoName.str().c_str()) != NULL)
673 map_H_vx_z_err2[String.str().c_str()] = (TH1F*)filename->Get(HistoName.str().c_str());
674
675 String.str("");
676 }
677
678 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
679 // @ Histograms containing ratios @
680 // @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
681 H_vx_pur_eff = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Purity Efficiency Vertices");
682 H_vx_x_diff_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex position difference along X vs Ntrks");
683 H_vx_y_diff_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex position difference along Y vs Ntrks");
684 H_vx_z_diff_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex position difference along Z vs Ntrks");
685 H_vx_x_diff_pull_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex position pull along X vs Ntrks");
686 H_vx_y_diff_pull_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex position pull along Y vs Ntrks");
687 H_vx_z_diff_pull_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex position pull along Z vs Ntrks");
688 H_vx_x_diff_orig_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Original vertex position difference along X vs Ntrks");
689 H_vx_y_diff_orig_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Original vertex position difference along Y vs Ntrks");
690 H_vx_z_diff_orig_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Original vertex position difference along Z vs Ntrks");
691 H_vx_x_diff_noerr_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex resolution along X vs Ntrks");
692 H_vx_y_diff_noerr_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex resolution along Y vs Ntrks");
693 H_vx_z_diff_noerr_vs_trk = (TH1F*)filename->Get("MyProcess/PixelVertices/CommonPxTkPlots/Vertex resolution along Z vs Ntrks");
694 H_vx_x_diff_vs_trk->Reset();
695 H_vx_y_diff_vs_trk->Reset();
696 H_vx_z_diff_vs_trk->Reset();
697 H_vx_x_diff_pull_vs_trk->Reset();
698 H_vx_y_diff_pull_vs_trk->Reset();
699 H_vx_z_diff_pull_vs_trk->Reset();
700 H_vx_x_diff_orig_vs_trk->Reset();
701 H_vx_y_diff_orig_vs_trk->Reset();
702 H_vx_z_diff_orig_vs_trk->Reset();
703 H_vx_x_diff_noerr_vs_trk->Reset();
704 H_vx_y_diff_noerr_vs_trk->Reset();
705 H_vx_z_diff_noerr_vs_trk->Reset();
706 H_vx_pur_eff->Reset();
707
708 return true;
709 }
710
711 return false;
712 }
713
714
715 void WriteHistos (string datatype, TFile* filename)
716 {
717 if (datatype.compare("Tk") == 0)
718 {
719 filename->cd("MyProcess/TrackFit");
720
721 for (myMap::iterator it1 = MapFit_d0_pull.begin(); it1 != MapFit_d0_pull.end(); it1++)
722 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
723 for (myMap::iterator it1 = MapFit_d0err.begin(); it1 != MapFit_d0err.end(); it1++)
724 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
725 for (myMap::iterator it1 = MapFit_d0.begin(); it1 != MapFit_d0.end(); it1++)
726 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
727
728 for (myMap::iterator it1 = MapFit_dz_pull.begin(); it1 != MapFit_dz_pull.end(); it1++)
729 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
730 for (myMap::iterator it1 = MapFit_dzerr.begin(); it1 != MapFit_dzerr.end(); it1++)
731 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
732 for (myMap::iterator it1 = MapFit_dz.begin(); it1 != MapFit_dz.end(); it1++)
733 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
734
735 for (myMap::iterator it1 = MapFit_pt_pull.begin(); it1 != MapFit_pt_pull.end(); it1++)
736 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
737 for (myMap::iterator it1 = MapFit_pterr.begin(); it1 != MapFit_pterr.end(); it1++)
738 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
739 for (myMap::iterator it1 = MapFit_pt.begin(); it1 != MapFit_pt.end(); it1++)
740 for (map<string, TH1F*>::iterator it2 = it1->second.begin(); it2 != it1->second.end(); it2++) it2->second->Write(it2->second->GetName(), TObject::kWriteDelete);
741
742 filename->cd("MyProcess/PixelTracks");
743
744 H_pix_trk_doubleCounter_eta_rel->Write(H_pix_trk_doubleCounter_eta_rel->GetName(), TObject::kWriteDelete);
745 H_pix_trk_doubleCounter_pt_rel->Write(H_pix_trk_doubleCounter_pt_rel->GetName(), TObject::kWriteDelete);
746
747 H_trk_trk_doubleCounter_eta_rel->Write(H_trk_trk_doubleCounter_eta_rel->GetName(), TObject::kWriteDelete);
748 H_trk_trk_doubleCounter_pt_rel->Write(H_trk_trk_doubleCounter_pt_rel->GetName(), TObject::kWriteDelete);
749
750 H_trk_eff_pt->Write(H_trk_eff_pt->GetName(), TObject::kWriteDelete);
751 H_trk_eff_eta->Write(H_trk_eff_eta->GetName(), TObject::kWriteDelete);
752 for (int i = 0; i < 3; i++) H_trk_eff_phi[i]->Write(H_trk_eff_phi[i]->GetName(), TObject::kWriteDelete);
753
754 H_trk_toteff_pt->Write(H_trk_eff_pt->GetName(), TObject::kWriteDelete);
755 H_trk_toteff_eta->Write(H_trk_eff_eta->GetName(), TObject::kWriteDelete);
756 for (int i = 0; i < 3; i++) H_trk_toteff_phi[i]->Write(H_trk_toteff_phi[i]->GetName(), TObject::kWriteDelete);
757
758 filename->cd("MyProcess/PixelTracks/CommonPxTkPlots");
759
760 for (int i = 0; i < 2; i++)
761 {
762 H_purity_trked_normchi2_pt_eta[i]->Write(H_purity_trked_normchi2_pt_eta[i]->GetName(), TObject::kWriteDelete);
763 H_purity_trked_pt_eta[i]->Write(H_purity_trked_pt_eta[i]->GetName(), TObject::kWriteDelete);
764
765 H_efficiency_trked_normchi2_pt_eta[i]->Write(H_efficiency_trked_normchi2_pt_eta[i]->GetName(), TObject::kWriteDelete);
766 H_efficiency_trked_pt_eta[i]->Write(H_efficiency_trked_pt_eta[i]->GetName(), TObject::kWriteDelete);
767 H_efficiency_trked_pt_eta_algo[i]->Write(H_efficiency_trked_pt_eta_algo[i]->GetName(), TObject::kWriteDelete);
768
769 H_ptres_pt_eta[i]->Write(H_ptres_pt_eta[i]->GetName(), TObject::kWriteDelete);
770 H_ptpull_pt_eta[i]->Write(H_ptpull_pt_eta[i]->GetName(), TObject::kWriteDelete);
771
772 H_dzres_pt_eta[i]->Write(H_dzres_pt_eta[i]->GetName(), TObject::kWriteDelete);
773 H_dzpull_pt_eta[i]->Write(H_dzpull_pt_eta[i]->GetName(), TObject::kWriteDelete);
774
775 H_d0res_pt_eta[i]->Write(H_d0res_pt_eta[i]->GetName(), TObject::kWriteDelete);
776 H_d0pull_pt_eta[i]->Write(H_d0pull_pt_eta[i]->GetName(), TObject::kWriteDelete);
777
778 H_skew_pt_eta[i]->Write(H_skew_pt_eta[i]->GetName(), TObject::kWriteDelete);
779 }
780
781 for (int i = 0; i < 5; i++)
782 {
783 H_ptres_eta_pt[i]->Write(H_ptres_eta_pt[i]->GetName(), TObject::kWriteDelete);
784 H_ptpull_eta_pt[i]->Write(H_ptpull_eta_pt[i]->GetName(), TObject::kWriteDelete);
785
786 H_dzres_eta_pt[i]->Write(H_dzres_eta_pt[i]->GetName(), TObject::kWriteDelete);
787 H_dzpull_eta_pt[i]->Write(H_dzpull_eta_pt[i]->GetName(), TObject::kWriteDelete);
788
789 H_d0res_eta_pt[i]->Write(H_d0res_eta_pt[i]->GetName(), TObject::kWriteDelete);
790 H_d0pull_eta_pt[i]->Write(H_d0pull_eta_pt[i]->GetName(), TObject::kWriteDelete);
791 }
792
793 for (int i = 0; i < 3; i++)
794 {
795 H_purity_trked_phi[i]->Write(H_purity_trked_phi[i]->GetName(), TObject::kWriteDelete);
796
797 H_efficiency_trked_phi[i]->Write(H_efficiency_trked_phi[i]->GetName(), TObject::kWriteDelete);
798 H_efficiency_trked_phi_algo[i]->Write(H_efficiency_trked_phi_algo[i]->GetName(), TObject::kWriteDelete);
799 }
800
801 H_purity_trked_eta->Write(H_purity_trked_eta->GetName(), TObject::kWriteDelete);
802 H_purity_trked_pt_whole_eta->Write(H_purity_trked_pt_whole_eta->GetName(), TObject::kWriteDelete);
803
804 H_efficiency_trked_eta->Write(H_efficiency_trked_eta->GetName(), TObject::kWriteDelete);
805 H_efficiency_trked_pt_whole_eta->Write(H_efficiency_trked_pt_whole_eta->GetName(), TObject::kWriteDelete);
806
807 H_efficiency_trked_eta_algo->Write(H_efficiency_trked_eta_algo->GetName(), TObject::kWriteDelete);
808 H_efficiency_trked_pt_whole_eta_algo->Write(H_efficiency_trked_pt_whole_eta_algo->GetName(), TObject::kWriteDelete);
809
810 H_tk_pur_eff->Write(H_tk_pur_eff->GetName(), TObject::kWriteDelete);
811 }
812 else if (datatype.compare("Vx") == 0)
813 {
814 filename->cd("MyProcess/VxNtrks");
815
816 for (map<string, TH1F*>::iterator it = map_H_vx_x_diff.begin(); it != map_H_vx_x_diff.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
817 for (map<string, TH1F*>::iterator it = map_H_vx_y_diff.begin(); it != map_H_vx_y_diff.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
818 for (map<string, TH1F*>::iterator it = map_H_vx_z_diff.begin(); it != map_H_vx_z_diff.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
819
820 for (map<string, TH1F*>::iterator it = map_H_vx_x_diff_pull.begin(); it != map_H_vx_x_diff_pull.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
821 for (map<string, TH1F*>::iterator it = map_H_vx_y_diff_pull.begin(); it != map_H_vx_y_diff_pull.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
822 for (map<string, TH1F*>::iterator it = map_H_vx_z_diff_pull.begin(); it != map_H_vx_z_diff_pull.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
823
824 for (map<string, TH1F*>::iterator it = map_H_vx_x_diff_orig.begin(); it != map_H_vx_x_diff_orig.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
825 for (map<string, TH1F*>::iterator it = map_H_vx_y_diff_orig.begin(); it != map_H_vx_y_diff_orig.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
826 for (map<string, TH1F*>::iterator it = map_H_vx_z_diff_orig.begin(); it != map_H_vx_z_diff_orig.end(); it++) it->second->Write(it->second->GetName(), TObject::kWriteDelete);
827
828 filename->cd("MyProcess/PixelVertices/CommonPxTkPlots");
829
830 H_vx_x_diff_vs_trk->Write(H_vx_x_diff_vs_trk->GetName(), TObject::kWriteDelete);
831 H_vx_y_diff_vs_trk->Write(H_vx_y_diff_vs_trk->GetName(), TObject::kWriteDelete);
832 H_vx_z_diff_vs_trk->Write(H_vx_z_diff_vs_trk->GetName(), TObject::kWriteDelete);
833 H_vx_x_diff_pull_vs_trk->Write(H_vx_x_diff_pull_vs_trk->GetName(), TObject::kWriteDelete);
834 H_vx_y_diff_pull_vs_trk->Write(H_vx_y_diff_pull_vs_trk->GetName(), TObject::kWriteDelete);
835 H_vx_z_diff_pull_vs_trk->Write(H_vx_z_diff_pull_vs_trk->GetName(), TObject::kWriteDelete);
836 H_vx_x_diff_orig_vs_trk->Write(H_vx_x_diff_orig_vs_trk->GetName(), TObject::kWriteDelete);
837 H_vx_y_diff_orig_vs_trk->Write(H_vx_y_diff_orig_vs_trk->GetName(), TObject::kWriteDelete);
838 H_vx_z_diff_orig_vs_trk->Write(H_vx_z_diff_orig_vs_trk->GetName(), TObject::kWriteDelete);
839 H_vx_x_diff_noerr_vs_trk->Write(H_vx_x_diff_noerr_vs_trk->GetName(), TObject::kWriteDelete);
840 H_vx_y_diff_noerr_vs_trk->Write(H_vx_y_diff_noerr_vs_trk->GetName(), TObject::kWriteDelete);
841 H_vx_z_diff_noerr_vs_trk->Write(H_vx_z_diff_noerr_vs_trk->GetName(), TObject::kWriteDelete);
842 H_vx_pur_eff->Write(H_vx_pur_eff->GetName(), TObject::kWriteDelete);
843 }
844 }
845
846
847 bool FitAlgorithm (TF1* FitFcn, TH1F* Histo)
848 {
849 double constant;
850 double mean;
851 double sigma;
852
853 constant = Histo->GetMaximum();
854 mean = Histo->GetMean();
855 sigma = Histo->GetRMS();
856
857 FitFcn->SetParameters(constant, mean, sigma);
858 Histo->Fit("Gauss", "QR");
859
860 if ((isnan(FitFcn->GetParameter(2)) == false) && (isnan(FitFcn->GetParError(2)) == false)) return true;
861 return false;
862 }
863
864
865 void FitAndFillMaps (string mapname,
866 double rangemin,
867 double rangemax)
868 {
869 double DeconvolvErrsq;
870 double ErrDeconvolvErrsq;
871 stringstream PtCutStr;
872 stringstream EtaCutStr;
873 PtCutStr.precision(3);
874 EtaCutStr.precision(3);
875 double Val, Err;
876
877 TF1* Gauss = new TF1("Gauss", "gaus", rangemin, rangemax);
878
879 if (mapname.compare("d0pulls") == 0) // @@@ d0pulls MAP @@@
880 {
881 TH1F* d0pulls[5];
882 d0pulls[0] = new TH1F("d0pulls_0", "d0pulls_0", 100, -pullsRange/2., pullsRange/2.);
883 d0pulls[1] = new TH1F("d0pulls_1", "d0pulls_1", 100, -pullsRange/2., pullsRange/2.);
884 d0pulls[2] = new TH1F("d0pulls_2", "d0pulls_2", 100, -pullsRange/2., pullsRange/2.);
885 d0pulls[3] = new TH1F("d0pulls_3", "d0pulls_3", 100, -pullsRange/2., pullsRange/2.);
886 d0pulls[4] = new TH1F("d0pulls_4", "d0pulls_4", 100, -pullsRange/2., pullsRange/2.);
887
888 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
889 {
890 PtCutStr.str("");
891 PtCutStr << PtCut;
892
893 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
894 {
895 EtaCutStr.str("");
896 EtaCutStr << EtaCut;
897
898 if ((MapFit_d0_pull.find(PtCutStr.str().c_str()) != MapFit_d0_pull.end()) &&
899 (MapFit_d0_pull[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_d0_pull[PtCutStr.str().c_str()].end()))
900 {
901 for (int i = 0; i < 2; i++)
902 if ((rint(EtaCut*100.) > rint(EtaThrBin[i]*100.)) && (rint(EtaCut*100.) <= rint(EtaThrBin[i+1]*100.)))
903 d0pulls[i]->Add(MapFit_d0_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
904 }
905 }
906
907 for (int i = 0; i < 2; i++)
908 {
909 if (FitAlgorithm (Gauss,d0pulls[i]) == true)
910 {
911 H_d0pull_pt_eta[i]->SetBinContent(H_d0pull_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Gauss->GetParameter(2));
912 H_d0pull_pt_eta[i]->SetBinError(H_d0pull_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Gauss->GetParError(2));
913 }
914
915 d0pulls[i]->Reset();
916 }
917 }
918
919 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
920 {
921 EtaCutStr.str("");
922 EtaCutStr << EtaCut;
923
924 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
925 {
926 PtCutStr.str("");
927 PtCutStr << PtCut;
928
929 if ((MapFit_d0_pull.find(PtCutStr.str().c_str()) != MapFit_d0_pull.end()) &&
930 (MapFit_d0_pull[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_d0_pull[PtCutStr.str().c_str()].end()))
931 {
932 for (int i = 0; i < 5; i++)
933 if ((rint(PtCut*100.) > rint(PtThr[i]*100.)) && (rint(PtCut*100.) <= rint(PtThr[i+1]*100.)))
934 d0pulls[i]->Add(MapFit_d0_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
935 }
936 }
937
938 for (int i = 0; i < 5; i++)
939 {
940 if (FitAlgorithm (Gauss,d0pulls[i]) == true)
941 {
942 H_d0pull_eta_pt[i]->SetBinContent(H_d0pull_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Gauss->GetParameter(2));
943 H_d0pull_eta_pt[i]->SetBinError(H_d0pull_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Gauss->GetParError(2));
944 }
945
946 d0pulls[i]->Reset();
947 }
948 }
949
950 for (int i = 0; i < 5; i++) delete d0pulls[i];
951 }
952 else if (mapname.compare("d0") == 0) // @@@ d0 MAP @@@
953 {
954 TH1F* d0[5];
955 d0[0] = new TH1F("d0_0", "d0_0", 100, -d0Range/2., d0Range/2.);
956 d0[1] = new TH1F("d0_1", "d0_1", 100, -d0Range/2., d0Range/2.);
957 d0[2] = new TH1F("d0_2", "d0_2", 100, -d0Range/2., d0Range/2.);
958 d0[3] = new TH1F("d0_3", "d0_3", 100, -d0Range/2., d0Range/2.);
959 d0[4] = new TH1F("d0_4", "d0_4", 100, -d0Range/2., d0Range/2.);
960
961 TH1F* d0err[5];
962 d0err[0] = new TH1F("d0err_0", "d0err_0", 200, 0., d0errRange);
963 d0err[1] = new TH1F("d0err_1", "d0err_1", 200, 0., d0errRange);
964 d0err[2] = new TH1F("d0err_2", "d0err_2", 200, 0., d0errRange);
965 d0err[3] = new TH1F("d0err_3", "d0err_3", 200, 0., d0errRange);
966 d0err[4] = new TH1F("d0err_4", "d0err_4", 200, 0., d0errRange);
967
968 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
969 {
970 PtCutStr.str("");
971 PtCutStr << PtCut;
972
973 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
974 {
975 EtaCutStr.str("");
976 EtaCutStr << EtaCut;
977
978 if ((MapFit_d0.find(PtCutStr.str().c_str()) != MapFit_d0.end()) &&
979 (MapFit_d0[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_d0[PtCutStr.str().c_str()].end()) &&
980 (MapFit_d0err.find(PtCutStr.str().c_str()) != MapFit_d0err.end()) &&
981 (MapFit_d0err[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_d0err[PtCutStr.str().c_str()].end()))
982 {
983 for (int i = 0; i < 2; i++)
984 if ((rint(EtaCut*100.) > rint(EtaThrBin[i]*100.)) && (rint(EtaCut*100.) <= rint(EtaThrBin[i+1]*100.)))
985 {
986 d0[i]->Add(MapFit_d0[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
987 d0err[i]->Add(MapFit_d0err[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
988 }
989 }
990 }
991
992 for (int i = 0; i < 2; i++)
993 {
994 DeconvolvErrsq = d0err[i]->GetMean() / 100000000.;
995 ErrDeconvolvErrsq = d0err[i]->GetRMS() / sqrt(d0err[i]->GetEntries()) / 100000000.;
996
997 if ((FitAlgorithm (Gauss,d0[i]) == true) && (Gauss->GetParameter(2)*Gauss->GetParameter(2) > DeconvolvErrsq))
998 {
999 Val = sqrt(Gauss->GetParameter(2)*Gauss->GetParameter(2) - DeconvolvErrsq);
1000 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) + powf(ErrDeconvolvErrsq,2.));
1001
1002 H_d0res_pt_eta[i]->SetBinContent(H_d0res_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Val);
1003 H_d0res_pt_eta[i]->SetBinError(H_d0res_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Err);
1004 }
1005
1006 d0[i]->Reset();
1007 d0err[i]->Reset();
1008 }
1009 }
1010
1011 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1012 {
1013 EtaCutStr.str("");
1014 EtaCutStr << EtaCut;
1015
1016 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1017 {
1018 PtCutStr.str("");
1019 PtCutStr << PtCut;
1020
1021 if ((MapFit_d0.find(PtCutStr.str().c_str()) != MapFit_d0.end()) &&
1022 (MapFit_d0[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_d0[PtCutStr.str().c_str()].end()) &&
1023 (MapFit_d0err.find(PtCutStr.str().c_str()) != MapFit_d0err.end()) &&
1024 (MapFit_d0err[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_d0err[PtCutStr.str().c_str()].end()))
1025 {
1026 for (int i = 0; i < 5; i++)
1027 if ((rint(PtCut*100.) > rint(PtThr[i]*100.)) && (rint(PtCut*100.) <= rint(PtThr[i+1]*100.)))
1028 {
1029 d0[i]->Add(MapFit_d0[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1030 d0err[i]->Add(MapFit_d0err[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1031 }
1032 }
1033 }
1034
1035 for (int i = 0; i < 5; i++)
1036 {
1037 DeconvolvErrsq = d0err[i]->GetMean() / 100000000.;
1038 ErrDeconvolvErrsq = d0err[i]->GetRMS() / sqrt(d0err[i]->GetEntries()) / 100000000.;
1039
1040 if ((FitAlgorithm (Gauss,d0[i]) == true) && (Gauss->GetParameter(2)*Gauss->GetParameter(2) > DeconvolvErrsq))
1041 {
1042 Val = sqrt(Gauss->GetParameter(2)*Gauss->GetParameter(2) - DeconvolvErrsq);
1043 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) + powf(ErrDeconvolvErrsq,2.));
1044
1045 H_d0res_eta_pt[i]->SetBinContent(H_d0res_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Val);
1046 H_d0res_eta_pt[i]->SetBinError(H_d0res_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Err);
1047 }
1048
1049 d0[i]->Reset();
1050 d0err[i]->Reset();
1051 }
1052 }
1053
1054 for (int i = 0; i < 5; i++) { delete d0[i]; delete d0err[i]; }
1055 }
1056 else if (mapname.compare("dzpulls") == 0) // @@@ dzpulls MAP @@@
1057 {
1058 TH1F* dzpulls[5];
1059 dzpulls[0] = new TH1F("dzpulls_0", "dzpulls_0", 100, -pullsRange/2., pullsRange/2.);
1060 dzpulls[1] = new TH1F("dzpulls_1", "dzpulls_1", 100, -pullsRange/2., pullsRange/2.);
1061 dzpulls[2] = new TH1F("dzpulls_2", "dzpulls_2", 100, -pullsRange/2., pullsRange/2.);
1062 dzpulls[3] = new TH1F("dzpulls_3", "dzpulls_3", 100, -pullsRange/2., pullsRange/2.);
1063 dzpulls[4] = new TH1F("dzpulls_4", "dzpulls_4", 100, -pullsRange/2., pullsRange/2.);
1064
1065 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1066 {
1067 PtCutStr.str("");
1068 PtCutStr << PtCut;
1069
1070 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1071 {
1072 EtaCutStr.str("");
1073 EtaCutStr << EtaCut;
1074
1075 if ((MapFit_dz_pull.find(PtCutStr.str().c_str()) != MapFit_dz_pull.end()) &&
1076 (MapFit_dz_pull[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_dz_pull[PtCutStr.str().c_str()].end()))
1077 {
1078 for (int i = 0; i < 2; i++)
1079 if ((rint(EtaCut*100.) > rint(EtaThrBin[i]*100.)) && (rint(EtaCut*100.) <= rint(EtaThrBin[i+1]*100.)))
1080 dzpulls[i]->Add(MapFit_dz_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1081 }
1082 }
1083
1084 for (int i = 0; i < 2; i++)
1085 {
1086 if (FitAlgorithm (Gauss,dzpulls[i]) == true)
1087 {
1088 H_dzpull_pt_eta[i]->SetBinContent(H_dzpull_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Gauss->GetParameter(2));
1089 H_dzpull_pt_eta[i]->SetBinError(H_dzpull_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Gauss->GetParError(2));
1090 }
1091
1092 dzpulls[i]->Reset();
1093 }
1094 }
1095
1096 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1097 {
1098 EtaCutStr.str("");
1099 EtaCutStr << EtaCut;
1100
1101 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1102 {
1103 PtCutStr.str("");
1104 PtCutStr << PtCut;
1105
1106 if ((MapFit_dz_pull.find(PtCutStr.str().c_str()) != MapFit_dz_pull.end()) &&
1107 (MapFit_dz_pull[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_dz_pull[PtCutStr.str().c_str()].end()))
1108 {
1109 for (int i = 0; i < 5; i++)
1110 if ((rint(PtCut*100.) > rint(PtThr[i]*100.)) && (rint(PtCut*100.) <= rint(PtThr[i+1]*100.)))
1111 dzpulls[i]->Add(MapFit_dz_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1112 }
1113 }
1114
1115 for (int i = 0; i < 5; i++)
1116 {
1117 if (FitAlgorithm (Gauss,dzpulls[i]) == true)
1118 {
1119 H_dzpull_eta_pt[i]->SetBinContent(H_dzpull_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Gauss->GetParameter(2));
1120 H_dzpull_eta_pt[i]->SetBinError(H_dzpull_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Gauss->GetParError(2));
1121 }
1122
1123 dzpulls[i]->Reset();
1124 }
1125 }
1126
1127 for (int i = 0; i < 5; i++) delete dzpulls[i];
1128 }
1129 else if (mapname.compare("dz") == 0) // @@@ dz MAP @@@
1130 {
1131 TH1F* dz[5];
1132 dz[0] = new TH1F("dz_0", "dz_0", 100, -dzRange/2., dzRange/2.);
1133 dz[1] = new TH1F("dz_1", "dz_1", 100, -dzRange/2., dzRange/2.);
1134 dz[2] = new TH1F("dz_2", "dz_2", 100, -dzRange/2., dzRange/2.);
1135 dz[3] = new TH1F("dz_3", "dz_3", 100, -dzRange/2., dzRange/2.);
1136 dz[4] = new TH1F("dz_4", "dz_4", 100, -dzRange/2., dzRange/2.);
1137
1138 TH1F* dzerr[5];
1139 dzerr[0] = new TH1F("dzerr_0", "dzerr_0", 200, 0., dzerrRange);
1140 dzerr[1] = new TH1F("dzerr_1", "dzerr_1", 200, 0., dzerrRange);
1141 dzerr[2] = new TH1F("dzerr_2", "dzerr_2", 200, 0., dzerrRange);
1142 dzerr[3] = new TH1F("dzerr_3", "dzerr_3", 200, 0., dzerrRange);
1143 dzerr[4] = new TH1F("dzerr_4", "dzerr_4", 200, 0., dzerrRange);
1144
1145 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1146 {
1147 PtCutStr.str("");
1148 PtCutStr << PtCut;
1149
1150 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1151 {
1152 EtaCutStr.str("");
1153 EtaCutStr << EtaCut;
1154
1155 if ((MapFit_dz.find(PtCutStr.str().c_str()) != MapFit_dz.end()) &&
1156 (MapFit_dz[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_dz[PtCutStr.str().c_str()].end()) &&
1157 (MapFit_dzerr.find(PtCutStr.str().c_str()) != MapFit_dzerr.end()) &&
1158 (MapFit_dzerr[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_dzerr[PtCutStr.str().c_str()].end()))
1159 {
1160 for (int i = 0; i < 2; i++)
1161 if ((rint(EtaCut*100.) > rint(EtaThrBin[i]*100.)) && (rint(EtaCut*100.) <= rint(EtaThrBin[i+1]*100.)))
1162 {
1163 dz[i]->Add(MapFit_dz[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1164 dzerr[i]->Add(MapFit_dzerr[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1165 }
1166 }
1167 }
1168
1169 for (int i = 0; i < 2; i++)
1170 {
1171 DeconvolvErrsq = dzerr[i]->GetMean() / 100000000.;
1172 ErrDeconvolvErrsq = dzerr[i]->GetRMS() / sqrt(dzerr[i]->GetEntries()) / 100000000.;
1173
1174 if ((FitAlgorithm (Gauss,dz[i]) == true) && (Gauss->GetParameter(2)*Gauss->GetParameter(2) > DeconvolvErrsq))
1175 {
1176 Val = sqrt(Gauss->GetParameter(2)*Gauss->GetParameter(2) - DeconvolvErrsq);
1177 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) + powf(ErrDeconvolvErrsq,2.));
1178
1179 H_dzres_pt_eta[i]->SetBinContent(H_dzres_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Val);
1180 H_dzres_pt_eta[i]->SetBinError(H_dzres_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Err);
1181 }
1182
1183 dz[i]->Reset();
1184 dzerr[i]->Reset();
1185 }
1186 }
1187
1188 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1189 {
1190 EtaCutStr.str("");
1191 EtaCutStr << EtaCut;
1192
1193 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1194 {
1195 PtCutStr.str("");
1196 PtCutStr << PtCut;
1197
1198 if ((MapFit_dz.find(PtCutStr.str().c_str()) != MapFit_dz.end()) &&
1199 (MapFit_dz[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_dz[PtCutStr.str().c_str()].end()) &&
1200 (MapFit_dzerr.find(PtCutStr.str().c_str()) != MapFit_dzerr.end()) &&
1201 (MapFit_dzerr[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_dzerr[PtCutStr.str().c_str()].end()))
1202 {
1203 for (int i = 0; i < 5; i++)
1204 if ((rint(PtCut*100.) > rint(PtThr[i]*100.)) && (rint(PtCut*100.) <= rint(PtThr[i+1]*100.)))
1205 {
1206 dz[i]->Add(MapFit_dz[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1207 dzerr[i]->Add(MapFit_dzerr[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1208 }
1209 }
1210 }
1211
1212 for (int i = 0; i < 5; i++)
1213 {
1214 DeconvolvErrsq = dzerr[i]->GetMean() / 100000000.;
1215 ErrDeconvolvErrsq = dzerr[i]->GetRMS() / sqrt(dzerr[i]->GetEntries()) / 100000000.;
1216
1217 if ((FitAlgorithm (Gauss,dz[i]) == true) && (Gauss->GetParameter(2)*Gauss->GetParameter(2) > DeconvolvErrsq))
1218 {
1219 Val = sqrt(Gauss->GetParameter(2)*Gauss->GetParameter(2) - DeconvolvErrsq);
1220 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) + powf(ErrDeconvolvErrsq,2.));
1221
1222 H_dzres_eta_pt[i]->SetBinContent(H_dzres_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Val);
1223 H_dzres_eta_pt[i]->SetBinError(H_dzres_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Err);
1224 }
1225
1226 dz[i]->Reset();
1227 dzerr[i]->Reset();
1228 }
1229 }
1230
1231 for (int i = 0; i < 5; i++) { delete dz[i]; delete dzerr[i]; }
1232 }
1233 else if (mapname.compare("ptpulls") == 0) // @@@ ptpulls MAP @@@
1234 {
1235 TH1F* ptpulls[5];
1236 ptpulls[0] = new TH1F("ptpulls_0", "ptpulls_0", 100, -pullsRange/2., pullsRange/2.);
1237 ptpulls[1] = new TH1F("ptpulls_1", "ptpulls_1", 100, -pullsRange/2., pullsRange/2.);
1238 ptpulls[2] = new TH1F("ptpulls_2", "ptpulls_2", 100, -pullsRange/2., pullsRange/2.);
1239 ptpulls[3] = new TH1F("ptpulls_3", "ptpulls_3", 100, -pullsRange/2., pullsRange/2.);
1240 ptpulls[4] = new TH1F("ptpulls_4", "ptpulls_4", 100, -pullsRange/2., pullsRange/2.);
1241
1242 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1243 {
1244 PtCutStr.str("");
1245 PtCutStr << PtCut;
1246
1247 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1248 {
1249 EtaCutStr.str("");
1250 EtaCutStr << EtaCut;
1251
1252 if ((MapFit_pt_pull.find(PtCutStr.str().c_str()) != MapFit_pt_pull.end()) &&
1253 (MapFit_pt_pull[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_pt_pull[PtCutStr.str().c_str()].end()))
1254 {
1255 for (int i = 0; i < 2; i++)
1256 if ((rint(EtaCut*100.) > rint(EtaThrBin[i]*100.)) && (rint(EtaCut*100.) <= rint(EtaThrBin[i+1]*100.)))
1257 ptpulls[i]->Add(MapFit_pt_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1258 }
1259 }
1260
1261 if ((rint(PtCut*100.) == rint(RangePt*100.)) || (H_ptpull_pt_eta[0]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.) != H_ptpull_pt_eta[0]->GetXaxis()->FindBin(rint(PtCut*100.)/100.)))
1262 {
1263 for (int i = 0; i < 2; i++)
1264 {
1265 if (FitAlgorithm (Gauss,ptpulls[i]) == true)
1266 {
1267 H_ptpull_pt_eta[i]->SetBinContent(H_ptpull_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Gauss->GetParameter(2));
1268 H_ptpull_pt_eta[i]->SetBinError(H_ptpull_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Gauss->GetParError(2));
1269 }
1270
1271 ptpulls[i]->Reset();
1272 }
1273 }
1274 }
1275
1276 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1277 {
1278 EtaCutStr.str("");
1279 EtaCutStr << EtaCut;
1280
1281 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1282 {
1283 PtCutStr.str("");
1284 PtCutStr << PtCut;
1285
1286 if ((MapFit_pt_pull.find(PtCutStr.str().c_str()) != MapFit_pt_pull.end()) &&
1287 (MapFit_pt_pull[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_pt_pull[PtCutStr.str().c_str()].end()))
1288 {
1289 for (int i = 0; i < 5; i++)
1290 if ((rint(PtCut*100.) > rint(PtThr[i]*100.)) && (rint(PtCut*100.) <= rint(PtThr[i+1]*100.)))
1291 ptpulls[i]->Add(MapFit_pt_pull[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1292 }
1293 }
1294
1295 for (int i = 0; i < 5; i++)
1296 {
1297 if (FitAlgorithm (Gauss,ptpulls[i]) == true)
1298 {
1299 H_ptpull_eta_pt[i]->SetBinContent(H_ptpull_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Gauss->GetParameter(2));
1300 H_ptpull_eta_pt[i]->SetBinError(H_ptpull_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Gauss->GetParError(2));
1301 }
1302
1303 ptpulls[i]->Reset();
1304 }
1305 }
1306
1307 for (int i = 0; i < 5; i++) delete ptpulls[i];
1308 }
1309 else if (mapname.compare("pt") == 0) // @@@ pt MAP @@@
1310 {
1311 TH1F* pt[5];
1312 pt[0] = new TH1F("pt_0", "pt_0", 600, -ptRange/2., ptRange/2.);
1313 pt[1] = new TH1F("pt_1", "pt_1", 600, -ptRange/2., ptRange/2.);
1314 pt[2] = new TH1F("pt_2", "pt_2", 600, -ptRange/2., ptRange/2.);
1315 pt[3] = new TH1F("pt_3", "pt_3", 600, -ptRange/2., ptRange/2.);
1316 pt[4] = new TH1F("pt_4", "pt_4", 600, -ptRange/2., ptRange/2.);
1317
1318 TH1F* pterr[5];
1319 pterr[0] = new TH1F("pterr_0", "pterr_0", 5000, 0., pterrRange);
1320 pterr[1] = new TH1F("pterr_1", "pterr_1", 5000, 0., pterrRange);
1321 pterr[2] = new TH1F("pterr_2", "pterr_2", 5000, 0., pterrRange);
1322 pterr[3] = new TH1F("pterr_3", "pterr_3", 5000, 0., pterrRange);
1323 pterr[4] = new TH1F("pterr_4", "pterr_4", 5000, 0., pterrRange);
1324
1325 double OldLowEdgeBin = MinPtTk;
1326
1327 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1328 {
1329 PtCutStr.str("");
1330 PtCutStr << PtCut;
1331
1332 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1333 {
1334 EtaCutStr.str("");
1335 EtaCutStr << EtaCut;
1336
1337 if ((MapFit_pt.find(PtCutStr.str().c_str()) != MapFit_pt.end()) &&
1338 (MapFit_pt[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_pt[PtCutStr.str().c_str()].end()) &&
1339 (MapFit_pterr.find(PtCutStr.str().c_str()) != MapFit_pterr.end()) &&
1340 (MapFit_pterr[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_pterr[PtCutStr.str().c_str()].end()))
1341 {
1342 for (int i = 0; i < 2; i++)
1343 if ((rint(EtaCut*100.) > rint(EtaThrBin[i]*100.)) && (rint(EtaCut*100.) <= rint(EtaThrBin[i+1]*100.)))
1344 {
1345 pt[i]->Add(MapFit_pt[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1346 pterr[i]->Add(MapFit_pterr[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1347 }
1348 }
1349 }
1350
1351 if ((rint(PtCut*100.) == rint(RangePt*100.)) || (H_ptres_pt_eta[0]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.) != H_ptres_pt_eta[0]->GetXaxis()->FindBin(rint(PtCut*100.)/100.)))
1352 {
1353 for (int i = 0; i < 2; i++)
1354 {
1355 DeconvolvErrsq = pterr[i]->GetMean() / 1000000.;
1356 ErrDeconvolvErrsq = pterr[i]->GetRMS() / sqrt(pterr[i]->GetEntries()) / 1000000.;
1357
1358 if ((FitAlgorithm (Gauss,pt[i]) == true) && (Gauss->GetParameter(2)*Gauss->GetParameter(2) > DeconvolvErrsq))
1359 {
1360 Val = sqrt(Gauss->GetParameter(2)*Gauss->GetParameter(2) - DeconvolvErrsq);
1361 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) + powf(ErrDeconvolvErrsq,2.));
1362
1363 H_ptres_pt_eta[i]->SetBinContent(H_ptres_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Val / ((PtCut + OldLowEdgeBin)/2.));
1364 H_ptres_pt_eta[i]->SetBinError(H_ptres_pt_eta[i]->GetXaxis()->FindBin(rint((PtCut-PtStep)*100.)/100.), Err / ((PtCut + OldLowEdgeBin)/2.));
1365
1366 H_skew_pt_eta[i]->Fill(rint((PtCut-PtStep)*100.)/100., pt[i]->GetSkewness());
1367 }
1368
1369 pt[i]->Reset();
1370 pterr[i]->Reset();
1371 }
1372
1373 OldLowEdgeBin = PtCut;
1374 }
1375 }
1376
1377 for (double EtaCut = EtaStep; rint(EtaCut*100.) <= rint(RangeEta*100.); EtaCut += EtaStep)
1378 {
1379 EtaCutStr.str("");
1380 EtaCutStr << EtaCut;
1381
1382 for (double PtCut = MinPtTk+PtStep; rint(PtCut*100.) <= rint(RangePt*100.); PtCut += PtStep)
1383 {
1384 PtCutStr.str("");
1385 PtCutStr << PtCut;
1386
1387 if ((MapFit_pt.find(PtCutStr.str().c_str()) != MapFit_pt.end()) &&
1388 (MapFit_pt[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_pt[PtCutStr.str().c_str()].end()) &&
1389 (MapFit_pterr.find(PtCutStr.str().c_str()) != MapFit_pterr.end()) &&
1390 (MapFit_pterr[PtCutStr.str().c_str()].find(EtaCutStr.str().c_str()) != MapFit_pterr[PtCutStr.str().c_str()].end()))
1391 {
1392 for (int i = 0; i < 5; i++)
1393 if ((rint(PtCut*100.) > rint(PtThr[i]*100.)) && (rint(PtCut*100.) <= rint(PtThr[i+1]*100.)))
1394 {
1395 pt[i]->Add(MapFit_pt[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1396 pterr[i]->Add(MapFit_pterr[PtCutStr.str().c_str()][EtaCutStr.str().c_str()]);
1397 }
1398 }
1399 }
1400
1401 for (int i = 0; i < 5; i++)
1402 {
1403 DeconvolvErrsq = pterr[i]->GetMean() / 1000000.;
1404 ErrDeconvolvErrsq = pterr[i]->GetRMS() / sqrt(pterr[i]->GetEntries()) / 1000000.;
1405
1406 if ((FitAlgorithm (Gauss,pt[i]) == true) && (Gauss->GetParameter(2)*Gauss->GetParameter(2) > DeconvolvErrsq))
1407 {
1408 Val = sqrt(Gauss->GetParameter(2)*Gauss->GetParameter(2) - DeconvolvErrsq);
1409 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) + powf(ErrDeconvolvErrsq,2.));
1410
1411 H_ptres_eta_pt[i]->SetBinContent(H_ptres_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Val / ((PtThr[i+1] + PtThr[i]) / 2.));
1412 H_ptres_eta_pt[i]->SetBinError(H_ptres_eta_pt[i]->GetXaxis()->FindBin(rint((EtaCut-EtaStep)*100.)/100.), Err / ((PtThr[i+1] + PtThr[i]) / 2.));
1413 }
1414
1415 pt[i]->Reset();
1416 pterr[i]->Reset();
1417 }
1418 }
1419
1420 for (int i = 0; i < 5; i++) { delete pt[i]; delete pterr[i]; }
1421 }
1422
1423 delete Gauss;
1424 }
1425
1426
1427 void MakePlots(string datatype)
1428 {
1429 double SumNum;
1430 double SumDen;
1431
1432 if (datatype.compare("Vx") == 0)
1433 {
1434 // #################
1435 // ### VERTEXING ###
1436 // #################
1437
1438 // ########################
1439 // # PURITY & EFFICIENCTY #
1440 // ########################
1441 SumNum = H_vx_counters->GetBinContent(H_vx_counters->GetXaxis()->FindBin("MatchEff"));
1442 SumDen = H_vx_counters->GetBinContent(H_vx_counters->GetXaxis()->FindBin("Strip"));
1443 H_vx_pur_eff->Fill("Efficiency", (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1444 H_vx_pur_eff->SetBinError(H_vx_pur_eff->GetXaxis()->FindBin("Efficiency"), (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1445
1446 SumNum = H_vx_counters->GetBinContent(H_vx_counters->GetXaxis()->FindBin("MatchPur"));
1447 SumDen = H_vx_counters->GetBinContent(H_vx_counters->GetXaxis()->FindBin("Pixel"));
1448 H_vx_pur_eff->Fill("Purity", (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1449 H_vx_pur_eff->SetBinError(H_vx_pur_eff->GetXaxis()->FindBin("Purity"), (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1450
1451 // #######################################
1452 // # Vertex Resolutions and Pulls vs Trk #
1453 // #######################################
1454 stringstream String;
1455 double Val,Err;
1456
1457 TF1* Gauss = new TF1("Gauss", "gaus", -zRange/2., zRange/2.);
1458
1459 for (map<string,TH1F*>::iterator it = map_H_vx_x_diff.begin(); it != map_H_vx_x_diff.end(); it++)
1460 {
1461 String.str("");
1462 String << atoi(it->first.c_str());
1463
1464 Gauss->SetRange(-zRange/2., zRange/2.);
1465
1466 // #############################
1467 // # Compute simple difference #
1468 // #############################
1469 if (map_H_vx_x_diff.find(String.str().c_str()) != map_H_vx_x_diff.end() && (FitAlgorithm (Gauss, map_H_vx_x_diff[String.str().c_str()]) == true))
1470 {
1471 H_vx_x_diff_vs_trk->SetBinContent(H_vx_x_diff_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1472 H_vx_x_diff_vs_trk->SetBinError(H_vx_x_diff_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1473 }
1474 if (map_H_vx_y_diff.find(String.str().c_str()) != map_H_vx_y_diff.end() && (FitAlgorithm (Gauss, map_H_vx_y_diff[String.str().c_str()]) == true))
1475 {
1476 H_vx_y_diff_vs_trk->SetBinContent(H_vx_y_diff_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1477 H_vx_y_diff_vs_trk->SetBinError(H_vx_y_diff_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1478 }
1479 if (map_H_vx_z_diff.find(String.str().c_str()) != map_H_vx_z_diff.end() && (FitAlgorithm (Gauss, map_H_vx_z_diff[String.str().c_str()]) == true))
1480 {
1481 H_vx_z_diff_vs_trk->SetBinContent(H_vx_z_diff_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1482 H_vx_z_diff_vs_trk->SetBinError(H_vx_z_diff_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1483 }
1484
1485 // #################################################
1486 // # Compute difference deconvolving vx resolution #
1487 // #################################################
1488 if ((map_H_vx_x_err2.find(String.str().c_str()) != map_H_vx_x_err2.end()) &&
1489 (map_H_vx_x_diff.find(String.str().c_str()) != map_H_vx_x_diff.end()) &&
1490 (FitAlgorithm (Gauss, map_H_vx_x_diff[String.str().c_str()]) == true) &&
1491 (powf(Gauss->GetParameter(2),2.) > map_H_vx_x_err2[String.str().c_str()]->GetMean()))
1492 {
1493 Val = sqrt(powf(Gauss->GetParameter(2),2.) - map_H_vx_x_err2[String.str().c_str()]->GetMean());
1494
1495 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) +
1496 powf(map_H_vx_x_err2[String.str().c_str()]->GetRMS()/sqrt(map_H_vx_x_err2[String.str().c_str()]->GetEntries()),2.));
1497
1498 H_vx_x_diff_noerr_vs_trk->SetBinContent(H_vx_x_diff_noerr_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Val);
1499 H_vx_x_diff_noerr_vs_trk->SetBinError(H_vx_x_diff_noerr_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Err);
1500 }
1501 if ((map_H_vx_y_err2.find(String.str().c_str()) != map_H_vx_y_err2.end()) &&
1502 (map_H_vx_y_diff.find(String.str().c_str()) != map_H_vx_y_diff.end()) &&
1503 (FitAlgorithm (Gauss, map_H_vx_y_diff[String.str().c_str()]) == true) &&
1504 (powf(Gauss->GetParameter(2),2.) > map_H_vx_y_err2[String.str().c_str()]->GetMean()))
1505 {
1506 Val = sqrt(powf(Gauss->GetParameter(2),2.) - map_H_vx_y_err2[String.str().c_str()]->GetMean());
1507
1508 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) +
1509 powf(map_H_vx_y_err2[String.str().c_str()]->GetRMS()/sqrt(map_H_vx_y_err2[String.str().c_str()]->GetEntries()),2.));
1510
1511 H_vx_y_diff_noerr_vs_trk->SetBinContent(H_vx_y_diff_noerr_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Val);
1512 H_vx_y_diff_noerr_vs_trk->SetBinError(H_vx_y_diff_noerr_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Err);
1513 }
1514 if ((map_H_vx_z_err2.find(String.str().c_str()) != map_H_vx_z_err2.end()) &&
1515 (map_H_vx_z_diff.find(String.str().c_str()) != map_H_vx_z_diff.end()) &&
1516 (FitAlgorithm (Gauss, map_H_vx_z_diff[String.str().c_str()]) == true) &&
1517 (powf(Gauss->GetParameter(2),2.) > map_H_vx_z_err2[String.str().c_str()]->GetMean()))
1518 {
1519 Val = sqrt(powf(Gauss->GetParameter(2),2.) - map_H_vx_z_err2[String.str().c_str()]->GetMean());
1520
1521 Err = 1./(2.*Val) * sqrt(powf(2.*Gauss->GetParameter(2)*Gauss->GetParError(2),2.) +
1522 powf(map_H_vx_z_err2[String.str().c_str()]->GetRMS()/sqrt(map_H_vx_z_err2[String.str().c_str()]->GetEntries()),2.));
1523
1524 H_vx_z_diff_noerr_vs_trk->SetBinContent(H_vx_z_diff_noerr_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Val);
1525 H_vx_z_diff_noerr_vs_trk->SetBinError(H_vx_z_diff_noerr_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Err);
1526 }
1527
1528 Gauss->SetRange(-pullsRange/2., pullsRange/2.);
1529
1530 // #################
1531 // # Compute pulls #
1532 // #################
1533 if (map_H_vx_x_diff_pull.find(String.str().c_str()) != map_H_vx_x_diff_pull.end() && (FitAlgorithm (Gauss, map_H_vx_x_diff_pull[String.str().c_str()]) == true))
1534 {
1535 H_vx_x_diff_pull_vs_trk->SetBinContent(H_vx_x_diff_pull_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1536 H_vx_x_diff_pull_vs_trk->SetBinError(H_vx_x_diff_pull_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1537 }
1538 if (map_H_vx_y_diff_pull.find(String.str().c_str()) != map_H_vx_y_diff_pull.end() && (FitAlgorithm (Gauss, map_H_vx_y_diff_pull[String.str().c_str()]) == true))
1539 {
1540 H_vx_y_diff_pull_vs_trk->SetBinContent(H_vx_y_diff_pull_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1541 H_vx_y_diff_pull_vs_trk->SetBinError(H_vx_y_diff_pull_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1542 }
1543 if (map_H_vx_z_diff_pull.find(String.str().c_str()) != map_H_vx_z_diff_pull.end() && (FitAlgorithm (Gauss, map_H_vx_z_diff_pull[String.str().c_str()]) == true))
1544 {
1545 H_vx_z_diff_pull_vs_trk->SetBinContent(H_vx_z_diff_pull_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1546 H_vx_z_diff_pull_vs_trk->SetBinError(H_vx_z_diff_pull_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1547 }
1548 }
1549
1550 Gauss->SetRange(-zRange/2., zRange/2.);
1551
1552 for (map<string,TH1F*>::iterator it = map_H_vx_x_diff_orig.begin(); it != map_H_vx_x_diff_orig.end(); it++)
1553 {
1554 String.str("");
1555 String << atoi(it->first.c_str());
1556
1557 // ####################################################
1558 // # Compute simple difference with original vertices #
1559 // ####################################################
1560 if (map_H_vx_x_diff_orig.find(String.str().c_str()) != map_H_vx_x_diff_orig.end() && (FitAlgorithm (Gauss, map_H_vx_x_diff_orig[String.str().c_str()]) == true))
1561 {
1562 H_vx_x_diff_orig_vs_trk->SetBinContent(H_vx_x_diff_orig_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1563 H_vx_x_diff_orig_vs_trk->SetBinError(H_vx_x_diff_orig_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1564 }
1565 if (map_H_vx_y_diff_orig.find(String.str().c_str()) != map_H_vx_y_diff_orig.end() && (FitAlgorithm (Gauss, map_H_vx_y_diff_orig[String.str().c_str()]) == true))
1566 {
1567 H_vx_y_diff_orig_vs_trk->SetBinContent(H_vx_y_diff_orig_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1568 H_vx_y_diff_orig_vs_trk->SetBinError(H_vx_y_diff_orig_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1569 }
1570 if (map_H_vx_z_diff_orig.find(String.str().c_str()) != map_H_vx_z_diff_orig.end() && (FitAlgorithm (Gauss, map_H_vx_z_diff_orig[String.str().c_str()]) == true))
1571 {
1572 H_vx_z_diff_orig_vs_trk->SetBinContent(H_vx_z_diff_orig_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParameter(2));
1573 H_vx_z_diff_orig_vs_trk->SetBinError(H_vx_z_diff_orig_vs_trk->GetXaxis()->FindBin(atoi(it->first.c_str())), Gauss->GetParError(2));
1574 }
1575 }
1576
1577 delete Gauss;
1578 }
1579 else if (datatype.compare("Tk") == 0)
1580 {
1581 // ################
1582 // ### TRACKING ###
1583 // ################
1584
1585 double TotalSumNum;
1586 double TotalSumDen;
1587
1588 // ################################
1589 // # 3-PIXEL-HIT EFFICIENCY VS PT #
1590 // ################################
1591 for (int i = 1; i <= H_trk_eff_pt->GetNbinsX(); i++)
1592 {
1593 H_trk_eff_pt->SetBinContent(i, (H_trk_eff_pt_den->GetBinContent(i) > H_trk_eff_pt_num->GetBinContent(i)) ?
1594 H_trk_eff_pt_num->GetBinContent(i) / H_trk_eff_pt_den->GetBinContent(i) : 0.);
1595 H_trk_eff_pt->SetBinError(i, (H_trk_eff_pt_den->GetBinContent(i) > H_trk_eff_pt_num->GetBinContent(i)) ?
1596 sqrt((H_trk_eff_pt_num->GetBinContent(i)*(H_trk_eff_pt_den->GetBinContent(i) - H_trk_eff_pt_num->GetBinContent(i))) /
1597 powf(H_trk_eff_pt_den->GetBinContent(i),3.)) : 0.);
1598 }
1599
1600 // #################################
1601 // # 3-PIXEL-HIT EFFICIENCY VS ETA #
1602 // #################################
1603 for (int i = 1; i <= H_trk_eff_eta->GetNbinsX(); i++)
1604 {
1605 H_trk_eff_eta->SetBinContent(i, (H_trk_eff_eta_den->GetBinContent(i) > H_trk_eff_eta_num->GetBinContent(i)) ?
1606 H_trk_eff_eta_num->GetBinContent(i) / H_trk_eff_eta_den->GetBinContent(i) : 0.);
1607 H_trk_eff_eta->SetBinError(i, (H_trk_eff_eta_den->GetBinContent(i) > H_trk_eff_eta_num->GetBinContent(i)) ?
1608 sqrt((H_trk_eff_eta_num->GetBinContent(i)*(H_trk_eff_eta_den->GetBinContent(i) - H_trk_eff_eta_num->GetBinContent(i))) /
1609 powf(H_trk_eff_eta_den->GetBinContent(i),3.)) : 0.);
1610 }
1611
1612 // #################################
1613 // # 3-PIXEL-HIT EFFICIENCY VS PHI #
1614 // #################################
1615 for (int i = 1; i <= H_trk_eff_phi[0]->GetNbinsX(); i++)
1616 {
1617 for (int j = 0; j < 3; j++)
1618 {
1619 H_trk_eff_phi[j]->SetBinContent(i, (H_trk_eff_phi_den[j]->GetBinContent(i) > H_trk_eff_phi_num[j]->GetBinContent(i)) ?
1620 H_trk_eff_phi_num[j]->GetBinContent(i) / H_trk_eff_phi_den[j]->GetBinContent(i) : 0.);
1621 H_trk_eff_phi[j]->SetBinError(i, (H_trk_eff_phi_den[j]->GetBinContent(i) > H_trk_eff_phi_num[j]->GetBinContent(i)) ?
1622 sqrt((H_trk_eff_phi_num[j]->GetBinContent(i)*(H_trk_eff_phi_den[j]->GetBinContent(i) - H_trk_eff_phi_num[j]->GetBinContent(i))) /
1623 powf(H_trk_eff_phi_den[j]->GetBinContent(i),3.)) : 0.);
1624 }
1625 }
1626
1627 // ########################
1628 // # ALL EFFICIENCY VS PT #
1629 // ########################
1630 for (int i = 1; i <= H_trk_toteff_pt->GetNbinsX(); i++)
1631 {
1632 H_trk_toteff_pt->SetBinContent(i, (H_trk_toteff_pt_den->GetBinContent(i) > H_trk_toteff_pt_num->GetBinContent(i)) ?
1633 H_trk_toteff_pt_num->GetBinContent(i) / H_trk_toteff_pt_den->GetBinContent(i) : 0.);
1634 H_trk_toteff_pt->SetBinError(i, (H_trk_toteff_pt_den->GetBinContent(i) > H_trk_toteff_pt_num->GetBinContent(i)) ?
1635 sqrt((H_trk_toteff_pt_num->GetBinContent(i)*(H_trk_toteff_pt_den->GetBinContent(i) - H_trk_toteff_pt_num->GetBinContent(i))) /
1636 powf(H_trk_toteff_pt_den->GetBinContent(i),3.)) : 0.);
1637 }
1638
1639 // #########################
1640 // # ALL EFFICIENCY VS ETA #
1641 // #########################
1642 for (int i = 1; i <= H_trk_toteff_eta->GetNbinsX(); i++)
1643 {
1644 H_trk_toteff_eta->SetBinContent(i, (H_trk_toteff_eta_den->GetBinContent(i) > H_trk_toteff_eta_num->GetBinContent(i)) ?
1645 H_trk_toteff_eta_num->GetBinContent(i) / H_trk_toteff_eta_den->GetBinContent(i) : 0.);
1646 H_trk_toteff_eta->SetBinError(i, (H_trk_toteff_eta_den->GetBinContent(i) > H_trk_toteff_eta_num->GetBinContent(i)) ?
1647 sqrt((H_trk_toteff_eta_num->GetBinContent(i)*(H_trk_toteff_eta_den->GetBinContent(i) - H_trk_toteff_eta_num->GetBinContent(i))) /
1648 powf(H_trk_toteff_eta_den->GetBinContent(i),3.)) : 0.);
1649 }
1650
1651 // #########################
1652 // # ALL EFFICIENCY VS PHI #
1653 // #########################
1654 for (int i = 1; i <= H_trk_toteff_phi[0]->GetNbinsX(); i++)
1655 {
1656 for (int j = 0; j < 3; j++)
1657 {
1658 H_trk_toteff_phi[j]->SetBinContent(i, (H_trk_toteff_phi_den[j]->GetBinContent(i) > H_trk_toteff_phi_num[j]->GetBinContent(i)) ?
1659 H_trk_toteff_phi_num[j]->GetBinContent(i) / H_trk_toteff_phi_den[j]->GetBinContent(i) : 0.);
1660 H_trk_toteff_phi[j]->SetBinError(i, (H_trk_toteff_phi_den[j]->GetBinContent(i) > H_trk_toteff_phi_num[j]->GetBinContent(i)) ?
1661 sqrt((H_trk_toteff_phi_num[j]->GetBinContent(i)*(H_trk_toteff_phi_den[j]->GetBinContent(i) - H_trk_toteff_phi_num[j]->GetBinContent(i))) /
1662 powf(H_trk_toteff_phi_den[j]->GetBinContent(i),3.)) : 0.);
1663 }
1664 }
1665
1666 // ################################
1667 // # PIXEL-TRACKS DOUBLE COUNTING #
1668 // ################################
1669 for (int i = 1; i <= H_pix_trk_doubleCounter_eta->GetNbinsX(); i++)
1670 {
1671 H_pix_trk_doubleCounter_eta_rel->SetBinContent(i, (H_pix_trk_eta->GetBinContent(i) > H_pix_trk_doubleCounter_eta->GetBinContent(i)) ?
1672 H_pix_trk_doubleCounter_eta->GetBinContent(i) / H_pix_trk_eta->GetBinContent(i) : 0.);
1673 H_pix_trk_doubleCounter_eta_rel->SetBinError(i, (H_pix_trk_eta->GetBinContent(i) > H_pix_trk_doubleCounter_eta->GetBinContent(i)) ?
1674 sqrt((H_pix_trk_doubleCounter_eta->GetBinContent(i)*(H_pix_trk_eta->GetBinContent(i) - H_pix_trk_doubleCounter_eta->GetBinContent(i))) /
1675 powf(H_pix_trk_eta->GetBinContent(i),3.)) : 0.);
1676 }
1677 for (int i = 1; i <= H_pix_trk_doubleCounter_pt->GetNbinsX(); i++)
1678 {
1679 H_pix_trk_doubleCounter_pt_rel->SetBinContent(i, (H_pix_trk_pt->GetBinContent(i) > H_pix_trk_doubleCounter_pt->GetBinContent(i)) ?
1680 H_pix_trk_doubleCounter_pt->GetBinContent(i) / H_pix_trk_pt->GetBinContent(i) : 0.);
1681 H_pix_trk_doubleCounter_pt_rel->SetBinError(i, (H_pix_trk_pt->GetBinContent(i) > H_pix_trk_doubleCounter_pt->GetBinContent(i)) ?
1682 sqrt((H_pix_trk_doubleCounter_pt->GetBinContent(i)*(H_pix_trk_pt->GetBinContent(i) - H_pix_trk_doubleCounter_pt->GetBinContent(i))) /
1683 powf(H_pix_trk_pt->GetBinContent(i),3.)) : 0.);
1684 }
1685
1686 // ##################################
1687 // # TRACKER-TRACKS DOUBLE COUNTING #
1688 // ##################################
1689 for (int i = 1; i <= H_trk_trk_doubleCounter_eta->GetNbinsX(); i++)
1690 {
1691 H_trk_trk_doubleCounter_eta_rel->SetBinContent(i, (H_trk_trk_eta->GetBinContent(i) > H_trk_trk_doubleCounter_eta->GetBinContent(i)) ?
1692 H_trk_trk_doubleCounter_eta->GetBinContent(i) / H_trk_trk_eta->GetBinContent(i) : 0.);
1693 H_trk_trk_doubleCounter_eta_rel->SetBinError(i, (H_trk_trk_eta->GetBinContent(i) > H_trk_trk_doubleCounter_eta->GetBinContent(i)) ?
1694 sqrt((H_trk_trk_doubleCounter_eta->GetBinContent(i)*(H_trk_trk_eta->GetBinContent(i) - H_trk_trk_doubleCounter_eta->GetBinContent(i))) /
1695 powf(H_trk_trk_eta->GetBinContent(i),3.)) : 0.);
1696 }
1697 for (int i = 1; i <= H_trk_trk_doubleCounter_pt->GetNbinsX(); i++)
1698 {
1699 H_trk_trk_doubleCounter_pt_rel->SetBinContent(i, (H_trk_trk_pt->GetBinContent(i) > H_trk_trk_doubleCounter_pt->GetBinContent(i)) ?
1700 H_trk_trk_doubleCounter_pt->GetBinContent(i) / H_trk_trk_pt->GetBinContent(i) : 0.);
1701 H_trk_trk_doubleCounter_pt_rel->SetBinError(i, (H_trk_trk_pt->GetBinContent(i) > H_trk_trk_doubleCounter_pt->GetBinContent(i)) ?
1702 sqrt((H_trk_trk_doubleCounter_pt->GetBinContent(i)*(H_trk_trk_pt->GetBinContent(i) - H_trk_trk_doubleCounter_pt->GetBinContent(i))) /
1703 powf(H_trk_trk_pt->GetBinContent(i),3.)) : 0.);
1704 }
1705
1706 // #####################
1707 // # PURITY vs CHI2-Pt #
1708 // #####################
1709 for (int k = 0; k < 2; k++)
1710 for (int i = 1; i <= H_purity_trked_normchi2_pt_eta[k]->GetNbinsX(); i++)
1711 for (int j = 1; j <= H_purity_trked_normchi2_pt_eta[k]->GetNbinsY(); j++)
1712 if (H_pix_trk_normchi2_pt_eta[k]->GetBinContent(i,j) >= MinEntriesEffPur)
1713 H_purity_trked_normchi2_pt_eta[k]->SetBinContent(i,j, H_trked_normchi2_pt_pur_eta[k]->GetBinContent(i,j) / H_pix_trk_normchi2_pt_eta[k]->GetBinContent(i,j));
1714
1715 // #################
1716 // # PURITY vs ETA #
1717 // #################
1718 TotalSumNum = 0.;
1719 TotalSumDen = 0.;
1720 for (int i = 1; i <= H_pix_trk_eta_phi->GetNbinsX(); i++)
1721 {
1722 SumNum = 0.;
1723 SumDen = 0.;
1724 for (int j = 1; j <= H_pix_trk_eta_phi->GetNbinsY(); j++)
1725 {
1726 SumNum = SumNum + H_trked_eta_phi_pur->GetBinContent(i,j);
1727 SumDen = SumDen + H_pix_trk_eta_phi->GetBinContent(i,j);
1728 }
1729
1730 TotalSumNum = TotalSumNum + SumNum;
1731 TotalSumDen = TotalSumDen + SumDen;
1732
1733 H_purity_trked_eta->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1734 H_purity_trked_eta->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1735 }
1736 H_tk_pur_eff->Fill("Purity", (TotalSumDen >= MinEntriesEffPur) ? TotalSumNum / TotalSumDen : 0.);
1737 H_tk_pur_eff->SetBinError(H_tk_pur_eff->GetXaxis()->FindBin("Purity"), (TotalSumDen >= MinEntriesEffPur) ? sqrt((TotalSumNum*(TotalSumDen-TotalSumNum))/powf(TotalSumDen,3.)) : 0.);
1738
1739 // ################
1740 // # PURITY vs Pt #
1741 // ################
1742 for (int j = H_pix_trk_eta_pt->GetYaxis()->FindBin(rint(MinPtTk*100.)/100.); j <= H_pix_trk_eta_pt->GetNbinsY(); j++)
1743 {
1744 SumNum = 0.;
1745 SumDen = 0.;
1746 for (int i = H_pix_trk_eta_pt->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); i <= H_pix_trk_eta_pt->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); i++)
1747 {
1748 SumNum = SumNum + H_trked_eta_pt_pur->GetBinContent(i,j);
1749 SumDen = SumDen + H_pix_trk_eta_pt->GetBinContent(i,j);
1750 }
1751 H_purity_trked_pt_eta[0]->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1752 H_purity_trked_pt_eta[0]->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1753
1754 SumNum = 0.;
1755 SumDen = 0.;
1756 for (int i = 1; i < H_pix_trk_eta_pt->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); i++)
1757 {
1758 SumNum = SumNum + H_trked_eta_pt_pur->GetBinContent(i,j);
1759 SumDen = SumDen + H_pix_trk_eta_pt->GetBinContent(i,j);
1760 }
1761 for (int i = H_pix_trk_eta_pt->GetXaxis()->FindBin(rint(EtaThr*100.)/100.) + 1; i <= H_pix_trk_eta_pt->GetNbinsX(); i++)
1762 {
1763 SumNum = SumNum + H_trked_eta_pt_pur->GetBinContent(i,j);
1764 SumDen = SumDen + H_pix_trk_eta_pt->GetBinContent(i,j);
1765 }
1766 H_purity_trked_pt_eta[1]->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1767 H_purity_trked_pt_eta[1]->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1768 }
1769
1770 for (int j = H_pix_trk_eta_pt->GetYaxis()->FindBin(rint(MinPtTk*100.)/100.); j <= H_pix_trk_eta_pt->GetNbinsY(); j++)
1771 {
1772 SumNum = 0.;
1773 SumDen = 0.;
1774 for (int i = H_pix_trk_eta_pt->GetXaxis()->FindBin(rint(-MaxEtaTkTrk*100.)/100.); i <= H_pix_trk_eta_pt->GetXaxis()->FindBin(rint(MaxEtaTkTrk*100.)/100.); i++)
1775 {
1776 SumNum = SumNum + H_trked_eta_pt_pur->GetBinContent(i,j);
1777 SumDen = SumDen + H_pix_trk_eta_pt->GetBinContent(i,j);
1778 }
1779 H_purity_trked_pt_whole_eta->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1780 H_purity_trked_pt_whole_eta->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1781 }
1782
1783 // #################
1784 // # PURITY vs PHI #
1785 // #################
1786 for (int i = 1; i <= H_pix_trk_eta_phi->GetNbinsY(); i++)
1787 {
1788 SumNum = 0.;
1789 SumDen = 0.;
1790 for (int j = H_pix_trk_eta_phi->GetXaxis()->FindBin(rint(-MaxEtaTkTrk*100.)/100.); j <= H_pix_trk_eta_phi->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); j++)
1791 {
1792 SumNum = SumNum + H_trked_eta_phi_pur->GetBinContent(j,i);
1793 SumDen = SumDen + H_pix_trk_eta_phi->GetBinContent(j,i);
1794 }
1795 H_purity_trked_phi[0]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1796 H_purity_trked_phi[0]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1797
1798 SumNum = 0.;
1799 SumDen = 0.;
1800 for (int j = H_pix_trk_eta_phi->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); j <= H_pix_trk_eta_phi->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); j++)
1801 {
1802 SumNum = SumNum + H_trked_eta_phi_pur->GetBinContent(j,i);
1803 SumDen = SumDen + H_pix_trk_eta_phi->GetBinContent(j,i);
1804 }
1805 H_purity_trked_phi[1]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1806 H_purity_trked_phi[1]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1807
1808 SumNum = 0.;
1809 SumDen = 0.;
1810 for (int j = H_pix_trk_eta_phi->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); j <= H_pix_trk_eta_phi->GetXaxis()->FindBin(rint(MaxEtaTkTrk*100.)/100.); j++)
1811 {
1812 SumNum = SumNum + H_trked_eta_phi_pur->GetBinContent(j,i);
1813 SumDen = SumDen + H_pix_trk_eta_phi->GetBinContent(j,i);
1814 }
1815 H_purity_trked_phi[2]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1816 H_purity_trked_phi[2]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1817 }
1818
1819 // #########################
1820 // # EFFICIENCY vs CHI2-Pt #
1821 // #########################
1822 for (int k = 0; k < 2; k++)
1823 for (int i = 1; i <= H_efficiency_trked_normchi2_pt_eta[k]->GetNbinsX(); i++)
1824 for (int j = 1; j <= H_efficiency_trked_normchi2_pt_eta[k]->GetNbinsY(); j++)
1825 if (H_trk_trk_normchi2_pt_eta[k]->GetBinContent(i,j) >= MinEntriesEffPur)
1826 H_efficiency_trked_normchi2_pt_eta[k]->SetBinContent(i,j, H_trked_normchi2_pt_eff_eta[k]->GetBinContent(i,j) / H_trk_trk_normchi2_pt_eta[k]->GetBinContent(i,j));
1827
1828 // #####################
1829 // # EFFICIENCY vs ETA #
1830 // #####################
1831 TotalSumNum = 0.;
1832 TotalSumDen = 0.;
1833 for (int i = 1; i <= H_trk_trk_eta_phi->GetNbinsX(); i++)
1834 {
1835 SumNum = 0.;
1836 SumDen = 0.;
1837 for (int j = 1; j <= H_trk_trk_eta_phi->GetNbinsY(); j++)
1838 {
1839 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(i,j);
1840 SumDen = SumDen + H_trk_trk_eta_phi->GetBinContent(i,j);
1841 }
1842
1843 TotalSumNum = TotalSumNum + SumNum;
1844 TotalSumDen = TotalSumDen + SumDen;
1845
1846 H_efficiency_trked_eta->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1847 H_efficiency_trked_eta->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1848 }
1849 H_tk_pur_eff->Fill("Efficiency", (TotalSumDen >= MinEntriesEffPur) ? TotalSumNum / TotalSumDen : 0.);
1850 H_tk_pur_eff->SetBinError(H_tk_pur_eff->GetXaxis()->FindBin("Efficiency"), (TotalSumDen >= MinEntriesEffPur) ? sqrt((TotalSumNum*(TotalSumDen-TotalSumNum))/powf(TotalSumDen,3.)) : 0.);
1851
1852 // ##########################
1853 // # ALGO-EFFICIENCY vs ETA #
1854 // ##########################
1855 TotalSumNum = 0.;
1856 TotalSumDen = 0.;
1857 for (int i = 1; i <= H_trk_trk_eta_phi_algo->GetNbinsX(); i++)
1858 {
1859 SumNum = 0.;
1860 SumDen = 0.;
1861 for (int j = 1; j <= H_trk_trk_eta_phi_algo->GetNbinsY(); j++)
1862 {
1863 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(i,j);
1864 SumDen = SumDen + H_trk_trk_eta_phi_algo->GetBinContent(i,j);
1865 }
1866
1867 TotalSumNum = TotalSumNum + SumNum;
1868 TotalSumDen = TotalSumDen + SumDen;
1869
1870 H_efficiency_trked_eta_algo->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1871 H_efficiency_trked_eta_algo->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1872 }
1873 H_tk_pur_eff->Fill("Algo-Efficiency", (TotalSumDen >= MinEntriesEffPur) ? TotalSumNum / TotalSumDen : 0.);
1874 H_tk_pur_eff->SetBinError(H_tk_pur_eff->GetXaxis()->FindBin("Algo-Efficiency"), (TotalSumDen >= MinEntriesEffPur) ? sqrt((TotalSumNum*(TotalSumDen-TotalSumNum))/powf(TotalSumDen,3.)) : 0.);
1875
1876 // ####################
1877 // # EFFICIENCY vs Pt #
1878 // ####################
1879 for (int j = H_trk_trk_eta_pt->GetYaxis()->FindBin(rint(MinPtTk*100.)/100.); j <= H_trk_trk_eta_pt->GetNbinsY(); j++)
1880 {
1881 SumNum = 0.;
1882 SumDen = 0.;
1883 for (int i = H_trk_trk_eta_pt->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); i <= H_trk_trk_eta_pt->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); i++)
1884 {
1885 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1886 SumDen = SumDen + H_trk_trk_eta_pt->GetBinContent(i,j);
1887 }
1888 H_efficiency_trked_pt_eta[0]->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1889 H_efficiency_trked_pt_eta[0]->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1890
1891 SumNum = 0.;
1892 SumDen = 0.;
1893 for (int i = 1; i < H_trk_trk_eta_pt->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); i++)
1894 {
1895 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1896 SumDen = SumDen + H_trk_trk_eta_pt->GetBinContent(i,j);
1897 }
1898 for (int i = H_trk_trk_eta_pt->GetXaxis()->FindBin(rint(EtaThr*100.)/100.) + 1; i <= H_trk_trk_eta_pt->GetNbinsX(); i++)
1899 {
1900 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1901 SumDen = SumDen + H_trk_trk_eta_pt->GetBinContent(i,j);
1902 }
1903 H_efficiency_trked_pt_eta[1]->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1904 H_efficiency_trked_pt_eta[1]->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1905 }
1906
1907 for (int j = H_trk_trk_eta_pt->GetYaxis()->FindBin(rint(MinPtTk*100.)/100.); j <= H_trk_trk_eta_pt->GetNbinsY(); j++)
1908 {
1909 SumNum = 0.;
1910 SumDen = 0.;
1911 for (int i = H_trk_trk_eta_pt->GetXaxis()->FindBin(rint(-MaxEtaTkTrk*100.)/100.); i <= H_trk_trk_eta_pt->GetXaxis()->FindBin(rint(MaxEtaTkTrk*100.)/100.); i++)
1912 {
1913 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1914 SumDen = SumDen + H_trk_trk_eta_pt->GetBinContent(i,j);
1915 }
1916 H_efficiency_trked_pt_whole_eta->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1917 H_efficiency_trked_pt_whole_eta->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1918 }
1919
1920 // #########################
1921 // # ALGO-EFFICIENCY vs Pt #
1922 // #########################
1923 for (int j = H_trk_trk_eta_pt_algo->GetYaxis()->FindBin(rint(MinPtTk*100.)/100.); j <= H_trk_trk_eta_pt_algo->GetNbinsY(); j++)
1924 {
1925 SumNum = 0.;
1926 SumDen = 0.;
1927 for (int i = H_trk_trk_eta_pt_algo->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); i <= H_trk_trk_eta_pt_algo->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); i++)
1928 {
1929 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1930 SumDen = SumDen + H_trk_trk_eta_pt_algo->GetBinContent(i,j);
1931 }
1932 H_efficiency_trked_pt_eta_algo[0]->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1933 H_efficiency_trked_pt_eta_algo[0]->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1934
1935 SumNum = 0.;
1936 SumDen = 0.;
1937 for (int i = 1; i < H_trk_trk_eta_pt_algo->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); i++)
1938 {
1939 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1940 SumDen = SumDen + H_trk_trk_eta_pt_algo->GetBinContent(i,j);
1941 }
1942 for (int i = H_trk_trk_eta_pt_algo->GetXaxis()->FindBin(rint(EtaThr*100.)/100.) + 1; i <= H_trk_trk_eta_pt_algo->GetNbinsX(); i++)
1943 {
1944 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1945 SumDen = SumDen + H_trk_trk_eta_pt_algo->GetBinContent(i,j);
1946 }
1947 H_efficiency_trked_pt_eta_algo[1]->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1948 H_efficiency_trked_pt_eta_algo[1]->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1949 }
1950
1951 for (int j = H_trk_trk_eta_pt_algo->GetYaxis()->FindBin(rint(MinPtTk*100.)/100.); j <= H_trk_trk_eta_pt_algo->GetNbinsY(); j++)
1952 {
1953 SumNum = 0.;
1954 SumDen = 0.;
1955 for (int i = H_trk_trk_eta_pt_algo->GetXaxis()->FindBin(rint(-MaxEtaTkTrk*100.)/100.); i <= H_trk_trk_eta_pt_algo->GetXaxis()->FindBin(rint(MaxEtaTkTrk*100.)/100.); i++)
1956 {
1957 SumNum = SumNum + H_trked_eta_pt_eff->GetBinContent(i,j);
1958 SumDen = SumDen + H_trk_trk_eta_pt_algo->GetBinContent(i,j);
1959 }
1960 H_efficiency_trked_pt_whole_eta_algo->SetBinContent(j, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1961 H_efficiency_trked_pt_whole_eta_algo->SetBinError(j, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1962 }
1963
1964 // #####################
1965 // # EFFICIENCY vs PHI #
1966 // #####################
1967 for (int i = 1; i <= H_trk_trk_eta_phi->GetNbinsY(); i++)
1968 {
1969 SumNum = 0.;
1970 SumDen = 0.;
1971 for (int j = H_trk_trk_eta_phi->GetXaxis()->FindBin(rint(-MaxEtaTkTrk*100.)/100.); j <= H_trk_trk_eta_phi->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); j++)
1972 {
1973 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(j,i);
1974 SumDen = SumDen + H_trk_trk_eta_phi->GetBinContent(j,i);
1975 }
1976 H_efficiency_trked_phi[0]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1977 H_efficiency_trked_phi[0]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1978
1979 SumNum = 0.;
1980 SumDen = 0.;
1981 for (int j = H_trk_trk_eta_phi->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); j <= H_trk_trk_eta_phi->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); j++)
1982 {
1983 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(j,i);
1984 SumDen = SumDen + H_trk_trk_eta_phi->GetBinContent(j,i);
1985 }
1986 H_efficiency_trked_phi[1]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1987 H_efficiency_trked_phi[1]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1988
1989 SumNum = 0.;
1990 SumDen = 0.;
1991 for (int j = H_trk_trk_eta_phi->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); j <= H_trk_trk_eta_phi->GetXaxis()->FindBin(rint(MaxEtaTkTrk*100.)/100.); j++)
1992 {
1993 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(j,i);
1994 SumDen = SumDen + H_trk_trk_eta_phi->GetBinContent(j,i);
1995 }
1996 H_efficiency_trked_phi[2]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
1997 H_efficiency_trked_phi[2]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
1998 }
1999
2000 // ##########################
2001 // # ALGO-EFFICIENCY vs PHI #
2002 // ##########################
2003 for (int i = 1; i <= H_trk_trk_eta_phi_algo->GetNbinsY(); i++)
2004 {
2005 SumNum = 0.;
2006 SumDen = 0.;
2007 for (int j = H_trk_trk_eta_phi_algo->GetXaxis()->FindBin(rint(-MaxEtaTkTrk*100.)/100.); j <= H_trk_trk_eta_phi_algo->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); j++)
2008 {
2009 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(j,i);
2010 SumDen = SumDen + H_trk_trk_eta_phi_algo->GetBinContent(j,i);
2011 }
2012 H_efficiency_trked_phi_algo[0]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
2013 H_efficiency_trked_phi_algo[0]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
2014
2015 SumNum = 0.;
2016 SumDen = 0.;
2017 for (int j = H_trk_trk_eta_phi_algo->GetXaxis()->FindBin(rint(-EtaThr*100.)/100.); j <= H_trk_trk_eta_phi_algo->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); j++)
2018 {
2019 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(j,i);
2020 SumDen = SumDen + H_trk_trk_eta_phi_algo->GetBinContent(j,i);
2021 }
2022 H_efficiency_trked_phi_algo[1]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
2023 H_efficiency_trked_phi_algo[1]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
2024
2025 SumNum = 0.;
2026 SumDen = 0.;
2027 for (int j = H_trk_trk_eta_phi_algo->GetXaxis()->FindBin(rint(EtaThr*100.)/100.); j <= H_trk_trk_eta_phi_algo->GetXaxis()->FindBin(rint(MaxEtaTkTrk*100.)/100.); j++)
2028 {
2029 SumNum = SumNum + H_trked_eta_phi_eff->GetBinContent(j,i);
2030 SumDen = SumDen + H_trk_trk_eta_phi_algo->GetBinContent(j,i);
2031 }
2032 H_efficiency_trked_phi_algo[2]->SetBinContent(i, (SumDen >= MinEntriesEffPur) ? SumNum / SumDen : 0.);
2033 H_efficiency_trked_phi_algo[2]->SetBinError(i, (SumDen >= MinEntriesEffPur) ? sqrt((SumNum*(SumDen-SumNum))/powf(SumDen,3.)) : 0.);
2034 }
2035
2036 FitAndFillMaps ("d0pulls", -pullsRange/2., pullsRange/2.);
2037 FitAndFillMaps ("d0", -d0Range/2., d0Range/2.);
2038 FitAndFillMaps ("dzpulls", -pullsRange/2., pullsRange/2.);
2039 FitAndFillMaps ("dz", -dzRange/2., dzRange/2.);
2040 FitAndFillMaps ("ptpulls", -pullsRange/2., pullsRange/2.);
2041 FitAndFillMaps ("pt", -ptRange/2., ptRange/2.);
2042 }
2043 }
2044
2045
2046 int main (int argc, char** argv)
2047 {
2048 stringstream DataType;
2049 stringstream FileName;
2050
2051 if (argc == 3)
2052 {
2053 FileName << argv[1];
2054 TFile* File = TFile::Open(FileName.str().c_str(),"UPDATE");
2055
2056 DataType << argv[2];
2057
2058 InitVariables();
2059 if (LoadHistos(DataType.str(), File) == true)
2060 {
2061 MakePlots(DataType.str());
2062 WriteHistos(DataType.str(), File);
2063 File->Close();
2064
2065 return 1;
2066 }
2067 else
2068 {
2069 cout << "Error while loading histograms from: " << FileName.str().c_str() << endl;
2070 return 0;
2071 }
2072 }
2073 else
2074 {
2075 cout << "@@@ Synopsis @@@" << endl;
2076 cout << "./MyPixelAnalysisRatios filename.root datatype[Tk or Vx]" << endl;
2077 return 0;
2078 }
2079 }