17#include "CLHEP/GenericFunctions/CumulativeChiSquare.hh"
31 declareInterface<ITrackSelectorTool>(
this);
39 InDet::BeamSpotData temp(evt->beamStatus(), evt->beamPosX(), evt->beamPosY(), evt->beamPosZ(),
40 evt->beamPosSigmaX(), evt->beamPosSigmaY(), evt->beamPosSigmaZ(),
41 evt->beamTiltXZ(), evt->beamTiltYZ(), evt->beamPosSigmaXY());
44 ATH_MSG_WARNING(
" Cannot get beamSpot center from xAOD::EventInfo. Using (0,0,0)... " );
52 ATH_MSG_WARNING(
" Cannot get beamSpot center from BeamSpotData. Using (0,0,0)... " );
66 ATH_MSG_DEBUG(
"No TrackSummaryTool set. OK if running on AOD.");
72 ATH_MSG_DEBUG(
"No TrackParticleCreatorTool set but shared hit selection used. OK if running on AOD.");
82 ATH_MSG_ERROR(
" Eta dependent cut on number of TRT hits requested but TrtDCCutTool not specified. ");
83 return StatusCode::FAILURE;
86 return StatusCode::FAILURE;
90 ATH_MSG_DEBUG(
"Using eta dependent cut on number of TRT hits.");
93 ATH_MSG_DEBUG(
"Using eta dependent cut on number of TRT hits + outliers.");
106 ATH_MSG_ERROR(
"Number of cuts DOES NOT match the number of intervals to apply. Please check jobOptions. ");
107 return StatusCode::FAILURE;
109 ATH_MSG_ERROR(
"Zero vectors for number of cuts and pt intervals. Please check jobOptions. ");
110 return StatusCode::FAILURE;
113 return StatusCode::SUCCESS;
120 return StatusCode::SUCCESS;
127 const EventContext& ctx = Gaudi::Hive::currentContext();
134 if (!preselectionDecision) {
135 ATH_MSG_DEBUG(
"Track rejected because of preselection decision!");
139 ATH_MSG_DEBUG(
" Preselection was requested but cannot be made since no Perigee in Track is available. This is not an error." );
143 if (myVertex==
nullptr) {
148 for (
const auto *i : *track.trackParameters()){
149 if ( i->covariance() && !
dynamic_cast<const Trk::Perigee*
>(i)) {
157 firstmeaspar=track.perigeeParameters();
159 ATH_MSG_WARNING(
" First measurment on track is missing. Using perigee Parameters, but they are missing: 0 pointer! Track selection failed " );
161 if (myVertex!=vertex) {
173 track.info().particleHypothesis() ).release();
174 const Trk::Perigee* extrapolatedPerigee = extrapolatedParameters ?
dynamic_cast<const Trk::Perigee*
>(extrapolatedParameters) :
nullptr;
175 if (!extrapolatedPerigee || !extrapolatedPerigee->covariance() ) {
177 if (extrapolatedParameters) {
178 ATH_MSG_WARNING(
"The return object of the extrapolator was not a perigee even if a perigeeSurface was used!" );
179 delete extrapolatedParameters;
180 extrapolatedParameters=
nullptr;
186 bool dec =
decision(extrapolatedPerigee, recVertex ? &recVertex->covariancePosition() :
nullptr );
187 if (myVertex!=vertex) {
191 bool isInTrtAcceptance=
true;
193 isInTrtAcceptance=
false;
195 if (extrapolatedPerigee!=track.perigeeParameters()) {
196 delete extrapolatedPerigee;
197 extrapolatedPerigee=
nullptr;
200 ATH_MSG_DEBUG(
"Track rejected because of perigee parameters!");
205 if (TrkQuality==
nullptr) {
206 ATH_MSG_WARNING(
"Requested cut on track quality was not possible. Track has no FitQuality object attached. Selection failed." );
216 std::unique_ptr<Trk::TrackSummary> summaryUniquePtr;
220 summary = summaryUniquePtr.get();
222 if (
nullptr==summary ) {
223 ATH_MSG_FATAL(
"Track preselection: cannot create a track summary (but useTrackSummary is true). Selection failed." );
230 ATH_MSG_FATAL(
"Track preselection: cannot create a track particle (but useSharedHitInfo is true). Selection failed." );
236 nHitTrt =
m_trtDCTool->minNumberDCs( (*track.trackParameters())[0] );
246 nHitTrtPlusOutliers =
m_trtDCTool->minNumberDCs( (*track.trackParameters())[0] );
255 nHitTrt, nHitTrtPlusOutliers)) {
266 const EventContext& ctx = Gaudi::Hive::currentContext();
274 if (!preselectionDecision) {
275 ATH_MSG_DEBUG(
"Track rejected because of preselection decision!");
279 ATH_MSG_WARNING(
" Preselection was requested but cannot be made since the Perigee is not the defining Parameter of the TrackParticle. This is not an error." );
281 bool isInTrtAcceptance=
true;
283 isInTrtAcceptance=
false;
287 if (TrkQuality==
nullptr) {
288 ATH_MSG_WARNING(
"Requested cut on track quality was not possible. TrackParticleBase has no FitQuality object attached. Selection failed." );
298 if (
nullptr==summary ) {
299 ATH_MSG_WARNING(
"Track preselection: cannot create a track summary (but useTrackSummary is true). Selection failed." );
305 ATH_MSG_ERROR(
"Use of InDetDetailedTrackSelectorTool with Trk::TrackParticleBase and useSharedHitInfo is not supported");
311 nHitTrt =
m_trtDCTool->minNumberDCs( (track.trackParameters())[0] );
320 nHitTrtPlusOutliers =
m_trtDCTool->minNumberDCs( (track.trackParameters())[0] );
328 if ((!perigeeBeforeExtrapolation) or
330 nHitTrt, nHitTrtPlusOutliers))) {
336 if (vertex==
nullptr) {
341 for (
const auto *i : track.trackParameters()) {
342 if (i->covariance() &&
349 if (!extrapolatedPerigee || !extrapolatedPerigee->covariance() ) {
350 ATH_MSG_DEBUG(
" Track Paraemters at first measurement not found. Perigee not found. Cannot do TrackSelection..." );
351 if (myVertex!=vertex) {
358 firstmeaspar=&(track.definingParameters());
370 extrapolatedPerigee = extrapolatedParameters ?
dynamic_cast<const Trk::Perigee*
>(extrapolatedParameters) :
nullptr;
371 if (extrapolatedPerigee==
nullptr || !extrapolatedPerigee->covariance()) {
373 if (extrapolatedParameters) {
374 ATH_MSG_WARNING(
"The return object of the extrapolator was not a perigee even if a perigeeSurface was used!" );
375 delete extrapolatedParameters;
376 extrapolatedParameters =
nullptr;
379 if (extrapolatedParameters)
ATH_MSG_VERBOSE (
"Result: " << *extrapolatedParameters);
381 bool dec =
decision(extrapolatedPerigee, recVertex ? &recVertex->covariancePosition() :
nullptr );
382 if (myVertex!=vertex) {
386 if (extrapolatedPerigee!=&(track.definingParameters())) {
387 delete extrapolatedPerigee;
388 extrapolatedPerigee=
nullptr;
391 ATH_MSG_DEBUG(
"Track rejected because of perigee parameters!");
399 if(vertex)
return vertex->position();
403 InDet::BeamSpotData temp(evt->beamStatus(), evt->beamPosX(), evt->beamPosY(), evt->beamPosZ(),
404 evt->beamPosSigmaX(), evt->beamPosSigmaY(), evt->beamPosSigmaZ(),
405 evt->beamTiltXZ(), evt->beamTiltYZ(), evt->beamPosSigmaXY());
406 return temp.beamVtx().position();
408 ATH_MSG_WARNING(
" Cannot get beamSpot center from xAOD::EventInfo. Using (0,0,0)... " );
413 if (beamSpotHandle.
isValid()) {
414 return beamSpotHandle->beamVtx().position();
416 ATH_MSG_WARNING(
" Cannot get beamSpot center from BeamSpotData. Using (0,0,0)... " );
426 const EventContext& ctx = Gaudi::Hive::currentContext();
433 ATH_MSG_DEBUG(
"Track rejected because of preselection decision!");
454 nHitTrtPlusOutliers =
m_trtDCTool->minNumberDCs( &perigee );
476 ATH_MSG_DEBUG(
"Track rejected because of Pt-Dependent SCT Hit cut (CAREFUL! Excludes dead modules)") ;
483 ATH_MSG_DEBUG(
"Track rejected because of Pt-Dependent SCT Hit cut (CAREFUL! Excludes dead modules)") ;
495 ATH_MSG_DEBUG(
"and track rejected because at least one hit is expected in the innermost pixel layer") ;
497 }
else ATH_MSG_DEBUG(
"recovered track as no b-layer expected") ;
550 ATH_MSG_DEBUG(
"Track rejected because of nHitTrt "<<nh<<
" < "<<nHitTrt);
555 if (nhh<nHitTrtPlusOutliers) {
556 ATH_MSG_DEBUG(
"Track rejected because of nHitTrtPlusOutliers "<<nhh<<
" < "<<nHitTrtPlusOutliers);
576 ATH_MSG_DEBUG(
"Track rejected because numberOfTRTHits is zero.");
590 if (nheh>1.) nheh=1.;
628 perigee,perigeeSurface,
630 const Trk::Perigee* extrapolatedPerigee = extrapolatedParameters ?
dynamic_cast<const Trk::Perigee*
>(extrapolatedParameters) :
nullptr;
631 if (extrapolatedPerigee==
nullptr) {
632 ATH_MSG_WARNING(
"Extrapolation to the vertex failed: " << perigeeSurface << std::endl << perigee );
633 if (extrapolatedParameters!=
nullptr) {
634 ATH_MSG_WARNING(
"The return object of the extrapolator was not a perigee even if a perigeeSurface was used!" );
635 delete extrapolatedParameters;
636 extrapolatedParameters=
nullptr;
643 const AmgSymMatrix(3)& vertexError = vertex->covariancePosition();
644 dec =
decision(extrapolatedPerigee,&vertexError);
646 dec =
decision(extrapolatedPerigee,
nullptr);
649 delete extrapolatedPerigee;
652 ATH_MSG_DEBUG(
"Track rejected because of perigee parameters!");
664 if(
nullptr==track || !track->covariance()) {
665 ATH_MSG_WARNING(
"Decision on measured perigee: Zero pointer to measured perigee passed. Selection failed." );
669 const AmgVector(5)& perigeeParms = track->parameters();
672 const EventContext& ctx = Gaudi::Hive::currentContext();
675 if (fieldCondObj ==
nullptr) {
687 double p = std::fabs(1./perigeeParms[
Trk::qOverP]);
692 double pt = p*std::sin(perigeeParms[
Trk::theta]);
726 double sinTheta = std::sin(perigeeParms[
Trk::theta]);
727 double cosTheta = std::cos(perigeeParms[
Trk::theta]);
728 double d0wrtPriVtx = perigeeParms[
Trk::d0];
729 double deltaZ = perigeeParms[
Trk::z0];
730 double z0wrtPriVtx = deltaZ*sinTheta;
731 double testtrackSigD0 = sqrt( (*track->covariance())(
Trk::d0,
Trk::d0) );
732 double testtrackSigZ0 = sqrt( (*track->covariance())(
Trk::z0,
Trk::z0) );
735 double trackPhi = perigeeParms[
Trk::phi];
736 double dIPdx = std::sin(trackPhi);
737 double dIPdy = -std::cos(trackPhi);
738 double DD0 = testtrackSigD0*testtrackSigD0;
740 if (covariancePosition) {
741 double DXX = dIPdx*dIPdx* (*covariancePosition)(0,0);
742 double DYY = dIPdy*dIPdy* (*covariancePosition)(1,1);
743 double DXY = 2.*dIPdx*dIPdy* (*covariancePosition)(0,1);
744 newD0Err = DD0 + DXX + DYY + DXY;
749 double d0ErrwrtPriVtx = (newD0Err>0 ? sqrt(newD0Err) : -10e-9);
751 if (d0ErrwrtPriVtx<0) {
752 ATH_MSG_WARNING(
" error on d0 is negative: numeric error... (not expected. please report!)" );
765 double dZIPdTheta = deltaZ*cosTheta;
766 double dZIPdz0 = sinTheta;
767 double dZIPdzV = -sinTheta;
768 double DTheta2 = dZIPdTheta*dZIPdTheta*testtrackSigTh*testtrackSigTh;
769 double DZ02 = dZIPdz0*dZIPdz0*testtrackSigZ0*testtrackSigZ0;
770 double DThetaZ0 = 2.*dZIPdTheta*dZIPdz0*(*track->covariance())(
Trk::theta,
Trk::z0);
772 if (covariancePosition) {
773 double DZV2 = dZIPdzV*dZIPdzV* (*covariancePosition)(2,2);
774 newZ0Err = DTheta2 + DZ02 + DZV2 + DThetaZ0;
776 newZ0Err = DTheta2 + DZ02 + DThetaZ0;
779 double z0ErrwrtPriVtx = (newZ0Err>0 ? sqrt(newZ0Err) : -10e-9);
781 if (z0ErrwrtPriVtx<0) {
782 ATH_MSG_WARNING(
" error on z0 is negative: numeric error... (not expected. please report!)" );
793 if (std::fabs(track->momentum().eta())>
m_etaMax) {
794 ATH_MSG_DEBUG(
"Track rejected because of fabs(eta) " << std::fabs(track->momentum().eta()) <<
" > " <<
m_etaMax);
804 if(
nullptr == trkQuality) {
805 ATH_MSG_WARNING(
"Null FitQuality pointer passed. No track Quality cut possible. Selection failed." );
815 if(ndf>0 &&
chi2>=0.) {
816 Genfun::CumulativeChiSquare myCumulativeChiSquare(ndf);
817 proba = 1.-myCumulativeChiSquare(
chi2);
845 bool useSharedHitInfo,
849 const int nHitTrtPlusOutliers)
const
851 if (summary==
nullptr) {
852 ATH_MSG_WARNING(
"Null TrackSummary pointer passed. Selection failed." );
870 if (nhp < 0) nhp = 0;
873 if (nhs < 0) nhs = 0;
876 if (ndhs < 0) ndhs = 0;
884 const AmgVector(5)& perigeeParms = track->parameters();
885 double p = std::fabs(1./perigeeParms[
Trk::qOverP]);
886 double pt = p*std::sin(perigeeParms[
Trk::theta]);
891 ATH_MSG_DEBUG(
"Track rejected because of Pt-Dependent SCT Hit cut (CAREFUL! Excludes dead modules)") ;
898 ATH_MSG_DEBUG(
"Track rejected because of Pt-Dependent SCT Hit cut (CAREFUL! Excludes dead modules)") ;
911 ATH_MSG_DEBUG(
"and no blayer tool configured, so will not try to recover track");
914 ATH_MSG_DEBUG(
"and track rejected because at least one hit is expected in the innermost pixel layer") ;
916 }
else ATH_MSG_DEBUG(
"recovered track as no b-layer expected") ;
982 ATH_MSG_DEBUG(
"Track rejected because of nHitTrt "<<nh<<
" < "<<nHitTrt);
988 if (nhh<nHitTrtPlusOutliers) {
989 ATH_MSG_DEBUG(
"Track rejected because of nHitTrtPlusOutliers "<<nhh<<
" < "<<nHitTrtPlusOutliers);
994 if (nhthits<0) nhthits=0;
1001 if (nhthitsWithOutliers<0) nhthitsWithOutliers=0;
1007 if (summary->get( Trk :: numberOfTRTHits )>0) {
1016 if ( summary->get( Trk :: numberOfTRTHits ) + summary->get( Trk :: numberOfTRTOutliers ) > 0 ) {
1019 if(nheh<0.) nheh=0.;
1020 if (nheh>1.) nheh=1.;
1028 if (useSharedHitInfo) {
1030 ATH_MSG_DEBUG(
"Track rejected because xAOD::TrackParticle not available");
1035 if(nbs < 0) nbs = 0;
1043 if(nps < 0) nps = 0;
1050 if(nss < 0) nss = 0;
1056 int nst = nps + nss;
1070 const AmgVector(5)& perigeeParms = myPerigee.parameters();
1073 const EventContext& ctx = Gaudi::Hive::currentContext();
1076 if (fieldCondObj ==
nullptr) {
1085 ATH_MSG_DEBUG(
"Track rejected because of perigee qOverP == 0.");
1088 double p = std::fabs(1./perigeeParms[
Trk::qOverP]);
1093 double pt = p*std::sin(perigeeParms[
Trk::theta]);
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_WARNING(x)
#define AmgSymMatrix(dim)
void getInitializedCache(MagField::AtlasFieldCache &cache) const
get B field cache for evaluation as a function of 2-d or 3-d position.
Local cache for magnetic field (based on MagFieldServices/AtlasFieldSvcTLS.h).
bool solenoidOn() const
status of the magnets
Class to represent and store fit qualities from track reconstruction in terms of and number of degre...
int numberDoF() const
returns the number of degrees of freedom of the overall track or vertex fit as integer
double chiSquared() const
returns the of the overall track fit
const Amg::Vector3D & momentum() const
Access method for the momentum.
Class describing the Line to which the Perigee refers to.
Trk::RecVertex inherits from Trk::Vertex.
A summary of the information contained by a track.
This class is a simplest representation of a vertex candidate.
const Amg::Vector3D & position() const
return position of vertex
float numberDoF() const
Returns the number of degrees of freedom of the overall track or vertex fit as float.
const Trk::Perigee & perigeeParameters() const
Returns the Trk::MeasuredPerigee track parameters.
virtual double pt() const override final
The transverse momentum ( ) of the particle.
float chiSquared() const
Returns the of the overall track fit.
virtual double eta() const override final
The pseudorapidity ( ) of the particle.
double chi2(TH1 *h0, TH1 *h1)
Eigen::Matrix< double, 3, 1 > Vector3D
ParametersT< TrackParametersDim, Charged, PerigeeSurface > Perigee
ParametersBase< TrackParametersDim, Charged > TrackParameters
@ numberOfSCTHits
number of SCT holes
@ numberOfPixelHits
number of pixel layers on track with absence of hits
@ numberOfTRTHighThresholdOutliers
number of dead TRT straws crossed
@ numberOfTRTOutliers
number of TRT holes
@ numberOfTRTHits
number of TRT outliers
@ numberOfInnermostPixelLayerHits
these are the hits in the 1st pixel layer
@ numberOfTRTHighThresholdHits
total number of TRT hits which pass the high threshold
@ numberOfSCTHoles
number of Holes in both sides of a SCT module
@ numberOfSCTDeadSensors
number of TRT hits
@ numberOfPixelHoles
number of pixels which have a ganged ambiguity.
@ numberOfPixelDeadSensors
number of pixel hits with broad errors (width/sqrt(12))
TrackParticle_v1 TrackParticle
Reference the current persistent version:
Vertex_v1 Vertex
Define the latest version of the vertex class.
@ expectInnermostPixelLayerHit
Do we expect a 0th-layer barrel hit for this track?
@ numberOfPixelHoles
number of pixel layers on track with absence of hits [unit8_t].
@ numberOfTRTHighThresholdOutliers
number of TRT high threshold outliers (only xenon counted) [unit8_t].
@ numberOfInnermostPixelLayerSharedHits
number of Pixel 0th layer barrel hits shared by several tracks.
@ numberOfTRTHits
number of TRT hits [unit8_t].
@ numberOfSCTDeadSensors
number of dead SCT sensors crossed [unit8_t].
@ numberOfSCTHits
number of hits in SCT [unit8_t].
@ numberOfSCTDoubleHoles
number of Holes in both sides of a SCT module [unit8_t].
@ numberOfInnermostPixelLayerHits
these are the hits in the 0th pixel barrel layer
@ numberOfPixelHits
these are the pixel hits, including the b-layer [unit8_t].
@ numberOfPixelSharedHits
number of Pixel all-layer hits shared by several tracks [unit8_t].
@ numberOfSCTSharedHits
number of SCT hits shared by several tracks [unit8_t].
@ numberOfTRTHighThresholdHits
number of TRT hits which pass the high threshold (only xenon counted) [unit8_t].
@ numberOfTRTOutliers
number of TRT outliers [unit8_t].
@ numberOfPixelDeadSensors
number of dead pixel sensors crossed [unit8_t].
@ numberOfSCTHoles
number of SCT holes [unit8_t].