GENFIT  Rev:NoNumberAvailable
GblFitter.cc
Go to the documentation of this file.
1 //-*-mode: C++; c-basic-offset: 2; -*-
2 /* Copyright 2013-2014
3  * Authors: Sergey Yashchenko and Tadeas Bilka
4  *
5  * This is an interface to General Broken Lines
6  *
7  * Version: 6 --------------------------------------------------------------
8  * - complete rewrite to GblFitter using GblFitterInfo and GblFitStatus
9  * - mathematics should be the same except for additional iterations
10  * (not possible before)
11  * - track is populated with scatterers + additional points with
12  * scatterers and no measurement (optional)
13  * - final track contains GblFitStatus and fitted states from GBL prediction
14  * - 1D/2D hits supported (pixel, single strip, combined strips(2D), wire)
15  * - At point: Only the very first raw measurement is used and from
16  * that, constructed MeasurementOnPlane with heighest weight
17  * Version: 5 --------------------------------------------------------------
18  * - several bug-fixes:
19  * - Scatterers at bad points
20  * - Jacobians at a point before they should be (code reorganized)
21  * - Change of sign of residuals
22  * Version: 4 --------------------------------------------------------------
23  * Fixed calculation of equvivalent scatterers (solution by C. Kleinwort)
24  * Now a scatterer is inserted at each measurement (except last) and
25  * between each two measurements.
26  * TrueHits/Clusters. Ghost (1D) hits ignored. With or
27  * without magnetic field.
28  * Version: 3 --------------------------------------------------------------
29  * This version now supports both TrueHits and Clusters for VXD.
30  * It can be used for arbitrary material distribution between
31  * measurements. Moments of scattering distribution are computed
32  * and translated into two equivalent thin GBL scatterers placed
33  * at computed positions between measurement points.
34  * Version: 2 ... never published -----------------------------------------
35  * Scatterer at each boundary (tooo many scatterers). TrueHits/Clusters.
36  * Without global der.&MP2 output.
37  * Version: 1 --------------------------------------------------------------
38  * Scatterers at measurement planes. TrueHits
39  * Version: 0 --------------------------------------------------------------
40  * Without scatterers. Genfit 1.
41  * -------------------------------------------------------------------------
42  *
43  * This file is part of GENFIT.
44  *
45  * GENFIT is free software: you can redistribute it and/or modify
46  * it under the terms of the GNU Lesser General Public License as published
47  * by the Free Software Foundation, either version 3 of the License, or
48  * (at your option) any later version.
49  *
50  * GENFIT is distributed in the hope that it will be useful,
51  * but WITHOUT ANY WARRANTY; without even the implied warranty of
52  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
53  * GNU Lesser General Public License for more details.
54  *
55  * You should have received a copy of the GNU Lesser General Public License
56  * along with GENFIT. If not, see <http://www.gnu.org/licenses/>.
57  */
58 
59 #include "GblFitter.h"
60 #include "../include/GblFitStatus.h"
61 #include "GblFitterInfo.h"
62 #include "GblTrajectory.h"
63 #include "GblPoint.h"
65 
66 #include "Track.h"
67 #include <TFile.h>
68 #include <TH1F.h>
69 #include <TTree.h>
70 #include <string>
71 #include <list>
72 #include <FieldManager.h>
73 #include <HMatrixU.h>
74 #include <HMatrixV.h>
75 #include <Math/SMatrix.h>
76 #include <TMatrixD.h>
77 #include <TVectorDfwd.h>
78 #include <TMatrixT.h>
79 #include <TVector3.h>
80 //#include "boost/algorithm/string.hpp"
81 
82 //#define DEBUG
83 
84 using namespace gbl;
85 using namespace std;
86 using namespace genfit;
87 
91 GblFitter::~GblFitter() {
92  if (m_segmentController) {
93  delete m_segmentController;
94  m_segmentController = NULL;
95  }
96 }
97 
98 void GblFitter::setTrackSegmentController(GblTrackSegmentController* controler)
99 {
100  if (m_segmentController) {
101  delete m_segmentController;
102  m_segmentController = NULL;
103  }
104  m_segmentController = controler;
105 }
106 
107 void GblFitter::processTrackWithRep(Track* trk, const AbsTrackRep* rep, bool resortHits)
108 {
109  cleanGblInfo(trk, rep);
110 
111  if (resortHits)
112  sortHits(trk, rep);
113 
114  // This flag enables/disables fitting of q/p parameter in GBL
115  // It is switched off automatically if no B-field at (0,0,0) is detected.
116  bool fitQoverP = true;
117  //TODO: Use clever way to determine zero B-field
118  double Bfield = genfit::FieldManager::getInstance()->getFieldVal(TVector3(0., 0., 0.)).Mag();
119  if (!(Bfield > 1.e-16))
120  fitQoverP = false;
121  // degrees of freedom after fit
122  int Ndf = 0;
123  // Chi2 after fit
124  double Chi2 = 0.;
125  //FIXME: d-w's not used so far...
126  double lostWeight = 0.;
127 
128  // Preparation of points (+add reference states) for GBL fit
129  // -----------------------------------------------------------------
131  trk->setFitStatus(gblfs, rep);
132  gblfs->setCurvature(fitQoverP);
133  //
134  // Propagate reference seed, create scattering points, calc Jacobians
135  // and store everything in fitter infos. (ready to collect points and fit)
136  //
137  //
138  gblfs->setTrackLen(
139  //
140  constructGblInfo(trk, rep)
141  //
142  );
143  //
144  //
145  gblfs->setIsFittedWithReferenceTrack(true);
146  gblfs->setNumIterations(0); //default value, still valid, No GBL iteration
147  if (m_externalIterations < 1)
148  return;
149  // -----------------------------------------------------------------
150 
151 
152  unsigned int nFailed = 0;
153  int fitRes = 0;
154  std::vector<std::string> gblIterations;
155  gblIterations.push_back(m_gblInternalIterations);
156  //boost::split(gblIterations, m_gblInternalIterations, boost::is_any_of(","), boost::token_compress_off);
157 
158  // Iterations and updates of fitter infos and fit status
159  // -------------------------------------------------------------------
160  for (unsigned int iIter = 0; iIter < m_externalIterations; iIter++) {
161  // GBL refit (1st of reference, then refit of GBL trajectory itself)
162  int nscat = 0, nmeas = 0, ndummy = 0;
163  std::vector<gbl::GblPoint> points = collectGblPoints(trk, rep);
164  for(unsigned int ip = 0;ip<points.size(); ip++) {
165  GblPoint & p = points.at(ip);
166  if (p.hasScatterer())
167  nscat++;
168  if (p.hasMeasurement())
169  nmeas++;
170  if(!p.hasMeasurement()&&!p.hasScatterer())
171  ndummy++;
172  }
173  gbl::GblTrajectory traj(points, gblfs->hasCurvature());
174 
175  fitRes = traj.fit(Chi2, Ndf, lostWeight, (iIter == m_externalIterations - 1) ? m_gblInternalIterations : "");
176 
177  // Update fit results in fitterinfos
178  updateGblInfo(traj, trk, rep);
179 
180  // This repropagates to get new Jacobians,
181  // if planes changed, predictions are extrapolated to new planes
182  if (m_recalcJacobians > iIter) {
183  GblFitterInfo* prevFitterInfo = 0;
184  GblFitterInfo* currFitterInfo = 0;
185  for (unsigned int ip = 0; ip < trk->getNumPoints(); ip++) {
186  if (trk->getPoint(ip)->hasFitterInfo(rep) && (currFitterInfo = dynamic_cast<GblFitterInfo*>(trk->getPoint(ip)->getFitterInfo(rep)))) {
187 
188  currFitterInfo->recalculateJacobian(prevFitterInfo);
189  prevFitterInfo = currFitterInfo;
190  }
191  }
192  }
193 
194  gblfs->setIsFitted(true);
195  gblfs->setIsFitConvergedPartially(fitRes == 0);
196  nFailed = trk->getNumPointsWithMeasurement() - nmeas;
197  gblfs->setNFailedPoints(nFailed);
198  gblfs->setIsFitConvergedFully(fitRes == 0 && nFailed == 0);
199  gblfs->setNumIterations(iIter + 1);
200  gblfs->setChi2(Chi2);
201  gblfs->setNdf(Ndf);
202  gblfs->setCharge(trk->getFittedState().getCharge());
203 
204  #ifdef DEBUG
205  int npoints_meas = trk->getNumPointsWithMeasurement();
206  int npoints_all = trk->getNumPoints();
207 
208  cout << "-------------------------------------------------------" << endl;
209  cout << " GBL processed genfit::Track " << endl;
210  cout << "-------------------------------------------------------" << endl;
211  cout << " # Track Points : " << npoints_all << endl;
212  cout << " # Meas. Points : " << npoints_meas << endl;
213  cout << " # GBL points all : " << traj.getNumPoints();
214  if (ndummy)
215  cout << " (" << ndummy << " dummy) ";
216  cout << endl;
217  cout << " # GBL points meas : " << nmeas << endl;
218  cout << " # GBL points scat : " << nscat << endl;
219  cout << "-------------- GBL Fit Results ----------- Iteration " << iIter+1 << " " << ((iIter < gblIterations.size()) ? gblIterations[iIter] : "") << endl;
220  cout << " Fit q/p parameter : " << (gblfs->hasCurvature() ? ("True") : ("False")) << endl;
221  cout << " Valid trajectory : " << ((traj.isValid()) ? ("True") : ("False")) << endl;
222  cout << " Fit result : " << fitRes << " (0 for success)" << endl;
223  cout << " GBL track NDF : " << Ndf << " (-1 for failure)" << endl;
224  cout << " GBL track Chi2 : " << Chi2 << endl;
225  cout << " GBL track P-value : " << TMath::Prob(Chi2, Ndf) << endl;
226  cout << "-------------------------------------------------------" << endl;
227  #endif
228 
229  }
230  // -------------------------------------------------------------------
231 
232 }
233 
234 void GblFitter::cleanGblInfo(Track* trk, const AbsTrackRep* rep) const {
235 
236  for (int ip = trk->getNumPoints() - 1; ip >=0; ip--) {
237  trk->getPoint(ip)->setScatterer(NULL);
238  trk->getPoint(ip)->deleteFitterInfo(rep);
239  //TODO
240  if (!trk->getPoint(ip)->hasRawMeasurements())
241  trk->deletePoint(ip);
242  }
243 }
244 
245 void GblFitter::sortHits(Track* trk, const AbsTrackRep* rep) const {
246  // All measurement points in ref. track
247  int npoints_meas = trk->getNumPointsWithMeasurement();
248  // Prepare state for extrapolation of track seed
249  StateOnPlane reference(rep);
250  rep->setTime(reference, trk->getTimeSeed());
251  rep->setPosMom(reference, trk->getStateSeed());
252  // Take the state to first plane
253  SharedPlanePtr firstPlane(trk->getPointWithMeasurement(0)->getRawMeasurement(0)->constructPlane(reference));
254  reference.extrapolateToPlane(firstPlane);
255  //1st point is at arc-len=0
256  double arcLenPos = 0;
257 
258  // Loop only between meas. points
259  for (int ipoint_meas = 0; ipoint_meas < npoints_meas - 1; ipoint_meas++) {
260  // current measurement point
261  TrackPoint* point_meas = trk->getPointWithMeasurement(ipoint_meas);
262  // Current detector plane
263  SharedPlanePtr plane = point_meas->getRawMeasurement(0)->constructPlane(reference);
264  // Get the next plane
265  SharedPlanePtr nextPlane(trk->getPointWithMeasurement(ipoint_meas + 1)->getRawMeasurement(0)->constructPlane(reference));
266 
267  point_meas->setSortingParameter(arcLenPos);
268  arcLenPos += reference.extrapolateToPlane(nextPlane);
269 
270  } // end of loop over track points with measurement
271  trk->getPointWithMeasurement(npoints_meas - 1)->setSortingParameter(arcLenPos);
272  trk->sort();
273 }
274 
275 std::vector<gbl::GblPoint> GblFitter::collectGblPoints(genfit::Track* trk, const genfit::AbsTrackRep* rep) {
276  //TODO store collected points in in fit status? need streamer for GblPoint (or something like that)
277  std::vector<gbl::GblPoint> thePoints;
278  thePoints.clear();
279 
280  // Collect points from track and fitterInfo(rep)
281  for (unsigned int ip = 0; ip < trk->getNumPoints(); ip++) {
282  GblFitterInfo * gblfi = dynamic_cast<GblFitterInfo*>(trk->getPoint(ip)->getFitterInfo(rep));
283  if (!gblfi)
284  continue;
285  thePoints.push_back(gblfi->constructGblPoint());
286  }
287  return thePoints;
288 }
289 
290 void GblFitter::updateGblInfo(gbl::GblTrajectory& traj, genfit::Track* trk, const genfit::AbsTrackRep* rep) {
291  //FIXME
292  if (!traj.isValid())
293  return;
294 
295  // Update points in track and fitterInfo(rep)
296  int igblfi = -1;
297  for (unsigned int ip = 0; ip < trk->getNumPoints(); ip++) {
298  GblFitterInfo * gblfi = dynamic_cast<GblFitterInfo*>(trk->getPoint(ip)->getFitterInfo(rep));
299  if (!gblfi)
300  continue;
301  igblfi++;
302 
303  // The point will calculate its position on the track
304  // (counting fitter infos) which hopefully
305  gblfi->updateFitResults(traj);
306 
307  // This is agains logic. User can do this if he wants and it is recommended usually
308  // so that following fit could reuse the updated seed
309  //if (igblfi == 0) {
310  // trk->setStateSeed( gblfi->getFittedState(true).getPos(), gblfi->getFittedState(true).getMom() );
311  // trk->setCovSeed( gblfi->getFittedState(true).get6DCov() );
312  //}
313  }
314 }
315 
317  double& theta, double& s, double& ds,
318  const double p, const double mass, const double charge,
319  const std::vector<genfit::MatStep>& steps) const {
320  theta = 0.; s = 0.; ds = 0.; length = 0;
321  if (steps.empty()) return;
322 
323  // sum of step lengths
324  double len = 0.;
325  // normalization
326  double sumxx = 0.;
327  // first moment (non-normalized)
328  double sumx2x2 = 0.;
329  // (part of) second moment / variance (non-normalized)
330  double sumx3x3 = 0.;
331 
332  double xmin = 0.;
333  double xmax = 0.;
334 
335  for (unsigned int i = 0; i < steps.size(); i++) {
336  const MatStep step = steps.at(i);
337  // inverse of material radiation length ... (in 1/cm) ... "density of scattering"
338  double rho = 1. / step.materialProperties_.getRadLen();
339  len += fabs(step.stepSize_);
340  xmin = xmax;
341  xmax = xmin + fabs(step.stepSize_);
342  // Compute integrals
343 
344  // integral of rho(x)
345  sumxx += rho * (xmax - xmin);
346  // integral of x*rho(x)
347  sumx2x2 += rho * (xmax * xmax - xmin * xmin) / 2.;
348  // integral of x*x*rho(x)
349  sumx3x3 += rho * (xmax * xmax * xmax - xmin * xmin * xmin) / 3.;
350  }
351  // This ensures PDG formula still gives positive results (but sumxx should be > 1e-4 for it to hold)
352  if (sumxx < 1.0e-10) return;
353  // Calculate theta from total sum of radiation length
354  // instead of summimg squares of individual deflection angle variances
355  // PDG formula:
356  double beta = p / sqrt(p * p + mass * mass);
357  theta = (0.0136 / p / beta) * fabs(charge) * sqrt(sumxx) * (1. + 0.038 * log(sumxx));
358  //theta = (0.015 / p / beta) * fabs(charge) * sqrt(sumxx);
359 
360  // track length
361  length = len;
362  // Normalization factor
363  double N = 1. / sumxx;
364  // First moment
365  s = N * sumx2x2;
366  // Square of second moment (variance)
367  // integral of (x - s)*(x - s)*rho(x)
368  double ds_2 = N * (sumx3x3 - 2. * sumx2x2 * s + sumxx * s * s);
369  ds = sqrt(ds_2);
370 
371  #ifdef DEBUG
372  #endif
374 }
375 
376 double GblFitter::constructGblInfo(Track* trk, const AbsTrackRep* rep)
377 {
378  // All measurement points in ref. track
379  int npoints_meas = trk->getNumPointsWithMeasurement();
380  // Dimesion of representation/state
381  int dim = rep->getDim();
382  // Jacobian for point with measurement = how to propagate from previous point (scat/meas)
383  TMatrixD jacPointToPoint(dim, dim);
384  jacPointToPoint.UnitMatrix();
385 
386  // Prepare state for extrapolation of track seed
387  // Take the state to first plane
388  StateOnPlane reference(rep);
389  rep->setTime(reference, trk->getTimeSeed());
390  rep->setPosMom(reference, trk->getStateSeed());
391 
392  SharedPlanePtr firstPlane(trk->getPointWithMeasurement(0)->getRawMeasurement(0)->constructPlane(reference));
393  reference.extrapolateToPlane(firstPlane);
394 
395  double sumTrackLen = 0;
396  // NOT used but useful
397  TMatrixDSym noise; TVectorD deltaState;
398 
399  // Loop only between meas. points
400  for (int ipoint_meas = 0; ipoint_meas < npoints_meas; ipoint_meas++) {
401  // current measurement point
402  TrackPoint* point_meas = trk->getPointWithMeasurement(ipoint_meas);
403  // Current detector plane
404  SharedPlanePtr plane = point_meas->getRawMeasurement(0)->constructPlane(reference);
405  // track direction at plane (in global coords)
406  TVector3 trackDir = rep->getDir(reference);
407  // track momentum direction vector at plane (in global coords)
408  double trackMomMag = rep->getMomMag(reference);
409  // charge of particle
410  double particleCharge = rep->getCharge(reference);
411  // mass of particle
412  double particleMass = rep->getMass(reference);
413  // Parameters of a thick scatterer between measurements
414  double trackLen = 0., scatTheta = 0., scatSMean = 0., scatDeltaS = 0.;
415  // Parameters of two equivalent thin scatterers
416  double theta1 = 0., theta2 = 0., s1 = 0., s2 = 0.;
417  // jacobian from s1=0 to s2
418  TMatrixD jacMeas2Scat(dim, dim);
419  jacMeas2Scat.UnitMatrix();
420 
421  // Stop here if we are at last point (do not add scatterers to last point),
422  // just the fitter info
423  if (ipoint_meas >= npoints_meas - 1) {
424 
425  // Construction last measurement (no scatterer)
426  // --------------------------------------------
427  // Just add the fitter info of last plane
428  GblFitterInfo* gblfimeas(new GblFitterInfo(point_meas, rep, reference));
429  gblfimeas->setJacobian(jacPointToPoint);
430  point_meas->setFitterInfo(gblfimeas);
431  // --------------------------------------------
432 
433  break;
434  }
435  // Extrapolate to next measurement to get material distribution
436  // Use a temp copy of the StateOnPlane to propage between measurements
437  StateOnPlane refCopy(reference);
438  // Get the next plane
439  SharedPlanePtr nextPlane(trk->getPointWithMeasurement(ipoint_meas + 1)->getRawMeasurement(0)->constructPlane(reference));
440 
441  // Extrapolation for multiple scattering calculation
442  // Extrap to point + 1, do NOT stop at boundary
443  TVector3 segmentEntry = refCopy.getPos();
444  rep->extrapolateToPlane(refCopy, nextPlane, false, false);
445  TVector3 segmentExit = refCopy.getPos();
446 
447  getScattererFromMatList(trackLen,
448  scatTheta,
449  scatSMean,
450  scatDeltaS,
451  trackMomMag,
452  particleMass,
453  particleCharge,
454  rep->getSteps());
455  // Now calculate positions and scattering variance for equivalent scatterers
456  // (Solution from Claus Kleinwort (DESY))
457  s1 = 0.; s2 = scatSMean + scatDeltaS * scatDeltaS / (scatSMean - s1);
458  theta1 = sqrt(scatTheta * scatTheta * scatDeltaS * scatDeltaS / (scatDeltaS * scatDeltaS + (scatSMean - s1) * (scatSMean - s1)));
459  theta2 = sqrt(scatTheta * scatTheta * (scatSMean - s1) * (scatSMean - s1) / (scatDeltaS * scatDeltaS + (scatSMean - s1) * (scatSMean - s1)));
460 
461  // Call segment controller to set MS options:
462  if (m_segmentController)
463  m_segmentController->controlTrackSegment(segmentEntry, segmentExit, this);
464 
465  // Scattering options: OFF / THIN / THICK
466  if (m_enableScatterers && !m_enableIntermediateScatterer) {
467  theta1 = scatTheta;
468  theta2 = 0;
469  } else if (!m_enableScatterers) {
470  theta1 = 0.;
471  theta2 = 0.;
472  }
473 
474  // Construction of measurement (with scatterer)
475  // --------------------------------------------
476 
477  if (theta1 > scatEpsilon)
478  point_meas->setScatterer(new ThinScatterer(plane, MaterialProperties(theta1, 0., 0., 0., 0.)));
479 
480  GblFitterInfo* gblfimeas(new GblFitterInfo(point_meas, rep, reference));
481  gblfimeas->setJacobian(jacPointToPoint);
482  point_meas->setFitterInfo(gblfimeas);
483  // --------------------------------------------
484 
485 
486  // If not last measurement, extrapolate and get jacobians for scattering points between this and next measurement
487  if (theta2 > scatEpsilon) {
488  // First scatterer has been placed at current measurement point (see above)
489  // theta2 > 0 ... we want second scatterer:
490  // Extrapolate to s2 (we have s1 = 0)
491  rep->extrapolateBy(reference, s2, false, true);
492  rep->getForwardJacobianAndNoise(jacMeas2Scat, noise, deltaState);
493 
494  // Construction of intermediate scatterer
495  // --------------------------------------
496  TrackPoint* scattp = new TrackPoint(trk);
497  scattp->setSortingParameter(point_meas->getSortingParameter() + s2);
498  scattp->setScatterer(new ThinScatterer(reference.getPlane(), MaterialProperties(theta2, 0., 0., 0., 0.)));
499  // Add point to track before next point
500  int pointIndex = 0;
501  //TODO Deduce this rather than looping over all points
502  for (unsigned int itp = 0; itp < trk->getNumPoints(); itp++) {
503  if (trk->getPoint(itp) == point_meas) {
504  pointIndex = itp;
505  break;
506  }
507  }
508  trk->insertPoint(scattp, pointIndex + 1);
509  // Create and store fitter info
510  GblFitterInfo * gblfiscat(new GblFitterInfo(scattp, rep, reference));
511  gblfiscat->setJacobian(jacMeas2Scat);
512  scattp->setFitterInfo(gblfiscat);
513  // ---------------------------------------
514 
515  // Finish extrapolation to next measurement
516  double nextStep = rep->extrapolateToPlane(reference, nextPlane, false, true);
517  rep->getForwardJacobianAndNoise(jacPointToPoint, noise, deltaState);
518 
519  if (0. > nextStep) {
520  cout << " ERROR: The extrapolation to measurement point " << (ipoint_meas + 2) << " stepped back by " << nextStep << "cm !!! Track will be cut before this point." << endl;
521  // stop trajectory construction here
522  break;
523  }
524 
525  } else {
526  // No scattering: extrapolate whole distance between measurements
527  double nextStep = rep->extrapolateToPlane(reference, nextPlane, false, true);
528  rep->getForwardJacobianAndNoise(jacPointToPoint, noise, deltaState);
529 
530  if (0. > nextStep) {
531  cout << " ERROR: The extrapolation to measurement point " << (ipoint_meas + 2) << " stepped back by " << nextStep << "cm !!! Track will be cut before this point." << endl;
532  // stop trajectory construction here
533  break;
534  }
535  }
536  // Track length up to next point
537  sumTrackLen += trackLen;
538 
539  } // end of loop over track points with measurement
540  return sumTrackLen;
541 }
genfit::AbsTrackRep::getMomMag
virtual double getMomMag(const StateOnPlane &state) const =0
get the magnitude of the momentum in GeV.
genfit::GblFitterInfo
Collects information needed and produced by a GblFitter/GBL and is specific to one AbsTrackRep of the...
Definition: GblFitterInfo.h:55
genfit::AbsTrackRep::getSteps
virtual std::vector< genfit::MatStep > getSteps() const =0
Get stepsizes and material properties of crossed materials of the last extrapolation.
genfit::SharedPlanePtr
boost::shared_ptr< genfit::DetPlane > SharedPlanePtr
Shared Pointer to a DetPlane.
Definition: SharedPlanePtr.h:43
genfit::AbsTrackRep::extrapolateToPlane
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,...
genfit::MatStep::stepSize_
double stepSize_
Definition: AbsTrackRep.h:44
genfit::Track
Collection of TrackPoint objects, AbsTrackRep objects and FitStatus objects.
Definition: Track.h:71
genfit::GblFitterInfo::recalculateJacobian
void recalculateJacobian(GblFitterInfo *prevFitterInfo)
Re-extrapolates between prevFitterInfo and this point using forward state to update the Jacobian (if ...
Definition: GblFitterInfo.cc:345
getScattererFromMatList
void getScattererFromMatList(double &length, double &theta, double &s, double &ds, const double p, const double mass, const double charge, const std::vector< MatStep > &steps)
Evaluates moments of radiation length distribution from list of material steps and computes parameter...
Definition: GFGbl.cc:263
genfit::GblTrackSegmentController
TrackSegmentController for use with GblFitter.
Definition: GblTrackSegmentController.h:40
genfit::Track::deletePoint
void deletePoint(int id)
Definition: Track.cc:473
genfit::TrackPoint::deleteFitterInfo
void deleteFitterInfo(const AbsTrackRep *rep)
Definition: TrackPoint.h:117
genfit::Track::getTimeSeed
double getTimeSeed() const
Definition: Track.h:161
genfit::GblFitterInfo::updateFitResults
void updateFitResults(gbl::GblTrajectory &traj)
Update fitter info from GBL fit results.
Definition: GblFitterInfo.cc:245
ICalibrationParametersDerivatives.h
genfit::StateOnPlane::getCharge
double getCharge() const
Definition: StateOnPlane.h:119
gbl::GblPoint::hasMeasurement
unsigned int hasMeasurement() const
Check for measurement at a point.
Definition: GblPoint.cc:150
genfit::Track::setFitStatus
void setFitStatus(FitStatus *fitStatus, const AbsTrackRep *rep)
Definition: Track.cc:341
genfit
Defines for I/O streams used for error and debug printing.
Definition: AbsFinitePlane.cc:22
genfit::GblFitStatus::setTrackLen
void setTrackLen(double trackLen)
Definition: GblFitStatus.h:63
GblPoint.h
genfit::FieldManager::getFieldVal
TVector3 getFieldVal(const TVector3 &position)
This does NOT use the cache!
Definition: FieldManager.h:63
genfit::Track::getPoint
TrackPoint * getPoint(int id) const
Definition: Track.cc:212
genfit::FitStatus::setNFailedPoints
void setNFailedPoints(int nFailedPoints)
Definition: FitStatus.h:131
genfit::AbsTrackRep::getMass
double getMass(const StateOnPlane &state) const
Get tha particle mass in GeV/c^2.
Definition: AbsTrackRep.cc:100
genfit::FitStatus::setNdf
void setNdf(const double &ndf)
Definition: FitStatus.h:138
genfit::AbsTrackRep::setTime
virtual void setTime(StateOnPlane &state, double time) const =0
Set time at which the state was defined.
genfit::FitStatus::setIsFitted
void setIsFitted(bool fitted=true)
Definition: FitStatus.h:128
genfit::AbsTrackRep
Abstract base class for a track representation.
Definition: AbsTrackRep.h:66
genfit::TrackPoint::setSortingParameter
void setSortingParameter(double sortingParameter)
Definition: TrackPoint.h:111
genfit::GblFitterInfo::setJacobian
void setJacobian(TMatrixD jacobian)
Set the Jacobian for further GblPoint construction.
Definition: GblFitterInfo.cc:91
genfit::GblFitStatus
FitStatus for use with GblFitter.
Definition: GblFitStatus.h:39
gbl::GblTrajectory::isValid
bool isValid() const
Retrieve validity of trajectory.
Definition: GblTrajectory.cc:219
gbl::GblPoint
Point on trajectory.
Definition: GblPoint.h:48
gbl
Namespace for the general broken lines package.
Definition: BorderedBandMatrix.cc:11
genfit::GblFitStatus::setNumIterations
void setNumIterations(unsigned int numIterations)
Definition: GblFitStatus.h:61
genfit::MatStep
Simple struct containing MaterialProperties and stepsize in the material.
Definition: AbsTrackRep.h:42
genfit::StateOnPlane::extrapolateToPlane
double extrapolateToPlane(const SharedPlanePtr &plane, bool stopAtBoundary=false, bool calcJacobianNoise=false)
Definition: StateOnPlane.h:75
genfit::TrackPoint::setFitterInfo
void setFitterInfo(genfit::AbsFitterInfo *fitterInfo)
Takes Ownership.
Definition: TrackPoint.cc:193
genfit::ThinScatterer
Thin or thick scatterer.
Definition: ThinScatterer.h:38
genfit::TrackPoint::hasFitterInfo
bool hasFitterInfo(const AbsTrackRep *rep) const
Definition: TrackPoint.h:103
genfit::Track::insertPoint
void insertPoint(TrackPoint *point, int id=-1)
Insert TrackPoint BEFORE TrackPoint with position id, if id >= 0.
Definition: Track.cc:363
genfit::Track::getPointWithMeasurement
TrackPoint * getPointWithMeasurement(int id) const
Definition: Track.cc:220
HMatrixV.h
genfit::TrackPoint
Object containing AbsMeasurement and AbsFitterInfo objects.
Definition: TrackPoint.h:50
genfit::StateOnPlane
A state with arbitrary dimension defined in a DetPlane.
Definition: StateOnPlane.h:45
genfit::GblFitStatus::setCurvature
void setCurvature(bool useCurvature)
Definition: GblFitStatus.h:51
genfit::TrackPoint::hasRawMeasurements
bool hasRawMeasurements() const
Definition: TrackPoint.h:96
GblFitter.h
genfit::AbsMeasurement::constructPlane
virtual SharedPlanePtr constructPlane(const StateOnPlane &state) const =0
HMatrixU.h
genfit::AbsTrackRep::getDir
TVector3 getDir(const StateOnPlane &state) const
Get the direction vector of a state.
Definition: AbsTrackRep.h:246
genfit::TrackPoint::getRawMeasurement
AbsMeasurement * getRawMeasurement(int i=0) const
Definition: TrackPoint.cc:147
genfit::AbsTrackRep::extrapolateBy
virtual double extrapolateBy(StateOnPlane &state, double step, bool stopAtBoundary=false, bool calcJacobianNoise=false) const =0
Extrapolates the state by step (cm) and returns the extrapolation length and, via reference,...
genfit::FitStatus::setCharge
void setCharge(double charge)
Definition: FitStatus.h:133
genfit::GblFitStatus::setIsFittedWithReferenceTrack
void setIsFittedWithReferenceTrack(bool fittedWithReferenceTrack=true)
Definition: GblFitStatus.h:62
genfit::GblFitStatus::hasCurvature
bool hasCurvature()
Definition: GblFitStatus.h:52
genfit::AbsTrackRep::getCharge
virtual double getCharge(const StateOnPlane &state) const =0
Get the (fitted) charge of a state. This is not always equal the pdg charge (e.g. if the charge sign ...
genfit::MaterialProperties::getRadLen
double getRadLen() const
Definition: MaterialProperties.h:54
genfit::FitStatus::setChi2
void setChi2(const double &chi2)
Definition: FitStatus.h:137
genfit::Track::getNumPointsWithMeasurement
unsigned int getNumPointsWithMeasurement() const
Definition: Track.h:114
scatEpsilon
const double scatEpsilon
Definition: GFGbl.cc:168
genfit::StateOnPlane::getPos
TVector3 getPos() const
Definition: StateOnPlane.h:111
genfit::GblFitterInfo::constructGblPoint
gbl::GblPoint constructGblPoint()
Collect all data and create a GblPoint.
Definition: GblFitterInfo.cc:111
genfit::StateOnPlane::getPlane
const SharedPlanePtr & getPlane() const
Definition: StateOnPlane.h:65
genfit::MatStep::materialProperties_
MaterialProperties materialProperties_
Definition: AbsTrackRep.h:43
genfit::AbsTrackRep::getForwardJacobianAndNoise
virtual void getForwardJacobianAndNoise(TMatrixD &jacobian, TMatrixDSym &noise, TVectorD &deltaState) const =0
Get the jacobian and noise matrix of the last extrapolation.
genfit::AbsTrackRep::setPosMom
virtual void setPosMom(StateOnPlane &state, const TVector3 &pos, const TVector3 &mom) const =0
Set position and momentum of state.
genfit::FieldManager::getInstance
static FieldManager * getInstance()
Get singleton instance.
Definition: FieldManager.h:112
genfit::FitStatus::setIsFitConvergedPartially
void setIsFitConvergedPartially(bool fitConverged=true)
Definition: FitStatus.h:130
Track.h
GblFitterInfo.h
genfit::Track::getStateSeed
const TVectorD & getStateSeed() const
Definition: Track.h:164
FieldManager.h
genfit::FitStatus::setIsFitConvergedFully
void setIsFitConvergedFully(bool fitConverged=true)
Definition: FitStatus.h:129
genfit::AbsTrackRep::getDim
virtual unsigned int getDim() const =0
Get the dimension of the state vector used by the track representation.
genfit::Track::getFittedState
const MeasuredStateOnPlane & getFittedState(int id=0, const AbsTrackRep *rep=NULL, bool biased=true) const
Shortcut to get FittedStates.
Definition: Track.cc:288
gbl::GblTrajectory
GBL trajectory.
Definition: GblTrajectory.h:26
gbl::GblTrajectory::fit
unsigned int fit(double &Chi2, int &Ndf, double &lostWeight, std::string optionList="")
Perform fit of trajectory.
Definition: GblTrajectory.cc:955
GblTrajectory.h
genfit::TrackPoint::setScatterer
void setScatterer(ThinScatterer *scatterer)
Definition: TrackPoint.h:119
genfit::Track::sort
bool sort()
Sort TrackPoint and according to their sorting parameters.
Definition: Track.cc:634
genfit::TrackPoint::getFitterInfo
AbsFitterInfo * getFitterInfo(const AbsTrackRep *rep=NULL) const
Get fitterInfo for rep. Per default, use cardinal rep.
Definition: TrackPoint.cc:169
genfit::Track::getNumPoints
unsigned int getNumPoints() const
Definition: Track.h:110
genfit::TrackPoint::getSortingParameter
double getSortingParameter() const
Definition: TrackPoint.h:88
genfit::MaterialProperties
Material properties needed e.g. for material effects calculation.
Definition: MaterialProperties.h:35
gbl::GblPoint::hasScatterer
bool hasScatterer() const
Check for scatterer at a point.
Definition: GblPoint.cc:232