File indexing completed on 2023-03-17 11:23:21
0001 #include "RecoVertex/KalmanVertexFit/interface/SingleTrackVertexConstraint.h"
0002 #include "DataFormats/GeometryCommonDetAlgo/interface/GlobalError.h"
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004 #include "MagneticField/Engine/interface/MagneticField.h"
0005
0006 #include <algorithm>
0007 using namespace std;
0008 using namespace reco;
0009
0010 namespace {
0011
0012
0013
0014 const float TrackerBoundsRadius = 112;
0015 const float TrackerBoundsHalfLength = 273.5;
0016 bool insideTrackerBounds(const GlobalPoint& point) {
0017 return ((point.transverse() < TrackerBoundsRadius) && (abs(point.z()) < TrackerBoundsHalfLength));
0018 }
0019 }
0020
0021 SingleTrackVertexConstraint::BTFtuple SingleTrackVertexConstraint::constrain(const TransientTrack& track,
0022 const GlobalPoint& priorPos,
0023 const GlobalError& priorError) const {
0024 VertexState priorVertexState(priorPos, priorError);
0025 return constrain(track, priorVertexState);
0026 }
0027
0028 SingleTrackVertexConstraint::BTFtuple SingleTrackVertexConstraint::constrain(const TransientTrack& track,
0029 const VertexState priorVertexState) const {
0030
0031
0032 typedef CachingVertex<5>::RefCountedVertexTrack RefCountedVertexTrack;
0033 typedef VertexTrack<5>::RefCountedLinearizedTrackState RefCountedLinearizedTrackState;
0034
0035 double field = track.field()->inInverseGeV(track.impactPointState().globalPosition()).z();
0036 int nominalBfield = track.field()->nominalValue();
0037 if ((fabs(field) < 1e-4) && (fabs(nominalBfield) != 0)) {
0038 LogDebug("RecoVertex/SingleTrackVertexConstraint")
0039 << "Initial state is very far, field is close to zero (<1e-4): " << field << "\n";
0040 return BTFtuple(false, TransientTrack(), 0.);
0041 }
0042
0043 RefCountedLinearizedTrackState lTrData = theLTrackFactory.linearizedTrackState(priorVertexState.position(), track);
0044 RefCountedVertexTrack vertexTrack = theVTrackFactory.vertexTrack(lTrData, priorVertexState);
0045
0046
0047
0048 std::vector<RefCountedVertexTrack> initialTracks;
0049 CachingVertex<5> vertex(priorVertexState, priorVertexState, initialTracks, 0);
0050 vertex = vertexUpdator.add(vertex, vertexTrack);
0051 if (!vertex.isValid()) {
0052 return BTFtuple(false, TransientTrack(), 0.);
0053 } else if (doTrackerBoundCheck_ && (!insideTrackerBounds(vertex.position()))) {
0054 LogDebug("RecoVertex/SingleTrackVertexConstraint") << "Fitted position is out of tracker bounds.\n";
0055 return BTFtuple(false, TransientTrack(), 0.);
0056 }
0057
0058 RefCountedVertexTrack nTrack = theVertexTrackUpdator.update(vertex, vertexTrack);
0059 return BTFtuple(true, nTrack->refittedState()->transientTrack(), nTrack->smoothedChi2());
0060 }
0061
0062 SingleTrackVertexConstraint::BTFtuple SingleTrackVertexConstraint::constrain(const FreeTrajectoryState& fts,
0063 const GlobalPoint& priorPos,
0064 const GlobalError& priorError) const {
0065 return constrain(ttFactory.build(fts), priorPos, priorError);
0066 }
0067
0068 SingleTrackVertexConstraint::BTFtuple SingleTrackVertexConstraint::constrain(const TransientTrack& track,
0069 const reco::BeamSpot& spot) const {
0070 VertexState priorVertexState(spot);
0071 return constrain(track, priorVertexState);
0072 }
0073
0074 SingleTrackVertexConstraint::BTFtuple SingleTrackVertexConstraint::constrain(const FreeTrajectoryState& fts,
0075 const reco::BeamSpot& spot) const {
0076 VertexState priorVertexState(spot);
0077 return constrain(ttFactory.build(fts), priorVertexState);
0078 }