T_AnalysisConfigMT_Tier0.h

Go to the documentation of this file.

/* emacs: this is -*- c++ -*- */
 
#ifndef TrigInDetAnalysisExample_T_AnalysisConfigMT_Tier0_H
#define TrigInDetAnalysisExample_T_AnalysisConfigMT_Tier0_H
 
#include "TrigInDetAnalysis/TIDAEvent.h"
#include "TrigInDetAnalysis/TIDAVertex.h"
#include "TrigInDetAnalysis/TrackSelector.h"     
#include "TrigInDetAnalysisUtils/T_AnalysisConfig.h"
#include "TrigInDetAnalysisUtils/TagNProbe.h"
 
#include "TrigInDetAnalysisExample/Analysis_Tier0.h"
#include "TrigInDetAnalysisExample/AnalysisR3_Tier0.h"
#include "TrigInDetAnalysisExample/VtxAnalysis.h"
#include "TrigInDetAnalysisExample/ChainString.h"
#include "TrigInDetAnalysisExample/TIDATools.h"
 
#include "TTree.h"
#include "TFile.h"
 
#include "GaudiKernel/ToolHandle.h"
#include "AthenaMonitoringKernel/GenericMonitoringTool.h"
 
 
// McParticleEvent includes
#include "McParticleEvent/TruthParticleContainer.h"
 
#include "GeneratorObjects/McEventCollection.h"
#include "AtlasHepMC/GenEvent.h"
#include "AtlasHepMC/GenVertex.h"
#include "AtlasHepMC/GenParticle.h"
 
#include "xAODEventInfo/EventInfo.h"
 
 
 
#include "TrigInDetAnalysis/TIDDirectory.h"
 
#include "TrigInDetAnalysis/Filter_AcceptAll.h"
 
#include "TrigInDetAnalysisUtils/TIDARoiDescriptorBuilder.h"
#include "TrigInDetAnalysisUtils/Filter_etaPT.h"
#include "TrigInDetAnalysisUtils/Filter_RoiSelector.h"
#include "TrigInDetAnalysisUtils/Associator_BestMatch.h"
#include "TrigInDetAnalysisUtils/Filters.h"
 
                            
 
#include "VxVertex/VxContainer.h"
 
#include "TrigSteeringEvent/HLTResult.h"
#include "TrigDecisionTool/ExpertMethods.h"
 
// #include "TrigSteeringEvent/TrigRoiDescriptorCollection.h"
 
#include "xAODTracking/TrackParticle.h"
#include "xAODTracking/TrackParticleContainer.h"
 
#include "TrigCompositeUtils/TrigCompositeUtils.h"
 
 
 
 
 
template<typename T, typename A=Analysis_Tier0>
class T_AnalysisConfigMT_Tier0 : public T_AnalysisConfig<T> {
 
public:
 
  // Full constructor: test/reference/selection
  // - analysisInstanceName: the name of the analysis chain being created
  // - xxxChainName: the name of the chain to be used as test/reference/selection; must be "StoreGate" in case of direct access to SG containers
  // - xxxType: the type of tracks to be retrieved from the test/reference/selection chain or container
  // - xxxKey:  the key for tracks to be retrieved from the test/reference/selection chain or container
  // - all standard operations are performed in loops over 0=test 1=reference 2=selection
  T_AnalysisConfigMT_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,
             TagNProbe*      TnP_tool = 0) :
    T_AnalysisConfig<T>( analysisInstanceName,
             testChainName,      testType,      testKey,
             referenceChainName, referenceType, referenceKey,
             testFilter, referenceFilter,
             associator,
             analysis),
    _analysis(0),
    m_useBeamCondSvc(false),
    m_doOffline(true),
    m_doMuons(false),
    m_doElectrons(false),
    m_doTaus(false),
    m_doBjets(false),
    m_hasTruthMap(false),
    m_doTauThreeProng(false),
    m_tauEtCutOffline(false),
    m_NRois(0),
    m_NRefTracks(0),
    m_NTestTracks(0),
    m_runPurity(false),
    m_shifter(false),
    m_pTthreshold(0),
    m_first(true),
    m_containTracks(false),
    m_TnP_tool(TnP_tool),
    m_invmass(0),
    m_invmass_obj(0),
    m_monTool(0)
  {
    m_event = new TIDA::Event();
    m_chainNames.push_back(testChainName);
 
#if 0
    ChainString& chain = m_chainNames.back(); 
 
    std::cout << "\nT_AnalysisConfigMT_Tier0::name:                " << name() << "\t" << this << std::endl;
    std::cout <<  "T_AnalysisConfigMT_Tier0::chain specification: " << testChainName << " -> " << chain << "\t" << chain.raw() << std::endl;
    std::cout << "\tchain: " << chain.head()    << std::endl;
    std::cout << "\tkey:   " << chain.tail()    << std::endl;
    std::cout << "\troi:   " << chain.roi()     << std::endl;
    std::cout << "\tvtx:   " << chain.vtx()     << std::endl;
    std::cout << "\tte:    " << chain.element() << std::endl;
 
    std::cout << "\tpost:  " << chain.post()          << std::endl; 
    std::cout << "\tpt:    " << chain.postvalue("pt") << std::endl;
 
    std::cout << "\tcontainTracks: " << m_containTracks << std::endl; 
#endif
    
    m_testType = testType;
  }
 
  virtual ~T_AnalysisConfigMT_Tier0() {
    delete m_event;
    if ( m_TnP_tool != 0 ) delete m_TnP_tool ;
  }
 
  void setRunPurity( bool b ) { m_runPurity=b; }
 
  void setShifter( bool b )    { m_shifter=b; }
 
  void useBeamCondSvc( bool b ) { m_useBeamCondSvc = b; }
 
  void containTracks( bool b ) { m_containTracks = b; }
 
public:
 
  A* _analysis;
 
protected:
 
  using T_AnalysisConfig<T>::m_provider;
  using T_AnalysisConfig<T>::m_tdt;
  using T_AnalysisConfig<T>::m_mcTruth;
 
  using T_AnalysisConfig<T>::name;
  using T_AnalysisConfig<T>::m_analysis;
 
  using T_AnalysisConfig<T>::m_selectorTest;
  using T_AnalysisConfig<T>::m_selectorRef;
  using T_AnalysisConfig<T>::m_associator;
  using T_AnalysisConfig<T>::m_filters;
  
  
  virtual void loop() {
    
    bool TnP_flag = (m_TnP_tool != 0) ; // flag for tag and probe analysis
    
    if( m_provider->msg().level() <= MSG::VERBOSE) {
      m_provider->msg(MSG::VERBOSE) <<  "AnalysisConfigMT_Tier0::loop() for " << T_AnalysisConfig<T>::m_analysisInstanceName <<  endmsg;
    }
 
    // get (offline) beam position
    double xbeam = 0;
    double ybeam = 0;
 
    if ( m_first ) {      
 
      m_first = false;
      
      if ( m_provider->msg().level() <= MSG::VERBOSE ) {
    m_provider->msg(MSG::VERBOSE) << " using beam position\tx=" << xbeam << "\ty=" << ybeam << endmsg;
    
    std::vector<std::string> configuredChains  = (*(m_tdt))->getListOfTriggers("L2_.*, EF_.*, HLT_.*");
    
        for ( unsigned i=0 ; i<configuredChains.size() ; i++ ) {
      //    std::cout << "Configured chain " << configuredChains[i]  << std::endl;
      m_provider->msg(MSG::VERBOSE)  << "Chain " << configuredChains[i]  << endmsg;
    }
      }
    
      //      std::cout << "\tloop() analyse chains " << m_chainNames.size() << std::endl;
    
      std::vector<ChainString>::iterator chainitr = m_chainNames.begin();
      
      std::vector<ChainString> chains;
      
      while ( chainitr!=m_chainNames.end() ) {
    
        ChainString& chainName = (*chainitr);
 
    if (m_provider->msg().level() <= MSG::VERBOSE) {
      m_provider->msg(MSG::VERBOSE) << "process chain " << chainName << endmsg;
    }   
 
    if ( chainName.head() == "" ) { 
      
      std::string selectChain = chainName.raw(); 
 
      chains.push_back( ChainString(selectChain) );
        
    }
    else { 

      std::vector<std::string> selectChains  = (*(m_tdt))->getListOfTriggers( chainName.head() );
          
      for ( unsigned iselected=0 ; iselected<selectChains.size() ; iselected++ ) {
 
        selectChains[iselected] = chainName.subs( selectChains[iselected] ); 
        
#if 0
        std::cout << "sorting:: chain specification: " << chainName << "\traw:" << chainName.raw() << std::endl;
        std::cout << "\tchain: " << chainName.head()    << std::endl;
        std::cout << "\tkey:   " << chainName.tail()    << std::endl;
        std::cout << "\tind:   " << chainName.extra()   << std::endl;
        std::cout << "\troi:   " << chainName.roi()     << std::endl;
        std::cout << "\tvtx:   " << chainName.vtx()     << std::endl;
        std::cout << "\tte:    " << chainName.element() << std::endl;
#endif
        
        chains.push_back( ChainString(selectChains[iselected]) );
        
        if(m_provider->msg().level() <= MSG::VERBOSE) {
          m_provider->msg(MSG::VERBOSE) << "Matching chain " << selectChains[iselected] << " (" << chainName.head() << ")" << endmsg;
        }
 
        //      std::cout<< "\tMatching chain " << sc << " (" << chainName.head() << ")" << std::endl;
        
      }
    }
 
        ++chainitr;
      }
 
      // m_chainNames.insert( m_chainNames.end(), chains.begin(), chains.end() );
      m_chainNames = chains;
      
      if(m_provider->msg().level() <= MSG::VERBOSE) {
    for ( size_t ic=m_chainNames.size() ; ic-- ;  ) m_provider->msg(MSG::VERBOSE) << "Analyse chain " << m_chainNames[ic] << endmsg;
      }      
 
    } 
    
    //    std::cout << "\tloop() event analysis ..." << std::endl;

    Filter_True filter;
    
    Filter_etaPT     filter_etaPT(5,200);
    Filter_Combined filter_truth( &filter_etaPT,   &filter_etaPT);
    
    // Filter_Combined filterRef (&filter_offline, &filter_vertex);
 
    int iRefFilter  = 1;
    int iTestFilter = 0;
 
    if ( m_runPurity ) { 
      iRefFilter  = 0;
      iTestFilter = 1;
    }
    
    Filter_Combined filterRef(   m_filters[iRefFilter][0],  &filter );
    Filter_Combined filterTest(  m_filters[iTestFilter][0], &filter );
 
 
    TrigTrackSelector selectorTruth( &filter_truth );
    TrigTrackSelector selectorRef( &filterRef );
    m_selectorRef = &selectorRef;
    TrigTrackSelector selectorTest( &filterTest );
    m_selectorTest = &selectorTest;
 
    if ( xbeam!=0 || ybeam!=0 ) { 
      m_selectorRef->setBeamline(  xbeam, ybeam );
    }  

 
    // clear the ntuple TIDA::Event class
    m_event->clear();

 
    const xAOD::EventInfo* pEventInfo;
 
    unsigned           run_number        = 0;
    unsigned long long event_number      = 0;
    unsigned           lumi_block        = 0;
    unsigned           bunch_crossing_id = 0;
    unsigned           time_stamp        = 0;
    double             mu_val            = 0;
 
    //    std::cout << "\tloop() get EventInfo" << std::endl;
 
    if ( this->template retrieve(pEventInfo, "EventInfo").isFailure() ) {
      m_provider->msg(MSG::WARNING) << "Failed to get EventInfo " << endmsg;
    } else {
 
      run_number        = pEventInfo->runNumber();
      event_number      = pEventInfo->eventNumber();
      lumi_block        = pEventInfo->lumiBlock();
      time_stamp        = pEventInfo->timeStamp();
      bunch_crossing_id = pEventInfo->bcid();
      mu_val            = pEventInfo->averageInteractionsPerCrossing();
    }
  
    if(m_provider->msg().level() <= MSG::VERBOSE){
      m_provider->msg(MSG::VERBOSE) << "run "     << run_number
                                    << "\tevent " << event_number
                                    << "\tlb "    << lumi_block << endmsg;
    }
    
    //      m_provider->msg(MSG::INFO) << "run "     << run_number
    //               << "\tevent " << event_number
    //               << "\tlb "    << lumi_block << endmsg;
 
    //    std::cout << "\trun "     << run_number  << "\tevent " << event_number  << "\tlb "    << lumi_block << std::endl;
 
 
    // clear the ntuple TIDA::Event class
    m_event->clear();
 
    m_event->run_number(run_number);
    m_event->event_number(event_number);
    m_event->lumi_block(lumi_block);
    m_event->time_stamp(time_stamp);
    m_event->bunch_crossing_id(bunch_crossing_id);
    m_event->mu(mu_val);

 
    bool analyse = false;
  
    // Check HLTResult
 
    //    std::cout << "\tloop() loop over trigger chains to determine whether to process this event ..." << std::endl;  
 
    for ( unsigned ichain=0 ; ichain<m_chainNames.size() ; ichain++ ) {
 
      const std::string& chainname = m_chainNames[ichain].head();
 
      if ( chainname == "" ) analyse = true;
      else { 
 
    //Only for trigger chains
    if ( chainname.find("L2")  == std::string::npos &&
         chainname.find("EF")  == std::string::npos &&
         chainname.find("HLT") == std::string::npos ) continue;
    
    if ( m_provider->msg().level() <= MSG::DEBUG ) {
      m_provider->msg(MSG::DEBUG) << "Chain "  << chainname
                      << "\tpass " << (*m_tdt)->isPassed(chainname)
                      << "\tpres " << (*m_tdt)->getPrescale(chainname) << endmsg;
    }
    
    //  std::cout << "\tChain "  << chainname << "\tpass " << (*m_tdt)->isPassed(chainname)
    //        << "\tpres " << (*m_tdt)->getPrescale(chainname) << std::endl;
    
    if ( (*(m_tdt))->isPassed(chainname) || (*(m_tdt))->getPrescale(chainname) ) analyse = true;
    
      }
    }
    
 
    //    if ( (*m_tdt)->ExperimentalAndExpertMethods().isHLTTruncated() ) {
    //    m_provider->msg(MSG::WARNING) << "HLTResult truncated, skipping event" << endmsg;
    //   return;
    //  }
    
    if ( !this->m_keepAllEvents && !analyse ) {
      //     m_provider->msg(MSG::VERBOSE) << "No chains passed unprescaled - not processing this event" << endmsg;
      if(m_provider->msg().level() <= MSG::VERBOSE) {
        m_provider->msg(MSG::VERBOSE) << "No chains passed unprescaled - not processing this event" << endmsg;
      }
      //    std::cout << "\tNo chains passed unprescaled - not processing this event" << std::endl;
      return;
    }
 
    

 
    selectorTruth.clear();
 
    
 
    if(m_provider->msg().level() <= MSG::VERBOSE)
      m_provider->msg(MSG::VERBOSE) << "MC Truth flag " << m_mcTruth << endmsg;
 
    if ( m_mcTruth && m_provider->msg().level() <= MSG::VERBOSE ) m_provider->msg(MSG::VERBOSE) << "getting Truth" << endmsg;
       

 
    std::vector<TIDA::Vertex> vertices;
    std::vector<TIDA::Vertex> vertices_rec;
 
    std::vector<double> refbeamspot;
    std::vector<double> testbeamspot;

 
    m_provider->msg(MSG::VERBOSE) << "fetching AOD Primary vertex container" << endmsg;
 
    if ( !this->select( vertices, "PrimaryVertices" ) ) { 
      m_provider->msg(MSG::VERBOSE) << "could not retrieve vertex collection " "PrimaryVertices" << std::endl;
    }

 
    if ( m_mcTruth ) {
      m_event->addChain( "Truth" );
      m_event->back().addRoi(TIDARoiDescriptor());
      m_event->back().back().addTracks(selectorTruth.tracks());
    }

 
    if ( m_doOffline ) {
      for ( unsigned i=0 ; i<vertices.size() ; i++ )  {
        if(m_provider->msg().level() <= MSG::VERBOSE)
          m_provider->msg(MSG::VERBOSE) << "vertex " << i << " " << vertices[i] << endmsg;
        m_event->addVertex(vertices[i]);
      }
    }

 
    // int Noff = 0;
    std::vector<TIDA::Track*> offline_tracks;
    std::vector<TIDA::Track*> electron_tracks;
    std::vector<TIDA::Track*> muon_tracks;
 
    std::vector<TIDA::Track*> ref_tracks;
    std::vector<TIDA::Track*> test_tracks;
 
    offline_tracks.clear();
    electron_tracks.clear();
    muon_tracks.clear();
 
    ref_tracks.clear();
    test_tracks.clear();
 
    // offline track retrieval now done once for each chain rather than each roi
    if ( m_provider->msg().level() <= MSG::VERBOSE )
      m_provider->msg(MSG::VERBOSE) << "MC Truth flag " << m_mcTruth << endmsg;
 
    bool foundTruth = false;
 
    if ( !m_doOffline && m_mcTruth ) {
 
      filter_truth.setRoi( 0 ); // don't filter on RoI yet (or until needed)  
 
      selectorTruth.clear();
 
      if ( m_provider->msg().level() <= MSG::VERBOSE )
    m_provider->msg(MSG::VERBOSE) << "getting Truth" << endmsg;
 
      if ( m_provider->evtStore()->template contains<TruthParticleContainer>("INav4MomTruthEvent") ) {
    //ESD
    this->template selectTracks<TruthParticleContainer>( &selectorTruth, "INav4MomTruthEvent" );
    foundTruth = true;
      }
      else if ( m_provider->evtStore()->template contains<TruthParticleContainer>("SpclMC") ) {
    this->template selectTracks<TruthParticleContainer>( &selectorTruth, "SpclMC");
    foundTruth = true;
      }
      else if ( m_provider->evtStore()->template contains<xAOD::TruthParticleContainer>("TruthParticles") ) {
    this->template selectTracks<xAOD::TruthParticleContainer>( &selectorTruth, "TruthParticles");
    foundTruth = true;
      }
      else if ( m_provider->evtStore()->template contains<TruthParticleContainer>("") ) {
    this->template selectTracks<TruthParticleContainer>( &selectorTruth, "");
    foundTruth = true;
      }
      else
    if ( m_provider->msg().level() <= MSG::VERBOSE ) {
      m_provider->msg(MSG::VERBOSE) << "Truth not found - none whatsoever!" << endmsg;
    }
    }
      
    if ( !m_doOffline && m_mcTruth && !foundTruth ) {
          
      if ( m_provider->msg().level() <= MSG::VERBOSE ) { 
    m_provider->msg(MSG::VERBOSE) << "getting Truth" << endmsg;
      }


      const McEventCollection* mcevent = 0;

 
      std::string keys[4] = { "GEN_AOD", "TruthEvent", "", "G4Truth" };
 
      std::string key = "";
 
      bool foundcollection = false;
 
      for ( int ik=0 ; ik<4 ; ik++ ) {
         
    if ( m_provider->msg().level() <= MSG::VERBOSE ) {
      m_provider->msg(MSG::VERBOSE) << "Try McEventCollection: " << keys[ik] << endmsg;
    }
 
    if ( !m_provider->evtStore()->template contains<McEventCollection>(keys[ik]) ) {
      if( m_provider->msg().level() <= MSG::VERBOSE )
        m_provider->msg(MSG::VERBOSE) << "No McEventCollection: " << keys[ik] << endmsg;
      continue;
    }
 
    if ( m_provider->msg().level() <= MSG::VERBOSE )
      m_provider->msg(MSG::VERBOSE) << "evtStore()->retrieve( mcevent, " << keys[ik] << " )" << endmsg;
 
    if ( this->template retrieve( mcevent, keys[ik] ).isFailure() ) {
      if ( m_provider->msg().level() <= MSG::VERBOSE )
        m_provider->msg(MSG::VERBOSE) << "Failed to get McEventCollection: " << keys[ik] << endmsg;
    }
    else {
      key = keys[ik];
      if(m_provider->msg().level() <= MSG::VERBOSE)
        m_provider->msg(MSG::VERBOSE) << "Found McEventCollection: " << key << endmsg;
      foundcollection = true;
      break;
    }
      }

      if ( !foundcollection ) {
    if(m_provider->msg().level() <= MSG::VERBOSE)
      m_provider->msg(MSG::WARNING) << "No MC Truth Collections of any sort, whatsoever!!!" << endmsg;
 
    //    m_tree->Fill();
    //    return StatusCode::FAILURE;
 
    return;
      }
 
      if ( m_provider->msg().level() <= MSG::VERBOSE ) {
    m_provider->msg(MSG::VERBOSE) << "Found McEventCollection: " << key << "\tNevents " << mcevent->size() << endmsg;
      }
 
      McEventCollection::const_iterator evitr = mcevent->begin();
      McEventCollection::const_iterator evend = mcevent->end();
 
      unsigned ie = 0; 
      unsigned ip = 0; 
 
      unsigned ie_ip = 0; 
 
      while ( evitr!=evend ) {
 
    int _ip = 0; 
 
    int pid = HepMC::signal_process_id((*evitr));
 
    //The logic should be clarified here
    if ( pid!=0 ) { 
 
      for (auto pitr: *(*evitr) ) { 
 
        selectorTruth.selectTrack( pitr );
 
        ++_ip;
                
      }
 
    }
    ++ie;
    ++evitr;
 
    if ( _ip>0 ) {
      //      m_provider->msg(MSG::VERBOSE) << "Found " << ie << "\tpid " << pid << "\t with " << ip << " TruthParticles (GenParticles)" << endmsg;
      ++ie_ip;
      ip += _ip;
    }
      }
 
      if(m_provider->msg().level() <= MSG::VERBOSE){
    m_provider->msg(MSG::VERBOSE) << "Found " << ip << " TruthParticles (GenParticles) in " << ie_ip << " GenEvents out of " << ie << endmsg;
    m_provider->msg(MSG::VERBOSE) << "selected " << selectorTruth.size() << " TruthParticles (GenParticles)" << endmsg;
      }
 
      if(selectorTruth.size() > 0) foundTruth = true;
 
      if ( !(ip>0) ) {
    if (m_provider->msg().level() <= MSG::VERBOSE) m_provider->msg(MSG::WARNING) << "NO TRUTH PARTICLES - returning" << endmsg;
    return; 
      }
      
    }
 
    // m_provider->msg(MSG::VERBOSE) << " Offline tracks " << endmsg;
    
    if ( m_doOffline ) {
      
      if ( m_provider->evtStore()->template contains<xAOD::TrackParticleContainer>("InDetTrackParticles") ) {
    this->template selectTracks<xAOD::TrackParticleContainer>( m_selectorRef, "InDetTrackParticles" );
    refbeamspot = this->template getBeamspot<xAOD::TrackParticleContainer>( "InDetTrackParticles" );
      }
      else if (m_provider->evtStore()->template contains<Rec::TrackParticleContainer>("TrackParticleCandidate") ) {
    this->template selectTracks<Rec::TrackParticleContainer>( m_selectorRef, "TrackParticleCandidate" );
      }
      else if ( m_provider->msg().level() <= MSG::WARNING ) {
    m_provider->msg(MSG::WARNING) << " Offline tracks not found " << endmsg;
      }
 
    }
 
 
    //    std::cout << "\tloop() loop over chains proper ..." << std::endl;

    for ( unsigned ichain=0 ; ichain<m_chainNames.size() ; ichain++ ) {

 
      // std::string& chainname = chains[ichain];
      const std::string& chainname = m_chainNames[ichain].head();
      const std::string&       key = m_chainNames[ichain].tail();
      const std::string&  vtx_name = m_chainNames[ichain].vtx();
      //Not used left just in case
      //const std::string&  roi_name = m_chainNames[ichain].roi();
      //const std::string&  te_name = m_chainNames[ichain].element();
      m_pTthreshold = 0;  
 
      if ( m_chainNames[ichain].postcount() ) { 
        std::string ptvalue = m_chainNames[ichain].postvalue("pt");
    if ( ptvalue!="" ) m_pTthreshold = std::stod(ptvalue);
      }
 
 
      unsigned _decisiontype = TrigDefs::Physics;
      unsigned  decisiontype;
     
      if ( this->requireDecision() ) _decisiontype =  TrigDefs::requireDecision;
      
      
      if ( m_chainNames[ichain].passed() ) decisiontype = _decisiontype;
      else                                 decisiontype = TrigDefs::alsoDeactivateTEs;

      //      if ( decisiontype==TrigDefs::requireDecision ) std::cout << "\tSUTT TrigDefs::requireDecision " << decisiontype << std::endl;
      //      if ( decisiontype==TrigDefs::Physics )         std::cout << "\tSUTT TrigDefs::Physics "         << decisiontype << std::endl;
 
      if ( chainname!="" && m_provider->msg().level() <= MSG::VERBOSE ) {
 
        m_provider->msg(MSG::VERBOSE) << "status for chain " << chainname
                                      << "\tpass "           << (*m_tdt)->isPassed(chainname)
                                      << "\tprescale "       << (*m_tdt)->getPrescale(chainname) << endmsg;
 
        m_provider->msg(MSG::VERBOSE) << "fetching features for chain " << chainname << endmsg;
    
        m_provider->msg(MSG::VERBOSE) << chainname << "\tpassed: " << (*m_tdt)->isPassed( chainname ) << endmsg;
      }

      //      std::cout << "\tstatus for chain " << chainname
      //        << "\tpass "           << (*m_tdt)->isPassed( chainname )
      //        << "\tpassdt "         << (*m_tdt)->isPassed( chainname, decisiontype )
      //        << "\tprescale "       << (*m_tdt)->getPrescale( chainname ) << std::endl;
    
 
      //   m_provider->msg(MSG::INFO) << chainname << "\tpassed: " << (*m_tdt)->isPassed( chainname ) << "\t" << m_chainNames[ichain] << "\trun " << run_number << "\tevent " << event_number << endmsg;
 
 
      if ( chainname!="" && !this->m_keepAllEvents && !(*m_tdt)->isPassed( chainname, decisiontype ) ) continue;
 

      // Trig::FeatureContainer f = (*m_tdt)->features( chainname, TrigDefs::alsoDeactivateTEs);
     
 
      // tag and probe analysis processes multiple chains passed in the tag and probe tool at once so loop over vector of chains
      std::vector<std::string> chainNames ;
            
      if ( !TnP_flag ) {
    chainNames.push_back(m_chainNames[ichain].raw()) ;
      }
      else {
    chainNames.push_back(m_TnP_tool->tag()) ;
    chainNames.push_back(m_TnP_tool->probe()) ;
      }
 
      // loop over new chainNames vector but doing the same stuff
      for ( unsigned i = 0 ; i < chainNames.size() ; ++i ) {
 
    ChainString chainConfig = chainNames[i] ;
    std::string chainName = chainConfig.head();
    
    m_event->addChain( chainNames[i] ); 
    
    TIDA::Chain& chain = m_event->back();
    
    if ( chainName == "" ) { 

 
      m_selectorTest->clear();
    
 
      TIDARoiDescriptor* roiInfo = new TIDARoiDescriptor(true);
        
      chain.addRoi( *roiInfo );
    
      if ( m_provider->evtStore()->template contains<xAOD::TrackParticleContainer>(key) ) {
        this->template selectTracks<xAOD::TrackParticleContainer>( m_selectorTest, key );
        refbeamspot = this->template getBeamspot<xAOD::TrackParticleContainer>( key );
      }
 
      const std::vector<TIDA::Track*>& testtracks = m_selectorTest->tracks();
 
      chain.back().addTracks(testtracks);
    
      if ( vtx_name!="" ) { 
      
 
        m_provider->msg(MSG::VERBOSE) << "\tFetch xAOD::VertexContainer with key " << vtx_name << endmsg;
 
        std::vector<TIDA::Vertex> tidavertices;
 
        if ( this->select( tidavertices, vtx_name ) ) chain.back().addVertices( tidavertices );  
      }
    
 
      delete roiInfo;
 
    }
    else {

      
      //std::string roi_key = m_chainNames[ichain].roi();
 
      std::string roi_key = chainConfig.roi();
    
      unsigned feature_type =TrigDefs::lastFeatureOfType;
 
      if ( roi_key!="" ) feature_type= TrigDefs::allFeaturesOfType;

      
      int leg = -1;
      
      if ( chainConfig.element()!="" ) { 
        leg = std::atoi(chainConfig.element().c_str());
      }
      
      std::vector< TrigCompositeUtils::LinkInfo<TrigRoiDescriptorCollection> > rois = 
        (*m_tdt)->template features<TrigRoiDescriptorCollection>( Trig::FeatureRequestDescriptor( chainName,  
                                                      decisiontype,
                                                      roi_key, 
                                                      feature_type,
                                                      "roi", 
                                                      leg ) );
      
      //    const unsigned int featureCollectionMode = const std::string& navElementLinkKey = "roi") const;
      
      int iroi = 0; 
      
      for ( const TrigCompositeUtils::LinkInfo<TrigRoiDescriptorCollection>& roi_info : rois ) {
      
        iroi++;

 
        if ( roi_key=="SuperRoi" && iroi>1 ) continue; 
    
        //    std::cout << "\troi: get link " << roi_key << " ..." << std::endl;
 
        const ElementLink<TrigRoiDescriptorCollection> roi_link = roi_info.link;

        if ( roi_key!="" && roi_link.dataID()!=roi_key ) continue;
 
        const TrigRoiDescriptor* const* roiptr = roi_link.cptr();
 
        if ( roiptr == 0 ) { 
          continue;
        }  
 
        //    std::cout << "\troi: link deref ..." << *roiptr << std::endl;
 
        if (m_provider->msg().level() <= MSG::VERBOSE) {
          m_provider->msg(MSG::VERBOSE) << " RoI descriptor for seeded chain " << chainname << " " << **roiptr << endmsg;
        }
     
        TIDARoiDescriptor* roiInfo = new TIDARoiDescriptor( TIDARoiDescriptorBuilder(**roiptr) );
     
        //    if ( dbg ) std::cout << "\troi " << iroi << " " << *roiInfo << std::endl;

 
        m_selectorTest->clear();
 
        if ( this->template selectTracks<xAOD::TrackParticleContainer>( m_selectorTest, roi_link,  key ) ) { } 
 
        // beamspot stuff not needed for xAOD::TrackParticles

 
        chain.addRoi( *roiInfo );
    
 
        const std::vector<TIDA::Track*>& testtracks = m_selectorTest->tracks();
 
        chain.back().addTracks(testtracks);
    

      
        if ( vtx_name!="" ) { 
 
          std::vector<TIDA::Vertex> tidavertices;
 
          this->select( tidavertices, roi_link, vtx_name );
 
          chain.back().addVertices( tidavertices );
        
        } 
      
 
#if 0
        if ( dbg ) { 
          std::cout << "\tTIDA analysis for chain: " << chainname << "\t key: " << key << "\t" << **roiptr << std::endl;
          std::cout << "\tcollections: " << chain.back() << std::endl; 
        }
#endif
      
        delete roiInfo;
      
      }
 
 
    } 
 
      } // end of loop chainNames vector loop
    
    
      if ( m_provider->msg().level() <= MSG::VERBOSE ) {
    m_provider->msg(MSG::VERBOSE) << "event: " << *m_event << endmsg;
      }
 
    }
 
    // close previous loop over chains and open new one
    
    for ( unsigned ichain=0 ; ichain<m_event->size() ; ichain++ ) {
      
      TIDA::Chain& chain = (*m_event)[ichain];
      ChainString chainConfig(chain.name());        
      const std::string&  vtx_name = chainConfig.vtx();
 
      // skip tag chains to avoid performing standard analysis on them (done for tnp at the same time as probes)
      if ( TnP_flag && chainConfig.extra().find("_tag")!=std::string::npos ) continue ;
      
      std::vector<TIDA::Roi*> rois ;
      
      if (TnP_flag) {
    rois = m_TnP_tool->GetRois( m_event->chains(), m_selectorRef, &filterRef, m_invmass, m_invmass_obj );
    // needs to be done AFTER retrieving offline tracks as m_selectorRef passed as arguement, hence restructuring
      }
      else {
        rois.reserve( chain.size() );
        for ( size_t ir=0 ; ir<chain.size() ; ir++ ) {
      rois.push_back( &(chain.rois()[ir]) );
    }
      }
 
      // now loop over the rois (again)
 
      for ( unsigned iroi=0 ; iroi<rois.size() ; iroi++ ) {
 
    if ( this->filterOnRoi() ) { 
      filterRef.setRoi( &(rois.at(iroi)->roi() ) ); 
      filterRef.containtracks( m_containTracks ); 
    }
    else filterRef.setRoi( 0 );
    
    test_tracks.clear();
    
    // this block is before the track retrieval in the original, is it working the same here?
    
 
        if ( m_provider->msg().level() <= MSG::VERBOSE )
          m_provider->msg(MSG::VERBOSE) << "MC Truth flag " << m_mcTruth << endmsg;
 
    if ( !m_doOffline && m_mcTruth ) {
      if ( this->filterOnRoi() )  filter_truth.setRoi( &(rois.at(iroi)->roi() ) );
      ref_tracks = m_selectorRef->tracks(&filter_truth);
    }
    else { // ie. if ( m_doOffline )
      ref_tracks = m_selectorRef->tracks(&filterRef) ; 
    }
    
    if ( m_provider->msg().level() <= MSG::VERBOSE ) {
      m_provider->msg(MSG::VERBOSE) << "ref tracks.size() " << m_selectorRef->tracks().size() << endmsg;
      for ( int ii=m_selectorRef->tracks().size() ; ii-- ; ) {
        m_provider->msg(MSG::VERBOSE) << "  ref track " << ii << " " << *m_selectorRef->tracks()[ii] << endmsg;
      }
    }
 
        test_tracks.clear();
 
 
        for ( unsigned itrk=0 ; itrk<rois.at(iroi)->tracks().size() ; itrk++ ) {
          test_tracks.push_back(&(rois.at(iroi)->tracks().at(itrk)));
        }   
    
    //  std::cout << "sutt track multiplicities: offline " << offline_tracks.size() << "\ttest " << test_tracks.size() << std::endl;
 
        _analysis->setvertices( vertices.size() );  
 
    if ( refbeamspot.size()>0 )  _analysis->setBeamRef(   refbeamspot ); 
    if ( testbeamspot.size()>0 ) _analysis->setBeamTest( testbeamspot ); 

 
    if ( m_runPurity ) {
 
      if ( this->getUseHighestPT() ) HighestPTOnly( test_tracks );
 
      if ( m_pTthreshold>0 ) FilterPT( test_tracks, m_pTthreshold );

      m_NRois++;
      m_NRefTracks  += test_tracks.size();
      m_NTestTracks += ref_tracks.size();
 
              
      m_associator->match( test_tracks, ref_tracks );
      
      _analysis->execute( test_tracks, ref_tracks, m_associator );
 
    }
    else { 

      if ( this->getUseHighestPT() ) HighestPTOnly( ref_tracks );

      
      if ( m_pTthreshold>0 )         FilterPT( ref_tracks, m_pTthreshold );

      m_NRois++;
      m_NRefTracks  += ref_tracks.size();
      m_NTestTracks += test_tracks.size();
      
      m_associator->match( ref_tracks, test_tracks );
    
      _analysis->setroi( &rois.at(iroi)->roi() );  
      _analysis->execute( ref_tracks, test_tracks, m_associator );
 
      if ( vtx_name!="" ) { 
        std::vector<TIDA::Vertex> vr = rois.at(iroi)->vertices();
        std::vector<TIDA::Vertex*> vtx_rec;    
        for ( unsigned iv=0 ; iv<vr.size() ; iv++ ) vtx_rec.push_back( &vr[iv] );
 
        std::vector<TIDA::Vertex*> vtx;
        if ( this->getVtxIndex()<0 ) { 
          for ( unsigned iv=0 ; iv<vertices.size() ; iv++ ) vtx.push_back( &vertices[iv] );
        }
        else { 
          if ( vertices.size()>unsigned(this->getVtxIndex()) ) vtx.push_back( &vertices[this->getVtxIndex()] );
        }
 
        _analysis->execute_vtx( vtx, vtx_rec, m_event );
      }
 
    }
 
    if ( _analysis->debug() ) { 
      m_provider->msg(MSG::INFO) << "Missing track for " << m_chainNames[ichain]  
                     << "\trun "             << run_number 
                     << "\tevent "           << event_number 
                     << "\tlb "              << lumi_block << endmsg;     
    }
 
      }
 
    }
    
    if ( m_provider->msg().level() <= MSG::VERBOSE ) {
      m_provider->msg(MSG::VERBOSE) << "\n\nEvent " << *m_event << endmsg;
    }
  }
 
 
 
  virtual void book() {
    
    if(m_provider->msg().level() <= MSG::VERBOSE)
      m_provider->msg(MSG::VERBOSE) << "AnalysisConfigMT_Tier0::book() " << name() << endmsg;
 
    m_provider->msg(MSG::ERROR) << "AnalysisConfigMT_Tier0::book() should no longer be called: " << name() << endmsg;
 
 
  }
 
 
 
  virtual void finalize() {
 
    if(m_provider->msg().level() <= MSG::VERBOSE){
      m_provider->msg(MSG::VERBOSE) << "AnalysisConfigMT_Tier0::finalise() " << m_provider->name() << endmsg;
    }
    
    m_analysis->finalise();
 
    // deleting instance of TnP_tool and setting pointer to null
    if ( m_TnP_tool != 0 ) {
      delete m_TnP_tool ;
      m_TnP_tool = 0 ;
    }
    
    m_provider->msg(MSG::INFO) << m_provider->name() << " " << m_chainNames[0] << "   \tNRois processed: " << m_NRois << "\tRef tracks: " << m_NRefTracks << "\tTestTracks: " << m_NTestTracks << endmsg;
 
    if(m_provider->msg().level() <= MSG::VERBOSE) {
      m_provider->msg(MSG::VERBOSE) << m_provider->name() << " finalised" << endmsg;
    }
  }
 
  void set_monTool( ToolHandle<GenericMonitoringTool>* m ) { m_monTool=m; }
 
  ToolHandle<GenericMonitoringTool>* monTool() { return m_monTool; }
 
protected:
 
  bool           m_useBeamCondSvc;
 
  TIDA::Event*  m_event;
 
  std::vector<ChainString>     m_chainNames;
  std::vector<A*> m_analyses;
  std::string m_testType;
 
  bool m_doOffline;
  bool m_doMuons;
  bool m_doElectrons;
  bool m_doTaus;
  bool m_doBjets;
  bool m_hasTruthMap;
  bool m_doTauThreeProng;
  bool m_tauEtCutOffline;
 
  std::string m_outputFileName;

  int m_NRois;
  int m_NRefTracks;
  int m_NTestTracks;
 
  bool m_runPurity;
 
  bool m_shifter;
 
  double m_pTthreshold;
 
  bool   m_first;
  
  bool   m_containTracks;
 
  TagNProbe* m_TnP_tool; 
 
  TH1F*      m_invmass;
  TH1F*      m_invmass_obj;
 
  ToolHandle<GenericMonitoringTool>* m_monTool;
 
};
 
 
 
#endif  // TrigInDetAnalysisExample_T_AnalysisConfigMT_Tier0_H