AsgLeptonTrackSelectionAlg.cxx

Go to the documentation of this file.

/*
  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
*/
 
 
 
 
//
// includes
//
 
#include <AsgAnalysisAlgorithms/AsgLeptonTrackSelectionAlg.h>
 
#include <xAODEgamma/Electron.h>
#include <xAODMuon/Muon.h>
#include <xAODTracking/TrackParticlexAODHelpers.h>
 
//
// method implementations
//
 
namespace CP
{
 
  StatusCode AsgLeptonTrackSelectionAlg ::
  initialize ()
  {
 
    if (m_maxD0Significance < 0 || !std::isfinite (m_maxD0Significance))
    {
      ATH_MSG_ERROR ("invalid value of maxD0Significance: " << m_maxD0Significance);
      return StatusCode::FAILURE;
    }
    if (m_maxDeltaZ0SinTheta < 0 || !std::isfinite (m_maxDeltaZ0SinTheta))
    {
      ATH_MSG_ERROR ("invalid value of maxDeltaZ0SinTheta: " << m_maxDeltaZ0SinTheta);
      return StatusCode::FAILURE;
    }
 
    m_accept.addCut ("trackRetrieval", "whether the track retrieval failed");
    if (m_maxD0Significance > 0)
      m_accept.addCut ("maxD0Significance", "maximum D0 significance cut");
    if (m_maxDeltaZ0SinTheta > 0)
      m_accept.addCut ("maxDeltaZ0SinTheta", "maximum Delta z0 sin theta cut");
    if (m_nMinPixelHits != -1 || m_nMaxPixelHits != -1)
      m_accept.addCut ("numPixelHits", "Minimum and/or maxiumum Pixel hits");
    if (m_nMinSCTHits != -1 || m_nMaxSCTHits != -1)
      m_accept.addCut ("numSCTHits", "Minimum and/or maxiumum SCT hits");
 
    if(m_decorateTTVAVars){
      if (m_d0sigDecoration.empty()){
        ANA_MSG_ERROR ("No d0significance decoration name set");
        return StatusCode::FAILURE;
      } else {
        m_d0sigDecorator = std::make_unique<SG::AuxElement::Decorator<float> > (m_d0sigDecoration);
      }
      if (m_z0sinthetaDecoration.empty()){
        ANA_MSG_ERROR ("No z0sintheta decoration name set");
        return StatusCode::FAILURE;
      } else {
        m_z0sinthetaDecorator = std::make_unique<SG::AuxElement::Decorator<float> > (m_z0sinthetaDecoration);
      }
    } else{
      if (!m_d0sigDecoration.empty()){
        ANA_MSG_WARNING ("d0significance decoration name set, please set decorateTTVVars to True if you want this to be written out");
      }
      if (!m_z0sinthetaDecoration.empty()){
        ANA_MSG_WARNING ("z0sintheta decoration name set, please set decorateTTVVars to True if you want this to be written out");
      }
    }
 
    ANA_CHECK (m_particlesHandle.initialize (m_systematicsList));
    ANA_CHECK (m_preselection.initialize (m_systematicsList, m_particlesHandle, SG::AllowEmpty));
    ANA_CHECK (m_selectionHandle.initialize (m_systematicsList, m_particlesHandle));
    ANA_CHECK (m_systematicsList.initialize());
 
    ANA_CHECK (m_eventInfoKey.initialize());
    ANA_CHECK (m_primaryVerticesKey.initialize());
 
    if (!m_nameSvc.empty())
    {
      ANA_CHECK (m_nameSvc.retrieve());
      ANA_CHECK (m_nameSvc->addAcceptInfo (m_particlesHandle.getNamePattern(), m_selectionHandle.getLabel(),
          m_accept));
    }
 
    return StatusCode::SUCCESS;
  }
 
 
 
  StatusCode AsgLeptonTrackSelectionAlg ::
  execute ()
  {
    SG::ReadHandle<xAOD::EventInfo> eventInfo(m_eventInfoKey);
    SG::ReadHandle<xAOD::VertexContainer> vertices(m_primaryVerticesKey);
    const xAOD::Vertex *primaryVertex {nullptr};
 
    for (const xAOD::Vertex *vertex : *vertices)
    {
      if (vertex->vertexType() == xAOD::VxType::PriVtx)
      {
        // The default "PrimaryVertex" container is ordered in
        // sum-pt, and the tracking group recommends to pick the one
        // with the maximum sum-pt, so this will do the right thing.
        // If the user needs a different primary vertex, they need to
        // provide a reordered primary vertex container and point
        // this algorithm to it.  Currently there is no central
        // algorithm to do that, so users will have to write their
        // own (15 Aug 18).
        if (primaryVertex == nullptr)
        {
          primaryVertex = vertex;
          break;
        }
      }
    }
 
    for (const auto& sys : m_systematicsList.systematicsVector())
    {
      const xAOD::IParticleContainer *particles = nullptr;
      ANA_CHECK (m_particlesHandle.retrieve (particles, sys));
      for (const xAOD::IParticle *particle : *particles)
      {
        asg::AcceptData acceptData (&m_accept);
        float d0sig = -999;
        float deltaZ0SinTheta = -999;
 
        if (m_preselection.getBool (*particle, sys))
        {
          std::size_t cutIndex {0};
          
          const xAOD::TrackParticle *track {nullptr};
          if (const xAOD::Muon *muon = dynamic_cast<const xAOD::Muon *>(particle)){
            track = muon->primaryTrackParticle();
          } else if (const xAOD::Electron *electron = dynamic_cast<const xAOD::Electron *>(particle)){
            track = electron->trackParticle();
          } else {
            ANA_MSG_ERROR ("failed to cast input to electron or muon");
            return StatusCode::FAILURE;
          }
 
          acceptData.setCutResult (cutIndex ++, track != nullptr);
          
          if (track != nullptr) {
            try {
              d0sig = xAOD::TrackingHelpers::d0significance(track, eventInfo->beamPosSigmaX(), eventInfo->beamPosSigmaY(), eventInfo->beamPosSigmaXY());
              if (m_maxD0Significance > 0) acceptData.setCutResult (cutIndex ++, fabs( d0sig ) < m_maxD0Significance);
            
            } catch (const std::runtime_error &) {
              acceptData.setCutResult (cutIndex ++, false);
            }
            
            const double vertex_z = primaryVertex ? primaryVertex->z() : 0;
            deltaZ0SinTheta = (track->z0() + track->vz() - vertex_z) * sin (particle->p4().Theta());
            if (m_maxDeltaZ0SinTheta > 0) acceptData.setCutResult (cutIndex ++, fabs (deltaZ0SinTheta) < m_maxDeltaZ0SinTheta);
 
            if (m_nMinPixelHits != -1 || m_nMaxPixelHits != -1) {
              uint8_t nPixelHits;
              track->summaryValue(nPixelHits, xAOD::numberOfPixelHits);
              bool accept = true;
              if(m_nMinPixelHits != -1) {
                accept &= nPixelHits >= m_nMinPixelHits;
              }
              if(m_nMaxPixelHits != -1) {
                accept &= nPixelHits <= m_nMaxPixelHits;
              }
              acceptData.setCutResult (cutIndex++, accept);
            }
            if (m_nMinSCTHits != -1 || m_nMaxSCTHits != -1) {
              uint8_t nSCTHits;
              track->summaryValue(nSCTHits, xAOD::numberOfSCTHits);
              bool accept = true;
              if(m_nMinSCTHits != -1) {
                accept &= nSCTHits >= m_nMinSCTHits;
              }
              if(m_nMaxSCTHits != -1) {
                accept &= nSCTHits <= m_nMaxSCTHits;
              }
              acceptData.setCutResult (cutIndex++, accept);
            }
          }
        }
        if(m_decorateTTVAVars){
          (*m_z0sinthetaDecorator)(*particle) = deltaZ0SinTheta;
          (*m_d0sigDecorator)(*particle) = d0sig;
        }
        m_selectionHandle.setBits(*particle, selectionFromAccept (acceptData), sys);
      }
    }
 
    return StatusCode::SUCCESS;
  }
}