8 #ifndef GBLTRAJECTORY_H_
9 #define GBLTRAJECTORY_H_
16 #include "TMatrixDSymEigen.h"
28 GblTrajectory(
const std::vector<GblPoint> &aPointList,
bool flagCurv =
true,
29 bool flagU1dir =
true,
bool flagU2dir =
true);
30 GblTrajectory(
const std::vector<GblPoint> &aPointList,
unsigned int aLabel,
31 const TMatrixDSym &aSeed,
bool flagCurv =
true,
bool flagU1dir =
32 true,
bool flagU2dir =
true);
34 const std::vector<std::pair<std::vector<GblPoint>, TMatrixD> > &aPointaAndTransList);
36 const std::vector<std::pair<std::vector<GblPoint>, TMatrixD> > &aPointaAndTransList,
37 const TMatrixD &extDerivatives,
const TVectorD &extMeasurements,
38 const TVectorD &extPrecisions);
42 unsigned int getResults(
int aSignedLabel, TVectorD &localPar,
43 TMatrixDSym &localCov)
const;
44 unsigned int getMeasResults(
unsigned int aLabel,
unsigned int &numRes,
45 TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors,
46 TVectorD &aDownWeights);
47 unsigned int getScatResults(
unsigned int aLabel,
unsigned int &numRes,
48 TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors,
49 TVectorD &aDownWeights);
50 void getLabels(std::vector<unsigned int> &aLabelList);
51 void getLabels(std::vector<std::vector< unsigned int> > &aLabelList);
52 unsigned int fit(
double &Chi2,
int &Ndf,
double &lostWeight,
53 std::string optionList =
"");
86 std::pair<std::vector<unsigned int>, TMatrixD>
getJacobian(
87 int aSignedLabel)
const;
90 unsigned int nJacobian = 1)
const;
100 void getResAndErr(
unsigned int aData,
double &aResidual,
101 double &aMeadsError,
double &aResError,
double &aDownWeight);
void construct()
Construct trajectory from list of points.
TVectorD externalMeasurements
ROOT::Math::SMatrix< double, 5, 5 > SMatrix55
unsigned int numInnerTrans
Number of inner transformations to external parameters.
TMatrixD externalDerivatives
(Symmetric) Bordered Band Matrix.
Millepede-II (binary) record.
double downWeight(unsigned int aMethod)
Down-weight all points.
TVectorD externalPrecisions
unsigned int getMeasResults(unsigned int aLabel, unsigned int &numRes, TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors, TVectorD &aDownWeights)
Get residuals at point from measurement.
std::vector< std::vector< GblPoint > > thePoints
(list of) List of points on trajectory
std::vector< unsigned int > theDimension
List of active dimensions (0=u1, 1=u2) in fit.
unsigned int numAllPoints
Number of all points on trajectory.
void prepare()
Prepare fit for simple or composed trajectory.
void defineOffsets()
Define offsets from list of points.
std::vector< unsigned int > numPoints
Number of points on (sub)trajectory.
ROOT::Math::SMatrix< double, 2, 7 > SMatrix27
void printData()
Print GblData blocks for trajectory.
TMatrixDSym externalSeed
Precision (inverse covariance matrix) of external seed.
void predict()
Calculate predictions for all points.
unsigned int numLocals
Total number of (additional) local parameters.
Simple Vector based on std::vector
unsigned int getNumPoints() const
Retrieve number of point from trajectory.
GblTrajectory(const std::vector< GblPoint > &aPointList, bool flagCurv=true, bool flagU1dir=true, bool flagU2dir=true)
Create new (simple) trajectory from list of points.
void printPoints(unsigned int level=0)
Print GblPoints on trajectory.
BorderedBandMatrix theMatrix
(Bordered band) matrix of linear equation system
std::vector< unsigned int > measDataIndex
mapping points to data blocks from measurements
unsigned int fit(double &Chi2, int &Ndf, double &lostWeight, std::string optionList="")
Perform fit of trajectory.
Namespace for the general broken lines package.
bool constructOK
Trajectory has been successfully constructed (ready for fit/output)
unsigned int getResults(int aSignedLabel, TVectorD &localPar, TMatrixDSym &localCov) const
Get fit results at point.
unsigned int numTrajectories
Number of trajectories (in composed trajectory)
void getFitToKinkJacobian(std::vector< unsigned int > &anIndex, SMatrix27 &aJacobian, const GblPoint &aPoint) const
Get jacobian for transformation from (trajectory) fit to kink parameters at point.
unsigned int numParameters
Number of fit parameters.
unsigned int externalPoint
Label of external point (or 0)
unsigned int getScatResults(unsigned int aLabel, unsigned int &numRes, TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors, TVectorD &aDownWeights)
Get (kink) residuals at point from scatterer.
bool fitOK
Trajectory has been successfully fitted (results are valid)
std::vector< unsigned int > scatDataIndex
mapping points to data blocks from scatterers
unsigned int numMeasurements
Total number of measurements.
void printTrajectory(unsigned int level=0)
Print GblTrajectory.
void getLabels(std::vector< unsigned int > &aLabelList)
Get (list of) labels of points on (simple) trajectory.
void buildLinearEquationSystem()
Build linear equation system from data (blocks).
void calcJacobians()
Calculate Jacobians to previous/next scatterer from point to point ones.
unsigned int numOffsets
Number of (points with) offsets on trajectory.
void getFitToLocalJacobian(std::vector< unsigned int > &anIndex, SMatrix55 &aJacobian, const GblPoint &aPoint, unsigned int measDim, unsigned int nJacobian=1) const
Get (part of) jacobian for transformation from (trajectory) fit to track parameters at point...
std::vector< GblData > theData
List of data blocks.
std::pair< std::vector< unsigned int >, TMatrixD > getJacobian(int aSignedLabel) const
Get jacobian for transformation from fit to track parameters at point.
VVector theVector
Vector of linear equation system.
unsigned int numCurvature
Number of curvature parameters (0 or 1) or external parameters.
bool isValid() const
Retrieve validity of trajectory.
void milleOut(MilleBinary &aMille)
Write valid trajectory to Millepede-II binary file.
std::vector< TMatrixD > innerTransformations
Transformations at innermost points of.
void getResAndErr(unsigned int aData, double &aResidual, double &aMeadsError, double &aResError, double &aDownWeight)
Get residual and errors from data block.