CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
KFTrajectoryFitter.cc
Go to the documentation of this file.
10 
12 
14 
15  if(aTraj.empty()) return Trajectory();
16 
17  const TM& firstTM = aTraj.firstMeasurement();
18  TSOS firstTsos = TrajectoryStateWithArbitraryError()(firstTM.updatedState());
19 
20  return fitOne(aTraj.seed(), aTraj.recHits(), firstTsos,type);
21 }
22 
24  const RecHitContainer&, fitType) const{
25 
26  throw cms::Exception("TrackFitters",
27  "KFTrajectoryFitter::fit(TrajectorySeed, <TransientTrackingRecHit>) not implemented");
28 
29  return Trajectory();
30 }
31 
33  const RecHitContainer& hits,
34  const TSOS& firstPredTsos,fitType) const
35 {
36  if(hits.empty()) return Trajectory();
37 
38 
39  if unlikely(aSeed.direction() == anyDirection)
40  throw cms::Exception("KFTrajectoryFitter","TrajectorySeed::direction() requested but not set");
41 
42  std::unique_ptr<Propagator> p_cloned = SetPropagationDirection(*thePropagator,
43  aSeed.direction());
44 
45 #ifdef EDM_ML_DEBUG
46  LogDebug("TrackFitters")
47  <<" +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
48  <<" KFTrajectoryFitter::fit starting with " << hits.size() <<" HITS";
49 
50  for (unsigned int j=0;j<hits.size();j++) {
51  if (hits[j]->det())
52  LogTrace("TrackFitters") << "hit #:" << j+1 << " rawId=" << hits[j]->det()->geographicalId().rawId()
53  << " validity=" << hits[j]->isValid();
54  else
55  LogTrace("TrackFitters") << "hit #:" << j+1 << " Hit with no Det information";
56  }
57  LogTrace("TrackFitters") << " INITIAL STATE "<< firstPredTsos;
58 #endif
59 
60  Trajectory ret(aSeed, p_cloned->propagationDirection());
61  Trajectory & myTraj = ret;
62  myTraj.reserve(hits.size());
63 
64  TSOS predTsos(firstPredTsos);
65  TSOS currTsos;
66 
67  int hitcounter = 0;
68  for(auto ihit : hits) {
69  ++hitcounter;
70 
71  const TransientTrackingRecHit & hit = (*ihit);
72 
73  // if unlikely(hit.det() == nullptr) continue;
74 
75  if unlikely( (!hit.isValid()) && hit.surface() == nullptr) {
76  LogDebug("TrackFitters")<< " Error: invalid hit with no GeomDet attached .... skipping";
77  continue;
78  }
79  // if (hit.det() && hit.geographicalId()<1000U) LogDebug("TrackFitters")<< "Problem 0 det id for " << typeid(hit).name() << ' ' << hit.det()->geographicalId() ;
80  // if (hit.isValid() && hit.geographicalId()<1000U) LogDebug("TrackFitters")<< "Problem 0 det id for " << typeid(hit).name() << ' ' << hit.det()->geographicalId();
81 
82  if ( hitcounter != 1) //no propagation needed for the first hit
83  predTsos = p_cloned->propagate( currTsos, *(hit.surface()) );
84 
85 
86  if unlikely(!predTsos.isValid()) {
87  LogDebug("TrackFitters")
88  << "SOMETHING WRONG !" << "\n"
89  << "KFTrajectoryFitter: predicted tsos not valid!\n"
90  << "current TSOS: " << currTsos << "\n";
91 
92  if(hit.surface()) LogTrace("TrackFitters") << "next Surface: " << hit.surface()->position() << "\n";
93 
94  if( myTraj.foundHits() >= minHits_ ) {
95  LogDebug("TrackFitters") << " breaking trajectory" << "\n";
96  break;
97  } else {
98  LogDebug("TrackFitters") << " killing trajectory" << "\n";
99  return Trajectory();
100  }
101  }
102 
103 
104  if likely(hit.isValid()) {
105  assert( (hit.geographicalId()!=0U) | !hit.canImproveWithTrack() ) ;
106  assert(hit.surface()!=nullptr);
107  //update
108  LogTrace("TrackFitters") << "THE HIT IS VALID: updating hit with predTsos";
109  assert( (!hit.canImproveWithTrack()) | (nullptr!=theHitCloner));
110  assert( (!hit.canImproveWithTrack()) | (nullptr!=dynamic_cast<BaseTrackerRecHit const*>((ihit).get())));
111  auto preciseHit = theHitCloner->makeShared(ihit,predTsos);
112  dump(*preciseHit,hitcounter,"TrackFitters");
113  assert(preciseHit->isValid());
114  assert( (preciseHit->geographicalId()!=0U) | (!preciseHit->canImproveWithTrack()) );
115  assert(preciseHit->surface()!=nullptr);
116 
117  if unlikely(!preciseHit->isValid()){
118  LogTrace("TrackFitters") << "THE Precise HIT IS NOT VALID: using currTsos = predTsos" << "\n";
119  currTsos = predTsos;
120  myTraj.push(TM(predTsos, ihit,0,theGeometry->idToLayer((ihit)->geographicalId()) ));
121 
122  }else{
123  LogTrace("TrackFitters") << "THE Precise HIT IS VALID: updating currTsos" << "\n";
124  currTsos = updator()->update(predTsos, *preciseHit);
125  //check for valid hits with no det (refitter with constraints)
126  bool badState = (!currTsos.isValid())
127  || (hit.geographicalId().det() == DetId::Tracker
128  &&
129  (std::abs(currTsos.localParameters().qbp())>100
130  || std::abs(currTsos.localParameters().position().y()) > 1000
131  || std::abs(currTsos.localParameters().position().x()) > 1000
132  ) ) || edm::isNotFinite(currTsos.localParameters().qbp());
133  if unlikely(badState){
134  if (!currTsos.isValid()) {
135  edm::LogError("FailedUpdate") <<"updating with the hit failed. Not updating the trajectory with the hit";
136 
137  }
138  else if (edm::isNotFinite(currTsos.localParameters().qbp())) {
139  edm::LogError("TrajectoryNaN")<<"Trajectory has NaN";
140 
141  }
142  else{
143  LogTrace("FailedUpdate")<<"updated state is valid but pretty bad, skipping. currTsos " <<currTsos<<"\n predTsos "<<predTsos;
144  }
145  myTraj.push(TM(predTsos, ihit,0,theGeometry->idToLayer((ihit)->geographicalId()) ));
146  //There is a no-fail policy here. So, it's time to give up
147  //Keep the traj with invalid TSOS so that it's clear what happened
148  if( myTraj.foundHits() >= minHits_ ) {
149  LogDebug("TrackFitters") << " breaking trajectory" << "\n";
150  break;
151  } else {
152  LogDebug("TrackFitters") << " killing trajectory" << "\n";
153  return Trajectory();
154  }
155  } else{
156  if (preciseHit->det()){
157  myTraj.push(TM(predTsos, currTsos, preciseHit,
158  estimator()->estimate(predTsos, *preciseHit).second,
159  theGeometry->idToLayer(preciseHit->geographicalId()) ));
160  }
161  else{
162  myTraj.push(TM(predTsos, currTsos, preciseHit,
163  estimator()->estimate(predTsos, *preciseHit).second));
164  }
165  }
166  }
167  } else { // invalid hit
168  dump(hit,hitcounter,"TrackFitters");
169  //no update
170  LogDebug("TrackFitters") << "THE HIT IS NOT VALID: using currTsos" << "\n";
171  currTsos = predTsos;
172  assert( ((ihit)->det()==nullptr) || (ihit)->geographicalId()!=0U);
173  if ((ihit)->det()) myTraj.push(TM(predTsos, ihit,0,theGeometry->idToLayer((ihit)->geographicalId()) ));
174  else myTraj.push(TM(predTsos, ihit,0));
175  }
176  LogTrace("TrackFitters")
177  << "predTsos !" << "\n"
178  << predTsos
179  <<" with local position " << predTsos.localPosition()
180  <<"currTsos !" << "\n"
181  << currTsos
182  <<" with local position " << currTsos.localPosition();
183  }
184 
185  LogDebug("TrackFitters") << "Found 1 trajectory with " << myTraj.foundHits() << " valid hits\n";
186 
187  return ret;
188 }
189 
#define LogDebug(id)
PropagationDirection direction() const
bool empty() const
True if trajectory has no measurements.
Definition: Trajectory.h:254
type
Definition: HCALResponse.h:21
virtual const Surface * surface() const
TrajectorySeed const & seed() const
Access to the seed used to reconstruct the Trajectory.
Definition: Trajectory.h:285
const LocalTrajectoryParameters & localParameters() const
const TrajectoryStateUpdator * updator() const
LocalPoint position() const
Local x and y position coordinates.
const MeasurementEstimator * estimator() const
T y() const
Definition: PV3DBase.h:63
ConstRecHitContainer recHits() const
Definition: Trajectory.h:204
TrajectoryMeasurement TM
TrackingRecHit::ConstRecHitPointer makeShared(TrackingRecHit::ConstRecHitPointer const &hit, TrajectoryStateOnSurface const &tsos) const
Definition: TkCloner.h:23
void reserve(unsigned int n)
Definition: Trajectory.h:150
virtual TrajectoryStateOnSurface update(const TrajectoryStateOnSurface &, const TrackingRecHit &) const =0
Trajectory fitOne(const Trajectory &aTraj, fitType) const override
#define unlikely(x)
U second(std::pair< T, U > const &p)
#define likely(x)
std::unique_ptr< Propagator > SetPropagationDirection(Propagator const &iprop, PropagationDirection dir)
bool isNotFinite(T x)
Definition: isFinite.h:10
static const DetLayerGeometry dummyGeometry
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
#define LogTrace(id)
const DetLayerGeometry * theGeometry
TrajectoryMeasurement const & firstMeasurement() const
Definition: Trajectory.h:187
virtual bool canImproveWithTrack() const
bool isValid() const
const Propagator * thePropagator
TrajectoryStateOnSurface const & updatedState() const
DetId geographicalId() const
virtual const DetLayer * idToLayer(const DetId &detId) const
Detector det() const
get the detector field from this detid
Definition: DetId.h:35
T x() const
Definition: PV3DBase.h:62
const PositionType & position() const
Trajectory::RecHitContainer RecHitContainer
TkCloner const * theHitCloner