Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Track.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file Track.cc
1 /* Copyright 2008-2010, Technische Universitaet Muenchen,
2  Authors: Christian Hoeppner & Sebastian Neubert & Johannes Rauch
3 
4  This file is part of GENFIT.
5 
6  GENFIT is free software: you can redistribute it and/or modify
7  it under the terms of the GNU Lesser General Public License as published
8  by the Free Software Foundation, either version 3 of the License, or
9  (at your option) any later version.
10 
11  GENFIT is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU Lesser General Public License for more details.
15 
16  You should have received a copy of the GNU Lesser General Public License
17  along with GENFIT. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #include "Track.h"
21 
22 #include "Exception.h"
23 #include "IO.h"
24 #include "KalmanFitterInfo.h"
25 #include "KalmanFitStatus.h"
26 #include "PlanarMeasurement.h"
27 #include "AbsMeasurement.h"
28 
29 #include "WireTrackCandHit.h"
30 
31 #include <algorithm>
32 #include <map>
33 
34 #include <TDatabasePDG.h>
35 #include <TMath.h>
36 #include <TBuffer.h>
37 
38 //#include <glog/logging.h>
39 
40 //#define DEBUG
41 
42 
43 namespace genfit {
44 
46  cardinalRep_(0), fitStatuses_(), mcTrackId_(-1), timeSeed_(0), stateSeed_(6), covSeed_(6)
47 {
48  ;
49 }
50 
51 
52 Track::Track(const TrackCand& trackCand, const MeasurementFactory<AbsMeasurement>& factory, AbsTrackRep* rep) :
53  cardinalRep_(0), fitStatuses_(), stateSeed_(6), covSeed_(6)
54 {
55 
56  if (rep != nullptr)
57  addTrackRep(rep);
58 
59  createMeasurements(trackCand, factory);
60 
61  // Copy seed information from candidate
62  timeSeed_ = trackCand.getTimeSeed();
63  stateSeed_ = trackCand.getStateSeed();
64  covSeed_ = trackCand.getCovSeed();
65 
66  mcTrackId_ = trackCand.getMcTrackId();
67 
68  // fill cache
70 
72 }
73 
74 void
76 {
77  // create the measurements using the factory.
78  std::vector <AbsMeasurement*> factoryHits = factory.createMany(trackCand);
79 
80  if (factoryHits.size() != trackCand.getNHits()) {
81  Exception exc("Track::Track ==> factoryHits.size() != trackCand->getNHits()",__LINE__,__FILE__);
82  exc.setFatal();
83  throw exc;
84  }
85 
86  // create TrackPoints
87  for (unsigned int i=0; i<factoryHits.size(); ++i){
88  TrackPoint* tp = new TrackPoint(factoryHits[i], this);
89  tp->setSortingParameter(trackCand.getHit(i)->getSortingParameter());
90  insertPoint(tp);
91  }
92 }
93 
94 
95 Track::Track(AbsTrackRep* trackRep, const TVectorD& stateSeed) :
96  cardinalRep_(0), fitStatuses_(), mcTrackId_(-1), timeSeed_(0), stateSeed_(stateSeed),
97  covSeed_(TMatrixDSym::kUnit, TMatrixDSym(6))
98 {
99  addTrackRep(trackRep);
100 }
101 
102 
103 Track::Track(AbsTrackRep* trackRep, const TVector3& posSeed, const TVector3& momSeed) :
104  cardinalRep_(0), fitStatuses_(), mcTrackId_(-1), timeSeed_(0), stateSeed_(6),
105  covSeed_(TMatrixDSym::kUnit, TMatrixDSym(6))
106 {
107  addTrackRep(trackRep);
108  setStateSeed(posSeed, momSeed);
109 }
110 
111 
112 Track::Track(AbsTrackRep* trackRep, const TVectorD& stateSeed, const TMatrixDSym& covSeed) :
113  cardinalRep_(0), fitStatuses_(), mcTrackId_(-1), timeSeed_(0), stateSeed_(stateSeed),
114  covSeed_(covSeed)
115 {
116  addTrackRep(trackRep);
117 }
118 
119 
121  TObject(rhs),
122  cardinalRep_(rhs.cardinalRep_), mcTrackId_(rhs.mcTrackId_), timeSeed_(rhs.timeSeed_),
123  stateSeed_(rhs.stateSeed_), covSeed_(rhs.covSeed_)
124 {
125  rhs.checkConsistency();
126 
127  std::map<const AbsTrackRep*, AbsTrackRep*> oldRepNewRep;
128 
129  for (std::vector<AbsTrackRep*>::const_iterator it=rhs.trackReps_.begin(); it!=rhs.trackReps_.end(); ++it) {
130  AbsTrackRep* newRep = (*it)->clone();
131  addTrackRep(newRep);
132  oldRepNewRep[(*it)] = newRep;
133  }
134 
135  trackPoints_.reserve(rhs.trackPoints_.size());
136  for (std::vector<TrackPoint*>::const_iterator it=rhs.trackPoints_.begin(); it!=rhs.trackPoints_.end(); ++it) {
137  trackPoints_.push_back(new TrackPoint(**it, oldRepNewRep));
138  trackPoints_.back()->setTrack(this);
139  }
140 
141  for (std::map< const AbsTrackRep*, FitStatus* >::const_iterator it=rhs.fitStatuses_.begin(); it!=rhs.fitStatuses_.end(); ++it) {
142  setFitStatus(it->second->clone(), oldRepNewRep[it->first]);
143  }
144 
146 
148 }
149 
151  swap(other);
152 
153  for (std::vector<TrackPoint*>::const_iterator it=trackPoints_.begin(); it!=trackPoints_.end(); ++it) {
154  trackPoints_.back()->setTrack(this);
155  }
156 
158 
160 
161  return *this;
162 }
163 
164 void Track::swap(Track& other) {
165  std::swap(this->trackReps_, other.trackReps_);
166  std::swap(this->cardinalRep_, other.cardinalRep_);
167  std::swap(this->trackPoints_, other.trackPoints_);
169  std::swap(this->fitStatuses_, other.fitStatuses_);
170  std::swap(this->mcTrackId_, other.mcTrackId_);
171  std::swap(this->timeSeed_, other.timeSeed_);
172  std::swap(this->stateSeed_, other.stateSeed_);
173  std::swap(this->covSeed_, other.covSeed_);
174 
175 }
176 
178  this->Clear();
179 }
180 
181 void Track::Clear(Option_t*)
182 {
183  // This function is needed for TClonesArray embedding.
184  // FIXME: smarter containers or pointers needed ...
185  for (size_t i = 0; i < trackPoints_.size(); ++i)
186  delete trackPoints_[i];
187 
188  trackPoints_.clear();
190 
191  for (std::map< const AbsTrackRep*, FitStatus* >::iterator it = fitStatuses_.begin(); it!= fitStatuses_.end(); ++it)
192  delete it->second;
193  fitStatuses_.clear();
194 
195  for (size_t i = 0; i < trackReps_.size(); ++i)
196  delete trackReps_[i];
197  trackReps_.clear();
198 
199  cardinalRep_ = 0;
200 
201  mcTrackId_ = -1;
202 
203  timeSeed_ = 0;
204  stateSeed_.Zero();
205  covSeed_.Zero();
206 }
207 
208 
209 TrackPoint* Track::getPoint(int id) const {
210  if (id < 0)
211  id += trackPoints_.size();
212 
213  return trackPoints_.at(id);
214 }
215 
216 
218  if (id < 0)
219  id += trackPointsWithMeasurement_.size();
220 
221  return trackPointsWithMeasurement_.at(id);
222 }
223 
224 
226  if (rep == nullptr)
227  rep = getCardinalRep();
228 
229  if (id >= 0) {
230  int i = 0;
231  for (std::vector<TrackPoint*>::const_iterator it = trackPointsWithMeasurement_.begin(); it != trackPointsWithMeasurement_.end(); ++it) {
232  if ((*it)->hasFitterInfo(rep)) {
233  if (id == i)
234  return (*it);
235  ++i;
236  }
237  }
238  } else {
239  // Search backwards.
240  int i = -1;
241  for (std::vector<TrackPoint*>::const_reverse_iterator it = trackPointsWithMeasurement_.rbegin(); it != trackPointsWithMeasurement_.rend(); ++it) {
242  if ((*it)->hasFitterInfo(rep)) {
243  if (id == i)
244  return (*it);
245  --i;
246  }
247  }
248  }
249 
250  // Not found, i.e. abs(id) > number of fitted TrackPoints
251  return 0;
252 }
253 
254 
256  if (rep == nullptr)
257  rep = getCardinalRep();
258 
259  if (id >= 0) {
260  int i = 0;
261  for (std::vector<TrackPoint*>::const_iterator it = trackPoints_.begin(); it != trackPoints_.end(); ++it) {
262  if ((*it)->hasFitterInfo(rep)) {
263  if (id == i)
264  return (*it);
265  ++i;
266  }
267  }
268  } else {
269  // Search backwards.
270  int i = -1;
271  for (std::vector<TrackPoint*>::const_reverse_iterator it = trackPoints_.rbegin(); it != trackPoints_.rend(); ++it) {
272  if ((*it)->hasFitterInfo(rep)) {
273  if (id == i)
274  return (*it);
275  --i;
276  }
277  }
278  }
279 
280  // Not found, i.e. abs(id) > number of fitted TrackPoints
281  return 0;
282 }
283 
284 
285 const MeasuredStateOnPlane& Track::getFittedState(int id, const AbsTrackRep* rep, bool biased) const {
286  if (rep == nullptr)
287  rep = getCardinalRep();
288 
290  if (point == nullptr) {
291  Exception exc("Track::getFittedState ==> no trackPoint with fitterInfo for rep",__LINE__,__FILE__);
292  exc.setFatal();
293  throw exc;
294  }
295  return point->getFitterInfo(rep)->getFittedState(biased);
296 }
297 
298 
299 int Track::getIdForRep(const AbsTrackRep* rep) const
300 {
301  for (size_t i = 0; i < trackReps_.size(); ++i)
302  if (trackReps_[i] == rep)
303  return i;
304 
305  Exception exc("Track::getIdForRep ==> cannot find TrackRep in Track",__LINE__,__FILE__);
306  exc.setFatal();
307  throw exc;
308 }
309 
310 
311 bool Track::hasFitStatus(const AbsTrackRep* rep) const {
312  if (rep == nullptr)
313  rep = getCardinalRep();
314 
315  if (fitStatuses_.find(rep) == fitStatuses_.end())
316  return false;
317 
318  return (fitStatuses_.at(rep) != nullptr);
319 }
320 
321 
322 bool Track::hasKalmanFitStatus(const AbsTrackRep* rep) const {
323  if (rep == nullptr)
324  rep = getCardinalRep();
325 
326  if (fitStatuses_.find(rep) == fitStatuses_.end())
327  return false;
328 
329  return (dynamic_cast<KalmanFitStatus*>(fitStatuses_.at(rep)) != nullptr);
330 }
331 
332 
334  return dynamic_cast<KalmanFitStatus*>(getFitStatus(rep));
335 }
336 
337 
338 void Track::setFitStatus(FitStatus* fitStatus, const AbsTrackRep* rep) {
339  if (fitStatuses_.find(rep) != fitStatuses_.end())
340  delete fitStatuses_.at(rep);
341 
342  fitStatuses_[rep] = fitStatus;
343 }
344 
345 
346 void Track::setStateSeed(const TVector3& pos, const TVector3& mom) {
347  stateSeed_.ResizeTo(6);
348 
349  stateSeed_(0) = pos.X();
350  stateSeed_(1) = pos.Y();
351  stateSeed_(2) = pos.Z();
352 
353  stateSeed_(3) = mom.X();
354  stateSeed_(4) = mom.Y();
355  stateSeed_(5) = mom.Z();
356 }
357 
358 
359 
361 
362  point->setTrack(this);
363 
364  #ifdef DEBUG
365  debugOut << "Track::insertPoint at position " << id << "\n";
366  #endif
367  assert(point!=nullptr);
368  trackHasChanged();
369 
370  point->setTrack(this);
371 
372  if (trackPoints_.size() == 0) {
373  trackPoints_.push_back(point);
374 
375  if (point->hasRawMeasurements())
376  trackPointsWithMeasurement_.push_back(point);
377 
378  return;
379  }
380 
381  if (id == -1 || id == (int)trackPoints_.size()) {
382  trackPoints_.push_back(point);
383 
384  if (point->hasRawMeasurements())
385  trackPointsWithMeasurement_.push_back(point);
386 
387  deleteReferenceInfo(std::max(0, (int)trackPoints_.size()-2), (int)trackPoints_.size()-1);
388 
389  // delete fitter infos if inserted point has a measurement
390  if (point->hasRawMeasurements()) {
391  deleteForwardInfo(-1, -1);
392  deleteBackwardInfo(0, -2);
393  }
394 
395  return;
396  }
397 
398  // [-size, size-1] is the allowed range
399  assert(id < (ssize_t)trackPoints_.size() || -id-1 <= (ssize_t)trackPoints_.size());
400 
401  if (id < 0)
402  id += trackPoints_.size() + 1;
403 
404  // insert
405  trackPoints_.insert(trackPoints_.begin() + id, point); // insert inserts BEFORE
406 
407  // delete fitter infos if inserted point has a measurement
408  if (point->hasRawMeasurements()) {
409  deleteForwardInfo(id, -1);
410  deleteBackwardInfo(0, id);
411  }
412 
413  // delete reference info of neighbouring points
414  deleteReferenceInfo(std::max(0, id-1), std::min((int)trackPoints_.size()-1, id+1));
415 
417 }
418 
419 
420 void Track::insertPoints(std::vector<TrackPoint*> points, int id) {
421 
422  int nBefore = getNumPoints();
423  int n = points.size();
424 
425  if (n == 0)
426  return;
427  if (n == 1) {
428  insertPoint(points[0], id);
429  return;
430  }
431 
432  for (std::vector<TrackPoint*>::iterator p = points.begin(); p != points.end(); ++p)
433  (*p)->setTrack(this);
434 
435  if (id == -1 || id == (int)trackPoints_.size()) {
436  trackPoints_.insert(trackPoints_.end(), points.begin(), points.end());
437 
438  deleteReferenceInfo(std::max(0, nBefore-1), nBefore);
439 
440  deleteForwardInfo(nBefore, -1);
441  deleteBackwardInfo(0, std::max(0, nBefore-1));
442 
444 
445  return;
446  }
447 
448 
449  assert(id < (ssize_t)trackPoints_.size() || -id-1 <= (ssize_t)trackPoints_.size());
450 
451  if (id < 0)
452  id += trackPoints_.size() + 1;
453 
454 
455  // insert
456  trackPoints_.insert(trackPoints_.begin() + id, points.begin(), points.end()); // insert inserts BEFORE
457 
458  // delete fitter infos if inserted point has a measurement
459  deleteForwardInfo(id, -1);
460  deleteBackwardInfo(0, id+n);
461 
462  // delete reference info of neighbouring points
463  deleteReferenceInfo(std::max(0, id-1), std::min((int)trackPoints_.size()-1, id));
464  deleteReferenceInfo(std::max(0, id+n-1), std::min((int)trackPoints_.size()-1, id+n));
465 
467 }
468 
469 
470 void Track::deletePoint(int id) {
471 
472  #ifdef DEBUG
473  debugOut << "Track::deletePoint at position " << id << "\n";
474  #endif
475 
476  trackHasChanged();
477 
478  if (id < 0)
479  id += trackPoints_.size();
480  assert(id>0);
481 
482 
483  // delete forwardInfo after point (backwardInfo before point) if deleted point has a measurement
484  if (trackPoints_[id]->hasRawMeasurements()) {
485  deleteForwardInfo(id, -1);
486  deleteBackwardInfo(0, id-1);
487  }
488 
489  // delete reference info of neighbouring points
490  deleteReferenceInfo(std::max(0, id-1), std::min((int)trackPoints_.size()-1, id+1));
491 
492  // delete point
493  std::vector<TrackPoint*>::iterator it = std::find(trackPointsWithMeasurement_.begin(), trackPointsWithMeasurement_.end(), trackPoints_[id]);
494  if (it != trackPointsWithMeasurement_.end())
495  trackPointsWithMeasurement_.erase(it);
496 
497  delete trackPoints_[id];
498  trackPoints_.erase(trackPoints_.begin()+id);
499 
501 
502 }
503 
504 
505 void Track::insertMeasurement(AbsMeasurement* measurement, int id) {
506  insertPoint(new TrackPoint(measurement, this), id);
507 }
508 
510  if(hasFitStatus(rep)) {
511  delete fitStatuses_.at(rep);
512  fitStatuses_.erase(rep);
513  }
514 
515  // delete FitterInfos related to the deleted TrackRep
516  for (const auto& trackPoint : trackPoints_) {
517  if(trackPoint->hasFitterInfo(rep)) {
518  trackPoint->deleteFitterInfo(rep);
519  }
520  }
521 }
522 
523 
524 void Track::mergeTrack(const Track* other, int id) {
525 
526  #ifdef DEBUG
527  debugOut << "Track::mergeTrack\n";
528  #endif
529 
530  if (other->getNumPoints() == 0)
531  return;
532 
533  std::map<const AbsTrackRep*, AbsTrackRep*> otherRepThisRep;
534  std::vector<const AbsTrackRep*> otherRepsToRemove;
535 
536  for (std::vector<AbsTrackRep*>::const_iterator otherRep=other->trackReps_.begin(); otherRep!=other->trackReps_.end(); ++otherRep) {
537  bool found(false);
538  for (std::vector<AbsTrackRep*>::const_iterator thisRep=trackReps_.begin(); thisRep!=trackReps_.end(); ++thisRep) {
539  if ((*thisRep)->isSame(*otherRep)) {
540  otherRepThisRep[*otherRep] = *thisRep;
541  #ifdef DEBUG
542  debugOut << " map other rep " << *otherRep << " to " << (*thisRep) << "\n";
543  #endif
544  if (found) {
545  Exception exc("Track::mergeTrack ==> more than one matching rep.",__LINE__,__FILE__);
546  exc.setFatal();
547  throw exc;
548  }
549  found = true;
550  break;
551  }
552  }
553  if (!found) {
554  otherRepsToRemove.push_back(*otherRep);
555  #ifdef DEBUG
556  debugOut << " remove other rep " << *otherRep << "\n";
557  #endif
558  }
559  }
560 
561 
562  std::vector<TrackPoint*> points;
563  points.reserve(other->getNumPoints());
564 
565  for (std::vector<TrackPoint*>::const_iterator otherTp=other->trackPoints_.begin(); otherTp!=other->trackPoints_.end(); ++otherTp) {
566  points.push_back(new TrackPoint(**otherTp, otherRepThisRep, &otherRepsToRemove));
567  }
568 
569  insertPoints(points, id);
570 }
571 
572 
574  trackReps_.push_back(trackRep);
575  fitStatuses_[trackRep] = new FitStatus();
576 }
577 
578 
579 void Track::deleteTrackRep(int id) {
580  if (id < 0)
581  id += trackReps_.size();
582 
583  AbsTrackRep* rep = trackReps_.at(id);
584 
585  // update cardinalRep_
586  if (int(cardinalRep_) == id)
587  cardinalRep_ = 0; // reset
588  else if (int(cardinalRep_) > id)
589  --cardinalRep_; // make cardinalRep_ point to the same TrackRep before and after deletion
590 
591  // delete FitterInfos related to the deleted TrackRep
592  for (std::vector<TrackPoint*>::const_iterator pointIt = trackPoints_.begin(); pointIt != trackPoints_.end(); ++pointIt) {
593  (*pointIt)->deleteFitterInfo(rep);
594  }
595 
596  // delete fitStatus
597  delete fitStatuses_.at(rep);
598  fitStatuses_.erase(rep);
599 
600  // delete rep
601  delete rep;
602  trackReps_.erase(trackReps_.begin()+id);
603 }
604 
605 
606 void Track::setCardinalRep(int id) {
607 
608  if (id < 0)
609  id += trackReps_.size();
610 
611  if (id >= 0 && (unsigned int)id < trackReps_.size())
612  cardinalRep_ = id;
613  else {
614  cardinalRep_ = 0;
615  errorOut << "Track::setCardinalRep: Attempted to set cardinalRep_ to a value out of bounds. Resetting cardinalRep_ to 0." << std::endl;
616  }
617 }
618 
619 
621 
622  // Todo: test
623 
624  if (trackReps_.size() <= 1)
625  return;
626 
627  double minChi2(9.E99);
628  const AbsTrackRep* bestRep(nullptr);
629 
630  for (std::map< const AbsTrackRep*, FitStatus* >::const_iterator it = fitStatuses_.begin(); it != fitStatuses_.end(); ++it) {
631  if (it->second->isFitConverged()) {
632  if (it->second->getChi2() < minChi2) {
633  minChi2 = it->second->getChi2();
634  bestRep = it->first;
635  }
636  }
637  }
638 
639  if (bestRep != nullptr) {
640  setCardinalRep(getIdForRep(bestRep));
641  }
642 }
643 
644 
645 bool Track::sort() {
646  #ifdef DEBUG
647  debugOut << "Track::sort \n";
648  #endif
649 
650  int nPoints(trackPoints_.size());
651  // original order
652  std::vector<TrackPoint*> pointsBefore(trackPoints_);
653 
654  // sort
655  std::stable_sort(trackPoints_.begin(), trackPoints_.end(), TrackPointComparator());
656 
657  // see where order changed
658  int equalUntil(-1), equalFrom(nPoints);
659  for (int i = 0; i<nPoints; ++i) {
660  if (pointsBefore[i] == trackPoints_[i])
661  equalUntil = i;
662  else
663  break;
664  }
665 
666  if (equalUntil == nPoints-1)
667  return false; // sorting did not change anything
668 
669 
670  trackHasChanged();
671 
672  for (int i = nPoints-1; i>equalUntil; --i) {
673  if (pointsBefore[i] == trackPoints_[i])
674  equalFrom = i;
675  else
676  break;
677  }
678 
679  #ifdef DEBUG
680  debugOut << "Track::sort. Equal up to (including) hit " << equalUntil << " and from (including) hit " << equalFrom << " \n";
681  #endif
682 
683  deleteForwardInfo(equalUntil+1, -1);
684  deleteBackwardInfo(0, equalFrom-1);
685 
686  deleteReferenceInfo(std::max(0, equalUntil+1), std::min((int)trackPoints_.size()-1, equalFrom-1));
687 
689 
690  return true;
691 }
692 
693 
694 bool Track::udpateSeed(int id, AbsTrackRep* rep, bool biased) {
695  try {
696  const MeasuredStateOnPlane& fittedState = getFittedState(id, rep, biased);
697  setTimeSeed(fittedState.getTime());
698  setStateSeed(fittedState.get6DState());
699  setCovSeed(fittedState.get6DCov());
700 
701  double fittedCharge = fittedState.getCharge();
702 
703  for (unsigned int i = 0; i<trackReps_.size(); ++i) {
704  if (trackReps_[i]->getPDGCharge() * fittedCharge < 0) {
705  trackReps_[i]->switchPDGSign();
706  }
707  }
708  }
709  catch (Exception& e) {
710  // in this case the original track seed will be used
711  return false;
712  }
713  return true;
714 }
715 
716 
718 
719  std::reverse(trackPoints_.begin(),trackPoints_.end());
720 
721  deleteForwardInfo(0, -1);
722  deleteBackwardInfo(0, -1);
723  deleteReferenceInfo(0, -1);
724 
726 }
727 
728 
730  if (rep != nullptr) {
731  rep->switchPDGSign();
732  return;
733  }
734 
735  for (unsigned int i = 0; i<trackReps_.size(); ++i) {
736  trackReps_[i]->switchPDGSign();
737  }
738 }
739 
740 
742  udpateSeed(-1); // set fitted state of last hit as new seed
743  reverseMomSeed(); // flip momentum direction
744  switchPDGSigns();
745  reverseTrackPoints(); // also deletes all fitterInfos
746 }
747 
748 
749 void Track::deleteForwardInfo(int startId, int endId, const AbsTrackRep* rep) {
750  #ifdef DEBUG
751  debugOut << "Track::deleteForwardInfo from position " << startId << " to " << endId << "\n";
752  #endif
753 
754  trackHasChanged();
755 
756  if (startId < 0)
757  startId += trackPoints_.size();
758  if (endId < 0)
759  endId += trackPoints_.size();
760  endId += 1;
761 
762  assert (endId >= startId);
763 
764  for (std::vector<TrackPoint*>::const_iterator pointIt = trackPoints_.begin() + startId; pointIt != trackPoints_.begin() + endId; ++pointIt) {
765  if (rep != nullptr) {
766  if ((*pointIt)->hasFitterInfo(rep))
767  (*pointIt)->getFitterInfo(rep)->deleteForwardInfo();
768  }
769  else {
770  const std::vector<AbsFitterInfo*> fitterInfos = (*pointIt)->getFitterInfos();
771  for (std::vector<AbsFitterInfo*>::const_iterator fitterInfoIt = fitterInfos.begin(); fitterInfoIt != fitterInfos.end(); ++fitterInfoIt) {
772  (*fitterInfoIt)->deleteForwardInfo();
773  }
774  }
775  }
776 }
777 
778 void Track::deleteBackwardInfo(int startId, int endId, const AbsTrackRep* rep) {
779 
780  #ifdef DEBUG
781  debugOut << "Track::deleteBackwardInfo from position " << startId << " to " << endId << "\n";
782  #endif
783 
784  trackHasChanged();
785 
786  if (startId < 0)
787  startId += trackPoints_.size();
788  if (endId < 0)
789  endId += trackPoints_.size();
790  endId += 1;
791 
792  assert (endId >= startId);
793 
794 
795  for (std::vector<TrackPoint*>::const_iterator pointIt = trackPoints_.begin() + startId; pointIt != trackPoints_.begin() + endId; ++pointIt) {
796  if (rep != nullptr) {
797  if ((*pointIt)->hasFitterInfo(rep))
798  (*pointIt)->getFitterInfo(rep)->deleteBackwardInfo();
799  }
800  else {
801  const std::vector<AbsFitterInfo*> fitterInfos = (*pointIt)->getFitterInfos();
802  for (std::vector<AbsFitterInfo*>::const_iterator fitterInfoIt = fitterInfos.begin(); fitterInfoIt != fitterInfos.end(); ++fitterInfoIt) {
803  (*fitterInfoIt)->deleteBackwardInfo();
804  }
805  }
806  }
807 }
808 
809 void Track::deleteReferenceInfo(int startId, int endId, const AbsTrackRep* rep) {
810 
811  #ifdef DEBUG
812  debugOut << "Track::deleteReferenceInfo from position " << startId << " to " << endId << "\n";
813  #endif
814 
815  trackHasChanged();
816 
817  if (startId < 0)
818  startId += trackPoints_.size();
819  if (endId < 0)
820  endId += trackPoints_.size();
821  endId += 1;
822 
823  assert (endId >= startId);
824 
825  for (std::vector<TrackPoint*>::const_iterator pointIt = trackPoints_.begin() + startId; pointIt != trackPoints_.begin() + endId; ++pointIt) {
826  if (rep != nullptr) {
827  if ((*pointIt)->hasFitterInfo(rep))
828  (*pointIt)->getFitterInfo(rep)->deleteReferenceInfo();
829  }
830  else {
831  std::vector<AbsFitterInfo*> fitterInfos = (*pointIt)->getFitterInfos();
832  for (std::vector<AbsFitterInfo*>::const_iterator fitterInfoIt = fitterInfos.begin(); fitterInfoIt != fitterInfos.end(); ++fitterInfoIt) {
833  (*fitterInfoIt)->deleteReferenceInfo();
834  }
835  }
836  }
837 }
838 
839 void Track::deleteMeasurementInfo(int startId, int endId, const AbsTrackRep* rep) {
840 
841  #ifdef DEBUG
842  debugOut << "Track::deleteMeasurementInfo from position " << startId << " to " << endId << "\n";
843  #endif
844 
845  trackHasChanged();
846 
847  if (startId < 0)
848  startId += trackPoints_.size();
849  if (endId < 0)
850  endId += trackPoints_.size();
851  endId += 1;
852 
853  assert (endId >= startId);
854 
855  for (std::vector<TrackPoint*>::const_iterator pointIt = trackPoints_.begin() + startId; pointIt != trackPoints_.begin() + endId; ++pointIt) {
856  if (rep != nullptr) {
857  if ((*pointIt)->hasFitterInfo(rep))
858  (*pointIt)->getFitterInfo(rep)->deleteMeasurementInfo();
859  }
860  else {
861  std::vector<AbsFitterInfo*> fitterInfos = (*pointIt)->getFitterInfos();
862  for (std::vector<AbsFitterInfo*>::const_iterator fitterInfoIt = fitterInfos.begin(); fitterInfoIt != fitterInfos.end(); ++fitterInfoIt) {
863  (*fitterInfoIt)->deleteMeasurementInfo();
864  }
865  }
866  }
867 }
868 
869 void Track::deleteFitterInfo(int startId, int endId, const AbsTrackRep* rep) {
870 
871  #ifdef DEBUG
872  debugOut << "Track::deleteFitterInfo from position " << startId << " to " << endId << "\n";
873  #endif
874 
875  trackHasChanged();
876 
877  if (startId < 0)
878  startId += trackPoints_.size();
879  if (endId < 0)
880  endId += trackPoints_.size();
881  endId += 1;
882 
883  assert (endId >= startId);
884 
885  for (std::vector<TrackPoint*>::const_iterator pointIt = trackPoints_.begin() + startId; pointIt != trackPoints_.begin() + endId; ++pointIt) {
886  if (rep != nullptr) {
887  if ((*pointIt)->hasFitterInfo(rep))
888  (*pointIt)->deleteFitterInfo(rep);
889  }
890  else {
891  for (std::vector<AbsTrackRep*>::const_iterator repIt = trackReps_.begin(); repIt != trackReps_.end(); ++repIt) {
892  if ((*pointIt)->hasFitterInfo(*repIt))
893  (*pointIt)->deleteFitterInfo(*repIt);
894  }
895  }
896  }
897 }
898 
899 
900 double Track::getTrackLen(AbsTrackRep* rep, int startId, int endId) const {
901 
902  if (startId < 0)
903  startId += trackPoints_.size();
904  if (endId < 0)
905  endId += trackPoints_.size();
906 
907  bool backwards(false);
908  if (startId > endId) {
909  double temp = startId;
910  startId = endId;
911  endId = temp;
912  backwards = true;
913  }
914 
915  endId += 1;
916 
917  if (rep == nullptr)
918  rep = getCardinalRep();
919 
920  double trackLen(0);
922 
923  for (std::vector<TrackPoint*>::const_iterator pointIt = trackPoints_.begin() + startId; pointIt != trackPoints_.begin() + endId; ++pointIt) {
924  if (! (*pointIt)->hasFitterInfo(rep)) {
925  Exception e("Track::getTracklength: trackPoint has no fitterInfo", __LINE__,__FILE__);
926  throw e;
927  }
928 
929  if (pointIt != trackPoints_.begin() + startId) {
930  trackLen += rep->extrapolateToPlane(state, (*pointIt)->getFitterInfo(rep)->getPlane());
931  }
932 
933  state = (*pointIt)->getFitterInfo(rep)->getFittedState();
934  }
935 
936  if (backwards)
937  trackLen *= -1.;
938 
939  return trackLen;
940 }
941 
942 
944  TrackCand* cand = new TrackCand();
945 
947  cand->setCovSeed(covSeed_);
948  cand->setMcTrackId(mcTrackId_);
949 
950  for (unsigned int i = 0; i < trackPointsWithMeasurement_.size(); ++i) {
952  const std::vector< AbsMeasurement* >& measurements = tp->getRawMeasurements();
953 
954  for (unsigned int j = 0; j < measurements.size(); ++j) {
955  const AbsMeasurement* m = measurements[j];
956  TrackCandHit* tch;
957 
958  int planeId = -1;
959  if (dynamic_cast<const PlanarMeasurement*>(m)) {
960  planeId = static_cast<const PlanarMeasurement*>(m)->getPlaneId();
961  }
962 
963  if (m->isLeftRightMeasurement()) {
964  tch = new WireTrackCandHit(m->getDetId(), m->getHitId(), planeId,
966  }
967  else {
968  tch = new TrackCandHit(m->getDetId(), m->getHitId(), planeId,
969  tp->getSortingParameter());
970  }
971  cand->addHit(tch);
972  }
973  }
974 
975  return cand;
976 }
977 
978 
979 double Track::getTOF(AbsTrackRep* rep, int startId, int endId) const {
980 
981  if (startId < 0)
982  startId += trackPoints_.size();
983  if (endId < 0)
984  endId += trackPoints_.size();
985 
986  if (startId > endId) {
987  std::swap(startId, endId);
988  }
989 
990  endId += 1;
991 
992  if (rep == nullptr)
993  rep = getCardinalRep();
994 
996 
997  const TrackPoint* startPoint(trackPoints_[startId]);
998  const TrackPoint* endPoint(trackPoints_[endId]);
999 
1000  if (!startPoint->hasFitterInfo(rep)
1001  || !endPoint->hasFitterInfo(rep)) {
1002  Exception e("Track::getTOF: trackPoint has no fitterInfo", __LINE__,__FILE__);
1003  throw e;
1004  }
1005 
1006  double tof = (rep->getTime(endPoint->getFitterInfo(rep)->getFittedState())
1007  - rep->getTime(startPoint->getFitterInfo(rep)->getFittedState()));
1008  return tof;
1009 }
1010 
1011 
1012 void Track::fixWeights(AbsTrackRep* rep, int startId, int endId) {
1013 
1014  if (startId < 0)
1015  startId += trackPoints_.size();
1016  if (endId < 0)
1017  endId += trackPoints_.size();
1018 
1019  assert(startId >= 0);
1020  assert(startId <= endId);
1021  assert(endId <= (int)trackPoints_.size());
1022 
1023  std::vector< AbsFitterInfo* > fis;
1024 
1025  for (std::vector<TrackPoint*>::iterator tp = trackPoints_.begin() + startId; tp != trackPoints_.begin() + endId; ++tp) {
1026  fis.clear();
1027  if (rep == nullptr) {
1028  fis = (*tp)->getFitterInfos();
1029  }
1030  else if ((*tp)->hasFitterInfo(rep)) {
1031  fis.push_back((*tp)->getFitterInfo(rep));
1032  }
1033 
1034  for (std::vector< AbsFitterInfo* >::iterator fi = fis.begin(); fi != fis.end(); ++fi) {
1035  KalmanFitterInfo* kfi = dynamic_cast<KalmanFitterInfo*>(*fi);
1036  if (kfi == nullptr)
1037  continue;
1038 
1039  kfi->fixWeights();
1040  }
1041  }
1042 }
1043 
1044 
1045 void Track::prune(const Option_t* option) {
1046 
1047  PruneFlags f;
1048  f.setFlags(option);
1049 
1050  for (std::map< const AbsTrackRep*, FitStatus* >::const_iterator it=fitStatuses_.begin(); it!=fitStatuses_.end(); ++it) {
1051  it->second->getPruneFlags().setFlags(option);
1052  }
1053 
1054  // prune trackPoints
1055  if (f.hasFlags("F") || f.hasFlags("L")) {
1056  TrackPoint* firstPoint = getPointWithFitterInfo(0);
1057  TrackPoint* lastPoint = getPointWithFitterInfo(-1);
1058  for (unsigned int i = 0; i<trackPoints_.size(); ++i) {
1059  if (trackPoints_[i] == firstPoint && f.hasFlags("F"))
1060  continue;
1061 
1062  if (trackPoints_[i] == lastPoint && f.hasFlags("L"))
1063  continue;
1064 
1065  delete trackPoints_[i];
1066  trackPoints_.erase(trackPoints_.begin()+i);
1067  --i;
1068  }
1069  }
1070 
1071  // prune TrackReps
1072  if (f.hasFlags("C")) {
1073  for (unsigned int i = 0; i < trackReps_.size(); ++i) {
1074  if (i != cardinalRep_) {
1075  deleteTrackRep(i);
1076  --i;
1077  }
1078  }
1079  }
1080 
1081 
1082  // from remaining trackPoints: prune measurementsOnPlane, unneeded fitterInfoStuff
1083  for (unsigned int i = 0; i<trackPoints_.size(); ++i) {
1084  if (f.hasFlags("W"))
1085  trackPoints_[i]->deleteRawMeasurements();
1086 
1087  std::vector< AbsFitterInfo* > fis = trackPoints_[i]->getFitterInfos();
1088  for (unsigned int j = 0; j<fis.size(); ++j) {
1089 
1090  if (i == 0 && f.hasFlags("FLI"))
1091  fis[j]->deleteForwardInfo();
1092  else if (i == trackPoints_.size()-1 && f.hasFlags("FLI"))
1093  fis[j]->deleteBackwardInfo();
1094  else if (f.hasFlags("FI"))
1095  fis[j]->deleteForwardInfo();
1096  else if (f.hasFlags("LI"))
1097  fis[j]->deleteBackwardInfo();
1098 
1099  if (f.hasFlags("U") && dynamic_cast<KalmanFitterInfo*>(fis[j]) != nullptr) {
1100  static_cast<KalmanFitterInfo*>(fis[j])->deletePredictions();
1101  }
1102 
1103  // also delete reference info if points have been removed since it is invalid then!
1104  if (f.hasFlags("R") or f.hasFlags("F") or f.hasFlags("L"))
1105  fis[j]->deleteReferenceInfo();
1106  if (f.hasFlags("M"))
1107  fis[j]->deleteMeasurementInfo();
1108  }
1109  }
1110 
1112 
1113  #ifdef DEBUG
1114  debugOut << "pruned Track: "; Print();
1115  #endif
1116 
1117 }
1118 
1119 
1120 void Track::Print(const Option_t* option) const {
1121  TString opt = option;
1122  opt.ToUpper();
1123  if (opt.Contains("C")) { // compact
1124 
1125  printOut << "\n ";
1126  for (unsigned int i=0; i<trackPoints_.size(); ++i) {
1127 
1128  int color = 32*(size_t)(trackPoints_[i]) % 15;
1129  switch (color) {
1130  case 0:
1131  printOut<<"\033[1;30m";
1132  break;
1133  case 1:
1134  printOut<<"\033[0;34m";
1135  break;
1136  case 2:
1137  printOut<<"\033[1;34m";
1138  break;
1139  case 3:
1140  printOut<<"\033[0;32m";
1141  break;
1142  case 4:
1143  printOut<<"\033[1;32m";
1144  break;
1145  case 5:
1146  printOut<<"\033[0;36m";
1147  break;
1148  case 6:
1149  printOut<<"\033[1;36m";
1150  break;
1151  case 7:
1152  printOut<<"\033[0;31m";
1153  break;
1154  case 8:
1155  printOut<<"\033[1;31m";
1156  break;
1157  case 9:
1158  printOut<<"\033[0;35m";
1159  break;
1160  case 10:
1161  printOut<<"\033[1;35m";
1162  break;
1163  case 11:
1164  printOut<<"\033[0;33m";
1165  break;
1166  case 12:
1167  printOut<<"\033[1;33m";
1168  break;
1169  case 13:
1170  printOut<<"\033[0;37m";
1171  break;
1172  default:
1173  ;
1174  }
1175  printOut << trackPoints_[i] << "\033[00m ";
1176  }
1177  printOut << "\n";
1178 
1179  printOut << " ";
1180  for (unsigned int i=0; i<trackPoints_.size(); ++i) {
1181  printf("% -9.3g ", trackPoints_[i]->getSortingParameter());
1182  }
1183 
1184  for (std::vector<AbsTrackRep*>::const_iterator rep = trackReps_.begin(); rep != trackReps_.end(); ++rep) {
1185  printOut << "\n" << getIdForRep(*rep) << " ";
1186  for (unsigned int i=0; i<trackPoints_.size(); ++i) {
1187  if (! trackPoints_[i]->hasFitterInfo(*rep)) {
1188  printOut << " ";
1189  continue;
1190  }
1191  AbsFitterInfo* fi = trackPoints_[i]->getFitterInfo(*rep);
1192  if (fi->hasMeasurements())
1193  printOut << "M";
1194  else
1195  printOut << " ";
1196 
1197  if (fi->hasReferenceState())
1198  printOut << "R";
1199  else
1200  printOut << " ";
1201 
1202  printOut << " ";
1203  }
1204  printOut << "\n";
1205 
1206  printOut << " -> ";
1207  for (unsigned int i=0; i<trackPoints_.size(); ++i) {
1208  if (! trackPoints_[i]->hasFitterInfo(*rep)) {
1209  printOut << " ";
1210  continue;
1211  }
1212  AbsFitterInfo* fi = trackPoints_[i]->getFitterInfo(*rep);
1213  if (fi->hasForwardPrediction())
1214  printOut << "P";
1215  else
1216  printOut << " ";
1217 
1218  if (fi->hasForwardUpdate())
1219  printOut << "U";
1220  else
1221  printOut << " ";
1222 
1223  printOut << " ";
1224  }
1225  printOut << "\n";
1226 
1227  printOut << " <- ";
1228  for (unsigned int i=0; i<trackPoints_.size(); ++i) {
1229  if (! trackPoints_[i]->hasFitterInfo(*rep)) {
1230  printOut << " ";
1231  continue;
1232  }
1233  AbsFitterInfo* fi = trackPoints_[i]->getFitterInfo(*rep);
1234  if (fi->hasBackwardPrediction())
1235  printOut << "P";
1236  else
1237  printOut << " ";
1238 
1239  if (fi->hasBackwardUpdate())
1240  printOut << "U";
1241  else
1242  printOut << " ";
1243 
1244  printOut << " ";
1245  }
1246 
1247  printOut << "\n";
1248 
1249  } //end loop over reps
1250 
1251  printOut << "\n";
1252  return;
1253  }
1254 
1255 
1256 
1257  printOut << "=======================================================================================\n";
1258  printOut << "genfit::Track, containing " << trackPoints_.size() << " TrackPoints and " << trackReps_.size() << " TrackReps.\n";
1259  printOut << " Seed state: "; stateSeed_.Print();
1260 
1261  for (unsigned int i=0; i<trackReps_.size(); ++i) {
1262  printOut << " TrackRep Nr. " << i;
1263  if (i == cardinalRep_)
1264  printOut << " (This is the cardinal rep)";
1265  printOut << "\n";
1266  trackReps_[i]->Print();
1267  }
1268 
1269  printOut << "---------------------------------------------------------------------------------------\n";
1270 
1271  for (unsigned int i=0; i<trackPoints_.size(); ++i) {
1272  printOut << "TrackPoint Nr. " << i << "\n";
1273  trackPoints_[i]->Print();
1274  printOut << "..........................................................................\n";
1275  }
1276 
1277  for (std::map< const AbsTrackRep*, FitStatus* >::const_iterator it=fitStatuses_.begin(); it!=fitStatuses_.end(); ++it) {
1278  it->second->Print();
1279  }
1280 
1281  printOut << "=======================================================================================\n";
1282 
1283 }
1284 
1285 
1287 
1288  // cppcheck-suppress unreadVariable
1289  bool consistent = true;
1290  std::stringstream failures;
1291 
1292  std::map<const AbsTrackRep*, const KalmanFitterInfo*> prevFis;
1293 
1294  // check if seed is 6D
1295  if (stateSeed_.GetNrows() != 6) {
1296  failures << "Track::checkConsistency(): stateSeed_ dimension != 6" << std::endl;
1297  // cppcheck-suppress unreadVariable
1298  consistent = false;
1299  }
1300 
1301  if (covSeed_.GetNrows() != 6) {
1302  failures << "Track::checkConsistency(): covSeed_ dimension != 6" << std::endl;
1303  // cppcheck-suppress unreadVariable
1304  consistent = false;
1305  }
1306 
1307  if (covSeed_.Max() == 0.) {
1308  // Nota bene: The consistency is not set to false when this occurs, because it does not break the consistency of
1309  // the track. However, when something else fails we keep this as additional error information.
1310  failures << "Track::checkConsistency(): Warning: covSeed_ is zero" << std::endl;
1311  }
1312 
1313  // check if correct number of fitStatuses
1314  if (fitStatuses_.size() != trackReps_.size()) {
1315  failures << "Track::checkConsistency(): Number of fitStatuses is != number of TrackReps " << std::endl;
1316  // cppcheck-suppress unreadVariable
1317  consistent = false;
1318  }
1319 
1320  // check if cardinalRep_ is in range of trackReps_
1321  if (trackReps_.size() && cardinalRep_ >= trackReps_.size()) {
1322  failures << "Track::checkConsistency(): cardinalRep id " << cardinalRep_ << " out of bounds" << std::endl;
1323  // cppcheck-suppress unreadVariable
1324  consistent = false;
1325  }
1326 
1327  for (std::vector<AbsTrackRep*>::const_iterator rep = trackReps_.begin(); rep != trackReps_.end(); ++rep) {
1328  // check for nullptr
1329  if ((*rep) == nullptr) {
1330  failures << "Track::checkConsistency(): TrackRep is nullptr" << std::endl;
1331  // cppcheck-suppress unreadVariable
1332  consistent = false;
1333  }
1334 
1335  // check for valid pdg code
1336  TParticlePDG* particle = TDatabasePDG::Instance()->GetParticle((*rep)->getPDG());
1337  if (particle == nullptr) {
1338  failures << "Track::checkConsistency(): TrackRep pdg ID " << (*rep)->getPDG() << " is not valid" << std::endl;
1339  // cppcheck-suppress unreadVariable
1340  consistent = false;
1341  }
1342 
1343  // check if corresponding FitStatus is there
1344  if (fitStatuses_.find(*rep) == fitStatuses_.end() and fitStatuses_.find(*rep)->second != nullptr) {
1345  failures << "Track::checkConsistency(): No FitStatus for Rep or FitStatus is nullptr" << std::endl;
1346  // cppcheck-suppress unreadVariable
1347  consistent = false;
1348  }
1349  }
1350 
1351  // check TrackPoints
1352  for (std::vector<TrackPoint*>::const_iterator tp = trackPoints_.begin(); tp != trackPoints_.end(); ++tp) {
1353  // check for nullptr
1354  if ((*tp) == nullptr) {
1355  failures << "Track::checkConsistency(): TrackPoint is nullptr" << std::endl;
1356  // cppcheck-suppress unreadVariable
1357  consistent = false;
1358  }
1359  // check if trackPoint points back to this track
1360  if ((*tp)->getTrack() != this) {
1361  failures << "Track::checkConsistency(): TrackPoint does not point back to this track" << std::endl;
1362  // cppcheck-suppress unreadVariable
1363  consistent = false;
1364  }
1365 
1366  // check rawMeasurements
1367  const std::vector<AbsMeasurement*>& rawMeasurements = (*tp)->getRawMeasurements();
1368  for (std::vector<AbsMeasurement*>::const_iterator m = rawMeasurements.begin(); m != rawMeasurements.end(); ++m) {
1369  // check for nullptr
1370  if ((*m) == nullptr) {
1371  failures << "Track::checkConsistency(): Measurement is nullptr" << std::endl;
1372  // cppcheck-suppress unreadVariable
1373  consistent = false;
1374  }
1375  // check if measurement points back to TrackPoint
1376  if ((*m)->getTrackPoint() != *tp) {
1377  failures << "Track::checkConsistency(): Measurement does not point back to correct TrackPoint" << std::endl;
1378  // cppcheck-suppress unreadVariable
1379  consistent = false;
1380  }
1381  }
1382 
1383  // check fitterInfos
1384  std::vector<AbsFitterInfo*> fitterInfos = (*tp)->getFitterInfos();
1385  for (std::vector<AbsFitterInfo*>::const_iterator fi = fitterInfos.begin(); fi != fitterInfos.end(); ++fi) {
1386  // check for nullptr
1387  if ((*fi) == nullptr) {
1388  failures << "Track::checkConsistency(): FitterInfo is nullptr. TrackPoint: " << *tp << std::endl;
1389  // cppcheck-suppress unreadVariable
1390  consistent = false;
1391  }
1392 
1393  // check if fitterInfos point to valid TrackReps in trackReps_
1394  int mycount (0);
1395  for (std::vector<AbsTrackRep*>::const_iterator rep = trackReps_.begin(); rep != trackReps_.end(); ++rep) {
1396  if ( (*rep) == (*fi)->getRep() ) {
1397  ++mycount;
1398  }
1399  }
1400  if (mycount == 0) {
1401  failures << "Track::checkConsistency(): fitterInfo points to TrackRep which is not in Track" << std::endl;
1402  // cppcheck-suppress unreadVariable
1403  consistent = false;
1404  }
1405 
1406  if (!( (*fi)->checkConsistency(&(this->getFitStatus((*fi)->getRep())->getPruneFlags())) ) ) {
1407  failures << "Track::checkConsistency(): FitterInfo not consistent. TrackPoint: " << *tp << std::endl;
1408  // cppcheck-suppress unreadVariable
1409  consistent = false;
1410  }
1411 
1412  if (dynamic_cast<KalmanFitterInfo*>(*fi) != nullptr) {
1413  if (prevFis[(*fi)->getRep()] != nullptr &&
1414  static_cast<KalmanFitterInfo*>(*fi)->hasReferenceState() &&
1415  prevFis[(*fi)->getRep()]->hasReferenceState() ) {
1416  double len = static_cast<KalmanFitterInfo*>(*fi)->getReferenceState()->getForwardSegmentLength();
1417  double prevLen = prevFis[(*fi)->getRep()]->getReferenceState()->getBackwardSegmentLength();
1418  if (fabs(prevLen + len) > 1E-10 ) {
1419  failures << "Track::checkConsistency(): segment lengths of reference states for rep " << (*fi)->getRep() << " (id " << getIdForRep((*fi)->getRep()) << ") at TrackPoint " << (*tp) << " don't match" << std::endl;
1420  failures << prevLen << " + " << len << " = " << prevLen + len << std::endl;
1421  failures << "TrackPoint " << *tp << ", FitterInfo " << *fi << ", rep " << getIdForRep((*fi)->getRep()) << std::endl;
1422  // cppcheck-suppress unreadVariable
1423  consistent = false;
1424  }
1425  }
1426 
1427  prevFis[(*fi)->getRep()] = static_cast<KalmanFitterInfo*>(*fi);
1428  }
1429  else
1430  prevFis[(*fi)->getRep()] = nullptr;
1431 
1432  } // end loop over FitterInfos
1433 
1434  } // end loop over TrackPoints
1435 
1436 
1437  // check trackPointsWithMeasurement_
1438  std::vector<TrackPoint*> trackPointsWithMeasurement;
1439 
1440  for (std::vector<TrackPoint*>::const_iterator it = trackPoints_.begin(); it != trackPoints_.end(); ++it) {
1441  if ((*it)->hasRawMeasurements()) {
1442  trackPointsWithMeasurement.push_back(*it);
1443  }
1444  }
1445 
1446  if (trackPointsWithMeasurement.size() != trackPointsWithMeasurement_.size()) {
1447  failures << "Track::checkConsistency(): trackPointsWithMeasurement_ has incorrect size" << std::endl;
1448  // cppcheck-suppress unreadVariable
1449  consistent = false;
1450  }
1451 
1452  for (unsigned int i = 0; i < trackPointsWithMeasurement.size(); ++i) {
1453  if (trackPointsWithMeasurement[i] != trackPointsWithMeasurement_[i]) {
1454  failures << "Track::checkConsistency(): trackPointsWithMeasurement_ is not correct" << std::endl;
1455  failures << "has id " << i << ", address " << trackPointsWithMeasurement_[i] << std::endl;
1456  failures << "should have id " << i << ", address " << trackPointsWithMeasurement[i] << std::endl;
1457  // cppcheck-suppress unreadVariable
1458  consistent = false;
1459  }
1460  }
1461 
1462  if (not consistent) {
1463  throw genfit::Exception(failures.str(), __LINE__, __FILE__);
1464  }
1465 }
1466 
1467 
1469 
1470  #ifdef DEBUG
1471  debugOut << "Track::trackHasChanged \n";
1472  #endif
1473 
1474  if (fitStatuses_.empty())
1475  return;
1476 
1477  for (std::map< const AbsTrackRep*, FitStatus* >::const_iterator it=fitStatuses_.begin(); it!=fitStatuses_.end(); ++it) {
1478  it->second->setHasTrackChanged();
1479  }
1480 }
1481 
1482 
1485  trackPointsWithMeasurement_.reserve(trackPoints_.size());
1486 
1487  for (std::vector<TrackPoint*>::const_iterator it = trackPoints_.begin(); it != trackPoints_.end(); ++it) {
1488  if ((*it)->hasRawMeasurements()) {
1489  trackPointsWithMeasurement_.push_back(*it);
1490  }
1491  }
1492 }
1493 
1494 
1495 void Track::Streamer(TBuffer &R__b)
1496 {
1497  // Stream an object of class genfit::Track.
1498  const bool streamTrackPoints = true; // debugging option
1499  //This works around a msvc bug and should be harmless on other platforms
1500  typedef ::genfit::Track thisClass;
1501  UInt_t R__s, R__c;
1502  if (R__b.IsReading()) {
1503  this->Clear();
1504  Version_t R__v = R__b.ReadVersion(&R__s, &R__c); if (R__v) { }
1505  TObject::Streamer(R__b);
1506  {
1507  std::vector<AbsTrackRep*> &R__stl = trackReps_;
1508  R__stl.clear();
1509  TClass *R__tcl1 = TBuffer::GetClass(typeid(genfit::AbsTrackRep));
1510  if (R__tcl1==0) {
1511  Error("trackReps_ streamer","Missing the TClass object for genfit::AbsTrackRep!");
1512  return;
1513  }
1514  int R__i, R__n;
1515  R__b >> R__n;
1516  R__stl.reserve(R__n);
1517  for (R__i = 0; R__i < R__n; R__i++) {
1518  genfit::AbsTrackRep* R__t;
1519  R__b >> R__t;
1520  R__stl.push_back(R__t);
1521  }
1522  }
1523  R__b >> cardinalRep_;
1524  if (streamTrackPoints)
1525  {
1526  std::vector<TrackPoint*> &R__stl = trackPoints_;
1527  R__stl.clear();
1528  TClass *R__tcl1 = TBuffer::GetClass(typeid(genfit::TrackPoint));
1529  if (R__tcl1==0) {
1530  Error("trackPoints_ streamer","Missing the TClass object for genfit::TrackPoint!");
1531  return;
1532  }
1533  int R__i, R__n;
1534  R__b >> R__n;
1535  R__stl.reserve(R__n);
1536  for (R__i = 0; R__i < R__n; R__i++) {
1537  genfit::TrackPoint* R__t;
1538  R__t = (genfit::TrackPoint*)R__b.ReadObjectAny(R__tcl1);
1539  R__t->setTrack(this);
1540  R__t->fixupRepsForReading();
1541  R__stl.push_back(R__t);
1542  }
1543  }
1544  {
1545  std::map<const AbsTrackRep*,FitStatus*> &R__stl = fitStatuses_;
1546  R__stl.clear();
1547  TClass *R__tcl1 = TBuffer::GetClass(typeid(genfit::AbsTrackRep));
1548  if (R__tcl1==0) {
1549  Error("fitStatuses_ streamer","Missing the TClass object for genfit::AbsTrackRep!");
1550  return;
1551  }
1552  TClass *R__tcl2 = TBuffer::GetClass(typeid(genfit::FitStatus));
1553  if (R__tcl2==0) {
1554  Error("fitStatuses_ streamer","Missing the TClass object for genfit::FitStatus!");
1555  return;
1556  }
1557  int R__i, R__n;
1558  R__b >> R__n;
1559  for (R__i = 0; R__i < R__n; R__i++) {
1560  int id;
1561  R__b >> id;
1562  genfit::FitStatus* R__t2;
1563  R__t2 = (genfit::FitStatus*)R__b.ReadObjectAny(R__tcl2);
1564 
1565  R__stl[this->getTrackRep(id)] = R__t2;
1566  }
1567  }
1568  // timeSeed_ was introduced in version 3. If reading an earlier
1569  // version, default to zero.
1570  if (R__v >= 3) {
1571  R__b >> timeSeed_;
1572  } else {
1573  timeSeed_ = 0;
1574  }
1575  stateSeed_.Streamer(R__b);
1576  covSeed_.Streamer(R__b);
1577  R__b.CheckByteCount(R__s, R__c, thisClass::IsA());
1578 
1580  } else {
1581  R__c = R__b.WriteVersion(thisClass::IsA(), kTRUE);
1582  TObject::Streamer(R__b);
1583  {
1584  std::vector<AbsTrackRep*> &R__stl = trackReps_;
1585  int R__n=int(R__stl.size());
1586  R__b << R__n;
1587  if(R__n) {
1588  TClass *R__tcl1 = TBuffer::GetClass(typeid(genfit::AbsTrackRep));
1589  if (R__tcl1==0) {
1590  Error("trackReps_ streamer","Missing the TClass object for genfit::AbsTrackRep!");
1591  return;
1592  }
1593  std::vector<AbsTrackRep*>::iterator R__k;
1594  for (R__k = R__stl.begin(); R__k != R__stl.end(); ++R__k) {
1595  R__b << *R__k;
1596  }
1597  }
1598  }
1599  R__b << cardinalRep_;
1600  if (streamTrackPoints)
1601  {
1602  std::vector<TrackPoint*> &R__stl = trackPoints_;
1603  int R__n=int(R__stl.size());
1604  R__b << R__n;
1605  if(R__n) {
1606  std::vector<TrackPoint*>::iterator R__k;
1607  for (R__k = R__stl.begin(); R__k != R__stl.end(); ++R__k) {
1608  R__b << (*R__k);
1609  }
1610  }
1611  }
1612  {
1613  std::map<const AbsTrackRep*,FitStatus*> &R__stl = fitStatuses_;
1614  int R__n=int(R__stl.size());
1615  R__b << R__n;
1616  if(R__n) {
1617  TClass *R__tcl1 = TBuffer::GetClass(typeid(genfit::AbsTrackRep));
1618  if (R__tcl1==0) {
1619  Error("fitStatuses_ streamer","Missing the TClass object for genfit::AbsTrackRep!");
1620  return;
1621  }
1622  TClass *R__tcl2 = TBuffer::GetClass(typeid(genfit::FitStatus));
1623  if (R__tcl2==0) {
1624  Error("fitStatuses_ streamer","Missing the TClass object for genfit::FitStatus!");
1625  return;
1626  }
1627  std::map<const AbsTrackRep*,FitStatus*>::iterator R__k;
1628  for (R__k = R__stl.begin(); R__k != R__stl.end(); ++R__k) {
1629  int id = this->getIdForRep((*R__k).first);
1630  R__b << id;
1631  R__b << ((*R__k).second);
1632  }
1633  }
1634  }
1635  R__b << timeSeed_;
1636  stateSeed_.Streamer(R__b);
1637  covSeed_.Streamer(R__b);
1638  R__b.SetByteCount(R__c, kTRUE);
1639  }
1640 }
1641 
1643  for (size_t i = 0; i < trackPoints_.size(); ++i)
1644  delete trackPoints_[i];
1645 
1646  trackPoints_.clear();
1648 
1649  for (std::map< const AbsTrackRep*, FitStatus* >::iterator it = fitStatuses_.begin(); it!= fitStatuses_.end(); ++it)
1650  delete it->second;
1651  fitStatuses_.clear();
1652 }
1653 
1654 } /* End of namespace genfit */