ATLAS Offline Software
Loading...
Searching...
No Matches
ConfAnalysis.cxx
Go to the documentation of this file.
1
9
10
11#include "ConfAnalysis.h"
12
13
14#include "ConfVtxAnalysis.h"
15
19
20#include <fstream>
21
22#include "BinConfig.h"
23
24#include "TF1.h"
25#include "TMath.h"
26
27#include "globals.h"
28
29bool PRINT_BRESIDUALS = false;
30
31std::ofstream dumpfile("dumpfile.log");
32
33void Normalise(TH1* h) {
34
35 for ( int i=1 ; i<=h->GetNbinsX() ; i++ ) {
36
37 double del = h->GetBinLowEdge(i+1)-h->GetBinLowEdge(i);
38
39 h->SetBinContent( i, h->GetBinContent(i)/del );
40 h->SetBinError( i, h->GetBinError(i)/del );
41 }
42
43}
44
45
46
47
49
55
56
58
60 if ( !m_initialiseFirstEvent ) {
61 std::cout << "ConfAnalysis::initialise() " << std::endl;
63 }
64}
65
67
68 if ( m_initialised ) return;
69
70 m_initialised = true;
71
72 // std::cout << "ConfAnalysis::initialiseInternal() " << name() << std::endl;
73
74 BinConfig& binConfig = g_binConfig;
75
76 if ( name().find("_e")!=std::string::npos ) binConfig = electronBinConfig;
77 else if ( name().find("_mu")!=std::string::npos ) binConfig = muonBinConfig;
78 else if ( name().find("_tau")!=std::string::npos ) binConfig = tauBinConfig;
79 else if ( name().find("_b")!=std::string::npos ) binConfig = bjetBinConfig;
80 else if ( name().find("cosmic")!=std::string::npos ) binConfig = cosmicBinConfig;
81
82
83 //+++ pT ranges
84 // double tmp_maxPt = 50000.;
85 // double tmp_absResPt = 0.0005;
86 //double tmp_maxPt = 50.;
87 double tmp_absResPt = 0.5;
88
89 const int pTResBins = int(100*binConfig.ptres_NScale);
90
91
92 //+++ Eta ranges
95 double tmp_maxEta = 5.;
96 double tmp_absResEta = 0.04; // 0.0005;
97
98 //+++ Phi ranges
99 double tmp_maxPhi = 3.142;
100 double tmp_absResPhi = 0.02; // 0.0001;
101
102
103 // std::cout << "ConfAnalysis::initialise() " << name() << " config: " << binConfig << std::endl;
104
107 int etaBins = int(50*binConfig.eta_NScale);
108 const int etaResBins = int(600*binConfig.eta_NScale);
109
110 const int phiBins = int(30*binConfig.phi_NScale);
111 const int phiResBins = int(100*binConfig.phires_NScale);
112
113 const int zBins = int(150*binConfig.z0_NScale);
114 const double zMax = binConfig.z0Max;
115
116 const int zresBins = 100;
117 const double zresMax = 10;
118
119 const int d0Bins = int(100*binConfig.d0_NScale);
120 const double d0Max = binConfig.d0Max;
121
122 const int d0resBins = 100;
123 const double d0resMax = 5;
124
125 // beamspot corrected position
126 const int a0Bins = int(300*binConfig.a0_NScale);
127 const double a0Max = binConfig.a0Max;
128
129 const int a0resBins = 100;
130 const double a0resMax = 5;
131
132 //+++ Book histograms
133
134 // calculate a logarithmic binning in pt
135
136 int Npt = 0;
137 double pt_a = 1;
138 double pt_b = 1;
139
140 // Npt = int(40*binConfig.pt_NScale);
141 // pt_a = 3.5;
142 Npt = int(45*binConfig.pt_NScale);
143 pt_a = 4;
144 pt_b = 2;
145 // etaBins = 12;
146 // }
147 // else {
148 // Npt = 40;
149 // pt_a = 4;
150 // pt_b = 1;
151 // }
152
153
154
155 const int ptnbins = Npt;
156 std::vector<double> ptbinlimsv(ptnbins+1);
157 double* ptbinlims = &ptbinlimsv[0];
158 // for ( int i=0 ; i<=ptnbins ; i++ ) { ptbinlims[i] = std::pow(10, 2.0*i/ptnbins+2)/1000; }
159 // for ( int i=0 ; i<=ptnbins ; i++ ) { ptbinlims[i] = std::pow(10, 2.3*i/ptnbins+2); }
160 for ( int i=0 ; i<=ptnbins ; i++ ) ptbinlims[i] = std::pow(10, pt_a*i/ptnbins+pt_b)/1000;
161
162
163 // ADDED BY JK - FOR SIGNED PT PLOTS
164 //-----
165 const int ptnbins2 = (2*ptnbins);
166 // std::cout << "ptnbins2 = " << ptnbins2 << std::endl;
167 std::vector<double> ptbinlims2v(ptnbins2 + 1);
168 double* ptbinlims2 = &ptbinlims2v[0];
169 // std::cout << "ptbinlims2v.size() = " << ptbinlims2v.size() << std::endl;
170 int ptnbin_counter = 0;
171 for ( int i=ptnbins; i>0 ; i-- ) {
172 ptbinlims2[ptnbin_counter] = std::pow(10, pt_a*i/ptnbins+pt_b)/(-2000);
173 // std::cout << "ptbinlims[" << i << "] = " << ptbinlims[i] << " , so ptbinlims2[" << ptnbin_counter << "] = " << ptbinlims2[ptnbin_counter] << std::endl;
174 ptnbin_counter++;
175 }
176
177 for ( int i=0 ; i<ptnbins+1 ; i++ ) {
178 ptbinlims2[ptnbin_counter] = std::pow(10, pt_a*i/ptnbins+pt_b)/2000;
179 // std::cout << "ptbinlims[" << i << "] = " << ptbinlims[i] << " , so ptbinlims2[" << ptnbin_counter << "] = " << ptbinlims2[ptnbin_counter] << std::endl;
180 ptnbin_counter++;
181 }
182 //-----
183
184 const int iptnbins = 20;
185 const double minmaxipt=0.5;
186 std::vector<double> iptbinlimsv(iptnbins+1);
187 double* iptbinlims = &iptbinlimsv[0];
188 for ( int i=0 ; i<=iptnbins ; i++ ) {
189 iptbinlims[i] = -minmaxipt+i*minmaxipt*2./iptnbins;
190 }
191
192
193
194 TDirectory* dir = gDirectory;
195
196 // std::cout << "ConfAnalysis::initialize() Directory " << gDirectory->GetName() << " " << name() << std::endl;
197
198 m_dir = new TIDDirectory(name());
199 m_dir->push();
200
201 // std::cout << "ConfAnalysis::initialize() Directory " << gDirectory->GetName() << " " << name() << std::endl;
202
203 if ( !m_initialiseFirstEvent ) std::cout << "ConfAnalysis::initialize() Directory " << name() << std::endl;
204
205 if ( name() != gDirectory->GetName() ) std::cerr << "ConfAnalysis::initialize() Directory: problem with directory " << gDirectory->GetName() << " " << name() << std::endl;
206
207 // TIDDirectory d("histos");
208 // d.push();
209
210 // std::cout << "ConfAnalysis::initialize() Directory " << gDirectory->GetName() << " package directory, " << name() << std::endl;
211
212
213 m_res.push_back( m_rnpix_eta = new Resplot( "npix_eta", 2*etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
214 m_res.push_back( m_rnsct_eta = new Resplot( "nsct_eta", 2*etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
215 m_res.push_back( m_rntrt_eta = new Resplot( "ntrt_eta", 2*etaBins, -tmp_maxEta, tmp_maxEta, 100, -0.5, 99.5 ) );
216 m_res.push_back( m_rnsihit_eta= new Resplot( "nsihit_eta",etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
217
218 m_res.push_back( m_rnpix_lb = new Resplot( "npix_lb", 250, 0, 2500, 22, -0.5, 21.5 ) );
219
220 m_res.push_back( m_rnpix_phi = new Resplot( "npix_phi", etaBins, -M_PI, M_PI, 22, -0.5, 21.5 ) );
221 m_res.push_back( m_rnsct_phi = new Resplot( "nsct_phi", etaBins, -M_PI, M_PI, 22, -0.5, 21.5 ) );
222 m_res.push_back( m_rntrt_phi = new Resplot( "ntrt_phi", etaBins, -M_PI, M_PI, 100, -0.5, 99.5 ) );
223
224 m_res.push_back( m_rnpix_pt = new Resplot( "npix_pt", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
225 m_res.push_back( m_rnsct_pt = new Resplot( "nsct_pt", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
226 m_res.push_back( m_rntrt_pt = new Resplot( "ntrt_pt", ptnbins, ptbinlims, 100, -0.5, 99.5 ) );
227
228 m_res.push_back( m_rnpixh_pt = new Resplot( "npixh_pt", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
229 m_res.push_back( m_rnscth_pt = new Resplot( "nscth_pt", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
230
231
232 m_res.push_back( m_rnpix_pt_bad = new Resplot( "npix_pt_bad", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
233 m_res.push_back( m_rnsct_pt_bad = new Resplot( "nsct_pt_bad", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
234 m_res.push_back( m_rntrt_pt_bad = new Resplot( "ntrt_pt_bad", ptnbins, ptbinlims, 100, -0.5, 99.5 ) );
235
236
237 m_res.push_back( m_rnpix_eta_rec = new Resplot( "npix_eta_rec", 2*etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
238 m_res.push_back( m_rnsct_eta_rec = new Resplot( "nsct_eta_rec", 2*etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
239 m_res.push_back( m_rntrt_eta_rec = new Resplot( "ntrt_eta_rec", 2*etaBins, -tmp_maxEta, tmp_maxEta, 100, -0.5, 99.5 ) );
240 m_res.push_back( m_rnsihit_eta_rec= new Resplot( "nsihit_eta_rec", etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
241
242 m_res.push_back( m_rnpix_phi_rec = new Resplot( "npix_phi_rec", etaBins, -M_PI, M_PI, 22, -0.5, 21.5 ) );
243 m_res.push_back( m_rnsct_phi_rec = new Resplot( "nsct_phi_rec", etaBins, -M_PI, M_PI, 22, -0.5, 21.5 ) );
244 m_res.push_back( m_rntrt_phi_rec = new Resplot( "ntrt_phi_rec", etaBins, -M_PI, M_PI, 100, -0.5, 99.5 ) );
245
246 m_res.push_back( m_rnpix_pt_rec = new Resplot( "npix_pt_rec", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
247 m_res.push_back( m_rnsct_pt_rec = new Resplot( "nsct_pt_rec", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
248 m_res.push_back( m_rntrt_pt_rec = new Resplot( "ntrt_pt_rec", ptnbins, ptbinlims, 100, -0.5, 99.5 ) );
249
250 m_res.push_back( m_rnpixh_pt_rec = new Resplot( "npixh_pt_rec", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
251 m_res.push_back( m_rnscth_pt_rec = new Resplot( "nscth_pt_rec", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
252
253
254 m_res.push_back( m_rnsct_vs_npix = new Resplot( "nsct_vs_npix", 12, -0.5, 11.5, 22, -0.5, 21.5 ) );
255 m_res.push_back( m_rnsct_vs_npix_rec = new Resplot( "nsct_vs_npix_rec", 12, -0.5, 11.5, 22, -0.5, 21.5 ) );
256
257 m_res.push_back( m_rChi2prob = new Resplot( "Chi2prob", ptnbins, ptbinlims, 20, 0, 1 ) );
258 m_res.push_back( m_rChi2 = new Resplot( "Chi2", ptnbins, ptbinlims, 200, 0, 100 ) );
259 m_res.push_back( m_rChi2dof = new Resplot( "Chi2dof", ptnbins, ptbinlims, 100, 0, 10 ) );
260
261 m_res.push_back( m_rChi2prob_bad = new Resplot( "Chi2prob_bad", ptnbins, ptbinlims, 20, 0, 1 ) );
262 m_res.push_back( m_rChi2_bad = new Resplot( "Chi2_bad", ptnbins, ptbinlims, 200, 0, 100 ) );
263 m_res.push_back( m_rChi2dof_bad = new Resplot( "Chi2dof_bad", ptnbins, ptbinlims, 100, 0, 10 ) );
264
265 m_res.push_back( m_rChi2prob_rec = new Resplot( "Chi2prob_rec", ptnbins, ptbinlims, 20, 0, 1 ) );
266 m_res.push_back( m_rChi2_rec = new Resplot( "Chi2_rec", ptnbins, ptbinlims, 200, 0, 100 ) );
267 m_res.push_back( m_rChi2dof_rec = new Resplot( "Chi2dof_rec", ptnbins, ptbinlims, 100, 0, 10 ) );
268
269 m_res.push_back( m_rChi2d_vs_Chi2d = new Resplot( "Chi2d_vs_Chi2d", 200, 0, 100, 200, 0, 100 ) );
270
271 m_res.push_back( m_rDChi2dof = new Resplot( "DChi2dof", 200, 0, 100, 200, -100, 100 ) );
272
274
275 double d0bins[40] = { -5.0, -4.0, -3.0, -2.5,
276 -2.0, -1.8, -1.6, -1.4, -1.2,
277 -1.05, -0.95, -0.85, -0.75, -0.65, -0.55, -0.45, -0.35, -0.25, -0.15, -0.05,
278 0.05, 0.15, 0.25, 0.35, 0.45, 0.55, 0.65, 0.75, 0.85, 0.95, 1.05,
279 1.2, 1.4, 1.6, 1.8, 2.0,
280 2.5, 3.0, 4.0, 5.0 };
281
282
283
284 m_res.push_back( m_rnpix_d0 = new Resplot( "npix_d0", 39, d0bins, 22, -0.5, 21.5 ) );
285 m_res.push_back( m_rnsct_d0 = new Resplot( "nsct_d0", 39, d0bins, 22, -0.5, 21.5 ) );
286 m_res.push_back( m_rntrt_d0 = new Resplot( "ntrt_d0", 39, d0bins, 22, -0.5, 21.5 ) );
287
288 m_res.push_back( m_rnpixh_d0 = new Resplot( "npixh_d0", 39, d0bins, 22, -0.5, 21.5 ) );
289 m_res.push_back( m_rnscth_d0 = new Resplot( "nscth_d0", 39, d0bins, 22, -0.5, 21.5 ) );
290
291 m_res.push_back( m_rnsi_pt = new Resplot( "nsi_pt", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
292 m_res.push_back( m_rnsih_pt = new Resplot( "nsih_pt", ptnbins, ptbinlims, 22, -0.5, 21.5 ) );
293
294 m_res.push_back( m_rnsi_d0 = new Resplot( "nsi_d0", 39, d0bins, 22, -0.5, 21.5 ) );
295 m_res.push_back( m_rnsih_d0 = new Resplot( "nsih_d0", 39, d0bins, 22, -0.5, 21.5 ) );
296
297 m_res.push_back( m_rnsi_eta = new Resplot( "nsi_eta", etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
298 m_res.push_back( m_rnsih_eta = new Resplot( "nsih_eta", etaBins, -tmp_maxEta, tmp_maxEta, 22, -0.5, 21.5 ) );
299
300 m_res.push_back( m_rnbl_d0 = new Resplot( "nbl_d0", 39, d0bins, 5, -0.5, 4.5 ) );
301 m_res.push_back( m_rnblh_d0 = new Resplot( "nblh_d0", 39, d0bins, 5, -0.5, 4.5 ) );
302
303
304 m_res.push_back( m_rnsct_lb = new Resplot( "nsct_lb", 250, 0, 2500, 22, -0.5, 21.5 ) );
305
306 m_res.push_back( m_rnpix_lb_rec = new Resplot( "npix_lb_rec", 250, 0, 2500, 22, -0.5, 21.5 ) );
307 m_res.push_back( m_rnsct_lb_rec = new Resplot( "nsct_lb_rec", 250, 0, 2500, 22, -0.5, 21.5 ) );
308
309 m_res.push_back( m_rnpix_d0_rec = new Resplot( "npix_d0_rec", 39, d0bins, 22, -0.5, 21.5 ) );
310 m_res.push_back( m_rnsct_d0_rec = new Resplot( "nsct_d0_rec", 39, d0bins, 22, -0.5, 21.5 ) );
311 m_res.push_back( m_rntrt_d0_rec = new Resplot( "ntrt_d0_rec", 39, d0bins, 22, -0.5, 21.5 ) );
312
313 // int Nptbins = 7;
314 // double _ptlims[8] = { 0, 500, 1000, 1500, 2000, 5000, 8000, 12000 };
315
316
317 // addHistogram( hchi2=new TH1F("chi2", "chi2", 100, 0, 20) );
318
319 // "reference" quantities
320 addHistogram( new TH1F( "pT", "pT", ptnbins, ptbinlims ) );
321 addHistogram( new TH1F( "eta", "eta", etaBins, -tmp_maxEta, tmp_maxEta ) );
322 addHistogram( new TH1F( "phi", "phi", phiBins, -tmp_maxPhi, tmp_maxPhi ) );
323 addHistogram( new TH1F( "z0", "z0", zBins, -zMax, zMax ) );
324 addHistogram( new TH1F( "d0", "d0", d0Bins, -d0Max, d0Max ) );
325 addHistogram( new TH1F( "a0", "a0", a0Bins, -a0Max, a0Max ) );
326
327 // error study histograms (reference)
328 addHistogram( new TH1F( "dpT", "dpT", 80, 0, 20 ) );
329 addHistogram( new TH1F( "deta", "deta", 50, 0, 1 ) );
330 addHistogram( new TH1F( "dphi", "dphi", 50, 0, 1 ) );
331 addHistogram( new TH1F( "dz0", "dz0", 100, 0, 2 ) );
332 addHistogram( new TH1F( "dd0", "dd0", 50, 0, 0.5 ) );
333 addHistogram( new TH1F( "da0", "da0", 50, 0, 0.5 ) );
334
335 addHistogram( new TH1F( "roi_deta", "roi_deta", 50, -1, 1 ) );
336 addHistogram( new TH1F( "roi_dphi", "roi_dphi", 50, -1, 1 ) );
337 addHistogram( new TH1F( "roi_dR", "roi_dR", 50, 0, 1 ) );
338
339 // tag and probe invariant mass histograms
340 if ( m_TnP_tool ) {
341 m_invmass = new TH1F( "invmass", "invariant mass;mass [GeV]", 320, 0, 200 );
342 m_invmassObj = new TH1F( "invmassObj", "invariant mass;mass [GeV]", 320, 0, 200 );
345 }
346
348 TH2F* heta_vs_pt = new TH2F( "eta_vs_pt", "p_{T} [GeV],#eta", ptnbins, ptbinlims, 10, -3, 3 );
349 TH2F* hd0_vs_pt = new TH2F( "d0_vs_pt", "p_{T} [GeV],d_{0} [mm]", ptnbins, ptbinlims, 30, -300, 300 );
350
351
352 m_eff_eta_vs_pt = new Efficiency2D( heta_vs_pt, "eff_eta_vs_pt" );
353 m_eff_d0_vs_pt = new Efficiency2D( hd0_vs_pt, "eff_d0_vs_pt" );
354
355 delete heta_vs_pt;
356 delete hd0_vs_pt;
357
358
359 // efficiencies and purities
360 m_eff_pt = new Efficiency1D( find("pT"), "pT_eff" );
361 m_eff_pt->Hist()->GetXaxis()->SetTitle("P_{T} [GeV]");
362 m_eff_pt->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
363
364 m_eff_eta = new Efficiency1D( find("eta"), "eta_eff" );
365 m_eff_eta->Hist()->GetXaxis()->SetTitle("#eta");
366 m_eff_eta->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
367
368 m_eff_phi = new Efficiency1D( find("phi"), "phi_eff" );
369 m_eff_phi->Hist()->GetXaxis()->SetTitle("#phi");
370 m_eff_phi->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
371
372 m_eff_z0 = new Efficiency1D( find("z0"), "z0_eff" );
373 m_eff_z0->Hist()->GetXaxis()->SetTitle("z0");
374 m_eff_z0->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
375
376 m_eff_d0 = new Efficiency1D( find("d0"), "d0_eff" );
377 m_eff_d0->Hist()->GetXaxis()->SetTitle("d0");
378 m_eff_d0->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
379
380 m_eff_a0 = new Efficiency1D( find("a0"), "a0_eff" );
381 m_eff_a0->Hist()->GetXaxis()->SetTitle("a0");
382 m_eff_a0->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
383
384 m_eff_ptm = new Efficiency1D( find("pT"), "pTm_eff" );
385 m_eff_ptm->Hist()->GetXaxis()->SetTitle("Negative P_{T} [GeV]");
386 m_eff_ptm->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
387
388 m_eff_ptp = new Efficiency1D( find("pT"), "pTp_eff" );
389 m_eff_ptp->Hist()->GetXaxis()->SetTitle("Positive P_{T} [GeV]");
390 m_eff_ptp->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
391
392 m_eff_roi_deta = new Efficiency1D( find("roi_deta"), "roi_deta_eff" );
393 m_eff_roi_deta->Hist()->GetXaxis()->SetTitle("RoI #Delta#eta");
394 m_eff_roi_deta->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
395
396 m_eff_roi_dphi = new Efficiency1D( find("roi_dphi"), "roi_dphi_eff" );
397 m_eff_roi_dphi->Hist()->GetXaxis()->SetTitle("RoI #Delta#phi");
398 m_eff_roi_dphi->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
399
400 m_eff_roi_dR = new Efficiency1D( find("roi_dR"), "roi_dR_eff" );
401 m_eff_roi_dR->Hist()->GetXaxis()->SetTitle("RoI #Delta R");
402 m_eff_roi_dR->Hist()->GetYaxis()->SetTitle("Efficiency [%]");
403
404 // addHistogram ( m_hDeltaR = new TH1F("DeltaR", "DeltaR", 100, 0, 0.1 ) );
405 addHistogram ( m_hDeltaR = new TH1F("DeltaR", "DeltaR", 100, 0, 0.2 ) );
406
407 m_purity_pt = new Efficiency1D( find("pT"), "pT_pur" );
408 m_purity_eta = new Efficiency1D( find("eta"), "eta_pur" );
409 m_purity_phi = new Efficiency1D( find("phi"), "phi_pur" );
410 m_purity_z0 = new Efficiency1D( find("z0"), "z0_pur" );
411 m_purity_d0 = new Efficiency1D( find("d0"), "d0_pur" );
412 m_purity_a0 = new Efficiency1D( find("a0"), "a0_pur" );
413
414 // "test" quantities
415 addHistogram( new TH1F( "pT_rec", "pT_rec", ptnbins, ptbinlims ) );
416 addHistogram( new TH1F( "eta_rec", "eta_rec", etaBins, -tmp_maxEta, tmp_maxEta ) );
417 addHistogram( new TH1F( "phi_rec", "phi_rec", phiBins, -tmp_maxPhi, tmp_maxPhi ) );
418 addHistogram( new TH1F( "z0_rec", "z0_rec", zBins, -zMax, zMax ) );
419 addHistogram( new TH1F( "d0_rec", "d0_rec", d0Bins, -d0Max, d0Max ) );
420 addHistogram( new TH1F( "a0_rec", "a0_rec", a0Bins, -a0Max, a0Max ) );
421
422 addHistogram2D( new TH2F( "eta_phi_rec", "eta_phi_rec", (tmp_maxEta+1)*30 , -tmp_maxEta-1, tmp_maxEta+1, (tmp_maxPhi+1)*30, -tmp_maxPhi-1, tmp_maxPhi+1 ) );
423 addHistogram2D( new TH2F( "phi_d0_rec", "phi_d0_rec", (2*tmp_maxPhi+2)*15, -tmp_maxPhi-1, tmp_maxPhi+1 ,d0Bins+20, -d0Max+7, d0Max-7 ));
424
425 // error study histograms (test)
426 addHistogram( new TH1F( "dpT_rec", "dpT_rec", 80, 0, 20.00 ) );
427 addHistogram( new TH1F( "deta_rec", "deta_rec", 50, 0, 0.02 ) );
428 addHistogram( new TH1F( "dphi_rec", "dphi_rec", 50, 0, 0.02 ) );
429 addHistogram( new TH1F( "dz0_rec", "dz0_rec", 100, 0, 1.5 ) );
430 addHistogram( new TH1F( "dd0_rec", "dd0_rec", 50, 0, 0.5 ) );
431 addHistogram( new TH1F( "da0_rec", "da0_rec", 50, 0, 0.5 ) );
432
433
434
435 // error study histograms (pull test-ref)
436 addHistogram( new TH1F("pT_pull", "pT_pull", 100, -10, 10) );
437 addHistogram( new TH1F("eta_pull", "eta_pull", 100, -10, 10) );
438 addHistogram( new TH1F("phi_pull", "phi_pull", 100, -10, 10) );
439 addHistogram( new TH1F("z0_pull", "z0_pull", 100, -10, 10) );
440 addHistogram( new TH1F("d0_pull", "d0_pull", 100, -10, 10) );
441 addHistogram( new TH1F("a0_pull", "a0_pull", 100, -10, 10) );
442
443 // error study histograms (pull test-ref) - SIMPLE VERSION
444 addHistogram( new TH1F("pT_pull_simple", "pT_pull_simple", 100, -10, 10) );
445 addHistogram( new TH1F("eta_pull_simple", "eta_pull_simple", 100, -10, 10) );
446 addHistogram( new TH1F("phi_pull_simple", "phi_pull_simple", 100, -10, 10) );
447 addHistogram( new TH1F("z0_pull_simple", "z0_pull_simple", 100, -10, 10) );
448 addHistogram( new TH1F("d0_pull_simple", "d0_pull_simple", 100, -10, 10) );
449 addHistogram( new TH1F("a0_pull_simple", "a0_pull_simple", 100, -10, 10) );
450
451
452 // resolutions
453 TH1F* pT_res = new TH1F( "pT_res", "pT_res", 4*pTResBins, -0.1, 0.1 );
454 pT_res->GetXaxis()->SetTitle("#Delta P_{T} [GeV]");
455 pT_res->GetYaxis()->SetTitle("Entries");
456 addHistogram( pT_res );
457
458
459 TH1F* spT_res = new TH1F( "spT_res", "spT_res", 4*pTResBins, -0.1, 0.1 );
460 spT_res->GetXaxis()->SetTitle("#Delta sP_{T} [GeV]");
461 spT_res->GetYaxis()->SetTitle("Entries");
462 addHistogram( spT_res );
463
464
465 TH1F* ipT_res = new TH1F( "ipT_res", "ipT_res", 4*pTResBins, -0.4, 0.4 );
466 ipT_res->GetXaxis()->SetTitle("#Delta 1/P_{T} [GeV^{-1}]");
467 ipT_res->GetYaxis()->SetTitle("Entries");
468 addHistogram( ipT_res );
469
470 TH1F* eta_res = new TH1F("eta_res", "eta_res", etaResBins, -2*tmp_absResEta, 2*tmp_absResEta );
471 eta_res->GetXaxis()->SetTitle("#Delta #eta");
472 eta_res->GetYaxis()->SetTitle("Entries");
473 addHistogram( eta_res );
474 addHistogram( new TH1F("etai_res", "etai_res", 1000, -0.04, 0.04 ) );
475
476 // TH1F* phi_res = new TH1F( "phi_res", "phi_res;#Delta #phi;Entries", 2*phiResBins, -2*tmp_absResPhi, 2*tmp_absResPhi );
477 // phi_res->GetXaxis()->SetTitle("#Delta #phi");
478 // phi_res->GetYaxis()->SetTitle("Entries");
479 // addHistogram( phi_res );
480
481 addHistogram( new TH1F( "phi_res", "phi_res;#Delta #phi;Entries", 2*phiResBins, -2*tmp_absResPhi, 2*tmp_absResPhi ) );
482 addHistogram( new TH1F( "z0_res", "z0_res;#Deltaz_{0};Entries", 16*zresBins, -8*zresMax, 8*zresMax ) );
483 addHistogram( new TH1F( "d0_res", "d0_res;#Deltad_{0};Entries", 4*d0resBins, -0.2*d0resMax, 0.2*d0resMax ) );
484 addHistogram( new TH1F( "a0_res", "a0_res;#Deltaa_{0};Entries", 4*a0resBins, -0.2*a0resMax, 0.2*a0resMax ) );
485
486
487 addHistogram( new TH1F( "dphi_res", "dphi_res;#Delta #phi;Entries", 2*phiResBins, -0.2*tmp_absResPhi, 0.2*tmp_absResPhi ) );
488 addHistogram( new TH1F( "dz0_res", "dz0_res;#Deltaz_{0};Entries", 8*zresBins, -0.8*zresMax, 0.8*zresMax ) );
489 addHistogram( new TH1F( "dd0_res", "dd0_res;#Deltad_{0};Entries", 4*d0resBins, -0.05*d0resMax, 0.05*d0resMax ) );
490 addHistogram( new TH1F( "da0_res", "da0_res;#Deltaa_{0};Entries", 4*a0resBins, -0.05*a0resMax, 0.05*a0resMax ) );
491
492
493 // std::cout << "booking resplots" << std::endl;
494
495 int factor = 10;
496
497 int wfactor = 2;
498 int zfactor = 16;
499
500
501 m_rDd0res.push_back( new Resplot("rDd0_vs_pt", ptnbins, ptbinlims, 1200, -0.1, 0.1 ) );
502 m_rDd0res.push_back( new Resplot("rDd0_vs_eta", 5*etaBins, -tmp_maxEta, tmp_maxEta, 1200, -0.1, 0.1 ) );
503 m_rDd0res.push_back( new Resplot("rDd0_vs_zed", zBins, -zMax, zMax, 1200, -0.1, 0.1 ) );
504 m_rDd0res.push_back( new Resplot("rDd0_vs_d0", 20, -1.5, 1.5, 1200, -0.1, 0.1 ) );
505 // m_rDd0res.push_back( new Resplot("rDd0_vs_phi", int(2*M_PI/0.1), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02), 1200, -0.01, 0.05 ) );
506 m_rDd0res.push_back( new Resplot("rDd0_vs_phi", 128, -M_PI, M_PI, 1200, -0.1, 0.1 ) );
507
508
509 m_rDa0res.push_back( new Resplot("rDa0_vs_pt", ptnbins, ptbinlims, 100, -0.1, 0.1 ) );
510 m_rDa0res.push_back( new Resplot("rDa0_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 100, -0.1, 0.1 ) );
511 m_rDa0res.push_back( new Resplot("rDa0_vs_zed", 0.2*zBins, -zMax, zMax, 100, -0.1, 0.1 ) );
512 m_rDa0res.push_back( new Resplot("rDa0_vs_da0", 20, -1.5, 1.5, 100, -0.1, 0.1 ) );
513 m_rDa0res.push_back( new Resplot("rDa0_vs_phi", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02), 100, -0.1, 0.1 ) );
514 m_rDa0res.push_back( new Resplot("rDa0_rec_vs_phi", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02), 100, -0.1, 0.1 ) );
515
516 m_rDa0res.push_back( new Resplot("rda0_vs_phi", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02), 100, -0.1, 0.1 ) );
517 m_rDa0res.push_back( new Resplot("rda0_rec_vs_phi", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02), 100, -0.1, 0.1 ) );
518
519 // m_rDd0res.push_back( new Resplot("rDd0_vs_ipt", iptnbins, iptbinlims, 100, -0.1, 0.1 ) );
520 // m_rDz0res.push_back( new Resplot("rDz0_vs_ipt", iptnbins, iptbinlims, 100, -1, 1 ) );
521
522
523 m_rDz0res.push_back( new Resplot("rDz0_vs_pt", ptnbins, ptbinlims, 501, -1, 1 ) );
524 m_rDz0res.push_back( new Resplot("rDz0_vs_eta", 5*etaBins, -tmp_maxEta, tmp_maxEta, 500, -1, 1 ) );
525 m_rDz0res.push_back( new Resplot("rDz0_vs_zed", zBins, -zMax, zMax, 500, -1, 1 ) );
526 m_rDz0res.push_back( new Resplot("rDz0_vs_phi", 128, -M_PI, M_PI, 500, -1, 1 ) );
527
528
529
532 m_retares.push_back( new Resplot("reta_vs_ipt", iptnbins, iptbinlims, 2*etaResBins, -wfactor*tmp_absResEta, wfactor*tmp_absResEta ) );
533 m_rphires.push_back( new Resplot("rphi_vs_ipt", iptnbins, iptbinlims, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
534 m_rzedres.push_back( new Resplot("rzed_vs_ipt", iptnbins, iptbinlims, 8*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
535 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_ipt", iptnbins, iptbinlims, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
536 m_riptres.push_back( new Resplot("ript_vs_ipt", iptnbins, iptbinlims, 16*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
537 m_rptres.push_back( new Resplot("rpt_vs_ipt", iptnbins, iptbinlims, 8*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
538 m_rd0res.push_back( new Resplot("rd0_vs_ipt", iptnbins, iptbinlims, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
539
540
541 m_retares.push_back( new Resplot("reta_vs_pt", ptnbins, ptbinlims, 8*etaResBins, -wfactor*tmp_absResEta, wfactor*tmp_absResEta ) );
542 m_rphires.push_back( new Resplot("rphi_vs_pt", ptnbins, ptbinlims, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
543 m_rzedres.push_back( new Resplot("rzed_vs_pt", ptnbins, ptbinlims, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
544 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_pt", ptnbins, ptbinlims, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
545 m_riptres.push_back( new Resplot("ript_vs_pt", ptnbins, ptbinlims, 16*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
546 m_rptres.push_back( new Resplot("rpt_vs_pt", ptnbins, ptbinlims, 8*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
547 m_rd0res.push_back( new Resplot("rd0_vs_pt", ptnbins, ptbinlims, factor*24*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
548
549 m_retares.push_back( new Resplot("reta_vs_ET", ptnbins, ptbinlims, 8*etaResBins, -wfactor*tmp_absResEta, wfactor*tmp_absResEta ) );
550 m_rphires.push_back( new Resplot("rphi_vs_ET", ptnbins, ptbinlims, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
551 m_rzedres.push_back( new Resplot("rzed_vs_ET", ptnbins, ptbinlims, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
552 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_ET", ptnbins, ptbinlims, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
553 m_riptres.push_back( new Resplot("ript_vs_ET", ptnbins, ptbinlims, 16*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
554 m_rptres.push_back( new Resplot("rpt_vs_ET", ptnbins, ptbinlims, 8*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
555 m_rd0res.push_back( new Resplot("rd0_vs_ET", ptnbins, ptbinlims, factor*24*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
556
557
558 // m_retares.push_back( new Resplot("reta_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 4*etaResBins, -tmp_absResEta, tmp_absResEta ) );
559 m_retares.push_back( new Resplot("reta_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 4*etaResBins, -wfactor*tmp_absResEta, wfactor*tmp_absResEta ) );
560 m_rphires.push_back( new Resplot("rphi_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
561 m_rzedres.push_back( new Resplot("rzed_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 12*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
562
563 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 24*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
564 //m_rzedres.push_back( new Resplot("rzed_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
565 m_riptres.push_back( new Resplot("ript_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 16*pTResBins, -tmp_absResPt, tmp_absResPt ) );
566 m_rptres.push_back( new Resplot("rpt_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 8*pTResBins, -tmp_absResPt, tmp_absResPt ) );
567 m_rd0res.push_back( new Resplot("rd0_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, factor*24*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
568
569
570 // rphivsDd0res = new Resplot( "rphi_vs_Dd0", 10, 0, 0.1, int(2*M_PI/0.02), -0.2*int(M_PI/0.02), 0.2*int(M_PI/0.02) );
571 m_hphivsDd0res[0] = new TH1F( "hphi_vs_Dd0_0", "phi for Dd0<0.1", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02) );
572 m_hphivsDd0res[1] = new TH1F( "hphi_vs_Dd0_1", "phi for Dd0>0.1", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02) );
573 m_hphivsDd0res[2] = new TH1F( "hphi_vs_Dd0_2", "all phi", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02) );
574
575 m_hphivsDa0res[0] = new TH1F( "hphi_vs_Da0_0", "phi for Da0<0.1", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02) );
576 m_hphivsDa0res[1] = new TH1F( "hphi_vs_Da0_1", "phi for Da0>0.1", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02) );
577 m_hphivsDa0res[2] = new TH1F( "hphi_vs_Da0_2", "all phi", int(2*M_PI/0.02), -0.02*int(M_PI/0.02), 0.02*int(M_PI/0.02) );
578
579
580 for ( unsigned ih=0 ; ih<3 ; ih++ ) {
581 m_hphivsDd0res[ih]->SetDirectory(0);
582 m_hphivsDa0res[ih]->SetDirectory(0);
583 }
584
585 m_retares.push_back( new Resplot("reta_vs_zed", 0.2*zBins, -zMax, zMax, 2*etaResBins, -tmp_absResEta, tmp_absResEta ) );
586 m_rphires.push_back( new Resplot("rphi_vs_zed", 0.2*zBins, -zMax, zMax, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
587 m_rzedres.push_back( new Resplot("rzed_vs_zed", 0.2*zBins, -zMax, zMax, 8*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
588 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_zed", 0.2*zBins, -zMax, zMax, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
589 //m_rzedres.push_back( new Resplot("rzed_vs_zed", 0.2*zBins, -zMax, zMax, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
590 //m_rzedres.push_back( new Resplot("rzed_vs_zed", zBins, -zMax, zMax, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
591 m_riptres.push_back( new Resplot("ript_vs_zed", 0.2*zBins, -zMax, zMax, 2*pTResBins, -2*tmp_absResPt, 2*tmp_absResPt ) );
592 m_rptres.push_back( new Resplot("rpt_vs_zed", 0.2*zBins, -zMax, zMax, 8*pTResBins, -tmp_absResPt, tmp_absResPt ) );
593 m_rd0res.push_back( new Resplot("rd0_vs_zed", 0.2*zBins, -zMax, zMax, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
594
595
596 m_retares.push_back( new Resplot("reta_vs_nvtx", 24, 0, 72, 4*etaResBins, -tmp_absResEta, tmp_absResEta ) );
597 m_rphires.push_back( new Resplot("rphi_vs_nvtx", 24, 0, 72, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
598 m_rzedres.push_back( new Resplot("rzed_vs_nvtx", 24, 0, 72, 4*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
599 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_nvtx", 24, 0, 72, 24*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
600 m_riptres.push_back( new Resplot("ript_vs_nvtx", 24, 0, 72, 4*pTResBins, -tmp_absResPt, tmp_absResPt ) );
601 m_rptres.push_back( new Resplot("rpt_vs_nvtx", 24, 0, 72, 8*pTResBins, -tmp_absResPt, tmp_absResPt ) );
602 m_rd0res.push_back( new Resplot("rd0_vs_nvtx", 24, 0, 72, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
603
604
605 m_retares.push_back( new Resplot("reta_vs_ntracks", 60, 0, 600, 4*etaResBins, -tmp_absResEta, tmp_absResEta ) );
606 m_rphires.push_back( new Resplot("rphi_vs_ntracks", 60, 0, 600, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
607 m_rzedres.push_back( new Resplot("rzed_vs_ntracks", 60, 0, 600, 4*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
608 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_ntracks", 60, 0, 600, 24*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
609 //m_rzedres.push_back( new Resplot("rzed_vs_ntracks", 60, 0, 600, 4*zfactor*zresBins, -zfactor*0.5*zresMax, zfactor*0.5*zresMax ) );
610 m_riptres.push_back( new Resplot("ript_vs_ntracks", 60, 0, 600, 4*pTResBins, -tmp_absResPt, tmp_absResPt ) );
611 m_rptres.push_back( new Resplot("rpt_vs_ntracks", 60, 0, 600, 8*pTResBins, -tmp_absResPt, tmp_absResPt ) );
612 m_rd0res.push_back( new Resplot("rd0_vs_ntracks", 60, 0, 600, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
613
614
615 m_retares.push_back( new Resplot("reta_vs_phi", 128, -M_PI, M_PI, 2*etaResBins, -wfactor*tmp_absResEta, wfactor*tmp_absResEta ) );
616 m_rphires.push_back( new Resplot("rphi_vs_phi", 128, -M_PI, M_PI, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
617 m_rzedres.push_back( new Resplot("rzed_vs_phi", 128, -M_PI, M_PI, 8*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
618 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_phi", 128, -M_PI, M_PI, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
619 m_riptres.push_back( new Resplot("ript_vs_phi", 128, -M_PI, M_PI, 16*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
620 m_rptres.push_back( new Resplot("rpt_vs_phi", 128, -M_PI, M_PI, 8*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
621 m_rd0res.push_back( new Resplot("rd0_vs_phi", 128, -M_PI, M_PI, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
622
623
624
625 m_retares.push_back( new Resplot("reta_vs_mu", 24, 0, 72, 4*etaResBins, -tmp_absResEta, tmp_absResEta ) );
626 m_rphires.push_back( new Resplot("rphi_vs_mu", 24, 0, 72, 8*phiResBins, -wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
627 m_rzedres.push_back( new Resplot("rzed_vs_mu", 24, 0, 72, 4*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
628 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_mu", 24, 0, 72, 24*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
629 m_riptres.push_back( new Resplot("ript_vs_mu", 24, 0, 72, 4*pTResBins, -tmp_absResPt, tmp_absResPt ) );
630 m_rptres.push_back( new Resplot("rpt_vs_mu", 24, 0, 72, 8*pTResBins, -tmp_absResPt, tmp_absResPt ) );
631 m_rd0res.push_back( new Resplot("rd0_vs_mu", 24, 0, 72, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
632
633
634 //ADDED BY JK
635 //-----
636 m_rzedres.push_back( new Resplot("rzed_vs_signed_pt", ptnbins2, ptbinlims2, 4*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
637 m_rzedres.push_back( new Resplot("rzed_vs_ABS_pt", ptnbins, ptbinlims, 4*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
638 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_signed_pt", ptnbins2, ptbinlims2, 24*zfactor*zresBins, -zfactor*zresMax, zfactor*zresMax ) );
639 m_rzedthetares.push_back( new Resplot("rzedtheta_vs_ABS_pt", ptnbins, ptbinlims, 24*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax ) );
640 //m_rzedres.push_back( new Resplot("rzed_vs_signed_pt", ptnbins2, ptbinlims2, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
641 //m_rzedres.push_back( new Resplot("rzed_vs_ABS_pt", ptnbins, ptbinlims, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
642 m_rd0res.push_back( new Resplot("rd0_vs_signed_pt", ptnbins2, ptbinlims2, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
643 m_rd0res.push_back( new Resplot("rd0_vs_ABS_pt", ptnbins, ptbinlims, factor*8*a0resBins, -wfactor*a0resMax, wfactor*a0resMax ) );
644 //-----
645
646
647 // std::cout << "booked" << std::endl;
648
649 m_retaresPull.push_back( new Resplot("retaPull_vs_ipt", iptnbins, iptbinlims, 2*etaResBins, -5,5));//-wfactor*tmp_absResEta, wfactor*tmp_absResEta ) );
650 m_rphiresPull.push_back( new Resplot("rphiPull_vs_ipt", iptnbins, iptbinlims, 8*phiResBins, -5,5));//-wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
651 m_rzedresPull.push_back( new Resplot("rzedPull_vs_ipt", iptnbins, iptbinlims, 4*zfactor*zresBins, -5,5));//-2*zfactor*zresMax, 2*zfactor*zresMax ) );
652 //m_riptresPull.push_back( new Resplot("ript_vs_ipt", iptnbins, iptbinlims, 16*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
653 m_rptresPull.push_back( new Resplot("rptPull_vs_ipt", iptnbins, iptbinlims, 8*pTResBins, -5,5));//-wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
654 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_ipt", iptnbins, iptbinlims, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ) ) ;
655
656 m_retaresPull.push_back( new Resplot("retaPull_vs_pt", ptnbins, ptbinlims, 2*etaResBins, -5,5));//-wfactor*tmp_absResEta, wfactor*tmp_absResEta ) );
657 m_rphiresPull.push_back( new Resplot("rphiPull_vs_pt", ptnbins, ptbinlims, 8*phiResBins, -5,5));//-wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
658 m_rzedresPull.push_back( new Resplot("rzedPull_vs_pt", ptnbins, ptbinlims, 4*zfactor*zresBins, -5,5));//-2*zfactor*zresMax, 2*zfactor*zresMax ) );
659 //rzedres.push_back( new Resplot("rzed_vs_pt", ptnbins, ptbinlims, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
660 //m_riptresPull.push_back( new Resplot("ript_vs_pt", ptnbins, ptbinlims, 16*pTResBins, -wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
661 m_rptresPull.push_back( new Resplot("rptPull_vs_pt", ptnbins, ptbinlims, 8*pTResBins, -5,5));//-wfactor*tmp_absResPt, wfactor*tmp_absResPt ) );
662 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_pt", ptnbins, ptbinlims, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ) );
663
664 m_retaresPull.push_back( new Resplot("retaPull_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 4*etaResBins, -5,5));//-tmp_absResEta, tmp_absResEta ) );
665 m_rphiresPull.push_back( new Resplot("rphiPull_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 8*phiResBins, -5,5));//-wfactor*tmp_absResPhi, wfactor*tmp_absResPhi) );
666 m_rzedresPull.push_back( new Resplot("rzedPull_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 4*zfactor*zresBins, -5,5));//-2*zfactor*zresMax, 2*zfactor*zresMax) );
667 //rzedres.push_back( new Resplot("rzed_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax) );
668 //riptresPull.push_back( new Resplot("ript_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 16*pTResBins, -tmp_absResPt, tmp_absResPt ) );
669 m_rptresPull.push_back( new Resplot("rptPull_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, 8*pTResBins, -5,5));//-tmp_absResPt, tmp_absResPt ) );
670 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_eta", etaBins, -tmp_maxEta, tmp_maxEta, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ));
671 m_retaresPull.push_back( new Resplot("retaPull_vs_zed", 0.2*zBins, -zMax, zMax, 2*etaResBins, -5,5));//-tmp_absResEta, tmp_absResEta ) );
672 m_rphiresPull.push_back( new Resplot("rphiPull_vs_zed", 0.2*zBins, -zMax, zMax, 8*phiResBins, -5,5));//-wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
673 m_rzedresPull.push_back( new Resplot("rzedPull_vs_zed", 0.2*zBins, -zMax, zMax, 4*zfactor*zresBins, -5,5));//-2*zfactor*zresMax, 2*zfactor*zresMax ) ) ;
674 //rzedres.push_back( new Resplot("rzed_vs_zed", 0.2*zBins, -zMax, zMax, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
675 //rzedres.push_back( new Resplot("rzed_vs_zed", zBins, -zMax, zMax, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
676 //riptresPull.push_back( new Resplot("ript_vs_zed", 0.2*zBins, -zMax, zMax, 2*pTResBins, -2*tmp_absResPt, 2*tmp_absResPt ) );
677 m_rptresPull.push_back( new Resplot("rptPull_vs_zed", 0.2*zBins, -zMax, zMax, 8*pTResBins, -5,5));//-tmp_absResPt, tmp_absResPt ) );
678 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_zed", 0.2*zBins, -zMax, zMax, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ) );
679
680 m_retaresPull.push_back( new Resplot("retaPull_vs_nvtx", 12, 0, 36, 4*etaResBins, -5,5));//-tmp_absResEta, tmp_absResEta ) );
681 m_rphiresPull.push_back( new Resplot("rphiPull_vs_nvtx", 12, 0, 36, 8*phiResBins, -5,5));//-wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
682 m_rzedresPull.push_back( new Resplot("rzedPull_vs_nvtx", 12, 0, 36, 4*zfactor*zresBins, -5,5));//-zfactor*zresMax, zfactor*zresMax ) );
683 //rzedres.push_back( new Resplot("rzed_vs_nvtx", 12, 0, 36, 4*zfactor*zresBins, -zfactor*0.5*zresMax, zfactor*0.5*zresMax ) );
684
685 //riptresPull.push_back( new Resplot("ript_vs_nvtx", 12, 0, 36, 4*pTResBins, -tmp_absResPt, tmp_absResPt ) );
686 m_rptresPull.push_back( new Resplot("rptPull_vs_nvtx", 12, 0, 36, 8*pTResBins, -5,5));//-tmp_absResPt, tmp_absResPt ) );
687 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_nvtx", 12, 0, 36, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ) );
688
689 m_retaresPull.push_back( new Resplot("retaPull_vs_ntracks", 60, 0, 600, 4*etaResBins, -5,5));//-tmp_absResEta, tmp_absResEta ) );
690 m_rphiresPull.push_back( new Resplot("rphiPull_vs_ntracks", 60, 0, 600, 8*phiResBins, -5,5));//-wfactor*tmp_absResPhi, wfactor*tmp_absResPhi ) );
691 m_rzedresPull.push_back( new Resplot("rzedPull_vs_ntracks", 60, 0, 600, 4*zfactor*zresBins, -5,5));//-zfactor*zresMax, zfactor*zresMax ) );
692 //rzedres.push_back( new Resplot("rzed_vs_ntracks", 60, 0, 600, 4*zfactor*zresBins, -zfactor*0.5*zresMax, zfactor*0.5*zresMax ) );
693 //riptresPull.push_back( new Resplot("ript_vs_ntracks", 60, 0, 600, 4*pTResBins, -tmp_absResPt, tmp_absResPt ) );
694 m_rptresPull.push_back( new Resplot("rptPull_vs_ntracks", 60, 0, 600, 8*pTResBins, -5,5));//-tmp_absResPt, tmp_absResPt ) );
695 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_ntracks", 60, 0, 600, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ) );
696
697 m_rzedresPull.push_back( new Resplot("rzedPull_vs_signed_pt", ptnbins2, ptbinlims2, 4*zfactor*zresBins, -5,5));//-2*zfactor*zresMax, 2*zfactor*zresMax ) );
698 m_rzedresPull.push_back( new Resplot("rzedPull_vs_ABS_pt", ptnbins, ptbinlims, 4*zfactor*zresBins, -5,5));//-2*zfactor*zresMax, 2*zfactor*zresMax ));
699 //rzedres.push_back( new Resplot("rzed_vs_signed_pt", ptnbins2, ptbinlims2, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ));
700 //rzedres.push_back( new Resplot("rzed_vs_ABS_pt", ptnbins, ptbinlims, 4*zfactor*zresBins, -2*zwidthfactor*zresMax, 2*zwidthfactor*zresMax ) );
701 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_signed_pt", ptnbins2, ptbinlims2, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ) );
702 m_rd0resPull.push_back( new Resplot("rd0Pull_vs_ABS_pt", ptnbins, ptbinlims, factor*8*a0resBins, -5,5));//-wfactor*a0resMax, wfactor*a0resMax ) );
703
704 m_rzedreslb = new Resplot("rzed_vs_lb", 301, -0.5, 3009.5, 8*zfactor*zresBins, -2*zfactor*zresMax, 2*zfactor*zresMax );
705
706 m_rzedlb = new Resplot("zed_vs_lb", 301, -0.5, 3009.5, 100, -300, 300 );
707 m_rzedlb_rec = new Resplot("zed_vs_lb_rec", 301, -0.5, 3009.5, 100, -300, 300 );
708
709
710 m_rd0_vs_phi = new Resplot( "d0_vs_phi", 20, -M_PI, M_PI, 200, -10, 10 );
711 m_rd0_vs_phi_rec = new Resplot( "d0_vs_phi_rec", 20, -M_PI, M_PI, 200, -10, 10 );
712
713
715
716 m_rRoi_deta_vs_eta = new Resplot( "rRoi_deta_vs_eta", 20, -2.5, 2.5, 101, -0.5, 0.5 );
717 m_rRoi_dphi_vs_eta = new Resplot( "rRoi_dphi_vs_eta", 20, -2.5, 2.5, 101, -0.5, 0.5 );
718 m_rRoi_dzed_vs_eta = new Resplot( "rRoi_dzed_vs_eta", 20, -2.5, 2.5, 401, -30.0, 30.0 );
719
720 // std::cout << "ROI RESPLOTS " << m_rRoi_deta_vs_eta->Name() << " " << m_rRoi_dphi_vs_eta->Name() << std::endl;
721
722
723 // hit occupancies
724
725 int NHits = 40;
726 int Ntracks = 10000;
727
728 addHistogram( new TH1F( "nsct", "nsct", NHits, -0.5, float(NHits-0.5) ) );
729 addHistogram( new TH1F( "nsct_rec", "nsct_rec", NHits, -0.5, float(NHits-0.5) ) );
730
731 addHistogram( new TH1F( "npix", "npix", NHits, -0.5, float(NHits-0.5) ) );
732 addHistogram( new TH1F( "npix_rec", "npix_rec", NHits, -0.5, float(NHits-0.5) ) );
733
734 addHistogram( new TH1F( "nsi", "nsi", NHits, -0.5, float(NHits-0.5) ) );
735 addHistogram( new TH1F( "nsi_rec", "nsi_rec", NHits, -0.5, float(NHits-0.5) ) );
736 addHistogram( new TH1F( "nsi_matched", "nsi_matched", NHits, -0.5, float(NHits-0.5) ) );
737
738
739 addHistogram( new TH1F( "ntrt", "ntrt", NHits, -0.5, float(NHits-0.5) ) );
740 addHistogram( new TH1F( "ntrt_rec", "ntrt_rec", NHits, -0.5, float(NHits-0.5) ) );
741
742 addHistogram( new TH1F( "nstraw", "nstraw", NHits*4, -0.5, float(4*NHits-0.5) ) );
743 addHistogram( new TH1F( "nstraw_rec", "nstraw_rec", NHits*4, -0.5, float(4*NHits-0.5) ) );
744
745 addHistogram( new TH1F( "ntracks", "ntracks", Ntracks+1, -0.5, float(Ntracks+0.5) ) );
746 addHistogram( new TH1F( "ntracks_rec", "ntracks_rec", Ntracks+1, -0.5, float(Ntracks+0.5) ) );
747
748
749 // beam offset fitting histos
750 m_h2 = new Resplot( "d0vphi", phiBins, -3.142, 3.142, d0Bins, -d0Max, d0Max );
751 m_h2r = new Resplot( "d0vphi_rec", phiBins, -3.142, 3.142, d0Bins, -d0Max, d0Max );
752 m_h2m = new Resplot( "d0vphi_match", phiBins, -3.142, 3.142, d0Bins, -d0Max, d0Max );
753
754 m_h2a0 = new Resplot( "a0vphi", phiBins, -3.142, 3.142, d0Bins, -d0Max, d0Max );
755 m_h2a0r = new Resplot( "a0vphi_rec", phiBins, -3.142, 3.142, d0Bins, -d0Max, d0Max );
756
758
759
760 TH1F heffvlb("eff vs lb", "eff vs lb", 301, -0.5, 3009.5 );
761
762 // 100,
763 // 1270515000, 1270560000
764 // 1272040000, 1272200000
765 // );
766 // 1270518944,
767 // 1270558762 );
768
769 // 1260470000,
770 // 1260680000 );
771 // 1260470000,
772 // 1260690000 );
773 // TH1F heffvlb("eff vs lb", "eff vs lb", 600, 0, 3000);
774
775 // TH1F heffvlb("eff vs lb", "eff vs lb", 600, 1260400000, 1260700000);
776
777
778 m_eff_vs_lb = new Efficiency1D( &heffvlb );
779
780 // m_z_vs_lb = new Resplot("z vs lb", 100, 1270515000, 1270560000, 100, -250, 250);
781 m_z_vs_lb = new Resplot("z vs lb", 301, -0.5, 3009.5, 100, -250, 250);
782
783 //m_rmap[142165] = 0;
784 //m_rmap[142166] = 500;
785 //m_rmap[142174] = 1000;
786 //m_rmap[142191] = 1500;
787 //m_rmap[142193] = 2000;
788 //m_rmap[142195] = 2500;
789
790
791 m_deltaR_v_eta = new Resplot("R v eta", 10, -2.5, 2.5, 100, 0, 0.1 );
792 m_deltaR_v_pt = new Resplot("R v pt", ptnbins, ptbinlims, 100, 0, 0.1 );
793
794
795 TH1F* eff_vs_mult = new TH1F( "eff_vs_mult", "eff_vs_mult", 25, 0, 25 );
796
797 m_eff_vs_mult = new Efficiency1D( eff_vs_mult, "eff_mult" );
798
799
800 m_n_vtx_tracks = new TH1F("nvtxtracks", "nvtxtracks", 150, 0, 600);
801 m_eff_vs_ntracks = new Efficiency1D( m_n_vtx_tracks, "eff_vs_ntracks");
802
803
804
805 //double oldnbins[21] = { 1, 5, 10, 15, 20,
806 // 25, 30, 35, 40, 45,
807 // 50, 60, 80, 100, 150,
808 // 200, 250, 300, 400, 500, 600 };
809
810
811 //double _nbins[23] = { 0, 29.4, 39.9179, 48.4358, 56.1813, 62.9524, 70.1377, 76.8907,
812 // 83.4667, 90.559, 97.4902, 105.737, 113.547, 121.281, 129.015, 139.824,
813 // 150.589, 164.093, 179.096, 206.867, 400, 500, 600 };
814
815 double nbins[23] = { 0, 29.5, 39.5, 48.5, 56.5, 63.5, 70.5, 77.5,
816 83.5, 91.5, 97.5, 106.5, 114.5, 121.5, 129.5, 140.5,
817 151.5, 164.5, 200.5, 250.5, 300.5, 400.5, 600 };
818
819 TH1F* n_vtx_tracks2 = new TH1F("nvtxtracks2", "nvtxtracks2", 22, nbins);
820 m_eff_vs_ntracks2 = new Efficiency1D( n_vtx_tracks2, "eff_vs_ntracks2");
821 delete n_vtx_tracks2;
822
823 m_n_vtx = new TH1F("nvtx", "nvtx", 81, -0.5, 80.5);
824 m_eff_vs_nvtx = new Efficiency1D( m_n_vtx, "eff_vs_nvtx");
825 //m_mu = new TH1F("mu", "mu", 3000, -0.5, 29.5);
826 m_mu = new TH1F("mu", "mu", 90, 0, 90);
827 m_eff_vs_mu = new Efficiency1D( m_mu, "eff_vs_mu");
828
829
831
832 double etovpt_bins[39] = {
833 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9,
834 1, 1.08571, 1.17877, 1.2798, 1.3895, 1.50859, 1.63789, 1.77828, 1.9307, 2.09618,
835 2.27585, 2.47091, 2.6827, 2.91263, 3.16228, 3.43332, 3.72759, 4.04709, 4.39397, 4.77058,
836 5.17947, 5.62341, 6.1054, 6.6287, 7.19686, 7.81371, 8.48343, 9.21055, 10
837 };
838
839 m_etovpt_raw = new TH1F("etovpt_raw", "ET / pT", 100, 0, 10 );
840
841 m_etovpt = new TH1F("etovpt", "ET / pT", 38, etovpt_bins );
842 m_eff_vs_etovpt = new Efficiency1D( m_etovpt, "eff_vs_etovpt");
843
844
845 m_et = new TH1F("ET", "ET; E_{T} [GeV]", ptnbins, ptbinlims );
846 m_eff_vs_et = new Efficiency1D( m_et, "eff_vs_ET" );
847
848
849 // std::cout << "initialize() Directory " << gDirectory->GetName() << " on leaving" << std::endl;
850
851 ConfVtxAnalysis* vtxanal = 0;
852 store().find( vtxanal, "rvtx" );
853 if ( vtxanal ) vtxanal->initialise();
854
855 // std::cout << "initialize() Directory " << gDirectory->GetName() << " on leaving" << std::endl;
856
857
858 m_dir->pop();
859
860 dir->cd();
861
862
863
864}
865
866
867
868
869
870// fit Gaussian to histogram
871
872TF1* FitFWGaussian(TH1D* s, double a, double b) {
873
874 // std::cout << "FitFWGaussian " << s->GetName() << std::endl;
875
876 // TF1* f1 = new TF1("gausr", "gaus");
877 TF1* f1 = new TF1("gausr", "[0]*exp(-(x-[1])*(x-[1])/([2]*[2]))");
878
879 f1->SetNpx(5000);
880 f1->SetLineWidth(1);
881
882 f1->SetParameter(0,s->GetBinContent(s->GetMaximumBin()) );
883 f1->SetParameter(1,s->GetBinCenter(s->GetMaximumBin()));
884 // f1->SetParameter(2,s->GetRMS());
885 f1->SetParameter(2,1.5);
886 f1->FixParameter(2,1.5);
887
888 f1->SetParName(0, "Norm");
889 f1->SetParName(1, "Mean");
890 f1->SetParName(2, "Sigma");
891
892 int nbins=0;
893 for (int j=1 ; j<=s->GetNbinsX() ; j++) if (s->GetBinContent(j)) nbins++;
894
895 if (nbins>2) {
896 if ( a==-999 || b==-999 ) s->Fit(f1,"Q");
897 else s->Fit(f1,"Q", "", a, b);
898 }
899 else for ( int j=0 ; j<3 ; j++ ) f1->SetParameter(j, 0);
900
901 // std::cout << "return" << std::endl;
902
903 return f1;
904}
905
906
907
908
909
910void fitSin( TH1D* h, const std::string& parent="" ) {
911
912 TF1* fsin = new TF1( "sinp", "sqrt([0]*[0])*sin([1]-x)+[2]" ); // , -M_PI, M_PI );
913
914 fsin->SetParameter(0,1);
915 fsin->SetParameter(1,0);
916 fsin->SetParameter(2,0);
917
918 fsin->SetLineWidth(1);
919
920 // h->SetTitle(";track #phi;d_{0} [mm]");
921
922 std::cout << __FUNCTION__ << " " << h->GetName() << std::endl;
923 h->Fit( fsin, "", "", -M_PI, M_PI );
924
925 double radius = fsin->GetParameter(0);
926 double phi = fsin->GetParameter(1);
927
928 double y = radius*sin(phi);
929 double x = radius*cos(phi);
930
931 std::cout << parent << "\t" << h->GetTitle() << "\tx = " << x << "\ty = " << y << std::endl;
932
933 delete fsin;
934}
935
936
937
938
940
942
943 // gDirectory->pwd();
944
945 if ( !m_initialised ) return;
946
947 std::cout << "ConfAnalysis::finalise() " << name();
948
949 if ( name().size()<19 ) std::cout << "\t";
950 if ( name().size()<30 ) std::cout << "\t";
951 if ( name().size()<41 ) std::cout << "\t";
952 if ( name().size()<52 ) std::cout << "\t";
953
954
955 std::cout << "\tNreco " << m_Nreco
956 << "\tNref " << m_Nref
957 << "\tNmatched " << m_Nmatched;
958
959 if (m_Nref) {
960 std::cout << " tracks approx " << (100*m_Nmatched)/m_Nref << "%" ;
961 }
962 std::cout << std::endl;
963
964 // if ( m_Nreco==0 ) return;
965
966 // TIDDirectory d( name() );
967 // d.push();
968
969 m_dir->push();
970
971 std::map<std::string, TH1F*>::iterator hitr=m_histos.begin();
972 std::map<std::string, TH1F*>::iterator hend=m_histos.end();
973 for ( ; hitr!=hend ; ++hitr ) hitr->second->Write();
974 // std::cout << "DBG >" << eff_pt->Hist()->GetName() << "< DBG" << std::endl;
975
976 // std::vector<Efficiency*> heff = { eff_pt,
977
978 const unsigned Neff = 11;
979 Efficiency1D* heff[Neff] = { m_eff_pt,
980 m_eff_eta,
981 m_eff_phi,
982 m_eff_z0,
983 m_eff_d0,
984 m_eff_a0,
985 m_eff_ptm,
986 m_eff_ptp,
989 m_eff_roi_dR };
990
991 for ( unsigned i=0 ; i<Neff ; i++ ) {
992 heff[i]->finalise();
993 heff[i]->Bayes()->Write( ( heff[i]->name()+"_tg" ).c_str() );
994 } // heff[i]->Hist()->Write(); }
995
996 // std::cout << "DBG >" << m_purity_pt->Hist()->GetName() << "< DBG" << std::endl;
997
998 m_eff_vs_mult->finalise();
999
1000 // Normalise(n_vtx_tracks);
1001
1002 m_eff_vs_ntracks->finalise();
1003 m_eff_vs_ntracks2->finalise();
1004
1005 m_eff_vs_nvtx->finalise();
1006 m_eff_vs_mu->finalise();
1007
1008 m_eff_vs_etovpt->finalise();
1009
1010 m_eff_vs_et->finalise();
1011
1012 m_eff_eta_vs_pt->finalise();
1013 m_eff_d0_vs_pt->finalise();
1014
1015
1016 const unsigned Npurity = 6;
1017 Efficiency1D* hpurity[Npurity] = {
1023 m_purity_a0 };
1024
1025 for ( unsigned i = 0 ; i<Npurity ; i++ ) hpurity[i]->finalise();
1026
1028 fitSin( m_rd0_vs_phi->Mean(), name()+"/rd0_vs_phi" );
1029 m_rd0_vs_phi->Write();
1030
1032 fitSin( m_rd0_vs_phi_rec->Mean(), name()+"/rd0_vs_phi_rec" );
1033 m_rd0_vs_phi_rec->Write();
1034
1035
1036 std::string spstr[5] = { "npix", "nsct", "nsi", "ntrt", "nbl" };
1037 for ( int i=m_res.size() ; i-- ; ) {
1038 TF1* (*resfit)(TH1D* s, double a, double b) = Resplot::FitNull95;
1039 for ( int ir=0 ; ir<5 ; ir++ ) if ( m_res[i]->Name().find(spstr[ir])!=std::string::npos ) { resfit = Resplot::FitNull; break; }
1040 m_res[i]->Finalise( resfit );
1041 m_res[i]->Write();
1042 delete m_res[i];
1043 }
1044
1045 m_deltaR_v_eta->Finalise(); m_deltaR_v_eta->Write();
1046 m_deltaR_v_pt->Finalise(); m_deltaR_v_pt->Write();
1047
1048 for ( unsigned i=m_rDd0res.size() ; i-- ; ) {
1049 m_rDd0res[i]->Finalise(Resplot::FitNull95);
1050 m_rDd0res[i]->Write();
1051 delete m_rDd0res[i];
1052 }
1053
1054 for ( unsigned i=m_rDa0res.size() ; i-- ; ) {
1055 m_rDa0res[i]->Finalise(Resplot::FitNull95);
1056 m_rDa0res[i]->Write();
1057 delete m_rDa0res[i];
1058 }
1059
1060 for ( unsigned i=m_rDz0res.size() ; i-- ; ) {
1061 m_rDz0res[i]->Finalise(Resplot::FitNull95);
1062 m_rDz0res[i]->Write();
1063 delete m_rDz0res[i];
1064 }
1065
1066
1067 for ( unsigned ih=0 ; ih<2 ; ih++ ) {
1068 m_hphivsDd0res[ih]->Divide( m_hphivsDd0res[2] );
1069 m_hphivsDa0res[ih]->Divide( m_hphivsDa0res[2] );
1070 }
1071
1072 for ( unsigned ih=0 ; ih<3 ; ih++ ) {
1073 m_hphivsDd0res[ih]->Write();
1074 m_hphivsDa0res[ih]->Write();
1075
1076 delete m_hphivsDd0res[ih];
1077 delete m_hphivsDa0res[ih];
1078 }
1079
1081
1085
1086 m_rRoi_deta_vs_eta->Write();
1087 m_rRoi_dphi_vs_eta->Write();
1088 m_rRoi_dzed_vs_eta->Write();
1089
1090 delete m_rRoi_deta_vs_eta;
1091 delete m_rRoi_dphi_vs_eta;
1092 delete m_rRoi_dzed_vs_eta;
1093
1095
1096
1097
1098 for ( unsigned i=m_retares.size() ; i-- ; ) {
1099
1100#if 1
1101
1102 m_retares[i]->Finalise(Resplot::FitNull95);
1103 m_rphires[i]->Finalise(Resplot::FitNull95);
1104 m_riptres[i]->Finalise(Resplot::FitNull95);
1105 m_rzedres[i]->Finalise(Resplot::FitNull95);
1106 m_rzedthetares[i]->Finalise(Resplot::FitNull95);
1107 // m_rptres[i]->Finalise(Resplot::FitBreit);
1108 //m_rptres[i]->Finalise(Resplot::FitNull95);
1109 m_rd0res[i]->Finalise(Resplot::FitNull95);
1110 // m_rd0res[i]->Finalise(Resplot::FitCentralGaussian);
1111 // m_rd0res_rms[i]->Finalise(Resplot::FitNull);
1112
1113#else
1114
1115 m_retares[i]->Finalise();
1116 m_rphires[i]->Finalise();
1117 m_rzedres[i]->Finalise();
1118 m_riptres[i]->Finalise();
1119 m_rptres[i]->Finalise();
1120 m_rd0res[i]->Finalise();
1121
1122#endif
1123
1124 m_retares[i]->Write();
1125 m_rphires[i]->Write();
1126 m_rzedres[i]->Write();
1127 m_rzedthetares[i]->Write();
1128 m_riptres[i]->Write();
1129 m_rptres[i]->Write();
1130 m_rd0res[i]->Write();
1131
1132 delete m_retares[i];
1133 delete m_rphires[i];
1134 delete m_rzedres[i];
1135 delete m_rzedthetares[i];
1136 delete m_riptres[i];
1137 delete m_rptres[i];
1138 delete m_rd0res[i];
1139
1140
1141 }
1142
1144
1145 m_rzedlb->Finalise(Resplot::FitNull95);
1147
1148 for ( unsigned i=m_retaresPull.size() ; i-- ; ) {
1149
1150 m_retaresPull[i]->Finalise(Resplot::FitNull);
1151 m_rphiresPull[i]->Finalise(Resplot::FitNull);
1152 m_rptresPull[i]->Finalise(Resplot::FitNull);
1153 m_rzedresPull[i]->Finalise(Resplot::FitNull);
1154 m_rd0resPull[i]->Finalise(Resplot::FitNull);
1155
1156 m_retaresPull[i]->Write();
1157 m_rphiresPull[i]->Write();
1158 m_rptresPull[i]->Write();
1159 m_rzedresPull[i]->Write();
1160 m_rd0resPull[i]->Write();
1161
1162 delete m_retaresPull[i];
1163 delete m_rphiresPull[i];
1164 delete m_rptresPull[i];
1165 delete m_rzedresPull[i];
1166 delete m_rd0resPull[i];
1167
1168 }
1169
1170
1171 m_rzedreslb->Write();
1172 delete m_rzedreslb;
1173
1174 m_rzedlb->Write();
1175 delete m_rzedlb;
1176
1177 m_rzedlb_rec->Write();
1178 delete m_rzedlb_rec;
1179
1180 // td->cd();
1181
1182 //ADDED BY JK
1183 //-----Only one more element in d0 and z0 vectors than eta now
1184 m_rzedres[m_rzedres.size()-1]->Finalise(Resplot::FitNull95);
1185 m_rzedres[m_rzedres.size()-1]->Write();
1186 m_rzedthetares.back()->Finalise(Resplot::FitNull95);
1187 m_rzedthetares.back()->Write();
1188 m_rd0res.back()->Finalise(Resplot::FitNull95);
1189 m_rd0res.back()->Write();
1190 //-----
1191
1192 m_eff_vs_lb->finalise();
1193
1194 m_z_vs_lb->Finalise(Resplot::FitNull95); m_z_vs_lb->Write();
1195 delete m_z_vs_lb;
1196
1197 // TH1F* hefflb = eff_vs_lb->Hist();
1198 // hefflb->Fit("pol0");
1199
1200 // ConfVtxAnalysis* vtxanal = 0;
1201 // store().find( vtxanal, "rvtx" );
1202 // if ( vtxanal ) vtxanal->finalise();
1203
1204 // d.pop();
1205
1206 ConfVtxAnalysis* vtxanal = 0;
1207 store().find( vtxanal, "rvtx" );
1208 if ( vtxanal ) vtxanal->finalise();
1209
1210 m_dir->pop();
1211
1212}
1213
1214extern int Nvtxtracks;
1215extern int NvtxCount;
1216
1217
1219
1220
1221double wrapphi( double phi ) {
1222 if ( phi<-M_PI ) phi += 2*M_PI;
1223 if ( phi> M_PI ) phi -= 2*M_PI;
1224 return phi;
1225}
1226
1227
1228void ConfAnalysis::execute( const std::vector<TIDA::Track*>& reftracks,
1229 const std::vector<TIDA::Track*>& testtracks,
1230 TrackAssociator* matcher,
1231 TrigObjectMatcher* objects ) {
1232
1233 // leave this commented code in for debug purposes ...
1234 // if ( objects ) std::cout << "TrigObjectMatcher: " << objects << std::endl;
1235
1237
1238 if ( m_print ) {
1239 std::cout << "ConfAnalysis::execute() \t " << name()
1240 << "\tref " << reftracks.size()
1241 << "\ttest " << testtracks.size();
1242 if ( groi ) std::cout << "\tgroi " << groi << " " << *groi;
1243 std::cout << std::endl;
1244 }
1245
1246 // std::cout << "ConfAnalysis (resolutions really) filling " << std::endl;
1247
1248 // should have these as a class variable
1249 static std::string varName[16] = { "pT", "eta", "phi", "z0", "d0", "a0",
1250 "nsct", "npix", "nsi", "ntrt", "nstraw",
1251 "dd0", "da0", "dz0", "deta", "dphi" };
1252
1253 std::map<std::string, TH1F*>::iterator hmitr = m_histos.find("ntracks");
1254 if ( hmitr!=m_histos.end() ) hmitr->second->Fill( reftracks.size() );
1255
1256 hmitr = m_histos.find("ntracks_rec");
1257 if ( hmitr!=m_histos.end() ) hmitr->second->Fill( testtracks.size() );
1258
1259 bool dump = false;
1260
1261 m_Nreco += testtracks.size();
1262 m_Nref += reftracks.size();
1263
1264
1265 // std::cout << "ConfAnalysis ref tracks " << std::endl;
1266
1267 // why don't we use this vertex position any more ???
1268 // m_Nref = 0;
1269 // for ( int i=reftracks.size() ; i-- ; ) {
1270 // double phit = reftracks[i]->phi();
1271 // // double a0t = reftracks[i]->a0() + sin(phit)*m_xBeamReference - cos(phit)*m_yBeamReference;
1272 // // if ( std::fabs(a0t)<a0 ) m_Nref++; ???
1273 // }
1274
1275 // if ( testtracks.size() ) std::cout << "NTRACKS " << testtracks.size() << std::endl;
1276
1277 for ( int i=reftracks.size() ; i-- ; ) {
1278
1280
1281 if ( groi!=0 ) {
1282 // std::cout << "ConfAnalysis::Fill() groi " << *groi << std::endl;
1283
1284 double deta = reftracks[i]->eta() - groi->eta();
1285 double dphi = wrapphi( reftracks[i]->phi() - groi->phi() );
1286 double dzed = reftracks[i]->z0() - groi->zed();
1287
1288 m_rRoi_deta_vs_eta->Fill( groi->eta(), deta );
1289 m_rRoi_dphi_vs_eta->Fill( groi->eta(), dphi );
1290 m_rRoi_dzed_vs_eta->Fill( groi->eta(), dzed );
1291 }
1292
1294 double ipTt = 1./(reftracks[i]->pT()/1000.);
1295 double pTt = reftracks[i]->pT()/1000;
1296
1297 double etat = reftracks[i]->eta();
1298 double phit = reftracks[i]->phi();
1299
1300 double thetat = 2*std::atan( exp( (-1)*etat ) );
1301
1306
1308
1309 double z0t = reftracks[i]->z0();
1310 double d0t = reftracks[i]->a0() - std::sin(phit)*m_xBeamReference + std::cos(phit)*m_yBeamReference;
1311 double a0t = reftracks[i]->a0();
1312
1313 if ( m_xBeamReference!=0 || m_yBeamReference!=0 ) {
1316 z0t = reftracks[i]->z0()+((std::cos(phit)*m_xBeamReference + std::sin(phit)*m_yBeamReference)/std::tan(thetat));
1317 d0t = reftracks[i]->a0();
1318 a0t = reftracks[i]->a0() + std::sin(phit)*m_xBeamReference - std::cos(phit)*m_yBeamReference;
1319 }
1320
1321 // if ( m_lfirst ) {
1322 // std::cout << "\na0t " << a0t << "(shifted " << d0t << ")" << std::endl;
1323 // std::cout << "\nz0t " << z0t << std::endl;
1324
1325 // std::cout << "\txBeamReference " << m_xBeamReference << "\tyBeamReference " << m_yBeamReference << std::endl;
1326 // }
1327
1328
1329
1330 // std::cout << "a0t " << a0t << std::endl;
1331
1332 m_rChi2prob->Fill( pTt, TMath::Prob(reftracks[i]->chi2(),reftracks[i]->dof()) );
1333 m_rChi2->Fill( pTt, reftracks[i]->chi2() );
1334 m_rChi2dof->Fill( pTt, reftracks[i]->chi2()/reftracks[i]->dof() );
1335
1336 // static double xbeamref = 0;
1337 // if ( m_lfirst || m_xBeamReference!=xbeamref) {
1338 // std::cout << __FUNCTION__ << "\tbeamline " << m_xBeamReference << " " << m_yBeamReference << " (ref)" << std::endl;
1339 // }
1340 // xbeamref = m_xBeamReference;
1341
1342
1344 double dpTt = reftracks[i]->dpT()/1000;
1345 double detat = reftracks[i]->deta();
1346 double dphit = reftracks[i]->dphi();
1347
1348 //RoI variables
1349 float droi_detat = groi->eta() - reftracks[i]->eta();
1350 float droi_dphit = groi->phi() - reftracks[i]->phi();
1351 if ( droi_dphit<-M_PI ) droi_dphit +=2*M_PI;
1352 if ( droi_dphit>M_PI ) droi_dphit -=2*M_PI;
1353 float droi_dRt = std::sqrt(droi_dphit*droi_dphit + droi_detat*droi_detat);
1354
1355
1356 // double dz0t = reftracks[i]->dz0()+((std::cos(phit)*m_xBeamReference + std::sin(phit)*m_yBeamReference)/std::tan(thetat));
1357 // double dd0t = reftracks[i]->da0();
1358 double dz0t = reftracks[i]->dz0();
1359 double dd0t = reftracks[i]->da0() - std::sin(phit)*m_xBeamReference + std::cos(phit)*m_yBeamReference;
1360
1361 // this will be changed when we know the beam spot position
1362 // double a0t = reftracks[i]->a0() + sin(phit)*m_xBeam - cos(phit)*m_yBeam;
1363 double da0t = reftracks[i]->da0();
1364
1365#if 0
1366 std::cout << "etat = " << etat << " +/- " << detat << std::endl;
1367 std::cout << "phit = " << phit << " +/- " << dphit << std::endl;
1368 std::cout << "z0t = " << z0t << " +/- " << dz0t << std::endl;
1369 std::cout << "d0t = " << d0t << " +/- " << dd0t << std::endl;
1370 std::cout << "a0t = " << a0t << " +/- " << da0t << std::endl;
1371 std::cout << "pTt = " << pTt << " +/- " << dpTt << std::endl;
1372#endif
1373
1374 if ( std::fabs(a0t)>a0 ) continue;
1375
1376 // double chi2t = reftracks[i]->chi2();
1377 // hchi2->Fill( chi2t );
1378
1379 double nsctt = reftracks[i]->sctHits();
1380 double npixt = reftracks[i]->pixelHits()*0.5;
1381 double nsit = reftracks[i]->pixelHits()*0.5 + reftracks[i]->sctHits();
1382
1383 double nsctht = reftracks[i]->sctHoles();
1384 double npixht = reftracks[i]->pixelHoles();
1385 double nsiht = reftracks[i]->pixelHoles() + reftracks[i]->sctHoles();
1386
1387 double nbl = reftracks[i]->bLayerHits();
1388 double nblh = ( ( reftracks[i]->expectBL() && reftracks[i]->bLayerHits()<1 ) ? 1 : 0 );
1389
1390 // double ntrtt = reftracks[i]->trHits();
1391 double nstrawt = reftracks[i]->strawHits();
1392
1393 // double ts_scale = (ts-1260400000)*3000.0/(1260700000-1260400000);
1394
1395 // std::cout << "Fill h2 " << " " << m_h2m << " " << *reftracks[i] << std::endl;
1396
1397 m_h2->Fill( phit, d0t );
1398 m_h2a0->Fill( phit, a0t );
1399
1400 m_rd0_vs_phi->Fill( phit, a0t );
1401
1402 double mu_val = gevent->mu();
1403
1404 m_mu->Fill( mu_val );
1405
1406
1407 const TIDA::Track* matchedreco = matcher->matched(reftracks[i]);
1408
1409 // std::cout << "\t\tConfAnalysis " << name() << "\t" << i << " " << *reftracks[i] << " -> ";
1410
1411 // raw reference track distributions
1412 double vpart[16] = { std::fabs(pTt), etat, phit, z0t, d0t, a0t, nsctt, npixt, nsctt, nsit, nstrawt, dd0t, da0t, dz0t, detat, dphit };
1413
1416 for ( int it=0 ; it<11 ; it++ ) {
1417 // std::string hname = varName[it];
1418 // std::map<std::string, TH1F*>::iterator hmitr = m_histos.find(hname);
1419 // if ( hmitr!=m_histos.end() ) hmitr->second->Fill( vpart[it] );
1420
1421 if ( TH1F* hptr = find( varName[it] ) ) hptr->Fill( vpart[it] );
1422 else std::cerr << "hmmm histo " << varName[it] << " not found" << std::endl;
1423 }
1424
1425
1426 m_rnpix_eta->Fill( etat, npixt*1.0 );
1427 m_rnsct_eta->Fill( etat, nsctt*1.0 );
1428 m_rntrt_eta->Fill( etat, nstrawt*1.0 );
1429 m_rnsihit_eta->Fill( etat, npixt + nsctt*1.);
1430
1431 m_rnpix_phi->Fill( phit, npixt*1.0 );
1432 m_rnsct_phi->Fill( phit, nsctt*1.0 );
1433 m_rntrt_phi->Fill( phit, nstrawt*1.0 );
1434
1435 m_rnpix_pt->Fill( std::fabs(pTt), npixt*1.0 );
1436 m_rnsct_pt->Fill( std::fabs(pTt), nsctt*1.0 );
1437 m_rntrt_pt->Fill( std::fabs(pTt), nstrawt*1.0 );
1438
1439
1440 m_rnpix_d0->Fill( a0t, npixt*1.0 );
1441 m_rnsct_d0->Fill( a0t, nsctt*1.0 );
1442 m_rntrt_d0->Fill( a0t, nstrawt*1.0 );
1443
1444 m_rnpixh_d0->Fill( a0t, npixht );
1445 m_rnscth_d0->Fill( a0t, nsctht );
1446
1447 m_rnsi_pt->Fill( std::fabs(pTt), nsit );
1448 m_rnsih_pt->Fill( std::fabs(pTt), nsiht );
1449
1450 m_rnsi_d0->Fill( a0t, nsit );
1451 m_rnsih_d0->Fill( a0t, nsiht );
1452
1453 m_rnsi_eta->Fill( etat, nsit );
1454 m_rnsih_eta->Fill(etat, nsiht );
1455
1456 m_rnbl_d0->Fill( a0t, nbl );
1457 m_rnblh_d0->Fill( a0t, nblh );
1458
1459
1460 m_rnpixh_pt->Fill( std::fabs(pTt), npixht );
1461 m_rnscth_pt->Fill( std::fabs(pTt), nsctht );
1462
1463 m_rnpix_lb->Fill( gevent->lumi_block(), npixt*1.0 );
1464 m_rnsct_lb->Fill( gevent->lumi_block(), nsctt*1.0 );
1465
1466
1467 m_rnsct_vs_npix->Fill( npixt, nsctt );
1468
1469 double etovpt_val = 0;
1470 const TrackTrigObject* tobj = 0;
1471
1472 double ET=0;
1473
1474 if ( objects ) {
1475 tobj = objects->object( reftracks[i]->id() );
1476 if ( tobj ) {
1479 etovpt_val = std::fabs( tobj->pt()/reftracks[i]->pT() );
1480 m_etovpt->Fill( etovpt_val );
1481 m_etovpt_raw->Fill( etovpt_val );
1482 m_et->Fill( tobj->pt()*0.001 );
1483 ET = std::fabs(tobj->pt()*0.001);
1484 }
1485 }
1486
1487
1488 if ( matchedreco ) {
1489
1490 m_Nmatched++;
1491
1492 // efficiency histos
1493 m_eff_pt->Fill(std::fabs(pTt));
1494 m_eff_z0->Fill(z0t);
1495 m_eff_eta->Fill(etat);
1496 m_eff_phi->Fill(phit);
1497 m_eff_d0->Fill(d0t);
1498 m_eff_a0->Fill(a0t);
1499 m_eff_roi_deta->Fill(droi_detat);
1500 m_eff_roi_dphi->Fill(droi_dphit);
1501 m_eff_roi_dR->Fill(droi_dRt);
1502
1503 // signed pT
1504 if ( pTt<0 ) m_eff_ptm->Fill(std::fabs(pTt));
1505 else m_eff_ptp->Fill(std::fabs(pTt));
1506
1507 m_eff_vs_mult->Fill( m_Nref );
1508
1509
1510 m_eff_eta_vs_pt->Fill( std::fabs(pTt), etat );
1511 m_eff_d0_vs_pt->Fill( std::fabs(pTt), a0t );
1512
1513
1514
1515 // m_eff_vs_lb->Fill( m_rmap[r]+lb );
1516 // m_eff_vs_lb->Fill( ts_scale );
1517 // m_eff_vs_lb->Fill( ts );
1518
1519 m_eff_vs_lb->Fill( gevent->lumi_block() );
1520
1521 if ( tobj ) {
1522 m_eff_vs_etovpt->Fill(etovpt_val);
1523 m_eff_vs_et->Fill( std::fabs(tobj->pt()*0.001) );
1524 }
1525
1526
1528
1530 double pTr = matchedreco->pT()/1000;
1531 double etar = matchedreco->eta();
1532 double phir = matchedreco->phi();
1533 //double z0r = matchedreco->z0() + std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest; ;
1534 double thetar = 2*std::atan( exp( (-1)*etar) );
1535
1536 // double z0r = matchedreco->z0()+((std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest)/std::tan(thetar));
1537 // double d0r = matchedreco->a0();
1538 // double a0r = matchedreco->a0() + sin(phir)*m_xBeamTest - cos(phir)*m_yBeamTest; // this will be changed when we know the beam spot position
1539
1540
1541
1542
1543
1544 // static bool tfirst = true;
1545 // static double xbeamtest = 0;
1546 // if ( m_lfirst || xbeamtest!=m_xBeamTest) {
1547 // std::cout << __FUNCTION__ << "\tbeamline " << m_xBeamTest << " " << m_yBeamTest << " (test)" << std::endl;
1548 // }
1549 // xbeamtest = m_xBeamTest;
1550
1551 double d0r = 0;
1552 double a0r = 0;
1553 double z0r = 0;
1554
1555 d0r = matchedreco->a0() - sin(phir)*m_xBeamTest + cos(phir)*m_yBeamTest; // this will be changed when we know the beam spot position
1556 a0r = matchedreco->a0();
1557 z0r = matchedreco->z0();
1558
1559 if ( m_xBeamTest!=0 || m_yBeamTest!=0 ) {
1560 d0r = matchedreco->a0();
1561 a0r = matchedreco->a0() + sin(phir)*m_xBeamTest - cos(phir)*m_yBeamTest; // this will be changed when we know the beam spot position
1562 z0r = matchedreco->z0()+((std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest)/std::tan(thetar));
1563 }
1564
1565
1566 // if ( m_lfirst ) {
1567 // std::cout << "\na0r " << a0r << "(shifted " << d0r << ")" << std::endl;
1568 // std::cout << "\nz0r " << z0r << std::endl;
1569 // std::cout << "\txBeamReference " << m_xBeamReference << "\tyBeamReference " << m_yBeamReference << std::endl;
1570 // }
1571
1572
1573 // static int it=0;
1574
1575 // std::cout << "it " << it++ << std::endl;
1576
1577 // m_lfirst = false;
1578
1579 double nsctr = matchedreco->sctHits();
1580 double npixr = matchedreco->pixelHits()*0.5;
1581 double nsir = matchedreco->pixelHits()*0.5 + matchedreco->sctHits();
1582
1583 double nscthr = matchedreco->sctHoles();
1584 double npixhr = matchedreco->pixelHoles();
1585
1586
1587 m_rnpix_lb_rec->Fill( gevent->lumi_block(), npixr*1.0 );
1588 m_rnsct_lb_rec->Fill( gevent->lumi_block(), nsctr*1.0 );
1589
1590
1591 //double ntrtr = matchedreco->trHits();
1592 double nstrawr = matchedreco->strawHits();
1593
1595
1596 double dpTr = matchedreco->dpT()/1000;
1597 double detar = matchedreco->deta();
1598 double dphir = matchedreco->dphi();
1599
1600 // double dz0r = matchedreco->dz0()+((std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest)/std::tan(thetar));
1601 // double dd0r = matchedreco->da0();
1602 // double da0r = matchedreco->da0() + sin(phir)*m_xBeamTest - cos(phir)*m_yBeamTest;
1603
1604 double dz0r = matchedreco->dz0();
1605 double dd0r = matchedreco->da0() - sin(phir)*m_xBeamTest + cos(phir)*m_yBeamTest;
1606 double da0r = matchedreco->da0();
1607
1608#if 0
1609 std::cout << "etar = " << etar << " +/- " << detar << std::endl;
1610 std::cout << "phir = " << phir << " +/- " << dphir << std::endl;
1611 std::cout << "pTr = " << pTr << " +/- " << dpTr << std::endl;
1612 std::cout << "a0r = " << a0r << " +/- " << da0r << std::endl;
1613 std::cout << "d0r = " << d0r << " +/- " << dd0r << std::endl;
1614 std::cout << "z0r = " << z0r << " +/- " << dz0r << std::endl;
1615#endif
1616
1617 if ( m_h2m ) m_h2m->Fill( phit, d0t );
1618
1619 m_rd0_vs_phi_rec->Fill( phir, a0r );
1620
1621
1623 double resfiller[9] = { std::fabs(ipTt), std::fabs(pTt), ET, etat, z0t, double(NvtxCount), double(Nvtxtracks), phit, mu_val };
1624
1625 double Delta_ipt = 1.0/pTr - 1.0/pTt;
1626 double Delta_pt = pTr - pTt;
1627
1628 if ( pTt<0 ) {
1629 Delta_ipt *= -1;
1630 Delta_pt *= -1;
1631 }
1632
1633 for ( int irfill=0 ; irfill<9 ; irfill++ ) {
1634 m_retares[irfill]->Fill( resfiller[irfill], etar-etat );
1635 m_rphires[irfill]->Fill( resfiller[irfill], phir-phit );
1636 m_rzedres[irfill]->Fill( resfiller[irfill], z0r-z0t );
1637 m_rzedthetares[irfill]->Fill( resfiller[irfill], z0r*std::sin(thetar)-z0t*std::sin(thetat) );
1638 m_riptres[irfill]->Fill( resfiller[irfill], Delta_ipt );
1639 m_rptres[irfill]->Fill( resfiller[irfill], Delta_pt );
1640 m_rd0res[irfill]->Fill( resfiller[irfill], a0r-a0t );
1641 }
1642
1643 double lb = gevent->lumi_block();
1644
1645 m_rzedreslb->Fill( lb, z0r-z0t );
1646 m_rzedlb->Fill( lb, z0t );
1647 m_rzedlb_rec->Fill( lb, z0r );
1648
1649 for ( int irfill=0 ; irfill<6 ; irfill++ ) {
1650 m_rphiresPull[irfill]->Fill( resfiller[irfill], (phir - phit) / sqrt( (dphit*dphit) + (dphir*dphir) ) );
1651 m_retaresPull[irfill]->Fill( resfiller[irfill], (etar - etat) / sqrt( (detat*detat) + (detar*detar) ) );
1652 m_rptresPull[irfill]->Fill( resfiller[irfill], Delta_pt / sqrt( (dpTt*dpTt) + (dpTr*dpTr) ) );
1653 m_rzedresPull[irfill]->Fill( resfiller[irfill], (z0r - z0t) / sqrt( (dz0t*dz0t) + (dz0r*dz0r) ) );
1654 m_rd0resPull[irfill]->Fill( resfiller[irfill], (a0r - a0t) / sqrt( (da0t*da0t) + (da0r*da0r) ) );
1655 }
1656
1657 m_rDz0res[0]->Fill( std::fabs(pTt), dz0r-dz0t );
1658 m_rDz0res[1]->Fill( etat, dz0r-dz0t );
1659 m_rDz0res[2]->Fill( z0t, dz0r-dz0t );
1660 m_rDz0res[3]->Fill( phit, dz0r-dz0t );
1661
1662 if ( dumpflag ) {
1663 std::ostream& dumpstream = dumpfile;
1664 if ( dz0t>0 && std::fabs( dz0r-dz0t )>0.04 ) {
1665 dump = true;
1666 dumpstream << "POOR sigma(z0) agreement \n\trefrack: " << *reftracks[i] << "\n\ttestrack: " << *matchedreco << std::endl;
1667 // std::cout << "dz0r dz0t" << dz0r << "\t" << dz0t << std::endl;
1668 }
1669 }
1670
1672 m_rDd0res[0]->Fill( std::fabs(pTt), dd0r-dd0t );
1673 m_rDd0res[1]->Fill( etat, dd0r-dd0t );
1674 m_rDd0res[2]->Fill( z0t, dd0r-dd0t );
1675 m_rDd0res[3]->Fill( d0t, dd0r-dd0t );
1676 m_rDd0res[4]->Fill( phit, dd0r-dd0t );
1677
1678 m_rDa0res[0]->Fill( std::fabs(pTt), da0r-da0t );
1679 m_rDa0res[1]->Fill( etat, da0r-da0t );
1680 m_rDa0res[2]->Fill( z0t, da0r-da0t );
1681 m_rDa0res[3]->Fill( da0t, da0r-da0t );
1682
1683 const double Deltaphi = 2*M_PI/NMod;
1684
1685 double phistart = 11.0819;
1686 if ( NMod==22 ) phistart = 7.05803;
1687
1688 double phit_wrap = phit - phistart*M_PI/180;
1689
1690 if ( phit_wrap<-M_PI ) phit_wrap += 2*M_PI;
1691
1692 double iphi = phit_wrap - (Deltaphi*int((phit_wrap+M_PI)/Deltaphi) - M_PI);
1693
1694 // double iphi = phit - M_PI*int(7*(phit+M_PI)/M_PI)/7 - 11.0819*M_PI/180 + M_PI;
1695 // double iphi = phit - M_PI*int(7*phit/M_PI)/7.0 - 11.0819*M_PI/180;
1696
1697 m_rDa0res[4]->Fill( phit, da0r-da0t );
1698 m_rDa0res[5]->Fill( iphi, da0r-da0t );
1699
1700 m_rDa0res[6]->Fill( iphi, da0t );
1701 m_rDa0res[7]->Fill( iphi, da0r );
1702
1703 //ADDED BY JK
1704 //-----
1707 m_rzedres[m_rphires.size()-1]->Fill( resfiller[1], z0r-z0t );
1708 m_rzedres[m_rphires.size()+1]->Fill( resfiller[1], z0r-z0t );
1709
1710 m_rzedthetares[m_rphires.size()-1]->Fill( resfiller[1], z0r*std::sin(thetar)-z0t*std::sin(thetat) );
1711 m_rzedthetares[m_rphires.size()]->Fill( resfiller[1], z0r*std::sin(thetar)-z0t*std::sin(thetat) );
1712
1713 m_rd0res[m_rphires.size()]->Fill( resfiller[1], a0r-a0t );
1714 m_rd0res[m_rphires.size()+1]->Fill( fabs(resfiller[1]), a0r-a0t ); //
1715
1716 m_rzedresPull[m_rphires.size()-3]->Fill( resfiller[1], (z0r - z0t) / std::sqrt( (dz0t*dz0t) + (dz0r*dz0r) ) );
1717 m_rzedresPull[m_rphires.size()-2]->Fill( fabs(resfiller[1]), (z0r - z0t) / std::sqrt( (dz0t*dz0t) + (dz0r*dz0r) ) );
1718 m_rd0resPull[m_rphires.size()-3]->Fill( resfiller[1], (a0r - a0t) / std::sqrt( (da0t*da0t) + (da0r*da0r) ) );
1719 m_rd0resPull[m_rphires.size()-2]->Fill( fabs(resfiller[1]), (a0r - a0t) / std::sqrt( (da0t*da0t) + (da0r*da0r) ) );
1720
1721 //-----
1722
1723 m_rnpix_eta_rec->Fill( etat, npixr*1.0 );
1724 m_rnsct_eta_rec->Fill( etat, nsctr*1.0 );
1725 m_rntrt_eta_rec->Fill( etat, nstrawr*1.0 );
1726 m_rnsihit_eta_rec->Fill( etat, npixr*0.5 + nsctr*1.0);
1727
1728 m_rnpix_phi_rec->Fill( phit, npixr*1.0 );
1729 m_rnsct_phi_rec->Fill( phit, nsctr*1.0 );
1730 m_rntrt_phi_rec->Fill( phit, nstrawr*1.0 );
1731
1732 m_rnpix_pt_rec->Fill( std::fabs(pTt), npixr*1.0 );
1733 m_rnsct_pt_rec->Fill( std::fabs(pTt), nsctr*1.0 );
1734 m_rntrt_pt_rec->Fill( std::fabs(pTt), nstrawr*1.0 );
1735
1736 m_rnpixh_pt_rec->Fill( std::fabs(pTt), npixhr*0.5 );
1737 m_rnscth_pt_rec->Fill( std::fabs(pTt), nscthr*1.0 );
1738
1739 m_rnpix_d0_rec->Fill( a0t, npixr*1.0 );
1740 m_rnsct_d0_rec->Fill( a0t, nsctr*1.0 );
1741 m_rntrt_d0_rec->Fill( a0t, nstrawr*1.0 );
1742
1745 m_n_vtx_tracks->Fill( Nvtxtracks );
1746
1747 m_eff_vs_nvtx->Fill( NvtxCount );
1748 m_n_vtx->Fill( NvtxCount );
1749
1750 m_eff_vs_mu->Fill( mu_val );
1751
1752 double vres[6] = { Delta_ipt, etar-etat, phir-phit, z0r-z0t, d0r-d0t, a0r-a0t };
1753 for ( int it=0 ; it<6 ; it++ ) {
1754 if ( it==0 ) {
1755 find("ipT_res")->Fill( vres[0] );
1756 find("spT_res")->Fill( 1.0/pTr-1.0/pTt );
1757 }
1758 if ( TH1F* hptr = find(varName[it]+"_res") ) hptr->Fill( vres[it] );
1759 else std::cerr << "hmmm histo " << varName[it]+"_res" << " not found" << std::endl;
1760 }
1761 //2D plot
1762 if ( TH2F* hptr = find2D("eta_phi_rec") ) {
1763 hptr->Fill( etar,phir );
1764 hptr->GetXaxis()->SetTitle("#eta");
1765 hptr->GetYaxis()->SetTitle("#phi");
1766 //hptr->SetFillStyle("COLZ");
1767 }
1768 if ( TH2F* hptr = find2D("phi_d0_Rec") ) {
1769 hptr->Fill( phir,d0r );
1770 hptr->GetXaxis()->SetTitle("#phi");
1771 hptr->GetYaxis()->SetTitle("d_{0} [mm]");
1772 //hptr->SetFillStyle("COLZ");
1773 }
1774
1775 // raw matched test track errors
1776 if ( TH1F* hptr = find("dpT_rec") ) hptr->Fill(dpTr);
1777 if ( TH1F* hptr =find("deta_rec")) hptr->Fill(detar);
1778 if ( TH1F* hptr =find("dphi_rec")) hptr->Fill(dphir);
1779 if ( TH1F* hptr =find("dz0_rec")) hptr->Fill(dz0r);
1780 if ( TH1F* hptr =find("dd0_rec")) hptr->Fill(dd0r);
1781 if ( TH1F* hptr =find("da0_rec")) hptr->Fill(da0r);
1782
1783 // raw matched reference track errors
1784 if ( TH1F* hptr = find("dpT") ) hptr->Fill(dpTt);
1785 if ( TH1F* hptr = find("deta")) hptr->Fill(detat);
1786 if ( TH1F* hptr = find("dphi")) hptr->Fill(dphit);
1787 if ( TH1F* hptr = find("dz0")) hptr->Fill(dz0t);
1788 if ( TH1F* hptr = find("dd0")) hptr->Fill(dd0t);
1789 if ( TH1F* hptr = find("da0")) hptr->Fill(da0t);
1790
1791
1792 if ( TH1F* hptr = find("dd0_res")) hptr->Fill(dd0r-dd0t);
1793 if ( TH1F* hptr = find("da0_res")) hptr->Fill(da0r-da0t);
1794 if ( TH1F* hptr = find("dz0_res")) hptr->Fill(dz0r-dz0t);
1795
1796 double Dd0 = dd0r-dd0t;
1797 double Da0 = da0r-da0t;
1798
1799 double Ddphi = dphir - dphit;
1800
1801 m_hphivsDd0res[2]->Fill( phit );
1802
1803 // if ( matchedreco->bLayerHits()<=3 ) std::cout << "\nov2\t" << Dd0 << " " << *reftracks[i] << std::endl;
1804 // else std::cout << "\nov4\t" << Dd0 << " " << *reftracks[i] << std::endl;
1805
1806
1807
1808 if ( PRINT_BRESIDUALS ) {
1809 if ( matchedreco->bLayerHits()<=3 ) std::cout << "\nov2\t" << Dd0 << " " << *matchedreco << std::endl;
1810 else std::cout << "\nov4\t" << Dd0 << " " << *matchedreco << std::endl;
1811 }
1812
1813 if ( std::fabs(Dd0)<0.01 ) {
1814 m_hphivsDd0res[0]->Fill( phit );
1815
1816 if ( PRINT_BRESIDUALS ) {
1817 std::cout << "close residual " << Dd0 << " " << Ddphi
1818 << " "<< reftracks[i]->bLayerHits()-matchedreco->bLayerHits()
1819 << " "<< reftracks[i]->pixelHits()-matchedreco->pixelHits();
1820 std::cout << "\nccr\t" << Dd0 << " " << Ddphi << " " << *reftracks[i];
1821 std::cout << "\ncct\t" << Dd0 << " " << Ddphi << " " << *matchedreco << std::endl;
1822 }
1823 }
1824 else {
1825 m_hphivsDd0res[1]->Fill( phit );
1826 if ( PRINT_BRESIDUALS ) {
1827 std::cout << "far residual " << Dd0 << " " << Ddphi
1828 << " "<< reftracks[i]->bLayerHits()-matchedreco->bLayerHits()
1829 << " "<< reftracks[i]->pixelHits()-matchedreco->pixelHits();
1830 std::cout << "\nffr\t" << Dd0 << " " << Ddphi << " " << *reftracks[i];
1831 std::cout << "\nfft\t" << Dd0 << " " << Ddphi << " " << *matchedreco << std::endl;
1832 }
1833 }
1834
1835
1836 m_hphivsDa0res[2]->Fill( iphi );
1837 if ( std::fabs(Da0)<0.01 ) m_hphivsDa0res[0]->Fill( iphi );
1838 else m_hphivsDa0res[1]->Fill( iphi );
1839
1840
1841 // pull stats
1842 double pull_pt = Delta_pt / std::sqrt( (dpTt*dpTt) + (dpTr*dpTr) );
1843 double pull_eta = (etar - etat) / std::sqrt( (detat*detat) + (detar*detar) );
1844 double pull_phi = (phir - phit) / std::sqrt( (dphit*dphit) + (dphir*dphir) );
1845 double pull_z0 = (z0r - z0t) / std::sqrt( (dz0t*dz0t) + (dz0r*dz0r) );
1846 double pull_d0 = (d0r - d0t) / std::sqrt( (dd0t*dd0t) + (dd0r*dd0r) );
1847 double pull_a0 = (a0r - a0t) / std::sqrt( (da0t*da0t) + (da0r*da0r) );
1848
1849 if ( TH1F* hptr = find("pT_pull") ) hptr->Fill(pull_pt);
1850 if ( TH1F* hptr = find("eta_pull")) hptr->Fill(pull_eta);
1851 if ( TH1F* hptr = find("phi_pull")) hptr->Fill(pull_phi);
1852 if ( TH1F* hptr = find("z0_pull")) hptr->Fill(pull_z0);
1853 if ( TH1F* hptr = find("d0_pull")) hptr->Fill(pull_d0);
1854 if ( TH1F* hptr = find("a0_pull")) hptr->Fill(pull_a0);
1855
1856 // pull stats - SIMPLE VERSION
1857 double pull_pt_simp = Delta_pt / sqrt( dpTr*dpTr );
1858 double pull_eta_simp = (etar - etat) / sqrt( detar*detar );
1859 double pull_phi_simp = (phir - phit) / sqrt( dphir*dphir );
1860 double pull_z0_simp = (z0r - z0t) / sqrt( dz0r*dz0r );
1861 double pull_d0_simp = (d0r - d0t) / sqrt( dd0r*dd0r );
1862 double pull_a0_simp = (a0r - a0t) / sqrt( da0r*da0r );
1863
1864 if ( TH1F* hptr = find("pT_pull_simple") ) hptr->Fill(pull_pt_simp);
1865 if ( TH1F* hptr = find("eta_pull_simple")) hptr->Fill(pull_eta_simp);
1866 if ( TH1F* hptr = find("phi_pull_simple")) hptr->Fill(pull_phi_simp);
1867 if ( TH1F* hptr = find("z0_pull_simple")) hptr->Fill(pull_z0_simp);
1868 if ( TH1F* hptr = find("d0_pull_simple")) hptr->Fill(pull_d0_simp);
1869 if ( TH1F* hptr = find("a0_pull_simple")) hptr->Fill(pull_a0_simp);
1870
1871
1872 if ( TH1F* hptr = find("etai_res") ) hptr->Fill( etat-etar );
1873
1874
1875 double Delphi = phit-phir;
1876 double Deleta = etat-etar;
1877
1878 if ( Delphi<-M_PI ) Delphi+=2*M_PI;
1879 if ( Delphi>M_PI ) Delphi -=2*M_PI;
1880
1881 double DeltaR = std::sqrt(Delphi*Delphi+Deleta*Deleta);
1882
1883 m_hDeltaR->Fill(DeltaR);
1884
1885 m_deltaR_v_eta->Fill(etat, DeltaR);
1886 m_deltaR_v_pt->Fill(std::fabs(pTt), DeltaR);
1887
1888 // in this loop over the reference tracks, could fill efficiency
1889 // histograms
1890
1891 // m_eff_vs_lb->Fill( m_rmap[r]+lb );
1892
1893 if ( TH1F* hptr = find("nsi_matched")) hptr->Fill(nsir);
1894
1896
1897
1898 m_rChi2prob_rec->Fill( std::fabs(pTr), TMath::Prob(matchedreco->chi2(),matchedreco->dof()) );
1899 m_rChi2_rec->Fill( std::fabs(pTr), matchedreco->chi2() );
1900 m_rChi2dof_rec->Fill( std::fabs(pTr), matchedreco->chi2()/matchedreco->dof() );
1901
1902 m_rChi2d_vs_Chi2d->Fill( reftracks[i]->chi2()/reftracks[i]->dof(),
1903 matchedreco->chi2()/matchedreco->dof() );
1904
1905 m_rDChi2dof->Fill( reftracks[i]->chi2()/reftracks[i]->dof(),
1906 (matchedreco->chi2()/matchedreco->dof())-(reftracks[i]->chi2()/reftracks[i]->dof()) );
1907
1908 }
1909 else {
1910
1912
1913 m_rnpix_pt_bad->Fill( std::fabs(pTt), npixt*0.5 );
1914 m_rnsct_pt_bad->Fill( std::fabs(pTt), nsctt*1.0 );
1915 m_rntrt_pt_bad->Fill( std::fabs(pTt), nstrawt*1.0 );
1916
1917 m_rChi2prob_bad->Fill( std::fabs(pTt), TMath::Prob(reftracks[i]->chi2(),reftracks[i]->dof()) );
1918 m_rChi2_bad->Fill( std::fabs(pTt), reftracks[i]->chi2() );
1919 m_rChi2dof_bad->Fill( std::fabs(pTt), reftracks[i]->chi2()/reftracks[i]->dof() );
1920
1921
1922 // fill efficiencies with unmatched histos
1923 // std::cout << "NULL" << std::endl;
1924 m_eff_pt->FillDenom(std::fabs(pTt));
1925 m_eff_z0->FillDenom(z0t);
1926 m_eff_eta->FillDenom(etat);
1927 m_eff_phi->FillDenom(phit);
1928 m_eff_d0->FillDenom(d0t);
1929 m_eff_a0->FillDenom(a0t);
1930
1931 // signed pT
1932 if ( pTt<0 ) m_eff_ptm->FillDenom(std::fabs(pTt));
1933 else m_eff_ptp->FillDenom(std::fabs(pTt));
1934
1935 m_eff_roi_deta->FillDenom(droi_detat);
1936 m_eff_roi_dphi->FillDenom(droi_dphit);
1937 m_eff_roi_dR->FillDenom(droi_dRt);
1938
1939 m_eff_vs_mult->FillDenom( m_Nref );
1940
1941 m_eff_eta_vs_pt->FillDenom( std::fabs(pTt), etat );
1942 m_eff_d0_vs_pt->FillDenom( std::fabs(pTt), a0t );
1943
1944 dump = false;
1945
1946 m_eff_vs_ntracks->FillDenom( Nvtxtracks );
1947 m_eff_vs_ntracks2->FillDenom( Nvtxtracks );
1948 m_n_vtx_tracks->Fill( Nvtxtracks );
1949
1950
1951 m_eff_vs_nvtx->FillDenom( NvtxCount );
1952 m_n_vtx->Fill( NvtxCount );
1953
1954 double mu_val = gevent->mu();
1955
1956 m_eff_vs_mu->FillDenom(mu_val);
1957
1958 if ( tobj ) {
1959 m_eff_vs_etovpt->FillDenom(etovpt_val);
1960 m_eff_vs_et->FillDenom( std::fabs(tobj->pt()*0.001) );
1961 }
1962
1963 if ( dumpflag ) {
1964 std::ostream& dumpstream = dumpfile;
1965
1966 if ( std::fabs(pTt)>1 ) {
1967 dump = true;
1968
1969 hipt = true;
1970 dumpstream << m_name << "\tMISSING TRACK run " << r << "\tevent " << ev
1971 << "\tlb " << lb << "\tN vertices " << NvtxCount << std::endl;
1972 dumpstream << m_name << "\tMISSING TRACK RoI " << *groi << std::endl;
1973 dumpstream << m_name << "\tMISSING TRACK Track " << *reftracks[i];
1974 if ( std::fabs(pTt)>=30 ) dumpstream << "\tvery high pt";
1975 if ( std::fabs(pTt)>4 &&
1976 std::fabs(pTt)<30 ) dumpstream << "\t high pt";
1977 dumpstream << std::endl;
1978
1979 if ( std::fabs(pTt)>=20 ){
1980 dumpstream << "Test tracks " << std::endl;
1981 for (unsigned int ii=0; ii<testtracks.size(); ii++){
1982 dumpstream << *testtracks[ii] << std::endl;
1983 }
1984 }
1985 }
1986 }
1987
1988
1989 // m_eff_vs_lb->FillDenom( ts );
1990 m_eff_vs_lb->FillDenom( gevent->lumi_block() );
1991 }
1992
1993 }
1994
1995 // return;
1996
1997
1998 // for fake/purity histograms, loop over the test tracks
1999 // and get the corresponding matched reference tracks from the
2000 // reverse map in the TrackAssociator class - revmatched()
2001
2002 static int icount = 0;
2003
2004 // if ( icount%1000 ) std::cout << "chain " << name() << "\t " << m_Nreco << " tracks" << std::endl;
2005 // if ( icount%1000 )
2006 if ( m_print ) std::cout << "ConfAnalysis::execute() \t " << name() << "\t " << icount << " events\t " << testtracks.size() << " tracks (" << m_Nreco << ")" << "\n---------------" << std::endl;
2007
2008 icount++;
2009
2010 for ( int i=testtracks.size() ; i-- ; ) {
2011
2012 // std::cout << "\t\tConfAnalysis purity " << name() << "\t" << i << " " << *testtracks[i] << " -> ";
2013
2014 // double pTr = std::fabs(testtracks[i]->pT());
2015 double pTr = testtracks[i]->pT()/1000;
2016 double etar = testtracks[i]->eta();
2017 double phir = testtracks[i]->phi();
2018 double thetar = 2*std::atan( exp( (-1)*etar) );
2019
2020 double z0r = testtracks[i]->z0(); // + ((std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest)/std::tan(thetar));
2021 double d0r = testtracks[i]->a0() - sin(phir)*m_xBeamTest + cos(phir)*m_yBeamTest; // this will be changed when we know the beam spot position
2022 double a0r = testtracks[i]->a0();
2023 // double a0rp = testtracks[i]->a0() - sin(phir)*m_xBeam - cos(phir)*m_yBeam; // this will be changed when we know the beam spot position
2024
2025 if ( m_xBeamTest!=0 && m_yBeamTest!=0 ) {
2026 d0r = testtracks[i]->a0();
2027 a0r = testtracks[i]->a0() + sin(phir)*m_xBeamTest - cos(phir)*m_yBeamTest; // this will be changed when we know the beam spot position
2028 z0r = testtracks[i]->z0()+((std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest)/std::tan(thetar));
2029 }
2030
2031
2032 // std::cout << "d0 " << d0r << "\tphi " << phir << "\tx " << m_xBeamTest << "\ty " << m_yBeamTest << std::endl;
2033
2034 double nsctr = testtracks[i]->sctHits();
2035 double npixr = testtracks[i]->pixelHits()*0.5;
2036 double nsir = testtracks[i]->pixelHits()*0.5 + testtracks[i]->sctHits();
2037
2038 double ntrtr = testtracks[i]->trHits();
2039 double nstrawr = testtracks[i]->strawHits();
2040
2041
2042 m_rnsct_vs_npix_rec->Fill( npixr, nsctr );
2043
2044
2045#if 0
2046 double dpTr_b = testtracks[i]->dpT()/1000;
2047 double detar_b = testtracks[i]->deta();
2048 double dphir_b = testtracks[i]->dphi();
2049 double dz0r_b = testtracks[i]->dz0(); // + ((std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest)/std::tan(thetar));
2050 double dd0r_b = testtracks[i]->da0() - sin(phir)*m_xBeamTest + cos(phir)*m_yBeamTest;
2051 double da0r_b = testtracks[i]->da0();
2052
2053
2054 std::cout << "pTr_b = " << pTr << " +/- " << dpTr_b << std::endl;
2055 std::cout << "phir_b = " << phir << " +/- " << dphir_b << std::endl;
2056 std::cout << "z0r_b = " << z0r << " +/- " << dz0r_b << std::endl;
2057 std::cout << "d0r_b = " << d0r << " +/- " << dd0r_b << std::endl;
2058 std::cout << "a0r_b = " << a0r << " +/- " << da0r_b << std::endl;
2059 std::cout << "etar_b = " << etar << " +/- " << detar_b << std::endl;
2060#endif
2061
2062 // double ts_scale = (ts-1260400000)*3000.0/(1260700000-1260400000);
2063
2064 // m_z_vs_lb->Fill( m_rmap[r]+lb, z0r );
2065 // m_z_vs_lb->Fill( ts, z0r );
2066 m_z_vs_lb->Fill( gevent->lumi_block(), z0r );
2067
2068 // hnpix_v_sct_rec->Fill( nsctr*0.5, npixr*1.0 );
2069
2070 if ( m_h2r ) m_h2r->Fill( phir, d0r );
2071 if ( m_h2a0r ) m_h2a0r->Fill( phir, a0r );
2072
2073 const TIDA::Track* matchedref = matcher->revmatched(testtracks[i]);
2074
2075 // if ( matchedref ) std::cout << *matchedref << std::endl;
2076 // else std::cout << "NULL" << std::endl;
2077
2078#if 1
2079 // raw test track distributions
2080 double vpart[11] = { std::fabs(pTr), etar, phir, z0r, d0r, a0r, nsctr, npixr, nsir, ntrtr, nstrawr };
2081 for ( int it=0 ; it<11 ; it++ ) {
2082 // std::string hname = name()+"_"+varName[it]+"_rec";
2083 // std::string hname = varName[it]+"_rec";
2084 // std::map<std::string, TH1F*>::iterator hmitr = m_histos.find(hname);
2085 // if ( hmitr!=m_histos.end() ) hmitr->second->Fill( vpar[it] );
2086 // else std::cerr << "hmmm histo " << hname << " not found" << std::endl;
2087 if ( TH1F* hptr = find(varName[it]+"_rec") ) hptr->Fill( vpart[it] );
2088 else std::cerr << "hmmm histo " << varName[it]+"_rec" << " not found" << std::endl;
2089 }
2090 //2D plot
2091 if ( TH2F* hptr = find2D("eta_phi_rec") ) {
2092 hptr->Fill( etar,phir );
2093 hptr->GetXaxis()->SetTitle("#eta");
2094 hptr->GetYaxis()->SetTitle("#phi");
2095 //hptr->SetFillStyle("COLZ");
2096 }
2097 if ( TH2F* hptr = find2D("phi_d0_rec") ) {
2098 hptr->Fill( phir,d0r );
2099 hptr->GetXaxis()->SetTitle("#phi");
2100 hptr->GetYaxis()->SetTitle("d_{0} [mm]");
2101 //hptr->SetFillStyle("COLZ");
2102 }
2103#endif
2104
2105
2106 // purities
2107 if ( matchedref ) {
2108
2109 // std::cout << *matchedref << std::endl;
2110
2111 m_purity_pt->Fill(std::fabs(pTr));
2112 m_purity_z0->Fill(z0r);
2113 m_purity_eta->Fill(etar);
2114 m_purity_phi->Fill(phir);
2115 m_purity_d0->Fill(d0r);
2116 m_purity_a0->Fill(a0r);
2117
2118 // hnpix_v_sct_match->Fill( nsctr*0.5, npixr*0.5 );
2119
2120 }
2121 else {
2122 // std::cout << "NULL" << std::endl;
2123 m_purity_pt->FillDenom(std::fabs(pTr));
2124 m_purity_z0->FillDenom(z0r);
2125 m_purity_eta->FillDenom(etar);
2126 m_purity_phi->FillDenom(phir);
2127 m_purity_d0->FillDenom(d0r);
2128 m_purity_a0->FillDenom(a0r);
2129 }
2130
2131 }
2132
2133 if ( dump && m_print ) {
2134
2135 std::cout << "ConfAnalysis::execute() missed a high pT track - dumping tracks" << std::endl;
2136
2137 for ( int i=reftracks.size() ; i-- ; ) {
2138
2139 if ( std::fabs( reftracks[i]->pT() ) > 1000 ) {
2140 std::cout << "\t dump " << *reftracks[i];
2141 const TIDA::Track* matchedreco = matcher->matched(reftracks[i]);
2142 if ( matchedreco ) std::cout << " <--> " << *matchedreco << std::endl;
2143 else std::cout << std::endl;
2144 }
2145
2146 }
2147
2148 for ( int i=testtracks.size() ; i-- ; ) {
2149 const TIDA::Track* matchedref = matcher->revmatched(testtracks[i]);
2150 if ( matchedref==0 ) std::cout << "\t\t\t\t\t " << *testtracks[i] << std::endl;
2151 }
2152
2153 }
2154
2155 if ( m_print ) std::cout << "ConfAnalysis::execute() exiting" << std::endl;
2156
2157}
2158
2159
2160
2161
2162
#define M_PI
Scalar phi() const
phi method
scales to modify the binning for histograms
int Nvtxtracks
Definition globals.cxx:17
void fitSin(TH1D *h, const std::string &parent="")
BinConfig cosmicBinConfig("cosmic")
TF1 * FitFWGaussian(TH1D *s, double a, double b)
BinConfig electronBinConfig("electron")
double wrapphi(double phi)
fill all the histograms - matched histograms, efficiencies etc
int NvtxCount
Definition globals.cxx:18
BinConfig muonBinConfig("muon")
bool PRINT_BRESIDUALS
stack trace headers
BinConfig g_binConfig("standard")
BinConfig bjetBinConfig("bjet")
std::ofstream dumpfile("dumpfile.log")
BinConfig tauBinConfig("tau")
void Normalise(TH1 *h)
TIDA::Event * gevent
Definition globals.cxx:13
JetDumper::Name Name
Definition JetDumper.cxx:19
static Double_t a
Basic event class to contain a vector of chains for trigger analysis.
TIDA::Associator< TIDA::Track > TrackAssociator
base class for a single track selection filter allowing parameter setting for complex track selection
#define y
#define x
Header file for AthHistogramAlgorithm.
Efficiency1D * m_eff_roi_deta
TH1F * m_n_vtx_tracks
Resplot * m_rnsih_pt
Resplot * m_rChi2prob_rec
Resplot * m_h2r
Efficiency1D * m_eff_ptp
Efficiency1D * m_eff_phi
Resplot * m_deltaR_v_eta
Efficiency1D * m_purity_phi
Resplot * m_rnpix_d0
Resplot * m_rnscth_pt
Resplot * m_rRoi_dzed_vs_eta
Resplot * m_rntrt_d0
Resplot * m_rnpix_eta_rec
Resplot * m_rnsih_eta
Resplot * m_rnsct_pt_rec
Efficiency1D * m_eff_vs_ntracks2
Resplot * m_rnsct_d0_rec
Resplot * m_rntrt_eta_rec
Efficiency1D * m_eff_vs_mu
Resplot * m_rnscth_pt_rec
Efficiency1D * m_eff_a0
Efficiency1D * m_purity_a0
Resplot * m_rnpixh_pt_rec
TH1F * m_etovpt_raw
electron specific ET/PT related stuff
Efficiency1D * m_eff_vs_lb
Resplot * m_rnsihit_eta
Resplot * m_rzedreslb
Resplot * m_h2a0
Efficiency1D * m_purity_pt
Resplot * m_rnsct_phi_rec
Efficiency1D * m_eff_vs_nvtx
Efficiency2D * m_eff_eta_vs_pt
std::vector< Resplot * > m_retaresPull
Resplot * m_h2m
int m_Nreco
number of reconstructed tracks
std::vector< Resplot * > m_rd0resPull
Resplot * m_deltaR_v_pt
TH1F * m_hphivsDd0res[3]
Resplot * m_rChi2_rec
void addHistogram(TH1F *h)
Resplot * m_rnpix_lb
std::vector< Resplot * > m_rptres
Resplot * m_rnpix_pt_bad
Resplot * m_rnpix_d0_rec
TH1F * find(const std::string &n)
Resplot * m_rntrt_pt_bad
Efficiency1D * m_eff_vs_ntracks
bool m_initialiseFirstEvent
Efficiency1D * m_eff_vs_etovpt
virtual void execute(const std::vector< TIDA::Track * > &reftracks, const std::vector< TIDA::Track * > &testtracks, TrackAssociator *matcher, TrigObjectMatcher *objects)
Resplot * m_rnblh_d0
Resplot * m_rnpix_eta
Efficiency1D * m_eff_eta
Resplot * m_rChi2dof
TagNProbe * m_TnP_tool
Resplot * m_rnsct_eta
Resplot * m_rChi2
std::vector< Resplot * > m_rDz0res
Resplot * m_rnbl_d0
TH1F * m_invmassObj
std::vector< Resplot * > m_rzedthetares
Resplot * m_rnpix_phi_rec
Efficiency1D * m_eff_vs_mult
virtual void initialise()
book all the histograms
Resplot * m_rChi2prob
Resplot * m_rnsi_eta
Resplot * m_rnsihit_eta_rec
Resplot * m_rnsct_lb_rec
Efficiency1D * m_eff_roi_dphi
Resplot * m_rnsct_phi
Efficiency1D * m_eff_roi_dR
Resplot * m_rntrt_phi_rec
std::vector< Resplot * > m_rptresPull
std::vector< Resplot * > m_res
Resplot * m_rnpix_lb_rec
Resplot * m_rd0_vs_phi
beam spot dependent
Resplot * m_rDChi2dof
Resplot * m_rntrt_eta
Resplot * m_rnsct_vs_npix_rec
void addHistogram2D(TH2F *h)
Resplot * m_rChi2dof_rec
TH1F * m_hphivsDa0res[3]
Efficiency1D * m_purity_d0
Efficiency1D * m_eff_pt
std::vector< Resplot * > m_rd0res
Resplot * m_rnsct_pt_bad
Efficiency1D * m_eff_vs_et
Resplot * m_rnsct_eta_rec
Resplot * m_rChi2prob_bad
std::vector< Resplot * > m_retares
Resplot * m_rnsct_lb
Resplot * m_rnpixh_d0
Resplot * m_rnsi_pt
Efficiency1D * m_eff_z0
Resplot * m_rnsih_d0
Resplot * m_rChi2d_vs_Chi2d
Resplot * m_rntrt_pt
Efficiency2D * m_eff_d0_vs_pt
Resplot * m_rd0_vs_phi_rec
virtual void finalise()
calculate the efficiencies and write them out with all the histograms
std::vector< Resplot * > m_riptres
Resplot * m_z_vs_lb
Resplot * m_rntrt_d0_rec
Resplot * m_rnscth_d0
Resplot * m_rnpixh_pt
Efficiency1D * m_purity_eta
std::vector< Resplot * > m_rzedresPull
Resplot * m_rnpix_phi
Resplot * m_rChi2_bad
Resplot * m_rnpix_pt_rec
Resplot * m_rnsct_d0
std::vector< Resplot * > m_rphires
Resplot * m_rChi2dof_bad
std::vector< Resplot * > m_rphiresPull
Resplot * m_h2a0r
virtual void initialiseInternal()
Resplot * m_h2
Resplot * m_rzedlb_rec
Resplot * m_rntrt_phi
Resplot * m_rnsi_d0
std::map< std::string, TH1F * > m_histos
bool m_print
flag to print out the matched tracks etc
Resplot * m_rzedlb
Resplot * m_rnpix_pt
Resplot * m_rntrt_pt_rec
TH2F * find2D(const std::string &n)
Resplot * m_rRoi_dphi_vs_eta
Resplot * m_rRoi_deta_vs_eta
Residuals.
TIDDirectory * m_dir
std::vector< Resplot * > m_rzedres
Efficiency1D * m_eff_d0
Resplot * m_rnsct_pt
std::vector< Resplot * > m_rDa0res
Efficiency1D * m_purity_z0
std::vector< Resplot * > m_rDd0res
Efficiency1D * m_eff_ptm
Resplot * m_rnsct_vs_npix
TGraphAsymmErrors * Bayes(double scale=100)
evaluate the uncertainties correctly ...
static TF1 * FitNull(TH1D *s, double a=-999, double b=-999)
Definition Resplot.cxx:1287
static TF1 * FitNull95(TH1D *s, double a=0, double b=0)
Definition Resplot.cxx:1672
virtual const T * revmatched(S *t)
virtual const S * matched(T *t)
int find(T *&t, const std::string &key)
void finalise(double scale=100)
actually calculate the efficiencies
TIDA::FeatureStore & store()
double m_yBeamReference
double m_xBeamReference
beamline positions reference sample
std::string m_name
identifier of the of the analysis - also used for the root directory into which the histograms are pu...
const std::string & name() const
return identifier
double m_xBeamTest
test sample
double pt() const
double chi2(TH1 *h0, TH1 *h1)
int ir
counter of the current depth
Definition fastadd.cxx:49
bool dumpflag
Definition globals.cxx:30
TIDARoiDescriptor * groi
all these externals initialised in globals.cxx
Definition globals.cxx:14
int NMod
Definition globals.cxx:20
int lb
Definition globals.cxx:23
bool hipt
Definition globals.cxx:29
int r
Definition globals.cxx:22
double a0
Definition globals.cxx:27
int ev
Definition globals.cxx:25
-event-from-file
double a0_NScale
Definition BinConfig.h:73
double phires_NScale
Definition BinConfig.h:79
double eta_NScale
Definition BinConfig.h:70
double a0Max
Definition BinConfig.h:85
double phi_NScale
Definition BinConfig.h:71
double d0Max
Definition BinConfig.h:84
double z0_NScale
Definition BinConfig.h:74
double z0Max
Definition BinConfig.h:86
double d0_NScale
Definition BinConfig.h:72
double ptres_NScale
scales for the residuals
Definition BinConfig.h:77
double pt_NScale
Definition BinConfig.h:69