ATLAS Offline Software
Loading...
Searching...
No Matches
ConfVtxAnalysis.cxx
Go to the documentation of this file.
1
9
10
11#include "ConfVtxAnalysis.h"
12
16
17
18ConfVtxAnalysis::ConfVtxAnalysis( const std::string& n, bool use_secVtx_limits ) :
19 VertexAnalysis( n ), m_initialised(false), m_finalised(false), m_use_secVtx_limits(use_secVtx_limits), m_dir(0) { }
20
21
22extern TIDA::Event* gevent;
23
24
26
27 // std::cerr << "ConfVtxAnalysis::initialise() " << name() << std::endl;
28
29 if ( m_initialised ) {
30 std::cerr << "ConfVtxAnalysis::initialise() " << name() << " already initialised" << std::endl;
31 return;
32 }
33
34 m_initialised = true;
35 m_finalised = false;
36
37 // std::cout << "ConfVtxAnalysis::initialise() " << name() << std::endl;
38
39 m_dir = new TIDDirectory(name());
40 m_dir->push();
41
42#if 0
43 double vnbins[41] = {
44 0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5,
45 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 16.5, 17.5, 18.5, 19.5,
46 21.5, 22.5, 24.5, 25.5, 27.5, 29.5,
47 31.5, 33.5, 35.5, 38.5,
48 40.5, 43.5, 46.5,
49 50.5, 53.5, 57.5,
50 61.5, 65.5, 69.5,
51 74.5,
52 80.5 };
53
54
55 double vnbins[81] = {
56 -0.5,
57 0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5,
58 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 16.5, 17.5, 18.5, 19.5,
59 20.5, 21.5, 22.5, 23.5, 24.5, 25.5, 26.5, 27.5, 28.5, 29.5,
60 31.5, 32.5, 33.5, 34.5, 36.5, 37.5, 39.5,
61 40.5, 42.5, 43.5, 45.5, 47.5, 49.5,
62 50.5, 52.5, 54.5, 57.5, 59.5,
63 61.5, 63.5, 66.5, 69.5,
64 71.5, 74.5, 77.5,
65 80.5, 83.5, 86.5,
66 90.5, 93.5, 97.5,
67 100.5, 104.5, 108.5,
68 113.5, 117.5,
69 122.5, 126.5,
70 131.5, 136.5,
71 142.5, 147.5,
72 153.5, 159.5,
73 165.5,
74 171.5, 178.5,
75 185.5,
76 192.5,
77 200.5 };
78
79#endif
80
81 double vnbins[101] = {
82 -0.5,
83 0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5, 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 17.5, 18.5, 19.5, 21.5,
84 23.5, 24.5, 26.5, 28.5, 30.5, 32.5, 35.5, 37.5, 40.5, 43.5, 46.5, 50.5, 53.5, 57.5, 61.5, 66.5, 71.5, 76.5, 81.5, 87.5,
85 93.5, 100.5, 107.5, 114.5, 123.5, 131.5, 141.5, 150.5, 161.5, 172.5, 185.5, 198.5, 211.5, 226.5, 242.5, 259.5, 277.5, 297.5, 317.5, 340.5,
86 363.5, 389.5, 416.5, 445.5, 476.5, 509.5,
87 544.5, 582.5, 623.5, 666.5, 713.5, 762.5, 815.5, 872.5, 933.5, 998.5, 1067.5,
88 1141.5, 1221.5, 1305.5, 1396.5, 1493.5, 1597.5,
89 1708.5, 1827.5, 1953.5, 2089.5,
90 2234.5, 2389.5, 2555.5,
91 2733.5, 2923.5, 3125.5,
92 3342.5, 3574.5,
93 3823.5, 4088.5,
94 4372.5, 4675.5,
95 5000.5
96 };
97
98 double vrbins[41] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9,
99 10, 11, 12, 13, 14, 15, 16, 17, 18, 19,
100 20, 21, 22, 23, 24, 25, 26, 27, 28, 29,
101 30, 33, 36, 39, 42, 45,
102 50, 55, 60,
103 70, 80 };
104
105 m_hnvtx = new TH1F( "nvtx", ";number of vertices", 101, -0.5, 100.5 );
106 m_hzed = new TH1F( "zed", ";vtx z [mm]", 200, -300, 300 );
107 m_hx = new TH1F( "x", ";vtx x [mm]", 800, -1, 1 );
108 m_hy = new TH1F( "y", ";vtx y [mm]", 800, -1, 1 );
109 m_hr = new TH1F( "r", ";vtx r [mm]", 40, vrbins );
110 // hntrax = new TH1F( "ntrax", ";number of tracks", 201, -0.5, 200.5 );
111 // AAAAARGH !! the mu values *all* lie on half integer values
112 // which means that you *can't* ever have "0" or "1" interactions, you can
113 // only ever have -0.5 or 0.5, or 1.5, which at low pile-up is completely
114 // ludicrous, leave the old code here as a reminder
115 // m_hmu = new TH1F( "mu", ";<mu>", 81, -0.5, 80.5 );
116 m_hmu = new TH1F( "mu", ";<mu>", 80, 0, 80 );
117 m_hlb = new TH1F( "lb", ";lumi block", 301, -0.5, 3009.5 );
118
119 m_hnvtx_rec = new TH1F( "nvtx_rec", ";number of vertices", 101, -0.5, 100.5 );
120 m_hzed_rec = new TH1F( "zed_rec", ";vtx z [mm]", 200, -300, 300 );
121 m_hx_rec = new TH1F( "x_rec", ";vtx x [mm]", 800, -1, 1 );
122 m_hy_rec = new TH1F( "y_rec", ";vtx y [mm]", 800, -1, 1 );
123 m_hr_rec = new TH1F( "r_rec", ";vtx r [mm]", 40, vrbins );
124
125 // different binning for primary and sec vtx analysis
127 m_hntrax = new TH1F( "ntrax", ";number of tracks", 14, -0.5, 13.5 );
128 m_hntrax_rec = new TH1F( "ntrax_rec", ";number of tracks", 14, -0.5, 13.5 );
129 m_h_dntrax = new TH1F( "dntrax", ";trigger - offline tracks", 61, -30.5, 30.5 );
130 }
131 else {
132 m_hntrax = new TH1F( "ntrax", ";number of tracks", 80, vnbins );
133 m_hntrax_rec = new TH1F( "ntrax_rec", ";number of tracks", 80, vnbins );
134 m_h_dntrax = new TH1F( "dntrax", ";trigger - offline tracks", 15, -7.5, 7.5 );
135 }
136
137
138 m_hzed_res = new TH1F( "zed_res", "Delta z [mm]", 400, -10, 10 );
139 m_hx_res = new TH1F( "x_res", "Delta x [mm]", 200, -0.1, 0.1 );
140 m_hy_res = new TH1F( "y_res", "Delta y [mm]", 200, -0.1, 0.1 );
141
142 m_rdntrax_vs_zed = new Resplot( "rdntrax_vs_zed", 100, -300, 300, 61, -30.5, 30.5);
143 m_rdntrax_vs_ntrax = new Resplot( "rdntrax_vs_ntrax", 14, -0.5, 13.5, 15, -30.5, 30.5);
144 m_rdntrax_vs_r = new Resplot( "rdntrax_vs_r", 40, vrbins, 15, -7.5, 7.5 );
145
146 m_rdz_vs_zed = new Resplot( "rdz_vs_zed", 100, -300, 300, 1600, -20, 20 );
147 m_rdz_vs_ntrax = new Resplot( "rdz_vs_ntrax", 201, -0.5, 200.5, 1600, -20, 20 );
148 m_rdz_vs_nvtx = new Resplot( "rdz_vs_nvtx", 81, -0.5, 80.5, 1600, -20, 20 );
149 m_rdz_vs_mu = new Resplot( "rdz_vs_mu", 30, 0, 30, 1600, -20, 20 );
150 m_rdz_vs_r = new Resplot( "rdz_vs_r", 40, vrbins, 1600, -20, 20 );
151
152 m_rdr_vs_zed = new Resplot( "rdr_vs_zed", 100, -300, 300, 800, -15, 15 );
153 m_rdr_vs_r = new Resplot( "rdr_vs_r", 40, vrbins, 800, -15, 15 );
154 m_rdr_vs_ntrax = new Resplot( "rdr_vs_ntrax", 14, -0.5, 13.5, 800, -15, 15 );
155
156 m_eff_zed = new Efficiency1D( m_hzed, "zed_eff" );
157 m_eff_x = new Efficiency1D( m_hx, "x_eff" );
158 m_eff_y = new Efficiency1D( m_hy, "y_eff" );
159 m_eff_ntrax = new Efficiency1D( m_hntrax, "ntrax_eff" );
160 m_eff_nvtx = new Efficiency1D( m_hnvtx, "nvtx_eff" );
161 m_eff_mu = new Efficiency1D( m_hmu, "mu_eff" );
162 m_eff_lb = new Efficiency1D( m_hlb, "lb_eff" );
163 m_eff_r = new Efficiency1D( m_hr, "r_eff" );
164
165 m_rnvtxrec_nvtx = new Resplot( "rnvtxrec_vs_nvtx", 81, -0.5, 80.5, 81, -0.5, 80.5 );
166
167 // double ntrax[10] = { 0, 2, 5, 10, 15, 20, 30, 15, 100 };
168
169 m_rdz_vs_lb = new Resplot( "rdz_vs_lb", 301, -0.5, 3009.5, 800, -20, 20 );
170 m_rdx_vs_lb = new Resplot( "rdx_vs_lb", 301, -0.5, 3009.5, 300, -3, 3 );
171 m_rdy_vs_lb = new Resplot( "rdy_vs_lb", 301, -0.5, 3009.5, 300, -3, 3 );
172
173 m_dir->pop();
174
175}
176
177template<typename T>
178std::ostream& operator<<( std::ostream& s, const std::vector<T*>& v) {
179 for ( size_t i=0 ; i<v.size() ; i++ ) {
180 if ( v[i] ) s << "\t" << *v[i] << "\n";
181 else s << "\t" << "----\n";
182 }
183 return s;
184}
185
186void ConfVtxAnalysis::execute( const std::vector<TIDA::Vertex*>& vtx0,
187 const std::vector<TIDA::Vertex*>& vtx1,
188 const TIDA::Event* tevt ) {
189
190
191 if ( !m_initialised ) return;
192
193 // if ( vtx1.size()<2 ) return;
194
195#if 0
196 std::cout << "ConfVtxAnalysis::execute() " << name()
197 << "\tvtx0.size() " << vtx0.size()
198 << "\tvtx1.size() " << vtx1.size()
199 << std::endl;
200
201 std::cout << "\tref: " << vtx0 << std::endl;
202 std::cout << "\ttest: " << vtx1 << std::endl;
203#endif
204
205 // use new matcher if tracks included
206 if ( vtx0[0]->tracks().size() > 0 ) {
207 VertexNewMatcher m( "vtx_matcher", 0.5 );
208 execute_internal<VertexNewMatcher>(vtx0, vtx1, m, tevt);
209 }
210 else {
211 VertexMatcher m( "vtx_matcher", 3 );
212 execute_internal<VertexMatcher>(vtx0, vtx1, m, tevt);
213 }
214}
215
216
217template<typename Matcher>
218void ConfVtxAnalysis::execute_internal( const std::vector<TIDA::Vertex*>& vtx0,
219 const std::vector<TIDA::Vertex*>& vtx1,
220 Matcher& m,
221 const TIDA::Event* tevt ) {
222
223 m.match( vtx0, vtx1 );
224
225 m_hnvtx->Fill( vtx0.size() );
226 m_hnvtx_rec->Fill( vtx1.size() );
227
228 m_rnvtxrec_nvtx->Fill( vtx0.size(), vtx1.size() );
229
230 // keep alternate functionality commented
231 // std::cout << "gevent " << gevent << std::endl;
232
234 // double mu = gevent->mu();
235 // double lb = gevent->lumi_block();
236 double mu = tevt->mu();
237 double lb = tevt->lumi_block();
238
239 for ( unsigned i=0 ; i<vtx0.size() ; i++ ) {
240
242 if ( vtx0[i]->Ntracks() == 0 ) continue;
243
244 // offline vertex radial position
245 double r = std::sqrt(vtx0[i]->x()*vtx0[i]->x() + vtx0[i]->y()*vtx0[i]->y());
246
247 m_hx->Fill( vtx0[i]->x() );
248 m_hy->Fill( vtx0[i]->y() );
249 m_hzed->Fill( vtx0[i]->z() );
250 m_hntrax->Fill( vtx0[i]->Ntracks() );
251 m_hr->Fill( r );
252
253 m_hlb->Fill( lb );
254 m_hmu->Fill( mu );
255
256 const TIDA::Vertex* mv = m.matched( vtx0[i] );
257
258 // std::cout << "\tvtx match: " << i << " " << mv << std::endl;
259 if ( mv ) {
260 // std::cout << "\ttest z " << mv->z() << " : delta z " << (mv->z()-vtx0[i]->z()) << std::endl;
261
264
265 // online vertex radial position
266 double r_rec = std::sqrt(mv->x()*mv->x() + mv->y()*mv->y());
267
268 int dntrax = mv->Ntracks() - vtx0[i]->Ntracks();
269
270 m_hx_rec->Fill( mv->x() );
271 m_hy_rec->Fill( mv->y() );
272 m_hzed_rec->Fill( mv->z() );
273 m_hntrax_rec->Fill( mv->Ntracks() );
274 m_hr_rec->Fill( r_rec );
275
276 m_hx_res->Fill( mv->x() - vtx0[i]->x() );
277 m_hy_res->Fill( mv->y() - vtx0[i]->y() );
278 m_hzed_res->Fill( mv->z() - vtx0[i]->z() );
279
280 m_h_dntrax->Fill( dntrax );
281
282 m_rdz_vs_zed->Fill( vtx0[i]->z(), mv->z() - vtx0[i]->z() );
283 m_rdz_vs_ntrax->Fill( vtx0[i]->Ntracks(), mv->z() - vtx0[i]->z() );
284 m_rdz_vs_nvtx->Fill( vtx0.size(), mv->z() - vtx0[i]->z() );
285 m_rdz_vs_mu->Fill( mu, mv->z() - vtx0[i]->z() );
286 m_rdz_vs_r->Fill( r, mv->z() - vtx0[i]->z() );
287
288 m_rdntrax_vs_zed->Fill( vtx0[i]->z(), dntrax );
289 m_rdntrax_vs_ntrax->Fill( vtx0[i]->Ntracks(), dntrax );
290 m_rdntrax_vs_r->Fill( r, dntrax );
291
292 m_rdr_vs_zed->Fill( vtx0[i]->z(), r_rec - r );
293 m_rdr_vs_r->Fill( r, r_rec - r );
294 m_rdr_vs_ntrax->Fill( vtx0[i]->Ntracks(), r_rec - r );
295
296 m_eff_zed->Fill( vtx0[i]->z() );
297 m_eff_x->Fill( vtx0[i]->x() );
298 m_eff_y->Fill( vtx0[i]->y() );
299 m_eff_r->Fill( r );
300
301 m_eff_ntrax->Fill( vtx0[i]->Ntracks() );
302 m_eff_nvtx->Fill( vtx0.size() );
303
304 m_eff_mu->Fill( mu );
305 m_eff_lb->Fill( lb );
306
307 // std::cout << "found vtx ref vertex size " << vtx0.size() << "\tonline " << vtx1.size() << std::endl;
308 // std::cout << "\tref: " << *vtx0[i] << std::endl;
309 // for ( unsigned iv=0 ; iv<vtx1.size() ; iv++ ) if ( vtx1[iv] ) std::cout << "\t" << iv << " : " << *vtx1[iv] << std::endl;
310
312 m_rdz_vs_lb->Fill( lb, mv->z() - vtx0[i]->z() );
313 m_rdx_vs_lb->Fill( lb, mv->x() - vtx0[i]->x() );
314 m_rdy_vs_lb->Fill( lb, mv->y() - vtx0[i]->y() );
315
316 } else {
317 // std::cout << "\t" << "------" << std::endl;
318
319 // std::cout << "missing vtx ref vertex size " << vtx0.size() << "\tonline " << vtx1.size() << std::endl;
320 // std::cout << "\tref: " << *vtx0[i] << std::endl;
321 // for ( unsigned iv=0 ; iv<vtx1.size() ; iv++ ) if ( vtx1[iv] ) std::cout << "\t" << iv << " : " << *vtx1[iv] << std::endl;
322
323#if 0
324 static int vtxcount=0;
325
326 if ( vtxcount<100 ) {
327 std::cout << "ConfVtxAnalysis::execute() " << name() << "\tnomatch\n"
328 << "\tvtx0.size() " << vtx0.size()
329 << "\tvtx1.size() " << vtx1.size()
330 << std::endl;
331
332 std::cout << "\tref: " << vtx0 << std::endl;
333 std::cout << "\ttest: " << vtx1 << std::endl;
334
335 vtxcount++;
336
337 // std::cout << *gevent << std::endl;
338 }
339
340
341#endif
342
343
344 m_eff_x->FillDenom( vtx0[i]->x() );
345 m_eff_y->FillDenom( vtx0[i]->y() );
346 m_eff_zed->FillDenom( vtx0[i]->z() );
347 m_eff_r->FillDenom( r );
348
349 m_eff_ntrax->FillDenom( vtx0[i]->Ntracks() );
350 m_eff_nvtx->FillDenom( vtx0.size() );
351
352 m_eff_mu->FillDenom( mu );
353 m_eff_lb->FillDenom( lb );
354 }
355 }
356}
357
359
360 if ( !m_initialised ) return;
361 if ( m_finalised ) return;
362
363 std::cout << "ConfVtxAnalysis::finalise() " << name() << std::endl;
364
365 m_finalised = true;
366
367 m_dir->push();
368
369 m_hnvtx->Write();
370 m_hzed->Write();
371 m_hntrax->Write();
372 m_hr->Write();
373
374 m_hnvtx_rec->Write();
375 m_hzed_rec->Write();
376 m_hntrax_rec->Write();
377 m_hr_rec->Write();
378
379 m_hmu->Write();
380 m_hlb->Write();
381 m_h_dntrax->Write();
382 m_hzed_res->Write();
383
384 std::cout << "finalising resplots" << std::endl;
385
386 m_rdz_vs_zed->Finalise( Resplot::FitNull95 ); m_rdz_vs_zed->Write();
387 m_rdz_vs_ntrax->Finalise( Resplot::FitNull95 ); m_rdz_vs_ntrax->Write();
388 m_rdz_vs_nvtx->Finalise( Resplot::FitNull95 ); m_rdz_vs_nvtx->Write();
389 m_rdz_vs_mu->Finalise( Resplot::FitNull95 ); m_rdz_vs_mu->Write();
390 m_rdz_vs_r->Finalise( Resplot::FitNull95 ); m_rdz_vs_r->Write();
391
394 m_rdntrax_vs_r->Finalise( Resplot::FitNull95 ); m_rdntrax_vs_r->Write();
395
396 m_rdr_vs_zed->Finalise( Resplot::FitNull95 ); m_rdr_vs_zed->Write();
397 m_rdr_vs_r->Finalise( Resplot::FitNull95 ); m_rdr_vs_r->Write();
398 m_rdr_vs_ntrax->Finalise( Resplot::FitNull95 ); m_rdr_vs_ntrax->Write();
399
401
402 m_eff_zed->finalise(); m_eff_zed->Bayes()->Write( (m_eff_zed->name()+"_tg").c_str() );
403 m_eff_x->finalise(); m_eff_x->Bayes()->Write( (m_eff_x->name()+"_tg").c_str() );
404 m_eff_y->finalise(); m_eff_y->Bayes()->Write( (m_eff_y->name()+"_tg").c_str() );
405 m_eff_r->finalise(); m_eff_r->Bayes()->Write( (m_eff_r->name()+"_tg").c_str());
406
407 m_eff_ntrax->finalise(); m_eff_ntrax->Bayes()->Write( (m_eff_ntrax->name()+"_tg").c_str() );
408 m_eff_nvtx->finalise(); m_eff_nvtx->Bayes()->Write( (m_eff_nvtx->name()+"_tg").c_str() );
409
410 m_eff_mu->finalise(); m_eff_mu->Bayes()->Write( (m_eff_mu->name()+"_tg").c_str() );
411 m_eff_lb->finalise(); m_eff_lb->Bayes()->Write( (m_eff_lb->name()+"_tg").c_str() );
412
413
414 m_rdx_vs_lb->Finalise( Resplot::FitNull95 ); m_rdx_vs_lb->Write();
415 m_rdy_vs_lb->Finalise( Resplot::FitNull95 ); m_rdy_vs_lb->Write();
416 m_rdz_vs_lb->Finalise( Resplot::FitNull95 ); m_rdz_vs_lb->Write();
417
418 m_dir->pop();
419}
420
std::ostream & operator<<(std::ostream &s, const std::vector< T * > &v)
TIDA::Event * gevent
Definition globals.cxx:13
Basic event class to contain a vector of chains for trigger analysis.
#define y
#define x
#define z
ConfVtxAnalysis(const std::string &n, bool use_secVtx_limits=false)
Resplot * m_rdr_vs_zed
Efficiency1D * m_eff_x
void execute_internal(const std::vector< TIDA::Vertex * > &vtx0, const std::vector< TIDA::Vertex * > &vtx1, Matcher &m, const TIDA::Event *tevt=0)
Efficiency1D * m_eff_zed
Resplot * m_rdntrax_vs_r
Efficiency1D * m_eff_mu
Resplot * m_rdz_vs_nvtx
Resplot * m_rdz_vs_ntrax
Efficiency1D * m_eff_y
Resplot * m_rnvtxrec_nvtx
void execute(const std::vector< TIDA::Vertex * > &vtx0, const std::vector< TIDA::Vertex * > &vtx1, const TIDA::Event *tevt=0)
Efficiency1D * m_eff_lb
Resplot * m_rdz_vs_zed
TIDDirectory * m_dir
Resplot * m_rdr_vs_ntrax
Efficiency1D * m_eff_nvtx
Efficiency1D * m_eff_ntrax
Resplot * m_rdntrax_vs_zed
Efficiency1D * m_eff_r
Resplot * m_rdz_vs_mu
Resplot * m_rdntrax_vs_ntrax
static TF1 * FitNull95(TH1D *s, double a=0, double b=0)
Definition Resplot.cxx:1672
void lumi_block(unsigned lb)
Definition TIDAEvent.h:44
void mu(double m)
Definition TIDAEvent.h:47
double y() const
Definition TIDAVertex.h:52
int Ntracks() const
Definition TIDAVertex.h:62
double z() const
Definition TIDAVertex.h:53
double x() const
Definition TIDAVertex.h:51
VertexAnalysis(const std::string &n)
const std::string & name() const
int lb
Definition globals.cxx:23
int r
Definition globals.cxx:22