fitters/src/KalmanFitterInfo.cc

0001 /* Copyright 2008-2013, Technische Universitaet Muenchen, Ludwig-Maximilians-Universität München
0002    Authors: Christian Hoeppner & Sebastian Neubert & Johannes Rauch & Tobias Schlüter
0003
0004    This file is part of GENFIT.
0005
0006    GENFIT is free software: you can redistribute it and/or modify
0007    it under the terms of the GNU Lesser General Public License as published
0008    by the Free Software Foundation, either version 3 of the License, or
0009    (at your option) any later version.
0010
0011    GENFIT is distributed in the hope that it will be useful,
0012    but WITHOUT ANY WARRANTY; without even the implied warranty of
0013    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0014    GNU Lesser General Public License for more details.
0015
0016    You should have received a copy of the GNU Lesser General Public License
0017    along with GENFIT.  If not, see <http://www.gnu.org/licenses/>.
0018 */
0019
0020 #include "KalmanFitterInfo.h"
0021
0022 #include <cassert>
0023 #include <TBuffer.h>
0024
0025 #include "IO.h"
0026 #include "Exception.h"
0027 #include "FitStatus.h"
0028 #include "Tools.h"
0029 #include "Track.h"
0030 #include "TrackPoint.h"
0031
0032 //#define DEBUG
0033
0034
0035 namespace genfit {
0036
0037 KalmanFitterInfo::KalmanFitterInfo() :
0038   AbsFitterInfo(), fixWeights_(false)
0039 {
0040   ;
0041 }
0042
0043 KalmanFitterInfo::KalmanFitterInfo(const TrackPoint* trackPoint, const AbsTrackRep* rep)  :
0044   AbsFitterInfo(trackPoint, rep), fixWeights_(false)
0045 {
0046   ;
0047 }
0048
0049 KalmanFitterInfo::~KalmanFitterInfo() {
0050   deleteMeasurementInfo();
0051 }
0052
0053
0054 KalmanFitterInfo* KalmanFitterInfo::clone() const {
0055   KalmanFitterInfo* retVal = new KalmanFitterInfo(this->getTrackPoint(), this->getRep());
0056   if (hasReferenceState())
0057     retVal->setReferenceState(new ReferenceStateOnPlane(*getReferenceState()));
0058   if (hasForwardPrediction())
0059     retVal->setForwardPrediction(new MeasuredStateOnPlane(*getForwardPrediction()));
0060   if (hasForwardUpdate())
0061     retVal->setForwardUpdate(new KalmanFittedStateOnPlane(*getForwardUpdate()));
0062   if (hasBackwardPrediction())
0063     retVal->setBackwardPrediction(new MeasuredStateOnPlane(*getBackwardPrediction()));
0064   if (hasBackwardUpdate())
0065     retVal->setBackwardUpdate(new KalmanFittedStateOnPlane(*getBackwardUpdate()));
0066
0067   retVal->measurementsOnPlane_.reserve(getNumMeasurements());
0068   for (std::vector<MeasurementOnPlane*>::const_iterator it = this->measurementsOnPlane_.begin(); it != this->measurementsOnPlane_.end(); ++it) {
0069     retVal->addMeasurementOnPlane(new MeasurementOnPlane(**it));
0070   }
0071
0072   retVal->fixWeights_ = this->fixWeights_;
0073
0074   return retVal;
0075 }
0076
0077
0078 MeasurementOnPlane KalmanFitterInfo::getAvgWeightedMeasurementOnPlane(bool ignoreWeights) const {
0079
0080   MeasurementOnPlane retVal(*(measurementsOnPlane_[0]));
0081
0082   if(measurementsOnPlane_.size() == 1) {
0083     if (ignoreWeights) {
0084       retVal.setWeight(1.);
0085     }
0086     else {
0087       double weight = (measurementsOnPlane_[0])->getWeight();
0088       if (weight != 1.) {
0089         retVal.getCov() *= 1. / weight;
0090       }
0091       retVal.setWeight(weight);
0092     }
0093     return retVal;
0094   }
0095
0096   // more than one hit
0097
0098   // cppcheck-suppress unreadVariable
0099   double sumOfWeights(0), weight(0);
0100
0101   retVal.getState().Zero();
0102   retVal.getCov().Zero();
0103
0104   TMatrixDSym covInv;
0105
0106   for(unsigned int i=0; i<measurementsOnPlane_.size(); ++i) {
0107
0108     if (i>0) {
0109       // make sure we have compatible measurement types
0110       // TODO: replace with Exceptions!
0111       assert(measurementsOnPlane_[i]->getPlane() == measurementsOnPlane_[0]->getPlane());
0112       assert(*(measurementsOnPlane_[i]->getHMatrix()) == *(measurementsOnPlane_[0]->getHMatrix()));
0113     }
0114
0115     tools::invertMatrix(measurementsOnPlane_[i]->getCov(), covInv); // invert cov
0116     if (ignoreWeights) {
0117       sumOfWeights += 1.;
0118     }
0119     else {
0120       weight = measurementsOnPlane_[i]->getWeight();
0121       sumOfWeights += weight;
0122       covInv *= weight; // weigh cov
0123     }
0124     retVal.getCov() += covInv; // covInv is already inverted and weighted
0125
0126     retVal.getState() += covInv * measurementsOnPlane_[i]->getState();
0127   }
0128
0129   // invert Cov
0130   tools::invertMatrix(retVal.getCov());
0131
0132   retVal.getState() *= retVal.getCov();
0133
0134   retVal.setWeight(sumOfWeights);
0135
0136   return retVal;
0137 }
0138
0139
0140 MeasurementOnPlane* KalmanFitterInfo::getClosestMeasurementOnPlane(const StateOnPlane* sop) const {
0141   if(measurementsOnPlane_.size() == 0)
0142     return nullptr;
0143
0144   if(measurementsOnPlane_.size() == 1)
0145     return getMeasurementOnPlane(0);
0146
0147   double normMin(9.99E99);
0148   unsigned int iMin(0);
0149   const AbsHMatrix* H = measurementsOnPlane_[0]->getHMatrix();
0150   for (unsigned int i=0; i<getNumMeasurements(); ++i) {
0151     if (*(measurementsOnPlane_[i]->getHMatrix()) != *H){
0152       Exception e("KalmanFitterInfo::getClosestMeasurementOnPlane: Cannot compare measurements with different H-Matrices. Maybe you have to select a different multipleMeasurementHandling.", __LINE__,__FILE__);
0153       e.setFatal();
0154       throw e;
0155     }
0156
0157     TVectorD res = measurementsOnPlane_[i]->getState() - H->Hv(sop->getState());
0158     double norm = sqrt(res.Norm2Sqr());
0159     if (norm < normMin) {
0160       normMin = norm;
0161       iMin = i;
0162     }
0163   }
0164
0165   return getMeasurementOnPlane(iMin);
0166 }
0167
0168
0169 std::vector<double> KalmanFitterInfo::getWeights() const {
0170   std::vector<double> retVal(getNumMeasurements());
0171
0172   for (unsigned int i=0; i<getNumMeasurements(); ++i) {
0173     retVal[i] = getMeasurementOnPlane(i)->getWeight();
0174   }
0175
0176   return retVal;
0177 }
0178
0179
0180 const MeasuredStateOnPlane& KalmanFitterInfo::getFittedState(bool biased) const {
0181
0182   // check if we can use cached states
0183   if (biased && fittedStateBiased_)
0184     return *fittedStateBiased_;
0185   if (!biased && fittedStateUnbiased_)
0186     return *fittedStateUnbiased_;
0187
0188
0189   const TrackPoint* tp = this->getTrackPoint();
0190   const Track* tr = tp->getTrack();
0191   const AbsTrackRep* rep =  this->getRep();
0192
0193   bool first(false), last(false);
0194   PruneFlags& flag = tr->getFitStatus(rep)->getPruneFlags();
0195   // if Track is pruned so that only one TrackPoint remains, see if it was the first or last one
0196   #ifdef DEBUG
0197   if (flag.isPruned()) {
0198     debugOut << "KalmanFitterInfo::getFittedState - Track is pruned and has " << tr->getNumPoints() << " TrackPoints \n";
0199     flag.Print();
0200   }
0201   #endif
0202   if (flag.isPruned() && tr->getNumPoints() == 1) {
0203     if (flag.hasFlags("F")) {
0204       first = true;
0205       #ifdef DEBUG
0206       debugOut << "KalmanFitterInfo::getFittedState - has flag F \n";
0207       #endif
0208     }
0209     else if (flag.hasFlags("L")) {
0210       last = true;
0211       #ifdef DEBUG
0212       debugOut << "KalmanFitterInfo::getFittedState - has flag L \n";
0213       #endif
0214     }
0215   }
0216   else { // otherwise check against TrackPoint order
0217     first = tr->getPointWithFitterInfo(0, rep) == tp;
0218     last = tr->getPointWithFitterInfo(-1, rep) == tp;
0219   }
0220
0221   #ifdef DEBUG
0222   debugOut << "KalmanFitterInfo::getFittedState first " << first << ", last " << last << "\n";
0223   debugOut << "KalmanFitterInfo::getFittedState forwardPrediction_ " << forwardPrediction_.get() << ", forwardUpdate_ " << forwardUpdate_.get() << "\n";
0224   debugOut << "KalmanFitterInfo::getFittedState backwardPrediction_ " << backwardPrediction_.get() << ", backwardUpdate_ " << backwardUpdate_.get() << "\n";
0225   #endif
0226
0227   /*
0228   if both needed forward prediction/update and backward prediction are available,
0229   use them to calculate smoothed state.
0230   Otherwise, if one is missing (i.e. has been pruned) and we are at first or last hit,
0231   use only backward or forward prediction (unbiased) of update (biased).
0232   */
0233
0234   if (biased) {
0235     // last measurement and no backward prediction -> use forward update
0236     if (last && !backwardPrediction_) {
0237       if(!forwardUpdate_) {
0238         Exception e("KalmanFitterInfo::getFittedState: Needed updates/predictions not available in this FitterInfo.", __LINE__,__FILE__);
0239         e.setFatal();
0240         throw e;
0241       }
0242       #ifdef DEBUG
0243       debugOut << "KalmanFitterInfo::getFittedState - biased at last measurement = forwardUpdate_ \n";
0244       #endif
0245       return *forwardUpdate_;
0246     }
0247
0248     // first measurement and no forward update -> use backward update
0249     if (first && (!forwardUpdate_ || (backwardUpdate_ && !forwardPrediction_) ) ) {
0250       if(!backwardUpdate_.get()) {
0251         Exception e("KalmanFitterInfo::getFittedState: Needed updates/predictions not available in this FitterInfo.", __LINE__,__FILE__);
0252         e.setFatal();
0253         throw e;
0254       }
0255       #ifdef DEBUG
0256       debugOut << "KalmanFitterInfo::getFittedState - biased at first measurement = backwardUpdate_ \n";
0257       //backwardUpdate_->Print();
0258       #endif
0259       return *backwardUpdate_;
0260     }
0261
0262     if(!forwardUpdate_ || !backwardPrediction_) {
0263       Exception e("KalmanFitterInfo::getFittedState: Needed updates/predictions not available in this FitterInfo.", __LINE__,__FILE__);
0264       e.setFatal();
0265       throw e;
0266     }
0267     #ifdef DEBUG
0268     debugOut << "KalmanFitterInfo::getFittedState - biased = mean(forwardUpdate_, backwardPrediction_) \n";
0269     #endif
0270     fittedStateBiased_.reset(new MeasuredStateOnPlane(calcAverageState(*forwardUpdate_, *backwardPrediction_)));
0271     return *fittedStateBiased_;
0272   }
0273
0274   // unbiased
0275
0276   // last measurement and no backward prediction -> use forward prediction
0277   if (last && !backwardPrediction_) {
0278     if(!forwardPrediction_) {
0279       Exception e("KalmanFitterInfo::getFittedState: Needed updates/predictions not available in this FitterInfo.", __LINE__,__FILE__);
0280       e.setFatal();
0281       throw e;
0282     }
0283     #ifdef DEBUG
0284     debugOut << "KalmanFitterInfo::getFittedState - unbiased at last measurement = forwardPrediction_ \n";
0285     #endif
0286     return *forwardPrediction_;
0287   }
0288
0289   // first measurement and no forward prediction -> use backward prediction
0290   if (first && !forwardPrediction_) {
0291     if(!backwardPrediction_) {
0292       Exception e("KalmanFitterInfo::getFittedState: Needed updates/predictions not available in this FitterInfo.", __LINE__,__FILE__);
0293       e.setFatal();
0294       throw e;
0295     }
0296     #ifdef DEBUG
0297     debugOut << "KalmanFitterInfo::getFittedState - unbiased at first measurement = backwardPrediction_ \n";
0298     #endif
0299     return *backwardPrediction_;
0300   }
0301
0302   if(!forwardPrediction_ || !backwardPrediction_) {
0303     Exception e("KalmanFitterInfo::getFittedState: Needed updates/predictions not available in this FitterInfo.", __LINE__,__FILE__);
0304     e.setFatal();
0305     throw e;
0306   }
0307   #ifdef DEBUG
0308   debugOut << "KalmanFitterInfo::getFittedState - unbiased = mean(forwardPrediction_, backwardPrediction_) \n";
0309   #endif
0310   fittedStateUnbiased_.reset(new MeasuredStateOnPlane(calcAverageState(*forwardPrediction_, *backwardPrediction_)));
0311   return *fittedStateUnbiased_;
0312 }
0313
0314
0315 MeasurementOnPlane KalmanFitterInfo::getResidual(unsigned int iMeasurement, bool biased, bool onlyMeasurementErrors) const {
0316
0317   const MeasuredStateOnPlane& smoothedState = getFittedState(biased);
0318   const MeasurementOnPlane* measurement = measurementsOnPlane_.at(iMeasurement);
0319   const SharedPlanePtr& plane = measurement->getPlane();
0320
0321   // check equality of planes and reps
0322   if(*(smoothedState.getPlane()) != *plane) {
0323     Exception e("KalmanFitterInfo::getResidual: smoothedState and measurement are not defined in the same plane.", __LINE__,__FILE__);
0324     throw e;
0325   }
0326   if(smoothedState.getRep() != measurement->getRep()) {
0327     Exception e("KalmanFitterInfo::getResidual: smoothedState and measurement are not defined wrt the same TrackRep.", __LINE__,__FILE__);
0328     throw e;
0329   }
0330
0331   const AbsHMatrix* H = measurement->getHMatrix();
0332
0333   //TODO: shouldn't the definition be (smoothed - measured) ?
0334   // res = -(H*smoothedState - measuredState)
0335   TVectorD res(H->Hv(smoothedState.getState()));
0336   res -= measurement->getState();
0337   res *= -1;
0338
0339   if (onlyMeasurementErrors) {
0340     return MeasurementOnPlane(res, measurement->getCov(), plane, smoothedState.getRep(), H->clone(), measurement->getWeight());
0341   }
0342
0343   TMatrixDSym cov(smoothedState.getCov());
0344   H->HMHt(cov);
0345   cov += measurement->getCov();
0346
0347   return MeasurementOnPlane(res, cov, plane, smoothedState.getRep(), H->clone(), measurement->getWeight());
0348 }
0349
0350
0351 double KalmanFitterInfo::getSmoothedChi2(unsigned int iMeasurement) const {
0352   const MeasurementOnPlane& res = getResidual(iMeasurement, true, false);
0353
0354   TMatrixDSym Rinv;
0355   tools::invertMatrix(res.getCov(), Rinv);
0356   return Rinv.Similarity(res.getState());
0357 }
0358
0359
0360 void KalmanFitterInfo::setReferenceState(ReferenceStateOnPlane* referenceState) {
0361   referenceState_.reset(referenceState);
0362   if (referenceState_)
0363     setPlane(referenceState_->getPlane());
0364
0365   // if plane has changed, delete outdated info
0366   /*if (forwardPrediction_ && forwardPrediction_->getPlane() != getPlane())
0367     setForwardPrediction(0);
0368
0369   if (forwardUpdate_ && forwardUpdate_->getPlane() != getPlane())
0370     setForwardUpdate(0);
0371
0372   if (backwardPrediction_ && backwardPrediction_->getPlane() != getPlane())
0373     setBackwardPrediction(0);
0374
0375   if (backwardUpdate_ && backwardUpdate_->getPlane() != getPlane())
0376     setBackwardUpdate(0);
0377
0378   if (measurementsOnPlane_.size() > 0 && measurementsOnPlane_[0]->getPlane() != getPlane())
0379     deleteMeasurementInfo();
0380   */
0381 }
0382
0383
0384 void KalmanFitterInfo::setForwardPrediction(MeasuredStateOnPlane* forwardPrediction) {
0385   forwardPrediction_.reset(forwardPrediction);
0386   fittedStateUnbiased_.reset();
0387   fittedStateBiased_.reset();
0388   if (forwardPrediction_)
0389     setPlane(forwardPrediction_->getPlane());
0390 }
0391
0392 void KalmanFitterInfo::setBackwardPrediction(MeasuredStateOnPlane* backwardPrediction) {
0393   backwardPrediction_.reset(backwardPrediction);
0394   fittedStateUnbiased_.reset();
0395   fittedStateBiased_.reset();
0396   if (backwardPrediction_)
0397     setPlane(backwardPrediction_->getPlane());
0398 }
0399
0400 void KalmanFitterInfo::setForwardUpdate(KalmanFittedStateOnPlane* forwardUpdate) {
0401   forwardUpdate_.reset(forwardUpdate);
0402   fittedStateUnbiased_.reset();
0403   fittedStateBiased_.reset();
0404   if (forwardUpdate_)
0405     setPlane(forwardUpdate_->getPlane());
0406 }
0407
0408 void KalmanFitterInfo::setBackwardUpdate(KalmanFittedStateOnPlane* backwardUpdate) {
0409   backwardUpdate_.reset(backwardUpdate);
0410   fittedStateUnbiased_.reset();
0411   fittedStateBiased_.reset();
0412   if (backwardUpdate_)
0413     setPlane(backwardUpdate_->getPlane());
0414 }
0415
0416
0417 void KalmanFitterInfo::setMeasurementsOnPlane(const std::vector< genfit::MeasurementOnPlane* >& measurementsOnPlane) {
0418   deleteMeasurementInfo();
0419
0420   for (std::vector<MeasurementOnPlane*>::const_iterator m = measurementsOnPlane.begin(), mend = measurementsOnPlane.end(); m < mend; ++m) {
0421     addMeasurementOnPlane(*m);
0422   }
0423 }
0424
0425
0426 void KalmanFitterInfo::addMeasurementOnPlane(MeasurementOnPlane* measurementOnPlane) {
0427   if (measurementsOnPlane_.size() == 0)
0428     setPlane(measurementOnPlane->getPlane());
0429
0430   measurementsOnPlane_.push_back(measurementOnPlane);
0431 }
0432
0433 void KalmanFitterInfo::addMeasurementsOnPlane(const std::vector< genfit::MeasurementOnPlane* >& measurementsOnPlane) {
0434   for (std::vector<MeasurementOnPlane*>::const_iterator m = measurementsOnPlane.begin(), mend = measurementsOnPlane.end(); m < mend; ++m) {
0435     addMeasurementOnPlane(*m);
0436   }
0437 }
0438
0439
0440 void KalmanFitterInfo::setRep(const AbsTrackRep* rep) {
0441   rep_ = rep;
0442
0443   if (referenceState_)
0444     referenceState_->setRep(rep);
0445
0446   if (forwardPrediction_)
0447     forwardPrediction_->setRep(rep);
0448
0449   if (forwardUpdate_)
0450     forwardUpdate_->setRep(rep);
0451
0452   if (backwardPrediction_)
0453     backwardPrediction_->setRep(rep);
0454
0455   if (backwardUpdate_)
0456     backwardUpdate_->setRep(rep);
0457
0458   for (std::vector<MeasurementOnPlane*>::iterator it = measurementsOnPlane_.begin(); it != measurementsOnPlane_.end(); ++it) {
0459     (*it)->setRep(rep);
0460   }
0461 }
0462
0463
0464 void KalmanFitterInfo::setWeights(const std::vector<double>& weights) {
0465
0466   if (weights.size() != getNumMeasurements()) {
0467     Exception e("KalmanFitterInfo::setWeights: weights do not have the same size as measurementsOnPlane", __LINE__,__FILE__);
0468     throw e;
0469   }
0470
0471   if (fixWeights_) {
0472     errorOut << "KalmanFitterInfo::setWeights - WARNING: setting weights even though weights are fixed!" << std::endl;
0473   }
0474
0475   for (unsigned int i=0; i<getNumMeasurements(); ++i) {
0476     getMeasurementOnPlane(i)->setWeight(weights[i]);
0477   }
0478 }
0479
0480
0481 void KalmanFitterInfo::deleteForwardInfo() {
0482   setForwardPrediction(nullptr);
0483   setForwardUpdate(nullptr);
0484   fittedStateUnbiased_.reset();
0485   fittedStateBiased_.reset();
0486 }
0487
0488 void KalmanFitterInfo::deleteBackwardInfo() {
0489   setBackwardPrediction(nullptr);
0490   setBackwardUpdate(nullptr);
0491   fittedStateUnbiased_.reset();
0492   fittedStateBiased_.reset();
0493 }
0494
0495 void KalmanFitterInfo::deletePredictions() {
0496   setForwardPrediction(nullptr);
0497   setBackwardPrediction(nullptr);
0498   fittedStateUnbiased_.reset();
0499   fittedStateBiased_.reset();
0500 }
0501
0502 void KalmanFitterInfo::deleteMeasurementInfo() {
0503   // FIXME: need smart pointers / smart containers here
0504   for (size_t i = 0; i < measurementsOnPlane_.size(); ++i)
0505     delete measurementsOnPlane_[i];
0506
0507   measurementsOnPlane_.clear();
0508 }
0509
0510
0511 void KalmanFitterInfo::Print(const Option_t*) const {
0512   printOut << "genfit::KalmanFitterInfo. Belongs to TrackPoint " << trackPoint_ << "; TrackRep " <<  rep_  << "\n";
0513
0514   if (fixWeights_)
0515     printOut << "Weights are fixed.\n";
0516
0517   for (unsigned int i=0; i<measurementsOnPlane_.size(); ++i) {
0518     printOut << "MeasurementOnPlane Nr " << i <<": "; measurementsOnPlane_[i]->Print();
0519   }
0520
0521   if (referenceState_) {
0522     printOut << "Reference state: "; referenceState_->Print();
0523   }
0524   if (forwardPrediction_) {
0525     printOut << "Forward prediction_: "; forwardPrediction_->Print();
0526   }
0527   if (forwardUpdate_) {
0528     printOut << "Forward update: "; forwardUpdate_->Print();
0529   }
0530   if (backwardPrediction_) {
0531     printOut << "Backward prediction_: "; backwardPrediction_->Print();
0532   }
0533   if (backwardUpdate_) {
0534     printOut << "Backward update: "; backwardUpdate_->Print();
0535   }
0536
0537 }
0538
0539
0540 bool KalmanFitterInfo::checkConsistency(const genfit::PruneFlags* flags) const {
0541
0542   bool retVal(true);
0543
0544   // check if in a TrackPoint
0545   if (!trackPoint_) {
0546     errorOut << "KalmanFitterInfo::checkConsistency(): trackPoint_ is nullptr" << std::endl;
0547     retVal = false;
0548   }
0549
0550   // check if there is a reference state
0551   /*if (!referenceState_) {
0552     errorOut << "KalmanFitterInfo::checkConsistency(): referenceState_ is nullptr" << std::endl;
0553     retVal = false;
0554   }*/
0555
0556   SharedPlanePtr plane = getPlane();
0557
0558   if (plane.get() == nullptr) {
0559     if (!(referenceState_ || forwardPrediction_ || forwardUpdate_ || backwardPrediction_ || backwardUpdate_ || measurementsOnPlane_.size() > 0))
0560       return true;
0561     errorOut << "KalmanFitterInfo::checkConsistency(): plane is nullptr" << std::endl;
0562     retVal = false;
0563   }
0564
0565   // cppcheck-suppress unreadVariable
0566   TVector3 oTest = plane->getO(); // see if the plane object is there
0567   // cppcheck-suppress unreadVariable
0568   oTest *= 47.11;
0569
0570   // if more than one measurement, check if they have the same dimensionality
0571   if (measurementsOnPlane_.size() > 1) {
0572     int dim = measurementsOnPlane_[0]->getState().GetNrows();
0573     for (unsigned int i = 1; i<measurementsOnPlane_.size(); ++i) {
0574       if(measurementsOnPlane_[i]->getState().GetNrows() != dim) {
0575         errorOut << "KalmanFitterInfo::checkConsistency(): measurementsOnPlane_ do not all have the same dimensionality" << std::endl;
0576         retVal = false;
0577       }
0578     }
0579     if (dim == 0) {
0580       errorOut << "KalmanFitterInfo::checkConsistency(): measurementsOnPlane_ have dimension 0" << std::endl;
0581       retVal = false;
0582     }
0583   }
0584
0585   // see if everything else is defined wrt this plane and rep_
0586   int dim = rep_->getDim(); // check dim
0587   if (referenceState_) {
0588     if(referenceState_->getPlane() != plane) {
0589       errorOut << "KalmanFitterInfo::checkConsistency(): referenceState_ is not defined with the correct plane " << referenceState_->getPlane().get() << " vs. " << plane.get() << std::endl;
0590       retVal = false;
0591     }
0592     if (referenceState_->getRep() != rep_) {
0593       errorOut << "KalmanFitterInfo::checkConsistency(): referenceState_ is not defined with the correct TrackRep" << std::endl;
0594       retVal = false;
0595     }
0596     if (referenceState_->getState().GetNrows() != dim) {
0597       errorOut << "KalmanFitterInfo::checkConsistency(): referenceState_ does not have the right dimension!" << std::endl;
0598       retVal = false;
0599     }
0600   }
0601
0602   if (forwardPrediction_) {
0603     if(forwardPrediction_->getPlane() != plane) {
0604       errorOut << "KalmanFitterInfo::checkConsistency(): forwardPrediction_ is not defined with the correct plane" << std::endl;
0605       retVal = false;
0606     }
0607     if(forwardPrediction_->getRep() != rep_) {
0608       errorOut << "KalmanFitterInfo::checkConsistency(): forwardPrediction_ is not defined with the correct TrackRep" << std::endl;
0609       retVal = false;
0610     }
0611     if (forwardPrediction_->getState().GetNrows() != dim || forwardPrediction_->getCov().GetNrows() != dim) {
0612       errorOut << "KalmanFitterInfo::checkConsistency(): forwardPrediction_ does not have the right dimension!" << std::endl;
0613       retVal = false;
0614     }
0615   }
0616   if (forwardUpdate_) {
0617     if(forwardUpdate_->getPlane() != plane) {
0618       errorOut << "KalmanFitterInfo::checkConsistency(): forwardUpdate_ is not defined with the correct plane" << std::endl;
0619       retVal = false;
0620     }
0621     if(forwardUpdate_->getRep() != rep_) {
0622       errorOut << "KalmanFitterInfo::checkConsistency(): forwardUpdate_ is not defined with the correct TrackRep" << std::endl;
0623       retVal = false;
0624     }
0625     if (forwardUpdate_->getState().GetNrows() != dim || forwardUpdate_->getCov().GetNrows() != dim) {
0626       errorOut << "KalmanFitterInfo::checkConsistency(): forwardUpdate_ does not have the right dimension!" << std::endl;
0627       retVal = false;
0628     }
0629   }
0630
0631   if (backwardPrediction_) {
0632     if(backwardPrediction_->getPlane() != plane) {
0633       errorOut << "KalmanFitterInfo::checkConsistency(): backwardPrediction_ is not defined with the correct plane" << std::endl;
0634       retVal = false;
0635     }
0636     if(backwardPrediction_->getRep() != rep_) {
0637       errorOut << "KalmanFitterInfo::checkConsistency(): backwardPrediction_ is not defined with the correct TrackRep" << std::endl;
0638       retVal = false;
0639     }
0640     if (backwardPrediction_->getState().GetNrows() != dim || backwardPrediction_->getCov().GetNrows() != dim) {
0641       errorOut << "KalmanFitterInfo::checkConsistency(): backwardPrediction_ does not have the right dimension!" << std::endl;
0642       retVal = false;
0643     }
0644   }
0645   if (backwardUpdate_) {
0646     if(backwardUpdate_->getPlane() != plane) {
0647       errorOut << "KalmanFitterInfo::checkConsistency(): backwardUpdate_ is not defined with the correct plane" << std::endl;
0648       retVal = false;
0649     }
0650     if(backwardUpdate_->getRep() != rep_) {
0651       errorOut << "KalmanFitterInfo::checkConsistency(): backwardUpdate_ is not defined with the correct TrackRep" << std::endl;
0652       retVal = false;
0653     }
0654     if (backwardUpdate_->getState().GetNrows() != dim || backwardUpdate_->getCov().GetNrows() != dim) {
0655       errorOut << "KalmanFitterInfo::checkConsistency(): backwardUpdate_ does not have the right dimension!" << std::endl;
0656       retVal = false;
0657     }
0658   }
0659
0660   for (std::vector<MeasurementOnPlane*>::const_iterator it = measurementsOnPlane_.begin(); it != measurementsOnPlane_.end(); ++it) {
0661     if((*it)->getPlane() != plane) {
0662       errorOut << "KalmanFitterInfo::checkConsistency(): measurement is not defined with the correct plane" << std::endl;
0663       retVal = false;
0664     }
0665     if((*it)->getRep() != rep_) {
0666       errorOut << "KalmanFitterInfo::checkConsistency(): measurement is not defined with the correct TrackRep" << std::endl;
0667       retVal = false;
0668     }
0669     if ((*it)->getState().GetNrows() == 0) {
0670       errorOut << "KalmanFitterInfo::checkConsistency(): measurement has dimension 0!" << std::endl;
0671       retVal = false;
0672     }
0673   }
0674
0675   if (flags == nullptr or !flags->hasFlags("U")) { // if predictions have not been pruned
0676     // see if there is an update w/o prediction or measurement
0677     if (forwardUpdate_ && !forwardPrediction_) {
0678       errorOut << "KalmanFitterInfo::checkConsistency(): forwardUpdate_ w/o forwardPrediction_" << std::endl;
0679       retVal = false;
0680     }
0681
0682
0683     if (backwardUpdate_ && !backwardPrediction_) {
0684       errorOut << "KalmanFitterInfo::checkConsistency(): backwardUpdate_ w/o backwardPrediction_" << std::endl;
0685       retVal = false;
0686     }
0687
0688     if (flags == nullptr or !flags->hasFlags("M")) {
0689       if (forwardUpdate_ && measurementsOnPlane_.size() == 0) {
0690         errorOut << "KalmanFitterInfo::checkConsistency(): forwardUpdate_ w/o measurement" << std::endl;
0691         retVal = false;
0692       }
0693
0694       if (backwardUpdate_ && measurementsOnPlane_.size() == 0) {
0695         errorOut << "KalmanFitterInfo::checkConsistency(): backwardUpdate_ w/o measurement" << std::endl;
0696         retVal = false;
0697       }
0698     }
0699   }
0700
0701
0702   return retVal;
0703 }
0704
0705
0706 // Modified from auto-generated Streamer to correctly deal with smart pointers.
0707 void KalmanFitterInfo::Streamer(TBuffer &R__b)
0708 {
0709    // Stream an object of class genfit::KalmanFitterInfo.
0710
0711    //This works around a msvc bug and should be harmless on other platforms
0712    typedef ::genfit::KalmanFitterInfo thisClass;
0713    UInt_t R__s, R__c;
0714    if (R__b.IsReading()) {
0715       Version_t R__v = R__b.ReadVersion(&R__s, &R__c); if (R__v) { }
0716       //This works around a msvc bug and should be harmless on other platforms
0717       typedef genfit::AbsFitterInfo baseClass0;
0718       baseClass0::Streamer(R__b);
0719       int flag;
0720       R__b >> flag;
0721       deleteForwardInfo();
0722       deleteBackwardInfo();
0723       deleteReferenceInfo();
0724       deleteMeasurementInfo();
0725       if (flag & 1) {
0726         referenceState_.reset(new ReferenceStateOnPlane());
0727         referenceState_->Streamer(R__b);
0728         referenceState_->setPlane(getPlane());
0729         // rep needs to be fixed up
0730       }
0731       if (flag & (1 << 1)) {
0732         forwardPrediction_.reset(new MeasuredStateOnPlane());
0733         forwardPrediction_->Streamer(R__b);
0734         forwardPrediction_->setPlane(getPlane());
0735         // rep needs to be fixed up
0736       }
0737       if (flag & (1 << 2)) {
0738         forwardUpdate_.reset(new KalmanFittedStateOnPlane());
0739         forwardUpdate_->Streamer(R__b);
0740         forwardUpdate_->setPlane(getPlane());
0741         // rep needs to be fixed up
0742       }
0743       if (flag & (1 << 3)) {
0744         backwardPrediction_.reset(new MeasuredStateOnPlane());
0745         backwardPrediction_->Streamer(R__b);
0746         backwardPrediction_->setPlane(getPlane());
0747         // rep needs to be fixed up
0748       }
0749       if (flag & (1 << 4)) {
0750         backwardUpdate_.reset(new KalmanFittedStateOnPlane());
0751         backwardUpdate_->Streamer(R__b);
0752         backwardUpdate_->setPlane(getPlane());
0753         // rep needs to be fixed up
0754       }
0755       {
0756         std::vector<genfit::MeasurementOnPlane*,std::allocator<genfit::MeasurementOnPlane*> > &R__stl =  measurementsOnPlane_;
0757         TClass *R__tcl1 = TBuffer::GetClass(typeid(genfit::MeasurementOnPlane));
0758         if (R__tcl1==0) {
0759           Error("measurementsOnPlane_ streamer","Missing the TClass object for genfit::MeasurementOnPlane!");
0760           return;
0761         }
0762         int R__i, R__n;
0763         R__b >> R__n;
0764         R__stl.reserve(R__n);
0765         for (R__i = 0; R__i < R__n; R__i++) {
0766           genfit::MeasurementOnPlane* R__t = new MeasurementOnPlane();
0767           R__t->Streamer(R__b);
0768           R__t->setPlane(getPlane());
0769           R__stl.push_back(R__t);
0770         }
0771       }
0772       R__b.CheckByteCount(R__s, R__c, thisClass::IsA());
0773    } else {
0774      R__c = R__b.WriteVersion(thisClass::IsA(), kTRUE);
0775      //This works around a msvc bug and should be harmless on other platforms
0776      typedef genfit::AbsFitterInfo baseClass0;
0777      baseClass0::Streamer(R__b);
0778      // "!!" forces the value to 1 or 0 (pointer != 0 or pointer == 0),
0779      // this value is then written as a bitfield.
0780      int flag = ((!!referenceState_)
0781          | (!!forwardPrediction_ << 1)
0782          | (!!forwardUpdate_ << 2)
0783          | (!!backwardPrediction_ << 3)
0784          | (!!backwardUpdate_ << 4));
0785      R__b << flag;
0786      if (flag & 1)
0787        referenceState_->Streamer(R__b);
0788      if (flag & (1 << 1))
0789        forwardPrediction_->Streamer(R__b);
0790      if (flag & (1 << 2))
0791        forwardUpdate_->Streamer(R__b);
0792      if (flag & (1 << 3))
0793        backwardPrediction_->Streamer(R__b);
0794      if (flag & (1 << 4))
0795        backwardUpdate_->Streamer(R__b);
0796      {
0797        std::vector<genfit::MeasurementOnPlane*,std::allocator<genfit::MeasurementOnPlane*> > &R__stl =  measurementsOnPlane_;
0798        int R__n=int(R__stl.size());
0799        R__b << R__n;
0800        if(R__n) {
0801          std::vector<genfit::MeasurementOnPlane*,std::allocator<genfit::MeasurementOnPlane*> >::iterator R__k;
0802          for (R__k = R__stl.begin(); R__k != R__stl.end(); ++R__k) {
0803            (*R__k)->Streamer(R__b);
0804          }
0805        }
0806      }
0807      R__b.SetByteCount(R__c, kTRUE);
0808   }
0809 }
0810
0811
0812 } /* End of namespace genfit */