ATLAS Offline Software
Loading...
Searching...
No Matches
T_AnalysisConfig_Tier0.h
Go to the documentation of this file.
1/* emacs: this is -*- c++ -*- */
22
23#ifndef TrigInDetAnalysisExample_T_AnalysisConfig_Tier0_H
24#define TrigInDetAnalysisExample_T_AnalysisConfig_Tier0_H
25
26
30
35
36#include "TTree.h"
37#include "TFile.h"
38
39
40// McParticleEvent includes
42
44#include "AtlasHepMC/GenEvent.h"
47
48#include "EventInfo/EventInfo.h"
49#include "EventInfo/EventID.h"
51
52
53
55
57
63
64
66
68
70
73
75
78
79
80
81
82
83template<typename T>
85
86public:
87
88 // Full constructor: test/reference/selection
89 // - analysisInstanceName: the name of the analysis chain being created
90 // - xxxChainName: the name of the chain to be used as test/reference/selection; must be "StoreGate" in case of direct access to SG containers
91 // - xxxType: the type of tracks to be retrieved from the test/reference/selection chain or container
92 // - xxxKey: the key for tracks to be retrieved from the test/reference/selection chain or container
93 // - all standard operations are performed in loops over 0=test 1=reference 2=selection
94 T_AnalysisConfig_Tier0(const std::string& analysisInstanceName,
95 const std::string& testChainName, const std::string& testType, const std::string& testKey,
96 const std::string& referenceChainName, const std::string& referenceType, const std::string& referenceKey,
97 TrackFilter* testFilter, TrackFilter* referenceFilter,
98 TrackAssociator* associator,
100 T_AnalysisConfig<T>( analysisInstanceName,
101 testChainName, testType, testKey,
102 referenceChainName, referenceType, referenceKey,
103 testFilter, referenceFilter,
104 associator,
105 analysis ),
106 _analysis(0),
107 m_useBeamCondSvc(false),
108 m_doOffline(true),
109 m_doMuons(false),
110 m_doElectrons(false),
111 m_doTaus(false),
112 m_doBjets(false),
113 m_hasTruthMap(false),
114 m_doTauThreeProng(false),
115 m_tauEtCutOffline(false),
116 m_NRois(0),
117 m_NRefTracks(0),
118 m_NTestTracks(0),
119 m_runPurity(false),
120 m_shifter(false),
121 m_pTthreshold(0),
122 m_first(true),
123 m_containTracks(false)
124 {
125 m_event = new TIDA::Event();
126 m_chainNames.push_back(testChainName);
127
128#if 0
129 ChainString& chain = m_chainNames.back();
130
131 std::cout << "\nT_AnalysisConfig_Tier0::name: " << name() << "\t" << this << std::endl;
132 std::cout << "T_AnalysisConfig_Tier0::chain specification: " << testChainName << " -> " << chain << "\t" << chain.raw() << std::endl;
133 std::cout << "\tchain: " << chain.head() << std::endl;
134 std::cout << "\tkey: " << chain.tail() << std::endl;
135 std::cout << "\troi: " << chain.roi() << std::endl;
136 std::cout << "\tvtx: " << chain.vtx() << std::endl;
137 std::cout << "\tte: " << chain.element() << std::endl;
138
139 std::cout << "\tpost: " << chain.post() << std::endl;
140 std::cout << "\tpt: " << chain.postvalue("pt") << std::endl;
141
142 std::cout << "\tcontainTracks: " << m_containTracks << std::endl;
143#endif
144
145 m_testType = testType;
146 }
147
148 virtual ~T_AnalysisConfig_Tier0() { delete m_event; }
149
150 void setRunPurity( bool b ) { m_runPurity=b; }
151
152 void setShifter( bool b ) { m_shifter=b; }
153
154 void useBeamCondSvc( bool b ) { m_useBeamCondSvc = b; }
155
156 void containTracks( bool b ) { m_containTracks = b; }
157
158public:
159
161
163 // void verbose( std::ostream& s ) {
164 // if( m_provider->msg().level() <= MSG::VERBOSE ) m_provider->msg(MSG::VERBOSE) << s << endmsg;
165 // }
166
167protected:
168
172
173 using T_AnalysisConfig<T>::name;
175
180
181 // using T_AnalysisConfig<T>::selectTracks;
182
183 virtual void loop() {
184
185 if( m_provider->msg().level() <= MSG::VERBOSE) {
186 m_provider->msg(MSG::VERBOSE) << "AnalysisConfig_Tier0::loop() for " << T_AnalysisConfig<T>::m_analysisInstanceName << endmsg;
187 }
188
189 // get (offline) beam position
190 double xbeam = 0;
191 double ybeam = 0;
192
193 if ( m_first ) {
194
195 m_first = false;
196
197 m_provider->msg(MSG::VERBOSE) << " using beam position\tx=" << xbeam << "\ty=" << ybeam << endmsg;
198
199 if (m_provider->msg().level() <= MSG::VERBOSE) {
200
201 std::vector<std::string> configuredChains = (*(m_tdt))->getListOfTriggers("L2_.*, EF_.*, HLT_.*");
202
203 for ( unsigned i=0 ; i<configuredChains.size() ; i++ ) {
204 m_provider->msg(MSG::VERBOSE) << "Chain " << configuredChains[i] << endmsg;
205 }
206
207 }
208
209
210 std::vector<ChainString>::iterator chainitr = m_chainNames.begin();
211
212 std::vector<ChainString> chains;
213
217 while ( chainitr!=m_chainNames.end() ) {
218
220 ChainString& chainName = (*chainitr);
221
222 m_provider->msg(MSG::INFO) << "process chain " << chainName << endmsg;
223
224 if ( chainName.head() == "" ) {
225
226 std::string selectChain = chainName.raw();
227
228 chains.push_back( ChainString(selectChain) );
229
230 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
231 m_provider->msg(MSG::VERBOSE) << "Matching chain " << selectChain << " (" << chainName.head() << ")" << endmsg;
232 }
233
234 }
235 else {
236
238 std::vector<std::string> selectChains = (*(m_tdt))->getListOfTriggers( chainName.head() );
239
240
241 for ( unsigned iselected=0 ; iselected<selectChains.size() ; iselected++ ) {
242
243 selectChains[iselected] = chainName.subs( selectChains[iselected] );
244
245#if 0
246 std::cout << "sorting:: chain specification: " << chainName << "\traw:" << chainName.raw() << std::endl;
247 std::cout << "\tchain: " << chainName.head() << std::endl;
248 std::cout << "\tkey: " << chainName.tail() << std::endl;
249 std::cout << "\troi: " << chainName.roi() << std::endl;
250 std::cout << "\tvtx: " << chainName.vtx() << std::endl;
251 std::cout << "\tte: " << chainName.element() << std::endl;
252 std::cout << "\tind: " << chainName.extra() << std::endl;
253#endif
254
256 chains.push_back( ChainString(selectChains[iselected]) );
257
258 if(m_provider->msg().level() <= MSG::VERBOSE) {
259 m_provider->msg(MSG::VERBOSE) << "Matching chain " << selectChains[iselected] << " (" << chainName.head() << ")" << endmsg;
260 }
261
262 }
263 }
264
265 ++chainitr;
266 }
267
268 // m_chainNames.insert( m_chainNames.end(), chains.begin(), chains.end() );
269 m_chainNames = chains;
270
271 for ( unsigned ic=0 ; ic<m_chainNames.size() ; ic++ ) m_provider->msg(MSG::VERBOSE) << "Analyse chain " << m_chainNames[ic] << endmsg;
272
273 }
274
275
277 Filter_True filter;
278
279 Filter_etaPT filter_etaPT(5,200);
280 Filter_Combined filter_truth( &filter_etaPT, &filter_etaPT);
281
283 // Filter_Combined filterRef (&filter_offline, &filter_vertex);
284
285 int iRefFilter = 1;
286 int iTestFilter = 0;
287
288 if ( m_runPurity ) {
289 iRefFilter = 0;
290 iTestFilter = 1;
291 }
292
293 Filter_Combined filterRef( m_filters[iRefFilter][0], &filter );
294 Filter_Combined filterTest( m_filters[iTestFilter][0], &filter );
295
296 TrigTrackSelector selectorTruth( &filter_truth );
297 TrigTrackSelector selectorRef( &filterRef );
298 m_selectorRef = &selectorRef;
299 TrigTrackSelector selectorTest( &filterTest );
300 m_selectorTest = &selectorTest;
301
302 if ( xbeam!=0 || ybeam!=0 ) {
303 m_selectorRef->setBeamline( xbeam, ybeam );
304 }
305
307
308 // clear the ntuple TIDA::Event class
309 m_event->clear();
310
312#ifndef XAODTRACKING_TRACKPARTICLE_H
313 const EventInfo* pEventInfo;
314#else
315 const xAOD::EventInfo* pEventInfo;
316#endif
317 unsigned run_number = 0;
318 uint64_t event_number = 0;
319 unsigned lumi_block = 0;
320 unsigned bunch_crossing_id = 0;
321 unsigned time_stamp = 0;
322 double mu_val = 0;
323
324 if ( this->template retrieve(pEventInfo, "EventInfo").isFailure() ) {
325 m_provider->msg(MSG::WARNING) << "Failed to get EventInfo " << endmsg;
326 } else {
327
328#ifndef XAODTRACKING_TRACKPARTICLE_H
329 run_number = pEventInfo->event_ID()->run_number();
330 event_number = pEventInfo->event_ID()->event_number();
331 lumi_block = pEventInfo->event_ID()->lumi_block();
332 time_stamp = pEventInfo->event_ID()->time_stamp();
333 bunch_crossing_id = pEventInfo->event_ID()->bunch_crossing_id();
334 mu_val = pEventInfo->averageInteractionsPerCrossing();
335#else
336 run_number = pEventInfo->runNumber();
337 event_number = pEventInfo->eventNumber();
338 lumi_block = pEventInfo->lumiBlock();
339 time_stamp = pEventInfo->timeStamp();
340 bunch_crossing_id = pEventInfo->bcid();
341 mu_val = pEventInfo->averageInteractionsPerCrossing();
342#endif
343 }
344
345 if(m_provider->msg().level() <= MSG::VERBOSE){
346 m_provider->msg(MSG::VERBOSE) << "run " << run_number
347 << "\tevent " << event_number
348 << "\tlb " << lumi_block << endmsg;
349 }
350
351 // m_provider->msg(MSG::INFO) << "run " << run_number
352 // << "\tevent " << event_number
353 // << "\tlb " << lumi_block << endmsg;
354
355 // std::cout << "run " << run_number << "\tevent " << event_number << "\tlb " << lumi_block << std::endl;
356
357
358 // clear the ntuple TIDA::Event class
359 m_event->clear();
360
361 m_event->run_number(run_number);
362 m_event->event_number(event_number);
363 m_event->lumi_block(lumi_block);
364 m_event->time_stamp(time_stamp);
365 m_event->bunch_crossing_id(bunch_crossing_id);
366 m_event->mu(mu_val);
367
370
371 bool analyse = false;
372
373 // Check HLTResult
374
375 for ( unsigned ichain=0 ; ichain<m_chainNames.size() ; ichain++ ) {
376
377 const std::string& chainname = m_chainNames[ichain].head();
378
379 if ( chainname == "" ) analyse = true;
380 else {
381
382 //Only for trigger chains
383 if ( chainname.find("L2") == std::string::npos &&
384 chainname.find("EF") == std::string::npos &&
385 chainname.find("HLT") == std::string::npos ) continue;
386
387 if ( m_provider->msg().level() <= MSG::DEBUG ) {
388 m_provider->msg(MSG::DEBUG) << "Chain " << chainname
389 << "\tpass " << (*m_tdt)->isPassed(chainname)
390 << "\tpres " << (*m_tdt)->getPrescale(chainname) << endmsg;
391 }
392
393 // std::cout << "Chain " << chainname << "\tpass " << (*m_tdt)->isPassed(chainname)
394 // << "\tpres " << (*m_tdt)->getPrescale(chainname) << std::endl;
395
396 if ( (*(m_tdt))->isPassed(chainname) || (*(m_tdt))->getPrescale(chainname) ) analyse = true;
397
398 }
399 }
400
401 // Remove this code to skip on truncated HLT results
402 // Need to leave the code here for the time being however, since we will
403 // still need in the future a more robust test to achieve this same
404 // functionality
405 //
406 // if ( (*m_tdt)->ExperimentalAndExpertMethods().isHLTTruncated() ) {
407 // m_provider->msg(MSG::WARNING) << "HLTResult truncated, skipping event" << endmsg;
408 // return;
409 // }
410
411 if ( !this->m_keepAllEvents && !analyse ) {
412 // m_provider->msg(MSG::VERBOSE) << "No chains passed unprescaled - not processing this event" << endmsg;
413 if(m_provider->msg().level() <= MSG::VERBOSE)
414 m_provider->msg(MSG::VERBOSE) << "No chains passed unprescaled - not processing this event" << endmsg;
415 return;
416 }
417
419
420 selectorTruth.clear();
421
422 if(m_provider->msg().level() <= MSG::VERBOSE)
423 m_provider->msg(MSG::VERBOSE) << "MC Truth flag " << m_mcTruth << endmsg;
424
425 if ( m_mcTruth && m_provider->msg().level() <= MSG::VERBOSE ) m_provider->msg(MSG::VERBOSE) << "getting Truth" << endmsg;
426
427
429
430 std::vector<TIDA::Vertex> vertices;
431 std::vector<TIDA::Vertex> vertices_rec;
432
433 std::vector<double> refbeamspot;
434 std::vector<double> testbeamspot;
435
436
437
438#ifndef XAODTRACKING_TRACKPARTICLE_H
439
440 const VxContainer* primaryVtxCollection;
441
442 if ( m_doOffline ) {
443 if ( m_provider->evtStore()->template contains<VxContainer>("VxPrimaryCandidate") ) {
444 if ( this->template retrieve(primaryVtxCollection, "VxPrimaryCandidate").isFailure()) {
445 if (m_provider->msg().level() <= MSG::WARNING) m_provider->msg(MSG::WARNING) << "Primary vertex container not found" << endmsg;
446 }
447 else {
448 VxContainer::const_iterator vtxitr = primaryVtxCollection->begin();
449 for ( ; vtxitr != primaryVtxCollection->end(); ++vtxitr) {
450 if ( (*vtxitr)->vxTrackAtVertex()->size()>0 ) {
451 vertices.push_back( TIDA::Vertex( (*vtxitr)->recVertex().position().x(),
452 (*vtxitr)->recVertex().position().y(),
453 (*vtxitr)->recVertex().position().z(),
454 0,0,0,
455 (*vtxitr)->vxTrackAtVertex()->size() ) );
456 }
457 }
458 }
459
460 // filter_vertex.setVertex(vertices);
461 }
462 }
463
464#else
465
466 //std::vector<TIDA::Vertex> vertices;
467
468 m_provider->msg(MSG::VERBOSE) << "fetching AOD Primary vertex container" << endmsg;
469
470 const xAOD::VertexContainer* xaodVtxCollection = 0;
471
472 if ( this->template retrieve( xaodVtxCollection, "PrimaryVertices" ).isFailure()) {
473 if (m_provider->msg().level() <= MSG::WARNING) m_provider->msg(MSG::WARNING) << "xAOD Primary vertex container not found with key " << "PrimaryVertices" << endmsg;
474 }
475
476 if ( xaodVtxCollection!=0 ) {
477
478 m_provider->msg(MSG::VERBOSE) << "xAOD Primary vertex container " << xaodVtxCollection->size() << " entries" << endmsg;
479
480 xAOD::VertexContainer::const_iterator vtxitr = xaodVtxCollection->begin();
481 for ( ; vtxitr != xaodVtxCollection->end(); ++vtxitr ) {
482 if ( (*vtxitr)->nTrackParticles()>0 && (*vtxitr)->vertexType()!=0 ) {
483 vertices.push_back( TIDA::Vertex( (*vtxitr)->x(),
484 (*vtxitr)->y(),
485 (*vtxitr)->z(),
487 (*vtxitr)->covariancePosition()(Trk::x,Trk::x),
488 (*vtxitr)->covariancePosition()(Trk::y,Trk::y),
489 (*vtxitr)->covariancePosition()(Trk::z,Trk::z),
490 (*vtxitr)->nTrackParticles(),
492 (*vtxitr)->chiSquared(),
493 (*vtxitr)->numberDoF() ) );
494 }
495 }
496 }
497
498
499#endif
500
501
503
504 if ( m_mcTruth ) {
505 m_event->addChain( "Truth" );
506 m_event->back().addRoi(TIDARoiDescriptor());
507 m_event->back().back().addTracks(selectorTruth.tracks());
508 }
509
511
512 if ( m_doOffline ) {
513 for ( unsigned i=0 ; i<vertices.size() ; i++ ) {
514 if(m_provider->msg().level() <= MSG::VERBOSE)
515 m_provider->msg(MSG::VERBOSE) << "vertex " << i << " " << vertices[i] << endmsg;
516 m_event->addVertex(vertices[i]);
517 }
518 }
519
521
522 // int Noff = 0;
523 std::vector<TIDA::Track*> offline_tracks;
524 std::vector<TIDA::Track*> electron_tracks;
525 std::vector<TIDA::Track*> muon_tracks;
526
527 std::vector<TIDA::Track*> ref_tracks;
528 std::vector<TIDA::Track*> test_tracks;
529
530 offline_tracks.clear();
531 electron_tracks.clear();
532 muon_tracks.clear();
533
534 ref_tracks.clear();
535 test_tracks.clear();
536
538 for ( unsigned ichain=0 ; ichain<m_chainNames.size() ; ichain++ ) {
539
540 test_tracks.clear();
541
543
544 // std::string& chainname = chains[ichain];
545 const std::string& chainname = m_chainNames[ichain].head();
546 const std::string& key = m_chainNames[ichain].tail();
547 const std::string& vtx_name = m_chainNames[ichain].vtx();
548 //no currently used but retained in case
549 //const std::string& roi_name = m_chainNames[ichain].roi();
550 //const std::string& te_name = m_chainNames[ichain].element();
551
552 m_pTthreshold = 0;
553
554 if ( m_chainNames[ichain].postcount() ) {
555 std::string ptvalue = m_chainNames[ichain].postvalue("pt");
556 if ( ptvalue!="" ) m_pTthreshold = std::stod(ptvalue);
557 }
558
559 // std::cout << "\tchain " << m_chainNames[ichain] << "\tchainname " << chainname << "\tvtx " << vtx_name << "\troi " << roi_name << std::endl;
560
561 unsigned _decisiontype = TrigDefs::Physics;
562 unsigned decisiontype;
563
564 if ( this->requireDecision() ) _decisiontype = TrigDefs::requireDecision;
565
566
567 if ( m_chainNames[ichain].passed() ) decisiontype = _decisiontype;
568 else decisiontype = TrigDefs::alsoDeactivateTEs;
569
570 if ( chainname!="" && m_provider->msg().level() <= MSG::VERBOSE ) {
571
572 m_provider->msg(MSG::VERBOSE) << "status for chain " << chainname
573 << "\tpass " << (*m_tdt)->isPassed(chainname)
574 << "\tprescale " << (*m_tdt)->getPrescale(chainname) << endmsg;
575
576 m_provider->msg(MSG::VERBOSE) << "fetching features for chain " << chainname << endmsg;
577
578 m_provider->msg(MSG::VERBOSE) << chainname << "\tpassed: " << (*m_tdt)->isPassed( chainname ) << endmsg;
579 }
580
581 // std::cout << "\tstatus for chain " << chainname
582 // << "\tpass " << (*m_tdt)->isPassed( chainname )
583 // << "\tpassdt " << (*m_tdt)->isPassed( chainname, decisiontype )
584 // << "\tprescale " << (*m_tdt)->getPrescale( chainname ) << std::endl;
585
586
587 // m_provider->msg(MSG::INFO) << chainname << "\tpassed: " << (*m_tdt)->isPassed( chainname ) << "\t" << m_chainNames[ichain] << "\trun " << run_number << "\tevent " << event_number << endmsg;
588
589
590 if ( chainname!="" && !this->m_keepAllEvents && !(*m_tdt)->isPassed( chainname, decisiontype ) ) continue;
591
594 // Trig::FeatureContainer f = (*m_tdt)->features( chainname, TrigDefs::alsoDeactivateTEs);
595
598
599 ChainString& chainConfig = m_chainNames[ichain];
600
601 std::string chainName = chainConfig.head();
602
603 m_event->addChain( chainConfig );
604
605 TIDA::Chain& chain = m_event->back();
606
607 if ( chainName == "" ) {
608
609 m_selectorTest->clear();
610
612
613 TIDARoiDescriptor* roiInfo = new TIDARoiDescriptor(true);
614
615 chain.addRoi( *roiInfo );
616
617
618# ifdef XAODTRACKING_TRACKPARTICLE_H
619 if ( m_provider->evtStore()->template contains<xAOD::TrackParticleContainer>(key) ) {
621 refbeamspot = this->template getBeamspot<xAOD::TrackParticleContainer>( key );
622 }
623# endif
624
625 const std::vector<TIDA::Track*>& testtracks = m_selectorTest->tracks();
626
627 chain.back().addTracks(testtracks);
628
629 delete roiInfo;
630
631 }
632 else {
633
634 Trig::FeatureContainer f = (*m_tdt)->features( chainname, decisiontype );
635 Trig::FeatureContainer::combination_const_iterator c(f.getCombinations().begin());
636 Trig::FeatureContainer::combination_const_iterator cEnd(f.getCombinations().end());
637
638
639 if ( c==cEnd ) {
640 if(m_provider->msg().level() <= MSG::VERBOSE){
641 m_provider->msg(MSG::VERBOSE) << "No combinations: skipping this chain " << chainname << endmsg;
642 }
643 continue;
644 }
645
646 if(m_provider->msg().level() <= MSG::VERBOSE) {
647 m_provider->msg(MSG::VERBOSE) << "combinations for chain " << chainname << " " << (cEnd-c) << endmsg;
648 }
649
650 unsigned icomb = 0;
651
652 for( ; c!=cEnd ; ++c ) {
653
654 icomb++;
655
656 // now add rois to this ntuple chain
657
658 // Get seeding RoI
659 // std::vector< Trig::Feature<TrigRoiDescriptor> > initRois = c->get<TrigRoiDescriptor>("initialRoI", TrigDefs::alsoDeactivateTEs);
660 // std::vector< Trig::Feature<TrigRoiDescriptor> > initRois = c->get<TrigRoiDescriptor>("forID", TrigDefs::alsoDeactivateTEs);
661
662 std::vector< Trig::Feature<TrigRoiDescriptor> > initRois;
663
664 std::string roi_key = m_chainNames[ichain].roi();
665
666 if ( roi_key=="SuperRoi" && icomb>1 ) continue;
667
668 if ( roi_key!="" ) {
669 initRois = c->get<TrigRoiDescriptor>(roi_key, decisiontype );
670 }
671 else {
672 initRois = c->get<TrigRoiDescriptor>("forID", decisiontype );
673 if ( initRois.empty() ) initRois = c->get<TrigRoiDescriptor>("", decisiontype );
674 if ( initRois.empty() ) initRois = c->get<TrigRoiDescriptor>("initialRoI", decisiontype );
675 }
676
677 // std::cout << "initRois.size() " << initRois.size() << std::endl;
678
679 if ( initRois.empty() ) continue;
680
681
682 // for ( unsigned ir=0 ; ir<initRois.size() ; ir++ ) {
683 // std::cout << "\t" << ir << "\t" << *initRois[ir].cptr() << std::endl;
684 // }
685
686 // Skip chains seeded by multiple RoIs: not yet implemented
687 if(initRois.size()>1 && roi_key!="SuperRoi" ) {
688 if(m_provider->msg().level() <= MSG::VERBOSE)
689 m_provider->msg(MSG::VERBOSE) << " More than one initial RoI found for seeded chain " << chainname << ": not yet supported" << endmsg;
690 continue;
691 }
692
693 TIDARoiDescriptor* roiInfo = 0;
694
695 if( !initRois.empty() ) {
696 const TrigRoiDescriptor* roid = initRois[0].cptr();
697
698 if(m_provider->msg().level() <= MSG::VERBOSE)
699 m_provider->msg(MSG::VERBOSE) << " RoI descriptor for seeded chain " << chainname << " " << *roid << endmsg;
700
701 roiInfo = new TIDARoiDescriptor(TIDARoiDescriptorBuilder(*roid));
702 // roiInfo->etaHalfWidth(m_roiInfo->etaHalfWidth());
703 // roiInfo->phiHalfWidth(m_roiInfo->phiHalfWidth());
704 // roiInfo->etaHalfWidth(roid->etaHalfWidth());
705 // roiInfo->phiHalfWidth(roid->phiHalfWidth());
706 // roiInfo->zedHalfWidth(roid->zedHalfWidth());
707
708 if(m_provider->msg().level() <= MSG::VERBOSE)
709 m_provider->msg(MSG::VERBOSE) << "using chain roi " << *roid << endmsg;
710
711 }
712 else {
713 roiInfo = new TIDARoiDescriptor();
714 // roiInfo->etaHalfWidth(5);
715 // roiInfo->phiHalfWidth(M_PI);
716 // roiInfo->zedHalfWidth(m_roiInfo->zedHalfWidth());
717
718 if(m_provider->msg().level() <= MSG::WARNING)
719 m_provider->msg(MSG::WARNING) << "roi not found" << endmsg;
720
721 }
722
723 if(m_provider->msg().level() <= MSG::VERBOSE) m_provider->msg(MSG::VERBOSE) << *roiInfo << endmsg;
724
725 m_selectorTest->clear();
726
727 m_provider->msg(MSG::VERBOSE) << "Searching for collection " << key << endmsg;
728 // std::cout << "Searching for collection " << key << std::endl;
729
731 if ( key.find("InDetTrigParticleCreation")!=std::string::npos ||
732 key.find("_IDTrig")!=std::string::npos ||
733 key.find("_EFID")!=std::string::npos ||
734 chainName.find("EF_")!=std::string::npos ||
735 chainName.find("HLT_")!=std::string::npos ) {
736# ifdef XAODTRACKING_TRACKPARTICLE_H
737 if ( this->template selectTracks<xAOD::TrackParticleContainer>( m_selectorTest, c, key ) ) testbeamspot = this->template getBeamspot<xAOD::TrackParticleContainer>( c, key );
738 else if ( this->template selectTracks<Rec::TrackParticleContainer>( m_selectorTest, c, key ) );
739# else
740 if ( this->template selectTracks<Rec::TrackParticleContainer>( m_selectorTest, c, key ) );
741# endif
742 else if ( this->template selectTracks<TrackCollection>( m_selectorTest, c, key ) );
743 else {
744 //m_provider->msg(MSG::WARNING) << "No track collection " << key << " found" << endmsg;
745 }
746 }
747 else {
749 if ( chainName.find("L2_")!=std::string::npos ) {
750 if ( this->template selectTracks<Rec::TrackParticleContainer>( m_selectorTest, c, key ) );
751 else if ( this->template selectTracks<TrackCollection>( m_selectorTest, c, key ) );
752# ifdef XAODTRACKING_TRACKPARTICLE_H
753 else if ( this->template selectTracks<xAOD::TrackParticleContainer>( m_selectorTest, c, key ) ) testbeamspot = this->template getBeamspot<xAOD::TrackParticleContainer>( c, key );
754# endif
755 else m_provider->msg(MSG::WARNING) << "No track collection " << key << " found" << endmsg;
756 }
757 }
758
759
760 const std::vector<TIDA::Track*>& testtracks = m_selectorTest->tracks();
761
762 m_provider->msg(MSG::VERBOSE) << "test tracks.size() " << testtracks.size() << endmsg;
763 // std::cout << "test tracks.size() " << testtracks.size() << std::endl;
764
765 // std::cout << "\ttest tracks.size() " << testtracks.size() << std::endl;
766
767 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
768 m_provider->msg(MSG::VERBOSE) << "test tracks.size() " << testtracks.size() << endmsg;
769 for ( int ii=testtracks.size() ; ii-- ; ) {
770 m_provider->msg(MSG::VERBOSE) << " test track " << ii << " " << *testtracks[ii] << endmsg;
771 //test_tracks.push_back(testtracks.at(ii));
772 }
773 }
774
775 chain.addRoi( *roiInfo );
776
777 chain.back().addTracks(testtracks);
778
779 delete roiInfo;
780
782
784
785 // std::cout << "vertex " << vtx_name << "\tchain " << chainName << "\tconfig " << chainConfig << std::endl;
786
787 if ( vtx_name!="" ) {
788
789 m_provider->msg(MSG::VERBOSE) << "\tFetch xAOD::VertexContainer for chain " << chainConfig << " with key " << vtx_name << endmsg;
790
791 std::vector< Trig::Feature<xAOD::VertexContainer> > xaodtrigvertices = c->get<xAOD::VertexContainer>(vtx_name);
792
793 if ( xaodtrigvertices.empty() ) {
794 if ( m_provider->msg().level() <= MSG::DEBUG ) {
795 m_provider->msg(MSG::WARNING) << "\tNo xAOD::VertexContainer for chain " << chainConfig << " for key " << vtx_name << endmsg;
796 }
797 }
798 else {
799
800 m_provider->msg(MSG::VERBOSE) << "\txAOD::VertexContainer found with size " << xaodtrigvertices.size() << "\t" << vtx_name << endmsg;
801
802 for ( unsigned iv=0 ; iv<xaodtrigvertices.size() ; iv++ ) {
803
804 const xAOD::VertexContainer* vert = xaodtrigvertices[iv].cptr();
805
806 m_provider->msg(MSG::VERBOSE) << "\t" << iv << " xAOD VxContainer for " << chainConfig << " " << vert << " key " << vtx_name << endmsg;
807
809
810 for ( ; vtxitr != vert->end(); ++vtxitr) {
813 // if ( ( (*vtxitr)->nTrackParticles()>0 && (*vtxitr)->vertexType()!=0 ) || vtx_name=="EFHistoPrmVtx" ) {
814 if ( (*vtxitr)->vertexType()!=0 || vtx_name=="EFHistoPrmVtx" ) {
815 chain.back().addVertex( TIDA::Vertex( (*vtxitr)->x(),
816 (*vtxitr)->y(),
817 (*vtxitr)->z(),
819 (*vtxitr)->covariancePosition()(Trk::x,Trk::x),
820 (*vtxitr)->covariancePosition()(Trk::y,Trk::y),
821 (*vtxitr)->covariancePosition()(Trk::z,Trk::z),
822 (*vtxitr)->nTrackParticles(),
824 (*vtxitr)->chiSquared(),
825 (*vtxitr)->numberDoF() ) );
826
827 }
828 }
829 }
830
831 }
832
833 }
834
835 }
836
837 }
838
839
840 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
841 m_provider->msg(MSG::VERBOSE) << "event: " << *m_event << endmsg;
842 }
843
845
846 for ( unsigned iroi=0 ; iroi<chain.size() ; iroi++ ) {
847
848 m_selectorRef->clear();
849
850 if ( this->filterOnRoi() ) {
851 filterRef.setRoi( &chain.rois().at(iroi).roi() );
852 filterRef.containtracks( m_containTracks );
853 }
854 else filterRef.setRoi( 0 );
855
856 test_tracks.clear();
857
858
872
873 if ( m_provider->msg().level() <= MSG::VERBOSE )
874 m_provider->msg(MSG::VERBOSE) << "MC Truth flag " << m_mcTruth << endmsg;
875
876 bool foundTruth = false;
877
878 if ( !m_doOffline && m_mcTruth ) {
879
880 if ( this->filterOnRoi() ) filter_truth.setRoi( &chain.rois().at(iroi).roi() );
881 else filter_truth.setRoi( 0 ); // don't filter on RoI unless needed
882
883 selectorTruth.clear();
884
885 if ( m_provider->msg().level() <= MSG::VERBOSE )
886 m_provider->msg(MSG::VERBOSE) << "getting Truth" << endmsg;
887
888 if ( m_provider->evtStore()->template contains<TruthParticleContainer>("INav4MomTruthEvent") ) {
889 //ESD
890 this->template selectTracks<TruthParticleContainer>( &selectorTruth, "INav4MomTruthEvent" );
891 foundTruth = true;
892 }
893 else if ( m_provider->evtStore()->template contains<TruthParticleContainer>("SpclMC") ) {
895 this->template selectTracks<TruthParticleContainer>( &selectorTruth, "SpclMC");
896 foundTruth = true;
897 }
898 else if ( m_provider->evtStore()->template contains<TruthParticleContainer>("") ) {
900 this->template selectTracks<TruthParticleContainer>( &selectorTruth, "");
901 foundTruth = true;
902 }
903 else
904 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
905 m_provider->msg(MSG::VERBOSE) << "Truth not found - none whatsoever!" << endmsg;
906 }
907 }
908
909
910 if ( !m_doOffline && m_mcTruth && !foundTruth ) {
911
912 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
913 m_provider->msg(MSG::VERBOSE) << "getting Truth" << endmsg;
914 }
915
917
919 const McEventCollection* mcevent = 0;
920
922
923 std::string keys[4] = { "GEN_AOD", "TruthEvent", "", "G4Truth" };
924
925 std::string key = "";
926
927 bool foundcollection = false;
928
929 for ( int ik=0 ; ik<4 ; ik++ ) {
930
931 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
932 m_provider->msg(MSG::VERBOSE) << "Try McEventCollection: " << keys[ik] << endmsg;
933 }
934
935 if ( !m_provider->evtStore()->template contains<McEventCollection>(keys[ik]) ) {
936 if( m_provider->msg().level() <= MSG::VERBOSE )
937 m_provider->msg(MSG::VERBOSE) << "No McEventCollection: " << keys[ik] << endmsg;
938 continue;
939 }
940
941 if ( m_provider->msg().level() <= MSG::VERBOSE )
942 m_provider->msg(MSG::VERBOSE) << "evtStore()->retrieve( mcevent, " << keys[ik] << " )" << endmsg;
943
944 if ( this->template retrieve( mcevent, keys[ik] ).isFailure() ) {
945 if ( m_provider->msg().level() <= MSG::VERBOSE )
946 m_provider->msg(MSG::VERBOSE) << "Failed to get McEventCollection: " << keys[ik] << endmsg;
947 }
948 else {
950 key = keys[ik];
951 if(m_provider->msg().level() <= MSG::VERBOSE)
952 m_provider->msg(MSG::VERBOSE) << "Found McEventCollection: " << key << endmsg;
953 foundcollection = true;
954 break;
955 }
956 }
957
959 if ( !foundcollection ) {
960 if(m_provider->msg().level() <= MSG::VERBOSE)
961 m_provider->msg(MSG::WARNING) << "No MC Truth Collections of any sort, whatsoever!!!" << endmsg;
962
963 // m_tree->Fill();
964 // return StatusCode::FAILURE;
965
966 return;
967 }
968
969 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
970 m_provider->msg(MSG::VERBOSE) << "Found McEventCollection: " << key << "\tNevents " << mcevent->size() << endmsg;
971 }
972
973 McEventCollection::const_iterator evitr = mcevent->begin();
974 McEventCollection::const_iterator evend = mcevent->end();
975
976 unsigned ie = 0;
977 unsigned ip = 0;
978
979 unsigned ie_ip = 0;
980
981 while ( evitr!=evend ) {
982
983 int _ip = 0;
984
985 int pid = HepMC::signal_process_id((*evitr));
986
987 //The logic should be clarified here
988 if ( pid!=0 ) {
989
990 for (auto pitr: *(*evitr)) {
991
992 selectorTruth.selectTrack( pitr );
993
994 ++_ip;
995
996 }
997
998 }
999 ++ie;
1000 ++evitr;
1001
1002 if ( _ip>0 ) {
1004 // m_provider->msg(MSG::VERBOSE) << "Found " << ie << "\tpid " << pid << "\t with " << ip << " TruthParticles (GenParticles)" << endmsg;
1005 ++ie_ip;
1006 ip += _ip;
1007 }
1008 }
1009
1010 if(m_provider->msg().level() <= MSG::VERBOSE){
1011 m_provider->msg(MSG::VERBOSE) << "Found " << ip << " TruthParticles (GenParticles) in " << ie_ip << " GenEvents out of " << ie << endmsg;
1012 m_provider->msg(MSG::VERBOSE) << "selected " << selectorTruth.size() << " TruthParticles (GenParticles)" << endmsg;
1013 }
1014
1015 if(selectorTruth.size() > 0) foundTruth = true;
1016
1017 if ( !(ip>0) ) {
1018 if (m_provider->msg().level() <= MSG::VERBOSE) m_provider->msg(MSG::WARNING) << "NO TRUTH PARTICLES - returning" << endmsg;
1019 return;
1020 }
1021
1022 }
1023
1025
1026 // m_provider->msg(MSG::VERBOSE) << " Offline tracks " << endmsg;
1027
1028 if ( m_doOffline ) {
1029
1030# ifdef XAODTRACKING_TRACKPARTICLE_H
1031 if ( m_provider->evtStore()->template contains<xAOD::TrackParticleContainer>("InDetTrackParticles") ) {
1032 this->template selectTracks<xAOD::TrackParticleContainer>( m_selectorRef, "InDetTrackParticles" );
1033 refbeamspot = this->template getBeamspot<xAOD::TrackParticleContainer>( "InDetTrackParticles" );
1034 }
1035 else if (m_provider->evtStore()->template contains<Rec::TrackParticleContainer>("TrackParticleCandidate") ) {
1036 this->template selectTracks<Rec::TrackParticleContainer>( m_selectorRef, "TrackParticleCandidate" );
1037 }
1038# else
1039 if (m_provider->evtStore()->template contains<Rec::TrackParticleContainer>("TrackParticleCandidate") ) {
1040 this->template selectTracks<Rec::TrackParticleContainer>( m_selectorRef, "TrackParticleCandidate" );
1041 }
1042# endif
1043 else if ( m_provider->msg().level() <= MSG::WARNING ) {
1044 m_provider->msg(MSG::WARNING) << " Offline tracks not found " << endmsg;
1045 }
1046
1047 ref_tracks = m_selectorRef->tracks();
1048
1049 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
1050 m_provider->msg(MSG::VERBOSE) << "ref tracks.size() " << m_selectorRef->tracks().size() << endmsg;
1051 for ( int ii=m_selectorRef->tracks().size() ; ii-- ; ) {
1052 m_provider->msg(MSG::VERBOSE) << " ref track " << ii << " " << *m_selectorRef->tracks()[ii] << endmsg;
1053 }
1054 }
1055
1056 }
1057 else {
1059 if ( m_mcTruth && foundTruth ){
1060 ref_tracks=selectorTruth.tracks();
1061 }
1062 }
1063
1064
1065
1066 test_tracks.clear();
1067
1068 for ( unsigned itrk=0 ; itrk<chain.rois().at(iroi).tracks().size() ; itrk++ ) {
1069 test_tracks.push_back(&(chain.rois().at(iroi).tracks().at(itrk)));
1070 }
1071
1072
1073 // std::cout << "sutt track multiplicities: offline " << offline_tracks.size() << "\ttest " << test_tracks.size() << std::endl;
1074
1075 _analysis->setvertices( vertices.size() );
1076
1077 if ( refbeamspot.size()>0 ) _analysis->setBeamRef( refbeamspot );
1078 if ( testbeamspot.size()>0 ) _analysis->setBeamTest( testbeamspot );
1079
1082
1083 if ( m_runPurity ) {
1084
1085 if ( this->getUseHighestPT() ) HighestPTOnly( test_tracks );
1086
1087 if ( m_pTthreshold>0 ) FilterPT( test_tracks, m_pTthreshold );
1088
1090 m_NRois++;
1091 m_NRefTracks += test_tracks.size();
1092 m_NTestTracks += ref_tracks.size();
1093
1095 m_associator->match( test_tracks, ref_tracks );
1096
1097 _analysis->execute( test_tracks, ref_tracks, m_associator );
1098
1099 }
1100 else {
1101
1103
1104 if ( this->getUseHighestPT() ) HighestPTOnly( ref_tracks );
1105
1107
1108 if ( m_pTthreshold>0 ) FilterPT( ref_tracks, m_pTthreshold );
1109
1111 m_NRois++;
1112 m_NRefTracks += ref_tracks.size();
1113 m_NTestTracks += test_tracks.size();
1114
1116 m_associator->match( ref_tracks, test_tracks );
1117
1118 // std::cout << "SUTT: execute : N tracks " << ref_tracks.size() << " " << test_tracks.size() << std::endl;
1119
1120 _analysis->setroi( &chain.rois().at(iroi).roi() );
1121 _analysis->execute( ref_tracks, test_tracks, m_associator );
1122
1123 // std::cout << "chain " << m_chainNames[ichain] << " " << "\tvtx name " << vtx_name << std::endl;
1124
1125 if ( vtx_name!="" ) {
1127 std::vector<TIDA::Vertex> vr = chain.rois().at(iroi).vertices();
1128 std::vector<TIDA::Vertex*> vtx_rec;
1129 for ( unsigned iv=0 ; iv<vr.size() ; iv++ ) vtx_rec.push_back( &vr[iv] );
1130
1131 std::vector<TIDA::Vertex*> vtx;
1132 if ( this->getVtxIndex()<0 ) {
1133 for ( unsigned iv=0 ; iv<vertices.size() ; iv++ ) vtx.push_back( &vertices[iv] );
1134 }
1135 else {
1136 if ( vertices.size()>unsigned(this->getVtxIndex()) ) vtx.push_back( &vertices[this->getVtxIndex()] );
1137 }
1138
1139 _analysis->execute_vtx( vtx, vtx_rec, m_event );
1140 }
1141
1142 }
1143
1144 if ( _analysis->debug() ) {
1145 m_provider->msg(MSG::INFO) << "Missing track for " << m_chainNames[ichain]
1146 << "\trun " << run_number
1147 << "\tevent " << event_number
1148 << "\tlb " << lumi_block << endmsg;
1149 }
1150
1151 }
1152
1153 }
1154
1155 if ( m_provider->msg().level() <= MSG::VERBOSE ) {
1156 m_provider->msg(MSG::VERBOSE) << "\n\nEvent " << *m_event << endmsg;
1157 }
1158 }
1159
1160
1161
1162 virtual void book() {
1163
1164 if(m_provider->msg().level() <= MSG::VERBOSE)
1165 m_provider->msg(MSG::VERBOSE) << "AnalysisConfig_Tier0::book() " << name() << endmsg;
1166
1167 m_provider->msg(MSG::ERROR) << "AnalysisConfig_Tier0::book() should no longer ever be called: " << name() << endmsg;
1168
1169 }
1170
1171
1172
1173 virtual void finalize() {
1174
1175 if(m_provider->msg().level() <= MSG::VERBOSE)
1176 m_provider->msg(MSG::VERBOSE) << "AnalysisConfig_Tier0::finalise() " << m_provider->name() << endmsg;
1177
1178 m_analysis->finalise();
1179
1180 m_provider->msg(MSG::INFO) << m_provider->name() << " " << m_chainNames[0] << " \tNRois processed: " << m_NRois << "\tRef tracks: " << m_NRefTracks << "\tTestTracks: " << m_NTestTracks << endmsg;
1181
1182 if(m_provider->msg().level() <= MSG::VERBOSE)
1183 m_provider->msg(MSG::VERBOSE) << m_provider->name() << " finalised" << endmsg;
1184
1185 }
1186
1187
1188protected:
1189
1191
1193
1194 std::vector<ChainString> m_chainNames;
1195 std::vector<Analysis_Tier0*> m_analyses;
1196 std::string m_testType;
1197
1206
1207 std::string m_outputFileName;
1208
1213
1215
1217
1219
1221
1223
1224};
1225
1226
1227
1228#endif // TrigInDetAnalysisExample_T_AnalysisConfig_Tier0_H
1229
#define endmsg
Erkcan's track matchers.
This class provides a unique identification for each event, in terms of run/event number and/or a tim...
generic track filter to filter on eta and PT only
bool passed(DecisionID id, const DecisionIDContainer &)
checks if required decision ID is in the set of IDs in the container
size_t size() const
Number of registered mappings.
Basic event class to contain a vector of chains for trigger analysis.
useful tool for the TrigInDetAnalysis class code
void HighestPTOnly(std::vector< T * > &tracks)
Definition TIDATools.h:20
void FilterPT(std::vector< T * > &tracks, double pt)
Definition TIDATools.h:41
TIDA::Associator< TIDA::Track > TrackAssociator
virtual void execute_vtx(const std::vector< TIDA::Vertex * > &vtx0, const std::vector< TIDA::Vertex * > &vtx1, const TIDA::Event *tevt=0)
const std::string & extra() const
Definition ChainString.h:38
const std::string & roi() const
Definition ChainString.h:35
std::string subs(std::string s) const
const std::string & raw() const
Definition ChainString.h:42
const std::string & head() const
Definition ChainString.h:33
const std::string & tail() const
Definition ChainString.h:34
const std::string & vtx() const
Definition ChainString.h:36
const std::string & element() const
Definition ChainString.h:37
DataModel_detail::const_iterator< DataVector > const_iterator
Definition DataVector.h:838
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
size_type size() const noexcept
Returns the number of elements in the collection.
EventID * event_ID()
the unique identification of the event.
float averageInteractionsPerCrossing() const
average interactions per crossing for all BCIDs - for out-of-time pile-up
void setRoi(TIDARoiDescriptor *r)
Definition Filters.h:236
void containtracks(bool b=true)
set / unset the flag to determine whether tracks should be fully contained in the RoI or not
Definition Filters.h:242
default simple filter which accepts all tracks
Definition Filters.h:26
This defines the McEventCollection, which is really just an ObjectVector of McEvent objectsFile: Gene...
Describes the Region of Interest geometry It has basically 8 parameters.
std::vector< ChainString > m_chainNames
T_AnalysisConfig_Tier0(const std::string &analysisInstanceName, const std::string &testChainName, const std::string &testType, const std::string &testKey, const std::string &referenceChainName, const std::string &referenceType, const std::string &referenceKey, TrackFilter *testFilter, TrackFilter *referenceFilter, TrackAssociator *associator, TrackAnalysis *analysis)
std::vector< Analysis_Tier0 * > m_analyses
const TrackAnalysis * analysis() const
TrackAnalysis * m_analysis
const std::string & name() const
TrigTrackSelector * m_selectorTest
std::string m_analysisInstanceName
bool selectTracks(TrigTrackSelector *selector, const ElementLink< TrigRoiDescriptorCollection > &roi_link, const std::string &key="")
std::vector< std::vector< TrackFilter * > > m_filters
TrigTrackSelector * m_selectorRef
StatusCode retrieve(Collection const *&collection, const std::string &key="")
std::vector< double > getBeamspot(const std::string &key)
ToolHandle< Trig::TrigDecisionTool > * m_tdt
T_AnalysisConfig(const std::string &analysisInstanceName, const std::string &testChainName, const std::string &testType, const std::string &testKey, const std::string &referenceChainName, const std::string &referenceType, const std::string &referenceKey, const std::string &selectionChainName, const std::string &selectionType, const std::string &selectionKey, TrackFilter *testFilter, TrackFilter *referenceFilter, TrackFilter *selectionFilter, TrackAssociator *associator, TrackAnalysis *analysis)
TrackAssociator * m_associator
nope - should be used for standalone also, perhaps need to protect the class def bits ifndef XAOD_ANA...
virtual void clear() override
bool selectTrack(const Rec::TrackParticle *track)
neater code to make use of vector function also for a single ancestor pdgid, instead of the full code...
std::vector< Combination >::const_iterator combination_const_iterator
bool contains(const std::string &s, const std::string &regx)
does a string contain the substring
Definition hcg.cxx:116
int signal_process_id(const GenEvent &evt)
Definition GenEvent.h:572
@ x
Definition ParamDefs.h:55
@ z
global position (cartesian)
Definition ParamDefs.h:57
@ y
Definition ParamDefs.h:56
EventInfo_v1 EventInfo
Definition of the latest event info version.
VertexContainer_v1 VertexContainer
Definition of the current "Vertex container version".