InDetSlidingWindowTrackTimeTool.cxx

Go to the documentation of this file.

/*
   Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
 */
 
// InDetSlidingWindowTrackTimeTool.cxx, (c) ATLAS Detector software
 
#include "InDetCosmicsEventPhase/InDetSlidingWindowTrackTimeTool.h"
 
#include "InDetRIO_OnTrack/TRT_DriftCircleOnTrack.h"
 
#include "TRT_ConditionsData/RtRelation.h"
#include "TRT_ConditionsData/BasicRtRelation.h"
 
#include "InDetIdentifier/TRT_ID.h"
#include "TrkTrack/Track.h"
#include "TrkSegment/TrackSegment.h"
 
#include "StoreGate/ReadCondHandle.h"
#include <cmath>
 
//================ Constructor =================================================
 
InDet::InDetSlidingWindowTrackTimeTool::InDetSlidingWindowTrackTimeTool(std::string const& t
                                                                        , std::string const& n
                                                                        , IInterface const* p
                                                                        ) : AthAlgTool(t, n, p)
  , m_caldbtool("TRT_CalDbTool", this) {
  declareInterface<IInDetCosmicsEventPhaseTool>(this);
  declareProperty("TRTCalDbTool", m_caldbtool);
  declareProperty("GlobalOffset", m_globalOffset = 10.);
  declareProperty("UseNewEP", m_useNewEP = true);
  declareProperty("NumberIterations", m_nIterations = 5);
  declareProperty("WindowSize", m_windowSize = 7.);
}
 
//================ Destructor =================================================
 
InDet::InDetSlidingWindowTrackTimeTool::~InDetSlidingWindowTrackTimeTool()
{}
 
 
//================ Initialisation =================================================
 
StatusCode InDet::InDetSlidingWindowTrackTimeTool::initialize() {
  ATH_CHECK( m_caldbtool.retrieve());
  // Read key
  ATH_CHECK( m_T0ReadKey.initialize() );
 
  return StatusCode::SUCCESS;
}
 
//================ Finalisation =================================================
 
StatusCode InDet::InDetSlidingWindowTrackTimeTool::finalize() {
 //
  return StatusCode::SUCCESS;
}
 
 
double InDet::InDetSlidingWindowTrackTimeTool::findPhase(Trk::Track const* track) const {
  ATH_MSG_DEBUG("Finding phase...");
  
  double aT0=0.;
  if (!m_useNewEP) {
    ATH_MSG_INFO(" Set averageT0 = 0 (m_useNewEP=false) ");
  } else {
    SG::ReadCondHandle<TRTCond::AverageT0> rhl{m_T0ReadKey};
    const TRTCond::AverageT0* avgT0{*rhl};
    aT0=avgT0->get();
  }
 
 
  TRTCond::RtRelation const* rtr = nullptr;
 
  double timeresidualsum = 0;
  size_t ntrthits = 0;
  double windowCenter = 0;
  for (int itr = 0; itr != m_nIterations; ++itr) {
    timeresidualsum = 0;
    ntrthits = 0;
    for (Trk::TrackStateOnSurface const* state : *track->trackStateOnSurfaces()) {
      Trk::MeasurementBase const* mesb = state->measurementOnTrack();
      if (!mesb
          || !state->type(Trk::TrackStateOnSurface::Measurement)
          ) continue;
      InDet::TRT_DriftCircleOnTrack const* trtcirc = dynamic_cast<InDet::TRT_DriftCircleOnTrack const*>(mesb);
      if (!trtcirc) continue;
      InDet::TRT_DriftCircle const* rawhit = trtcirc->prepRawData();
      if (!rawhit) continue;
      if (!rawhit->driftTimeValid()
          || rawhit->firstBinHigh()
          ) continue;
      Identifier const& ident = trtcirc->identify();
      double rawdrifttime = rawhit->rawDriftTime();
      double t0 = m_caldbtool->getT0(ident);
 
      ATH_MSG_DEBUG("T0 : " << t0);
      rtr = m_caldbtool->getRtRelation(ident);
 
      Trk::TrackParameters const* tparp = state->trackParameters();
      if (!tparp) continue;
      double trkdistance = tparp->parameters()[Trk::driftRadius];
      double trkdrifttime = rtr->drifttime(std::abs(trkdistance));
      double timeresidual = rawdrifttime - t0 + aT0 - trkdrifttime;
      ATH_MSG_DEBUG("trkdistance=" << trkdistance
                                   << "  trkdrifttime=" << trkdrifttime
                                   << "  timeresidual=" << timeresidual
                                   << " rawdrifttime=" << rawdrifttime);
      if (timeresidual < 2000
          && std::abs(trkdistance) < 2.8) {
        if (itr == 0
            || std::abs(timeresidual - windowCenter) < m_windowSize) {
          timeresidualsum += timeresidual;
          ++ntrthits;
        }
      }
    }
    if (ntrthits > 1) {
      windowCenter = timeresidualsum / ntrthits;
    } else break;
  }
 
  ATH_MSG_DEBUG("timeresidualsum = " << timeresidualsum);
  ATH_MSG_DEBUG("ntrtrhits = " << ntrthits);
 
  if (ntrthits > 1)
    return windowCenter + m_globalOffset;
  return 0.;
}
 
double InDet::InDetSlidingWindowTrackTimeTool::findPhase(Trk::Segment const* segment) const {
 
  double sum_tr = 0.;
  double sum_goodhits = 0.;
  int nhits = segment->numberOfMeasurementBases();
 
  std::vector<double> data{0., 0.1, -0.00087,0.};
  TRTCond::RtRelation const* rtr = nullptr;
 
  for (int i = 0; i < nhits; ++i) {
    Trk::RIO_OnTrack const* rio = dynamic_cast<Trk::RIO_OnTrack const*>(segment->measurement(i));
    InDet::TRT_DriftCircleOnTrack const* trtcirc = dynamic_cast<InDet::TRT_DriftCircleOnTrack const*>(rio);
    if (!trtcirc) continue;
    InDet::TRT_DriftCircle const* rawhit = trtcirc->prepRawData();
    if (!rawhit) continue;
    if (!rawhit->lastBinHigh() && !rawhit->isNoise()) {
      Identifier const& ident = trtcirc->identify();
      rtr = m_caldbtool->getRtRelation(ident);
      if (not rtr) {
        ATH_MSG_WARNING("Rt relation pointer is null!");
        return 0.;
      }
      double full_drifttime = rtr->drifttime(2.0);
      sum_tr += rawhit->trailingEdge() * 3.125 - full_drifttime;
      ATH_MSG_VERBOSE("Hit " << sum_goodhits
                             << " : " << rawhit->trailingEdge() * 3.125
                             << "  fulldrifttime=" << full_drifttime);
      sum_goodhits += 1;
    } else {
      ATH_MSG_VERBOSE("Hit has lastbin high");
    }
  }
  if (sum_goodhits > 1) return sum_tr / sum_goodhits + m_globalOffset;
  return 0;
}
 
double InDet::InDetSlidingWindowTrackTimeTool::findPhaseFromTE(Trk::Track const* track) const {
 
  ATH_MSG_DEBUG("Finding phase...");
 
  double aT0=0.;
  if (!m_useNewEP) {
    ATH_MSG_INFO(" Set averageT0 = 0 (m_useNewEP=false) ");
  } else {
    SG::ReadCondHandle<TRTCond::AverageT0> rhl{m_T0ReadKey};
    const TRTCond::AverageT0* avgT0{*rhl};
    aT0=avgT0->get();
  }
 
  double timeresidualsum = 0;
  size_t ntrthits = 0;
  for (Trk::TrackStateOnSurface const* state : *track->trackStateOnSurfaces()) {
    Trk::MeasurementBase const* mesb = state->measurementOnTrack();
    if (!mesb
        || !state->type(Trk::TrackStateOnSurface::Measurement)
        ) continue;
    InDet::TRT_DriftCircleOnTrack const* trtcirc = dynamic_cast<InDet::TRT_DriftCircleOnTrack const*>(mesb);
    if (!trtcirc) continue;
    InDet::TRT_DriftCircle const* rawhit = trtcirc->prepRawData();
    if (!rawhit) continue;
    Identifier const& ident = trtcirc->identify();
    double rawtrailingedge = rawhit->trailingEdge() * 3.125;
    double t0 = m_caldbtool->getT0(ident);
    ATH_MSG_DEBUG("T0 : " << t0);
    double timeresidual = rawtrailingedge - t0 + aT0;
    ATH_MSG_DEBUG("timeresidual=" << timeresidual);
    if (timeresidual < 2000) {
      timeresidualsum += timeresidual;
      ++ntrthits;
    }
  }
  if (ntrthits > 1) return timeresidualsum / ntrthits + m_globalOffset;
  return 0.;
}