GENFIT  Rev:NoNumberAvailable
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules
GblFitterInfo.cc
Go to the documentation of this file.
1 #include "GblFitterInfo.h"
2 #include "MeasurementOnPlane.h"
3 #include "HMatrixUV.h"
4 #include "HMatrixV.h"
5 #include "HMatrixU.h"
6 
8 namespace genfit {
9 
10  GblFitterInfo::GblFitterInfo() : AbsFitterInfo(), jacobian_(TMatrixD(5, 5)) {
11  reset();
12  }
13 
14  GblFitterInfo::GblFitterInfo(const TrackPoint* trackPoint, const AbsTrackRep* rep) : AbsFitterInfo(trackPoint, rep), jacobian_(TMatrixD(5, 5)) {
15  reset();
16 
17  // Retrieve jacobian from rep by default (may not be used ... at 1st point)
18  unsigned int dim = rep->getDim();
19  if (dim != 5)
20  throw new genfit::Exception("GblFitterInfo: Representation state is not 5D", __LINE__, __FILE__);
21  TMatrixDSym noise(dim, dim);
22  TVectorD dState(dim);
23  rep->getForwardJacobianAndNoise(jacobian_, noise, dState);
24 
25  }
26 
27  GblFitterInfo::GblFitterInfo(const TrackPoint* trackPoint, const AbsTrackRep* rep, StateOnPlane& referenceState) : AbsFitterInfo(trackPoint, rep), jacobian_(TMatrixD(5, 5)) {
28  reset();
29 
30  // Retrieve jacobian from rep by default (may not be used ... at 1st point)
31  unsigned int dim = rep->getDim();
32  if (dim != 5)
33  throw new genfit::Exception("GblFitterInfo: Representation state is not 5D", __LINE__, __FILE__);
34  TMatrixDSym noise(dim, dim);
35  TVectorD dState(dim);
36  rep->getForwardJacobianAndNoise(jacobian_, noise, dState);
37 
38  setReferenceState(referenceState);
39  }
40 
41  void GblFitterInfo::reset(unsigned int measurementDim, unsigned int repDim) {
42  refPrediction_.ResizeTo(repDim);
43  refPrediction_.Zero();
44 
45  measResiduals_.ResizeTo(measurementDim);
46  measResiduals_.Zero();
47  kinkResiduals_.ResizeTo(measurementDim);
48  kinkResiduals_.Zero();
49  measResidualErrors_.ResizeTo(measurementDim);
50  measResidualErrors_.Zero();
51  kinkResidualErrors_.ResizeTo(measurementDim);
52  kinkResidualErrors_.Zero();
53  measDownWeights_.ResizeTo(measurementDim);
54  measDownWeights_.Zero();
55  kinkDownWeights_.ResizeTo(measurementDim);
56  kinkDownWeights_.Zero();
57 
58  TMatrixDSym zero(repDim);
59  zero.Zero();
60 
61  fwdStateCorrection_.ResizeTo(repDim);
62  fwdCov_.ResizeTo(zero);
63  bwdStateCorrection_.ResizeTo(repDim);
64  bwdCov_.ResizeTo(zero);
65  fwdPrediction_.ResizeTo(repDim);
66  bwdPrediction_.ResizeTo(repDim);
67  fwdCov_.SetMatrixArray(zero.GetMatrixArray());
68  bwdCov_.SetMatrixArray(zero.GetMatrixArray());
69 
70  measurement_.ResizeTo(2);
71  measurement_.Zero();
72  measCov_.ResizeTo(TMatrixDSym(2));
73  measCov_.Zero();
74  hMatrix_.ResizeTo(HMatrixUV().getMatrix());
76  }
77 
79  updateMeasurementAndPlane(referenceState);
80 
81  refPrediction_ = referenceState.getState();
82  fwdPrediction_ = referenceState.getState();
83  bwdPrediction_ = referenceState.getState();
84 
85  //TODO: reset even if already fitted?
86  fittedStateBwd_.reset(new MeasuredStateOnPlane(referenceState, trackPoint_->getTrack()->getCovSeed()));
87  fittedStateFwd_.reset(new MeasuredStateOnPlane(referenceState, trackPoint_->getTrack()->getCovSeed()));
88 
89  }
90 
91  void GblFitterInfo::setJacobian(TMatrixD jacobian) {
92  jacobian_.ResizeTo(jacobian);
93  jacobian_ = jacobian;
94  }
95 
96  TMatrixDSym GblFitterInfo::getCovariance(double variance, TVector3 trackDirection, SharedPlanePtr measurementPlane) const {
97 
98  double c1 = trackDirection.Dot(measurementPlane->getU());
99  double c2 = trackDirection.Dot(measurementPlane->getV());
100 
101  TMatrixDSym scatCov(2);
102  scatCov(0, 0) = 1. - c2 * c2;
103  scatCov(1, 1) = 1. - c1 * c1;
104  scatCov(0, 1) = c1 * c2;
105  scatCov(1, 0) = c1 * c2;
106  scatCov *= variance * variance / (1. - c1 * c1 - c2 * c2) / (1. - c1 * c1 - c2 * c2) ;
107 
108  return scatCov;
109  }
110 
112  // All meas/scat info is added from genfit data again (to cope with possible RecoHit update)
113 
114  // Create GBL point with current jacobian
115  gbl::GblPoint thePoint(jacobian_);
116 
117  //NOTE: 3rd update and update anytime GblPoint is requested
118  // mostly likely will update with reference as on 2nd update
119  StateOnPlane sop(getFittedState().getState(), sharedPlane_, rep_);
120 
121 
122  // Scatterer
123  if (trackPoint_->hasThinScatterer()) {
124  if (!hasMeasurements()) {
125  //double scatVariance = trackPoint_->getMaterialInfo()->getMaterial().getDensity();
126  //TVectorD kinkPrec(2);
127  //kinkPrec(0) = 1./scatVariance/scatVariance; kinkPrec(1) = 1./scatVariance/scatVariance;
128  //thePoint.addScatterer(getKinks(), kinkPrec);
129  //TODO: if state at scatterer is updated, the direction of track might not be perpendicular anymore
130  // plane does not change at pure scatterer
131  TMatrixDSym kinkCov = getCovariance(trackPoint_->getMaterialInfo()->getMaterial().getDensity(), sop.getDir(), sop.getPlane());
132  thePoint.addScatterer(getKinks(), kinkCov.Invert());
133  } else {
134  TMatrixDSym kinkCov = getCovariance(trackPoint_->getMaterialInfo()->getMaterial().getDensity(), sop.getDir(), sop.getPlane());
135  thePoint.addScatterer(getKinks(), kinkCov.Invert());
136  }
137  }
138 
139 
140 
141  // Measurement
142  if (hasMeasurements()) {
143  MeasuredStateOnPlane measurement = getResidual(0, true, true);
144  TVectorD aResiduals(measurement.getState());
145  TMatrixDSym aPrecision(measurement.getCov().Invert());
146  if (HMatrixU().getMatrix() == hMatrix_) {
147  double res = aResiduals(0);
148  double prec = aPrecision(0, 0);
149  aResiduals.ResizeTo(2);
150  aPrecision.ResizeTo(TMatrixDSym(2));
151  aResiduals.Zero();
152  aResiduals(0) = res;
153  aPrecision.Zero();
154  aPrecision(0, 0) = prec;
155  }
156  if (HMatrixV().getMatrix() == hMatrix_) {
157  double res = aResiduals(0);
158  double prec = aPrecision(0, 0);
159  aResiduals.ResizeTo(2);
160  aPrecision.ResizeTo(TMatrixDSym(2));
161  aResiduals.Zero();
162  aResiduals(1) = res;
163  aPrecision.Zero();
164  aPrecision(1, 1) = prec;
165  }
166  // always 2D, U/V set by precision (=0 to disable)
167  thePoint.addMeasurement(aResiduals, aPrecision);
168  }
169 
170  // Derivatives
171  ICalibrationParametersDerivatives* globals = NULL;
172  if (hasMeasurements() && (globals = dynamic_cast<ICalibrationParametersDerivatives*>(trackPoint_->getRawMeasurement(0)) )) {
173  TMatrixD derivs = globals->derivatives(&sop);
174  std::vector<int> labels = globals->labels();
175 
176  if (derivs.GetNcols() > 0 && !labels.empty() && (unsigned int)derivs.GetNcols() == labels.size()) {
177  thePoint.addGlobals(labels, derivs);
178  }
179  TMatrixD locals = globals->localDerivatives(&sop);
180  if (locals.GetNcols() > 0) {
181  thePoint.addLocals(locals);
182  GblFitStatus* gblfs = dynamic_cast<GblFitStatus*>(trackPoint_->getTrack()->getFitStatus(rep_));
183  if (gblfs) {
184  if (gblfs->getMaxLocalFitParams() < locals.GetNcols())
185  gblfs->setMaxLocalFitParams(locals.GetNcols());
186  }
187  }
188  }
189 
190  return thePoint;
191 
192  }
193 
195  if (!trackPoint_)
196  return;
198  //no measurement, plane is supposed to come from state (is perpendicular)
199  setPlane(sop.getPlane());
200  return;
201  }
202  std::vector<MeasurementOnPlane*> allMeas = trackPoint_->getRawMeasurement(0)->constructMeasurementsOnPlane(sop);
203 
204  /*
205  double normMin(9.99E99);
206  unsigned int imop(0);
207  const AbsHMatrix* H = allMeas[0]->getHMatrix();
208  for (unsigned int i=0; i<allMeas.size(); ++i) {
209  if (*(allMeas[i]->getHMatrix()) != *H){
210  Exception e("GblFitterInfo::updateMeasurementAndPlane: Cannot compare measurements with different H-Matrices.", __LINE__,__FILE__);
211  e.setFatal();
212  throw e;
213  }
214 
215  TVectorD res = allMeas[i]->getState() - H->Hv(sop.getState());
216  double norm = sqrt(res.Norm2Sqr());
217  if (norm < normMin) {
218  normMin = norm;
219  imop = i;
220  }
221  }
222  */
223  unsigned int imop = 0;
224  double maxWeight = allMeas.at(0)->getWeight();
225  for (unsigned int i = 0; i < allMeas.size(); i++)
226  if (allMeas.at(i)->getWeight() > maxWeight)
227  imop = i;
228 
229  measurement_.ResizeTo(allMeas.at(imop)->getState());
230  measurement_ = allMeas.at(imop)->getState();
231  measCov_.ResizeTo(allMeas.at(imop)->getCov());
232  measCov_ = allMeas.at(imop)->getCov();
233  hMatrix_.ResizeTo(allMeas.at(imop)->getHMatrix()->getMatrix());
234  hMatrix_ = allMeas.at(imop)->getHMatrix()->getMatrix();
235 
236  setPlane(allMeas.at(imop)->getPlane());
237 
238  for (unsigned int imeas = 0; imeas < allMeas.size(); imeas++)
239  delete allMeas[imeas];
240  allMeas.clear();
241 
242 
243  }
244 
246  if (!traj.isValid())
247  return;
248 
249  // Deduce our position in the trajectory
250  //----------------------------------------------
251  unsigned int label = 0;
252  // Loop over track and find index of this fitter info
254  for (unsigned int ip = 0; ip < trk->getNumPoints(); ip++) {
255  // Indexing of GblFitter info in track (each gives one point to trajectory at trackpoint position)
256  if (dynamic_cast<GblFitterInfo*>( trk->getPoint(ip)->getFitterInfo(rep_) )) {
257  // First fitter info has label 1 (for point in GBL traj)
258  label++;
259  // We found itself = we have our label
260  if (trk->getPoint(ip)->getFitterInfo(rep_) == this)
261  break;
262  }
263  }
264  if (label == 0)
265  throw genfit::Exception("GblFitterInfo: fitter info did not found itself in track to update", __LINE__, __FILE__);
266  if (label > traj.getNumPoints())
267  throw genfit::Exception("GblFitterInfo: Deduced point label not valid", __LINE__, __FILE__);
268  //----------------------------------------------
269 
270  // Residuals (in plane)
271  unsigned int numMRes = 2;
272  TVectorD mResiduals(2), mMeasErrors(2), mResErrors(2), mDownWeights(2);
273  if (0 != traj.getMeasResults(label, numMRes, mResiduals, mMeasErrors, mResErrors, mDownWeights))
274  throw genfit::Exception(" NO measurement results ", __LINE__, __FILE__);
275 
276  // Kinks
277  unsigned int numKRes = 2;
278  TVectorD kResiduals(2), kMeasErrors(2), kResErrors(2), kDownWeights(2);
279  if (0 != traj.getScatResults(label, numKRes, kResiduals, kMeasErrors, kResErrors, kDownWeights))
280  throw genfit::Exception(" NO scattering results ", __LINE__, __FILE__);
281 
282  // Check for local derivatives
283  int nLocals = 0;
284  GblFitStatus* gblfs = dynamic_cast<GblFitStatus*>(trackPoint_->getTrack()->getFitStatus(rep_));
285  if (gblfs)
286  nLocals = gblfs->getMaxLocalFitParams();
287 
288  // Predictions (different at scatterer)
289  TVectorD bwdUpdate(5 + nLocals), fwdUpdate(5 + nLocals);
290  TMatrixDSym bwdCov(5 + nLocals), fwdCov(5 + nLocals);
291 
292  // forward prediction
293  if (0 != traj.getResults(label, fwdUpdate, fwdCov))
294  throw genfit::Exception(" NO forward results ", __LINE__, __FILE__);
295 
296  // backward prediction
297  if (0 != traj.getResults(-1 * label, bwdUpdate, bwdCov))
298  throw genfit::Exception(" NO backward results ", __LINE__, __FILE__);
299 
300  if (nLocals > 0) {
301  TVectorD _bwdUpdate(5 + nLocals), _fwdUpdate(5 + nLocals);
302  TMatrixDSym _bwdCov(5 + nLocals), _fwdCov(5 + nLocals);
303  _bwdUpdate = bwdUpdate;
304  _fwdUpdate = fwdUpdate;
305  _bwdCov = bwdCov;
306  _fwdCov = fwdCov;
307  bwdUpdate.ResizeTo(5);
308  fwdUpdate.ResizeTo(5);
309  bwdCov.ResizeTo(TMatrixDSym(5));
310  fwdCov.ResizeTo(TMatrixDSym(5));
311  for (int i = 0; i < 5; i++) {
312  bwdUpdate(i) = _bwdUpdate(i);
313  fwdUpdate(i) = _fwdUpdate(i);
314  for (int j = 0; j < 5; j++) {
315  bwdCov(i, j) = _bwdCov(i, j);
316  fwdCov(i, j) = _fwdCov(i, j);
317  }
318  }
319  }
320  // Now update the the fitter info
321  //
322  //-------------------------------------------------
323  // Backward/forward prediction (residual) (differs at scatterers) AFTER GBL fit
324  bwdStateCorrection_ = bwdUpdate;
325  fwdStateCorrection_ = fwdUpdate;
326  bwdCov_ = bwdCov;
327  fwdCov_ = fwdCov;
328  fwdPrediction_ += fwdStateCorrection_; // This is the update!
329  bwdPrediction_ += bwdStateCorrection_; // This is the update!
330 
333 
334  // Set scattering/measurement residual data
335  kinkResiduals_ = kResiduals;
336  measResiduals_ = mResiduals;
337  kinkResidualErrors_ = kResErrors;
338  measResidualErrors_ = mResErrors;
339  measDownWeights_ = mDownWeights;
340  kinkDownWeights_ = kDownWeights;
341 
342  //-------------------------------------------------
343  }
344 
346  {
347  // Invalidates errors and corrections from last iteration
348  // (will be defined in different plane). But fitted state and residual is ok.
349 
350  if (!prevFitterInfo) {
351  jacobian_.UnitMatrix();
352 
353  //TODO: For 1st plane this has no sense
354  //return;
355  prevFitterInfo = this;
356  }
357 
358  //TODO
360 
361  //updateMeasurementAndPlane(StateOnPlane(fwdPrediction_, sharedPlane_, rep_));
362  //
363 
364  TMatrixDSym noise;
365  TVectorD dstate;
366  // Take forward state from previous fitter info,
367  // its (maybe updated) plane
368  // and our rep
369  StateOnPlane prevState(prevFitterInfo->getFittedState(true).getState(), prevFitterInfo->getPlane(), rep_);
370 
371  if (hasMeasurements()) {
372  SharedPlanePtr newPlane = trackPoint_->getRawMeasurement(0)->constructPlane(prevState);
373  rep_->extrapolateToPlane(prevState, newPlane, false, true);
374  } else {
375  rep_->extrapolateToPlane(prevState, sharedPlane_, false, true);
376  }
377 
378  rep_->getForwardJacobianAndNoise(jacobian_, noise, dstate);
379  // Now update meas data
380  updateMeasurementAndPlane(prevState);
381 
382  //
383  // Extrap predictions to new plane
384  //
385  StateOnPlane oldFwdState(fwdPrediction_, oldPlane, rep_);
386  StateOnPlane oldBwdState(bwdPrediction_, oldPlane, rep_);
387  rep_->extrapolateToPlane(oldFwdState, sharedPlane_);
388  rep_->extrapolateToPlane(oldBwdState, sharedPlane_);
389  fwdPrediction_ = oldFwdState.getState();
390  bwdPrediction_ = oldBwdState.getState();
391  fittedStateBwd_.reset();
392  fittedStateFwd_.reset();
393  //
394  }
395 
396 
397  const MeasuredStateOnPlane& GblFitterInfo::getFittedState(bool afterKink) const {
398  // ALways biased from GBL (global fit!)
399 
400  if (!fittedStateFwd_ || !fittedStateBwd_) {
403  }
404 
405  if (afterKink) {
406  return *fittedStateFwd_;
407  }
408  else {
409  return *fittedStateBwd_;
410  }
411 
412  }
413 
414  MeasurementOnPlane GblFitterInfo::getResidual(unsigned int, bool, bool onlyMeasurementErrors) const {
415  // Return error of resduals (=0, none, for reference track, no errors known before fit (seed covariance might not be correct)
416  TMatrixDSym localCovariance(2);
417  localCovariance.Zero();
418  // combined meas+fit errors:
419  //TODO: 1D covariance + more dimensions of residuals
420 
421  localCovariance(0, 0) = measResidualErrors_(0) * measResidualErrors_(0);
422  localCovariance(1, 1) = measResidualErrors_(1) * measResidualErrors_(1);
423 
424  if (HMatrixU().getMatrix() == hMatrix_) {
425  localCovariance.ResizeTo(TMatrixDSym(1));
426  localCovariance(0, 0) = measResidualErrors_(0) * measResidualErrors_(0);
427  }
428  if (HMatrixV().getMatrix() == hMatrix_) {
429  localCovariance.ResizeTo(TMatrixDSym(1));
430  localCovariance(0, 0) = measResidualErrors_(1) * measResidualErrors_(1);
431  }
432 
433  if (hasMeasurements()){
434  // this means only for reference state before gbl fit, this way will be used
435  TVectorD res( measurement_ - hMatrix_ * getFittedState(true).getState() );
436 
437  if (onlyMeasurementErrors) {
438  localCovariance.ResizeTo(measCov_);
439  localCovariance = measCov_;
440  }
442  }
443  TVectorD zeroRes(2);
444  zeroRes.Zero();
445  // Else return 0's or whatever
446  //TODO: or throw?
447  return MeasurementOnPlane(zeroRes, localCovariance, sharedPlane_, rep_, new HMatrixUV());
448 
449  } // calculate residual, track and measurement errors are added if onlyMeasurementErrors is false
450 
452  TMatrixDSym localCovariance(2);
453  localCovariance.Zero();
454  localCovariance(0, 0) = kinkResidualErrors_(0) * kinkResidualErrors_(0);
455  localCovariance(1, 1) = kinkResidualErrors_(1) * kinkResidualErrors_(1);
456  return MeasurementOnPlane(getKinks(), localCovariance, sharedPlane_, rep_, new genfit::HMatrixUV());
457 
458  }
459 
460  TVectorD GblFitterInfo::getKinks() const {
461  TVectorD kinks(2);
462  kinks.Zero();
463 
464  TVectorD stateDiff(getFittedState(true).getState() - getFittedState(false).getState());
465  kinks(0) = -stateDiff(1);
466  kinks(1) = -stateDiff(2);
467 
468  return kinks;
469  }
470 
473  }
474 
476  //TODO
477  return true;
478  }
479 
481 
482  GblFitterInfo* retVal = new GblFitterInfo(this->getTrackPoint(), this->getRep());
483 
484  retVal->setPlane(sharedPlane_);
485 
486  retVal->jacobian_ = jacobian_;
487  retVal->measResiduals_ = measResiduals_;
489  retVal->kinkResiduals_ = kinkResiduals_;
495  retVal->bwdCov_ = bwdCov_;
496  retVal->fwdCov_ = fwdCov_;
497  retVal->fwdPrediction_ = fwdPrediction_;
498  retVal->bwdPrediction_ = bwdPrediction_;
499  retVal->refPrediction_ = refPrediction_;
500  retVal->measurement_.ResizeTo(measurement_);
501  retVal->measCov_.ResizeTo(measCov_);
502  retVal->hMatrix_.ResizeTo(hMatrix_);
503  retVal->measurement_ = measurement_;
504  retVal->measCov_ = measCov_;
505  retVal->hMatrix_ = hMatrix_;
506 
507  return retVal;
508  }
509 
510  void GblFitterInfo::Print(const Option_t*) const {
511  //TODO
512  std::cout << "=============================================================================================" << std::endl;
513  std::cout << " >>> GblFitterInfo " << std::endl;
514  std::cout << " ************* " << std::endl;
515 
516  std::cout << " rep: " << rep_ << ", trackpoint: " << trackPoint_ << ", plane: " << sharedPlane_.get() << std::endl;
517  sharedPlane_->Print();
518  std::cout << std::endl;
519 
520  std::cout << "=============================================================================================" << std::endl;
521  std::cout << " | PREDICTIONS | REFERENCE | Corrections from last iteration |" << std::endl;
522  std::cout << " | (+)prediction | (-)prediction | state | (+)correction | (-) correction |" << std::endl;
523  std::cout << "---------------------------------------------------------------------------------------------" << std::endl;
524 
525  for (int i = 0; i <5; i++) {
526  std::cout << std::left;
527  std::cout << " ";
528  if (i==0)
529  std::cout << "q/p";
530  if (i==1)
531  std::cout << "u' ";
532  if (i==2)
533  std::cout << "v' ";
534  if (i==3)
535  std::cout << "u ";
536  if (i==4)
537  std::cout << "v ";
538  std::cout << " | "
539  << std::setw(12) << fwdPrediction_(i) << " | "
540  << std::setw(12) << bwdPrediction_(i) << " | "
541  << std::setw(12) << refPrediction_(i) << " | "
542  << std::setw(12) << fwdStateCorrection_(i) << " | "
543  << std::setw(12) << bwdStateCorrection_(i) << std::endl;
544  }
545  std::cout << "=============================================================================================" << std::endl;
546 
547  if (hasMeasurements()) {
548  std::cout << " | Meas. residual | measurement - prediction | Down-weight | Fit+meas Err. |" << std::endl;
549  std::cout << " | | | | -diagonaliz. |" << std::endl;
550  std::cout << "---------------------------------------------------------------------------------------------" << std::endl;
551 
552  TVectorD residual = getResidual().getState();
553  if (residual.GetNoElements()<2) {
554  residual.ResizeTo(2);
555  residual.Zero();
557  residual(0) = getResidual().getState()(0);
558  else
559  residual(1) = getResidual().getState()(0);
560 
561  }
562  std::cout << " u | "
563  << std::setw(12) << measResiduals_(0) << " | "
564  << std::setw(12) << residual(0) << " | "
565  << std::setw(12) << measDownWeights_(0) << " | "
566  << std::setw(12) << measResidualErrors_(0)
567  << std::endl;
568 
569  std::cout << " v | "
570  << std::setw(12) << measResiduals_(1) << " | "
571  << std::setw(12) << residual(1) << " | "
572  << std::setw(12) << measDownWeights_(1) << " | "
573  << std::setw(12) << measResidualErrors_(1)
574  << std::endl;
575 
576  std::cout << "---------------------------------------------------------------------------------------------" << std::endl;
577  }
578 
579  std::cout << " | Kink residual | Residual of slope difference | Down-weight | Fit Kink Err. |" << std::endl;
580  std::cout << " | -diagonalized | - ( (+)pred - (-)pred ) | | -diagonaliz. |" << std::endl;
581  std::cout << "---------------------------------------------------------------------------------------------" << std::endl;
582 
583  std::cout << " u' | "
584  << std::setw(12) << kinkResiduals_(0) << " | "
585  << std::setw(12) << getKinks()(0) << " | "
586  << std::setw(12) << kinkDownWeights_(0) << " | "
587  << std::setw(12) << kinkResidualErrors_(0)
588  << std::endl;
589 
590  std::cout << " v' | "
591  << std::setw(12) << kinkResiduals_(1) << " | "
592  << std::setw(12) << getKinks()(1) << " | "
593  << std::setw(12) << kinkDownWeights_(1) << " | "
594  << std::setw(12) << kinkResidualErrors_(1)
595  << std::endl;
596  std::cout << "=============================================================================================" << std::endl;
597  std::cout << "Measurement: "; measurement_.Print();
598 
599  std::cout << "H Matrix: "; hMatrix_.Print();
600 
601  std::cout << "Measurement covariance: "; measCov_.Print();
602 
603  std::cout << "Jacobian: "; jacobian_.Print();
604  std::cout << "Backward covariance: "; bwdCov_.Print();
605  std::cout << "Forward covariance : "; fwdCov_.Print();
606 
607  std::cout << "=============================================================================================" << std::endl;
608 
609  }
610 
611 
612 } // end of namespace genfit
void addLocals(const TMatrixD &aDerivatives)
gbl::GblPoint constructGblPoint()
Collect all data and create a GblPoint.
TVector3 getDir() const
Definition: StateOnPlane.h:113
virtual void getForwardJacobianAndNoise(TMatrixD &jacobian, TMatrixDSym &noise, TVectorD &deltaState) const =0
Get the jacobian and noise matrix of the last extrapolation.
AbsHMatrix implementation for one-dimensional MeasurementOnPlane and RKTrackRep parameterization.
Definition: HMatrixU.h:37
boost::shared_ptr< genfit::DetPlane > SharedPlanePtr
Shared Pointer to a DetPlane.
GBL trajectory.
Definition: GblTrajectory.h:26
virtual SharedPlanePtr constructPlane(const StateOnPlane &state) const =0
void setMaxLocalFitParams(unsigned maxLocalFitParams)
Definition: GblFitStatus.h:53
FitStatus for use with GblFitter.
Definition: GblFitStatus.h:39
virtual std::vector< int > labels()=0
Vector of integer labels for calibration/alignment parameters available (must match #columns of deriv...
unsigned int getMeasResults(unsigned int aLabel, unsigned int &numRes, TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors, TVectorD &aDownWeights)
Get residuals at point from measurement.
void addMeasurement(const TMatrixD &aProjection, const TVectorD &aResiduals, const TVectorD &aPrecision, double minPrecision=0.)
Add a measurement to a point.
Definition: GblPoint.cc:47
const TVectorD & getState() const
Definition: StateOnPlane.h:61
bool hasThinScatterer() const
Definition: TrackPoint.h:108
virtual void Print(const Option_t *="") const
TMatrixDSym getCovariance(double variance, TVector3 trackDirection, SharedPlanePtr measurementPlane) const
Get scattering covariance projected into (measurement) plane.
Info which information has been pruned from the Track.
Definition: FitStatus.h:47
virtual double extrapolateToPlane(StateOnPlane &state, const genfit::SharedPlanePtr &plane, bool stopAtBoundary=false, bool calcJacobianNoise=false) const =0
Extrapolates the state to plane, and returns the extrapolation length and, via reference, the extrapolated state.
virtual bool checkConsistency(const genfit::PruneFlags *=NULL) const
unsigned int getNumPoints() const
Definition: Track.h:110
void addGlobals(const std::vector< int > &aLabels, const TMatrixD &aDerivatives)
Point on trajectory.
Definition: GblPoint.h:48
virtual std::vector< genfit::MeasurementOnPlane * > constructMeasurementsOnPlane(const StateOnPlane &state) const =0
Collects information needed and produced by a GblFitter/GBL and is specific to one AbsTrackRep of the...
Definition: GblFitterInfo.h:55
SharedPlanePtr sharedPlane_
No ownership.
Definition: AbsFitterInfo.h:94
GblFitterInfo()
Constructor for ROOT I/O.
virtual GblFitterInfo * clone() const
Deep copy ctor for polymorphic class.
const SharedPlanePtr & getPlane() const
Definition: AbsFitterInfo.h:74
boost::scoped_ptr< MeasuredStateOnPlane > fittedStateFwd_
cache
const TMatrixDSym & getCovSeed() const
Definition: Track.h:168
AbsMeasurement * getRawMeasurement(int i=0) const
Definition: TrackPoint.cc:147
TVectorD getKinks() const
Get kink (residual) (2D) = 0 - ( (+)pred - (-)pred )
virtual bool isEqual(const AbsHMatrix &other) const =0
boost::scoped_ptr< MeasuredStateOnPlane > fittedStateBwd_
StateOnPlane with additional covariance matrix.
unsigned int getNumPoints() const
Retrieve number of point from trajectory.
virtual TMatrixD derivatives(const genfit::StateOnPlane *sop)=0
Derivatives of residuals (local measurement coordinates) w.r.t. alignment/calibration parameters Matr...
AbsHMatrix implementation for two-dimensional MeasurementOnPlane and RKTrackRep parameterization.
Definition: HMatrixUV.h:39
TrackPoint * getPoint(int id) const
Definition: Track.cc:212
bool hasMeasurements() const
unsigned int getResults(int aSignedLabel, TVectorD &localPar, TMatrixDSym &localCov) const
Get fit results at point.
void updateMeasurementAndPlane(const StateOnPlane &sop)
SHOULD BE USED ONLY INTERNALY! Update the plane from measurement constructed with state or take plane...
const TMatrixD & getMatrix() const
Get the actual matrix representation.
Definition: HMatrixUV.cc:33
void reset(unsigned int measurementDim=2, unsigned int repDim=5)
(Initial) reset of fitter info
Abstract base class for a track representation.
Definition: AbsTrackRep.h:66
ThinScatterer * getMaterialInfo() const
Definition: TrackPoint.h:107
void updateFitResults(gbl::GblTrajectory &traj)
Update fitter info from GBL fit results.
const TMatrixD & getMatrix() const
Get the actual matrix representation.
Definition: HMatrixU.cc:33
unsigned int getScatResults(unsigned int aLabel, unsigned int &numRes, TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors, TVectorD &aDownWeights)
Get (kink) residuals at point from scatterer.
Exception class for error handling in GENFIT (provides storage for diagnostic information) ...
Definition: Exception.h:48
Object containing AbsMeasurement and AbsFitterInfo objects.
Definition: TrackPoint.h:50
MeasurementOnPlane getResidual(unsigned int=0, bool=false, bool onlyMeasurementErrors=true) const
Get the residual.
const MaterialProperties & getMaterial() const
Definition: ThinScatterer.h:48
void setReferenceState(StateOnPlane &referenceState)
Set the prediction and plane (from measurement if any) You should use the user constructor instead...
virtual TMatrixD localDerivatives(const genfit::StateOnPlane *)
Derivatives for additional local parameters to be fitted in global calibration algorithms together wi...
void addScatterer(const TVectorD &aResiduals, const TVectorD &aPrecision)
Add a (thin) scatterer to a point.
Definition: GblPoint.cc:188
MeasurementOnPlane getKink() const
Get kink (residual) with diagonalized covariance (2D) Covariance may be zero if not yet fitted or no ...
const TrackPoint * trackPoint_
Definition: AbsFitterInfo.h:88
const SharedPlanePtr & getPlane() const
Definition: StateOnPlane.h:65
Detector plane.
Definition: DetPlane.h:61
A state with arbitrary dimension defined in a DetPlane.
Definition: StateOnPlane.h:45
const MeasuredStateOnPlane & getFittedState(bool afterKink=true) const
Get the prediction at this point Always biased in GBL (global fit) There are 2 states, before and after kink (if there is a scatterer at this point). Per default the state after kink for forward propagation is returned.
Track * getTrack() const
Definition: TrackPoint.h:90
AbsHMatrix implementation for one-dimensional MeasurementOnPlane and RKTrackRep parameterization.
Definition: HMatrixV.h:37
bool hasRawMeasurements() const
Definition: TrackPoint.h:96
void recalculateJacobian(GblFitterInfo *prevFitterInfo)
Re-extrapolates between prevFitterInfo and this point using forward state to update the Jacobian (if ...
const AbsTrackRep * getRep() const
Definition: AbsFitterInfo.h:55
MeasurementOnPlane getMeasurement() const
Get the measurement on plane from stored measurement data (from last construction/update) ...
Collection of TrackPoint objects, AbsTrackRep objects and FitStatus objects.
Definition: Track.h:71
const TrackPoint * getTrackPoint() const
Definition: AbsFitterInfo.h:54
bool isValid() const
Retrieve validity of trajectory.
void setPlane(const SharedPlanePtr &plane)
Definition: AbsFitterInfo.h:78
Measured coordinates on a plane.
virtual const AbsHMatrix * constructHMatrix(const AbsTrackRep *) const =0
const TMatrixDSym & getCov() const
Abstract base class to establish an interface between physical representation of the detector for ali...
AbsFitterInfo * getFitterInfo(const AbsTrackRep *rep=NULL) const
Get fitterInfo for rep. Per default, use cardinal rep.
Definition: TrackPoint.cc:169
const AbsTrackRep * rep_
No ownership.
Definition: AbsFitterInfo.h:92
Defines for I/O streams used for error and debug printing.
This class collects all information needed and produced by a specific AbsFitter and is specific to on...
Definition: AbsFitterInfo.h:42
virtual unsigned int getDim() const =0
Get the dimension of the state vector used by the track representation.
FitStatus * getFitStatus(const AbsTrackRep *rep=NULL) const
Get FitStatus for a AbsTrackRep. Per default, return FitStatus for cardinalRep.
Definition: Track.h:152
void setJacobian(TMatrixD jacobian)
Set the Jacobian for further GblPoint construction.