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