ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/PixelOnlyTracking/MyPixelAnalysisRatios.cc
Revision: 1.6
Committed: Wed Jun 9 15:22:11 2010 UTC (14 years, 10 months ago) by dinardo
Content type: text/plain
Branch: MAIN
Changes since 1.5: +13 -14 lines
Log Message:
*** empty log message ***

File Contents

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