gbl::GblTrajectory Class Reference

GBL trajectory. More...

#include <GblTrajectory.h>

Public Member Functions
unsigned int	fit (double &Chi2, int &Ndf, double &lostWeight, std::string optionList="")
	Perform fit of trajectory. More...
	GblTrajectory (const std::vector< GblPoint > &aPointList, bool flagCurv=true, bool flagU1dir=true, bool flagU2dir=true)
	Create new (simple) trajectory from list of points. More...
	GblTrajectory (const std::vector< GblPoint > &aPointList, unsigned int aLabel, const TMatrixDSym &aSeed, bool flagCurv=true, bool flagU1dir=true, bool flagU2dir=true)
	Create new (simple) trajectory from list of points with external seed. More...
	GblTrajectory (const std::vector< std::pair< std::vector< GblPoint >, TMatrixD > > &aPointaAndTransList)
	Create new composed trajectory from list of points and transformations. More...
	GblTrajectory (const std::vector< std::pair< std::vector< GblPoint >, TMatrixD > > &aPointaAndTransList, const TMatrixD &extDerivatives, const TVectorD &extMeasurements, const TVectorD &extPrecisions)
	Create new composed trajectory from list of points and transformations with (independent) external measurements. More...
void	getLabels (std::vector< unsigned int > &aLabelList)
	Get (list of) labels of points on (simple) trajectory. More...
void	getLabels (std::vector< std::vector< unsigned int > > &aLabelList)
	Get (list of lists of) labels of points on (composed) trajectory. More...
unsigned int	getMeasResults (unsigned int aLabel, unsigned int &numRes, TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors, TVectorD &aDownWeights)
	Get residuals at point from measurement. More...
unsigned int	getNumPoints () const
	Retrieve number of point from trajectory. More...
unsigned int	getResults (int aSignedLabel, TVectorD &localPar, TMatrixDSym &localCov) const
	Get fit results at point. More...
unsigned int	getScatResults (unsigned int aLabel, unsigned int &numRes, TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors, TVectorD &aDownWeights)
	Get (kink) residuals at point from scatterer. More...
bool	isValid () const
	Retrieve validity of trajectory. More...
void	milleOut (MilleBinary &aMille)
	Write valid trajectory to Millepede-II binary file. More...
void	printData ()
	Print GblData blocks for trajectory. More...
void	printPoints (unsigned int level=0)
	Print GblPoints on trajectory. More...
void	printTrajectory (unsigned int level=0)
	Print GblTrajectory. More...
virtual	~GblTrajectory ()

Private Member Functions
void	buildLinearEquationSystem ()
	Build linear equation system from data (blocks). More...
void	calcJacobians ()
	Calculate Jacobians to previous/next scatterer from point to point ones. More...
void	construct ()
	Construct trajectory from list of points. More...
void	defineOffsets ()
	Define offsets from list of points. More...
double	downWeight (unsigned int aMethod)
	Down-weight all points. More...
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. More...
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. More...
std::pair< std::vector< unsigned int >, TMatrixD >	getJacobian (int aSignedLabel) const
	Get jacobian for transformation from fit to track parameters at point. More...
void	getResAndErr (unsigned int aData, double &aResidual, double &aMeadsError, double &aResError, double &aDownWeight)
	Get residual and errors from data block. More...
void	predict ()
	Calculate predictions for all points. More...
void	prepare ()
	Prepare fit for simple or composed trajectory. More...

Private Attributes
bool	constructOK
	Trajectory has been successfully constructed (ready for fit/output) More...
TMatrixD	externalDerivatives
TVectorD	externalMeasurements
unsigned int	externalPoint
	Label of external point (or 0) More...
TVectorD	externalPrecisions
TMatrixDSym	externalSeed
	Precision (inverse covariance matrix) of external seed. More...
bool	fitOK
	Trajectory has been successfully fitted (results are valid) More...
std::vector< TMatrixD >	innerTransformations
	Transformations at innermost points of. More...
std::vector< unsigned int >	measDataIndex
	mapping points to data blocks from measurements More...
unsigned int	numAllPoints
	Number of all points on trajectory. More...
unsigned int	numCurvature
	Number of curvature parameters (0 or 1) or external parameters. More...
unsigned int	numInnerTrans
	Number of inner transformations to external parameters. More...
unsigned int	numLocals
	Total number of (additional) local parameters. More...
unsigned int	numMeasurements
	Total number of measurements. More...
unsigned int	numOffsets
	Number of (points with) offsets on trajectory. More...
unsigned int	numParameters
	Number of fit parameters. More...
std::vector< unsigned int >	numPoints
	Number of points on (sub)trajectory. More...
unsigned int	numTrajectories
	Number of trajectories (in composed trajectory) More...
std::vector< unsigned int >	scatDataIndex
	mapping points to data blocks from scatterers More...
std::vector< GblData >	theData
	List of data blocks. More...
std::vector< unsigned int >	theDimension
	List of active dimensions (0=u1, 1=u2) in fit. More...
BorderedBandMatrix	theMatrix
	(Bordered band) matrix of linear equation system More...
std::vector< std::vector< GblPoint > >	thePoints
	(list of) List of points on trajectory More...
VVector	theVector
	Vector of linear equation system. More...

Detailed Description

GBL trajectory.

List of GblPoints ordered by arc length. Can be fitted and optionally written to MP-II binary file.

Definition at line 26 of file GblTrajectory.h.

Constructor & Destructor Documentation

gbl::GblTrajectory::GblTrajectory	(	const std::vector< GblPoint > &	aPointList,
		bool	flagCurv = true,
		bool	flagU1dir = true,
		bool	flagU2dir = true
	)

Create new (simple) trajectory from list of points.

Curved trajectory in space (default) or without curvature (q/p) or in one plane (u-direction) only.

Parameters

[in]	aPointList	List of points
[in]	flagCurv	Use q/p
[in]	flagU1dir	Use in u1 direction
[in]	flagU2dir	Use in u2 direction

Definition at line 117 of file GblTrajectory.cc.

References construct(), numAllPoints, numPoints, theDimension, and thePoints.

                                                       :
        numAllPoints(aPointList.size()), numPoints(), numOffsets(0), numInnerTrans(
                0), numCurvature(flagCurv ? 1 : 0), numParameters(0), numLocals(
                0), numMeasurements(0), externalPoint(0), theDimension(0), thePoints(), theData(), measDataIndex(), scatDataIndex(), externalSeed(), innerTransformations(), externalDerivatives(), externalMeasurements(), externalPrecisions() {

    if (flagU1dir)
        theDimension.push_back(0);
    if (flagU2dir)
        theDimension.push_back(1);
    // simple (single) trajectory
    thePoints.push_back(aPointList);
    numPoints.push_back(numAllPoints);
    construct(); // construct trajectory
}

gbl::GblTrajectory::GblTrajectory	(	const std::vector< GblPoint > &	aPointList,
		unsigned int	aLabel,
		const TMatrixDSym &	aSeed,
		bool	flagCurv = true,
		bool	flagU1dir = true,
		bool	flagU2dir = true
	)

Create new (simple) trajectory from list of points with external seed.

Curved trajectory in space (default) or without curvature (q/p) or in one plane (u-direction) only.

Parameters

[in]	aPointList	List of points
[in]	aLabel	(Signed) label of point for external seed (<0: in front, >0: after point, slope changes at scatterer!)
[in]	aSeed	Precision matrix of external seed
[in]	flagCurv	Use q/p
[in]	flagU1dir	Use in u1 direction
[in]	flagU2dir	Use in u2 direction

Definition at line 145 of file GblTrajectory.cc.

References construct(), numAllPoints, numPoints, theDimension, and thePoints.

                                        :
        numAllPoints(aPointList.size()), numPoints(), numOffsets(0), numInnerTrans(
                0), numCurvature(flagCurv ? 1 : 0), numParameters(0), numLocals(
                0), numMeasurements(0), externalPoint(aLabel), theDimension(0), thePoints(), theData(), measDataIndex(), scatDataIndex(), externalSeed(
                aSeed), innerTransformations(), externalDerivatives(), externalMeasurements(), externalPrecisions() {

    if (flagU1dir)
        theDimension.push_back(0);
    if (flagU2dir)
        theDimension.push_back(1);
    // simple (single) trajectory
    thePoints.push_back(aPointList);
    numPoints.push_back(numAllPoints);
    construct(); // construct trajectory
}

gbl::GblTrajectory::GblTrajectory

(

const std::vector< std::pair< std::vector< GblPoint >, TMatrixD > > &

aPointsAndTransList

)

Create new composed trajectory from list of points and transformations.

Composed of curved trajectory in space.

Parameters

[in]

aPointsAndTransList

List containing pairs with list of points and transformation (at inner (first) point)

Definition at line 168 of file GblTrajectory.cc.

References construct(), plotBeamSpotDB::first, innerTransformations, numAllPoints, numCurvature, numPoints, edm::second(), theDimension, and thePoints.

                                                                                       :
        numAllPoints(), numPoints(), numOffsets(0), numInnerTrans(
                aPointsAndTransList.size()), numParameters(0), numLocals(0), numMeasurements(
                0), externalPoint(0), theDimension(0), thePoints(), theData(), measDataIndex(), scatDataIndex(), externalSeed(), innerTransformations(), externalDerivatives(), externalMeasurements(), externalPrecisions() {

    for (unsigned int iTraj = 0; iTraj < aPointsAndTransList.size(); ++iTraj) {
        thePoints.push_back(aPointsAndTransList[iTraj].first);
        numPoints.push_back(thePoints.back().size());
        numAllPoints += numPoints.back();
        innerTransformations.push_back(aPointsAndTransList[iTraj].second);
    }
    theDimension.push_back(0);
    theDimension.push_back(1);
    numCurvature = innerTransformations[0].GetNcols();
    construct(); // construct (composed) trajectory
}

gbl::GblTrajectory::GblTrajectory	(	const std::vector< std::pair< std::vector< GblPoint >, TMatrixD > > &	aPointsAndTransList,
		const TMatrixD &	extDerivatives,
		const TVectorD &	extMeasurements,
		const TVectorD &	extPrecisions
	)

Create new composed trajectory from list of points and transformations with (independent) external measurements.

Composed of curved trajectory in space.

Parameters

[in]	aPointsAndTransList	List containing pairs with list of points and transformation (at inner (first) point)
[in]	extDerivatives	Derivatives of external measurements vs external parameters
[in]	extMeasurements	External measurements (residuals)
[in]	extPrecisions	Precision of external measurements

Definition at line 194 of file GblTrajectory.cc.

References construct(), plotBeamSpotDB::first, innerTransformations, numAllPoints, numCurvature, numPoints, edm::second(), theDimension, and thePoints.

                                       :
        numAllPoints(), numPoints(), numOffsets(0), numInnerTrans(
                aPointsAndTransList.size()), numParameters(0), numLocals(0), numMeasurements(
                0), externalPoint(0), theDimension(0), thePoints(), theData(), measDataIndex(), scatDataIndex(), externalSeed(), innerTransformations(), externalDerivatives(
                extDerivatives), externalMeasurements(extMeasurements), externalPrecisions(
                extPrecisions) {

    for (unsigned int iTraj = 0; iTraj < aPointsAndTransList.size(); ++iTraj) {
        thePoints.push_back(aPointsAndTransList[iTraj].first);
        numPoints.push_back(thePoints.back().size());
        numAllPoints += numPoints.back();
        innerTransformations.push_back(aPointsAndTransList[iTraj].second);
    }
    theDimension.push_back(0);
    theDimension.push_back(1);
    numCurvature = innerTransformations[0].GetNcols();
    construct(); // construct (composed) trajectory
}

gbl::GblTrajectory::~GblTrajectory ( )

virtual

Definition at line 216 of file GblTrajectory.cc.

                              {
}

Member Function Documentation

void gbl::GblTrajectory::buildLinearEquationSystem ( )

private

Build linear equation system from data (blocks).

Definition at line 692 of file GblTrajectory.cc.

References gbl::BorderedBandMatrix::addBlockMatrix(), j, numCurvature, numLocals, numParameters, gbl::VVector::resize(), gbl::BorderedBandMatrix::resize(), theData, theMatrix, and theVector.

Referenced by fit().

                                              {
    unsigned int nBorder = numCurvature + numLocals;
    theVector.resize(numParameters);
    theMatrix.resize(numParameters, nBorder);
    double aValue, aWeight;
    std::vector<unsigned int>* indLocal;
    std::vector<double>* derLocal;
    std::vector<GblData>::iterator itData;
    for (itData = theData.begin(); itData < theData.end(); ++itData) {
        itData->getLocalData(aValue, aWeight, indLocal, derLocal);
        for (unsigned int j = 0; j < indLocal->size(); ++j) {
            theVector((*indLocal)[j] - 1) += (*derLocal)[j] * aWeight * aValue;
        }
        theMatrix.addBlockMatrix(aWeight, indLocal, derLocal);
    }
}

void gbl::GblTrajectory::calcJacobians ( )

private

Calculate Jacobians to previous/next scatterer from point to point ones.

Definition at line 297 of file GblTrajectory.cc.

References gbl::GblPoint::addNextJacobian(), begin, end, gbl::GblPoint::getP2pJacobian(), numTrajectories, and thePoints.

Referenced by construct().

                                  {

    SMatrix55 scatJacobian;
    // loop over trajectories
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        // forward propagation (all)
        GblPoint* previousPoint = &thePoints[iTraj].front();
        unsigned int numStep = 0;
        std::vector<GblPoint>::iterator itPoint;
        for (itPoint = thePoints[iTraj].begin() + 1;
                itPoint < thePoints[iTraj].end(); ++itPoint) {
            if (numStep == 0) {
                scatJacobian = itPoint->getP2pJacobian();
            } else {
                scatJacobian = itPoint->getP2pJacobian() * scatJacobian;
            }
            numStep++;
            itPoint->addPrevJacobian(scatJacobian); // iPoint -> previous scatterer
            if (itPoint->getOffset() >= 0) {
                previousPoint->addNextJacobian(scatJacobian); // lastPoint -> next scatterer
                numStep = 0;
                previousPoint = &(*itPoint);
            }
        }
        // backward propagation (without scatterers)
        for (itPoint = thePoints[iTraj].end() - 1;
                itPoint > thePoints[iTraj].begin(); --itPoint) {
            if (itPoint->getOffset() >= 0) {
                scatJacobian = itPoint->getP2pJacobian();
                continue; // skip offsets
            }
            itPoint->addNextJacobian(scatJacobian); // iPoint -> next scatterer
            scatJacobian = scatJacobian * itPoint->getP2pJacobian();
        }
    }
}

void gbl::GblTrajectory::construct ( )

private

Construct trajectory from list of points.

Trajectory is prepared for fit or output to binary file, may consists of sub-trajectories.

Definition at line 233 of file GblTrajectory.cc.

References begin, calcJacobians(), constructOK, gather_cfg::cout, defineOffsets(), MillePedeFileConverter_cfg::e, fitOK, hpstanc_transforms::max, numAllPoints, numCurvature, numInnerTrans, numLocals, numMeasurements, numOffsets, numParameters, numTrajectories, prepare(), theDimension, and thePoints.

Referenced by GblTrajectory().

                              {

    constructOK = false;
    fitOK = false;
    unsigned int aLabel = 0;
    if (numAllPoints < 2) {
        std::cout << " GblTrajectory construction failed: too few GblPoints "
                << std::endl;
        return;
    }
    // loop over trajectories
    numTrajectories = thePoints.size();
    //std::cout << " numTrajectories: " << numTrajectories << ", " << innerTransformations.size() << std::endl;
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        std::vector<GblPoint>::iterator itPoint;
        for (itPoint = thePoints[iTraj].begin();
                itPoint < thePoints[iTraj].end(); ++itPoint) {
            numLocals = std::max(numLocals, itPoint->getNumLocals());
            numMeasurements += itPoint->hasMeasurement();
            itPoint->setLabel(++aLabel);
        }
    }
    defineOffsets();
    calcJacobians();
    try {
        prepare();
    } catch (std::overflow_error &e) {
        std::cout << " GblTrajectory construction failed: " << e.what()
                << std::endl;
        return;
    }
    constructOK = true;
    // number of fit parameters
    numParameters = (numOffsets - 2 * numInnerTrans) * theDimension.size()
            + numCurvature + numLocals;
}

void gbl::GblTrajectory::defineOffsets ( )

private

Define offsets from list of points.

Define offsets at points with scatterers and first and last point. All other points need interpolation from adjacent points with offsets.

Definition at line 275 of file GblTrajectory.cc.

References begin, numOffsets, numTrajectories, and thePoints.

Referenced by construct().

                                  {

    // loop over trajectories
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        // first point is offset
        thePoints[iTraj].front().setOffset(numOffsets++);
        // intermediate scatterers are offsets
        std::vector<GblPoint>::iterator itPoint;
        for (itPoint = thePoints[iTraj].begin() + 1;
                itPoint < thePoints[iTraj].end() - 1; ++itPoint) {
            if (itPoint->hasScatterer()) {
                itPoint->setOffset(numOffsets++);
            } else {
                itPoint->setOffset(-numOffsets);
            }
        }
        // last point is offset
        thePoints[iTraj].back().setOffset(numOffsets++);
    }
}

double gbl::GblTrajectory::downWeight ( unsigned int  aMethod )

private

Down-weight all points.

Parameters

[in]

aMethod

M-estimator (1: Tukey, 2:Huber, 3:Cauchy)

Definition at line 937 of file GblTrajectory.cc.

References theData.

Referenced by fit().

                                                     {
    double aLoss = 0.;
    std::vector<GblData>::iterator itData;
    for (itData = theData.begin(); itData < theData.end(); ++itData) {
        aLoss += (1. - itData->setDownWeighting(aMethod));
    }
    return aLoss;
}

unsigned int gbl::GblTrajectory::fit	(	double &	Chi2,
		int &	Ndf,
		double &	lostWeight,
		std::string	optionList = ""
	)

Perform fit of trajectory.

Optionally iterate for outlier down-weighting.

Parameters

[out]	Chi2	Chi2 sum (corrected for down-weighting)
[out]	Ndf	Number of degrees of freedom
[out]	lostWeight	Sum of weights lost due to down-weighting
[in]	optionList	Iterations for down-weighting (One character per iteration: t,h,c (or T,H,C) for Tukey, Huber or Cauchy function)

Returns

Error code (non zero value indicates failure of fit)

Definition at line 956 of file GblTrajectory.cc.

References buildLinearEquationSystem(), constructOK, gather_cfg::cout, downWeight(), MillePedeFileConverter_cfg::e, fitOK, i, numParameters, predict(), gbl::BorderedBandMatrix::solveAndInvertBorderedBand(), AlCaHLTBitMon_QueryRunRegistry::string, theData, theMatrix, and theVector.

Referenced by trackingPlots.Iteration::modules().

                              {
    const double normChi2[4] = { 1.0, 0.8737, 0.9326, 0.8228 };
    const std::string methodList = "TtHhCc";

    Chi2 = 0.;
    Ndf = -1;
    lostWeight = 0.;
    if (not constructOK)
        return 10;

    unsigned int aMethod = 0;

    buildLinearEquationSystem();
    lostWeight = 0.;
    unsigned int ierr = 0;
    try {

        theMatrix.solveAndInvertBorderedBand(theVector, theVector);
        predict();

        for (unsigned int i = 0; i < optionList.size(); ++i) // down weighting iterations
                {
            size_t aPosition = methodList.find(optionList[i]);
            if (aPosition != std::string::npos) {
                aMethod = aPosition / 2 + 1;
                lostWeight = downWeight(aMethod);
                buildLinearEquationSystem();
                theMatrix.solveAndInvertBorderedBand(theVector, theVector);
                predict();
            }
        }
        Ndf = theData.size() - numParameters;
        Chi2 = 0.;
        for (unsigned int i = 0; i < theData.size(); ++i) {
            Chi2 += theData[i].getChi2();
        }
        Chi2 /= normChi2[aMethod];
        fitOK = true;

    } catch (int e) {
        std::cout << " GblTrajectory::fit exception " << e << std::endl;
        ierr = e;
    }
    return ierr;
}

void gbl::GblTrajectory::getFitToKinkJacobian ( std::vector< unsigned int > &  anIndex, SMatrix27 &  aJacobian, const GblPoint &  aPoint  ) const

private

Get jacobian for transformation from (trajectory) fit to kink parameters at point.

Jacobian broken lines (q/p,..,u_i-1,u_i,u_i+1..) to kink (du') parameters.

Parameters

[out]	anIndex	List of fit parameters with non zero derivatives
[out]	aJacobian	Corresponding transformation matrix
[in]	aPoint	Point to use

Definition at line 514 of file GblTrajectory.cc.

References gbl::GblPoint::getDerivatives(), gbl::GblPoint::getOffset(), i, numCurvature, numLocals, and theDimension.

Referenced by prepare().

                                                            {

    unsigned int nDim = theDimension.size();
    unsigned int nCurv = numCurvature;
    unsigned int nLocals = numLocals;

    int nOffset = aPoint.getOffset();

    SMatrix22 prevW, prevWJ, nextW, nextWJ;
    SVector2 prevWd, nextWd;
    aPoint.getDerivatives(0, prevW, prevWJ, prevWd); // W-, W- * J-, W- * d-
    aPoint.getDerivatives(1, nextW, nextWJ, nextWd); // W-, W- * J-, W- * d-
    const SMatrix22 sumWJ(prevWJ + nextWJ); // W- * J- + W+ * J+
    const SVector2 sumWd(prevWd + nextWd); // W+ * d+ + W- * d-

    unsigned int iOff = (nOffset - 1) * nDim + nCurv + nLocals + 1; // first offset ('i' in u_i)

    // local offset
    if (nCurv > 0) {
        aJacobian.Place_in_col(-sumWd, 0, 0); // from curvature
        anIndex[0] = nLocals + 1;
    }
    aJacobian.Place_at(prevW, 0, 1); // from 1st Offset
    aJacobian.Place_at(-sumWJ, 0, 3); // from 2nd Offset
    aJacobian.Place_at(nextW, 0, 5); // from 1st Offset
    for (unsigned int i = 0; i < nDim; ++i) {
        anIndex[1 + theDimension[i]] = iOff + i;
        anIndex[3 + theDimension[i]] = iOff + nDim + i;
        anIndex[5 + theDimension[i]] = iOff + nDim * 2 + i;
    }
}

void gbl::GblTrajectory::getFitToLocalJacobian ( std::vector< unsigned int > &  anIndex, SMatrix55 &  aJacobian, const GblPoint &  aPoint, unsigned int  measDim, unsigned int  nJacobian = 1  ) const

private

Get (part of) jacobian for transformation from (trajectory) fit to track parameters at point.

Jacobian broken lines (q/p,..,u_i,u_i+1..) to local (q/p,u',u) parameters.

Parameters

[out]	anIndex	List of fit parameters with non zero derivatives
[out]	aJacobian	Corresponding transformation matrix
[in]	aPoint	Point to use
[in]	measDim	Dimension of 'measurement' (<=2: calculate only offset part, >2: complete matrix)
[in]	nJacobian	Direction (0: to previous offset, 1: to next offset)

Definition at line 419 of file GblTrajectory.cc.

References gbl::GblPoint::getDerivatives(), gbl::GblPoint::getOffset(), i, numCurvature, numLocals, Validation_hcalonly_cfi::sign, and theDimension.

Referenced by getJacobian(), and prepare().

                                      {

    unsigned int nDim = theDimension.size();
    unsigned int nCurv = numCurvature;
    unsigned int nLocals = numLocals;

    int nOffset = aPoint.getOffset();

    if (nOffset < 0) // need interpolation
            {
        SMatrix22 prevW, prevWJ, nextW, nextWJ, matN;
        SVector2 prevWd, nextWd;
        int ierr;
        aPoint.getDerivatives(0, prevW, prevWJ, prevWd); // W-, W- * J-, W- * d-
        aPoint.getDerivatives(1, nextW, nextWJ, nextWd); // W-, W- * J-, W- * d-
        const SMatrix22 sumWJ(prevWJ + nextWJ);
        matN = sumWJ.Inverse(ierr); // N = (W- * J- + W+ * J+)^-1
        // derivatives for u_int
        const SMatrix22 prevNW(matN * prevW); // N * W-
        const SMatrix22 nextNW(matN * nextW); // N * W+
        const SVector2 prevNd(matN * prevWd); // N * W- * d-
        const SVector2 nextNd(matN * nextWd); // N * W+ * d+

        unsigned int iOff = nDim * (-nOffset - 1) + nLocals + nCurv + 1; // first offset ('i' in u_i)

        // local offset
        if (nCurv > 0) {
            aJacobian.Place_in_col(-prevNd - nextNd, 3, 0); // from curvature
            anIndex[0] = nLocals + 1;
        }
        aJacobian.Place_at(prevNW, 3, 1); // from 1st Offset
        aJacobian.Place_at(nextNW, 3, 3); // from 2nd Offset
        for (unsigned int i = 0; i < nDim; ++i) {
            anIndex[1 + theDimension[i]] = iOff + i;
            anIndex[3 + theDimension[i]] = iOff + nDim + i;
        }

        // local slope and curvature
        if (measDim > 2) {
            // derivatives for u'_int
            const SMatrix22 prevWPN(nextWJ * prevNW); // W+ * J+ * N * W-
            const SMatrix22 nextWPN(prevWJ * nextNW); // W- * J- * N * W+
            const SVector2 prevWNd(nextWJ * prevNd); // W+ * J+ * N * W- * d-
            const SVector2 nextWNd(prevWJ * nextNd); // W- * J- * N * W+ * d+
            if (nCurv > 0) {
                aJacobian(0, 0) = 1.0;
                aJacobian.Place_in_col(prevWNd - nextWNd, 1, 0); // from curvature
            }
            aJacobian.Place_at(-prevWPN, 1, 1); // from 1st Offset
            aJacobian.Place_at(nextWPN, 1, 3); // from 2nd Offset
        }
    } else { // at point
        // anIndex must be sorted
        // forward : iOff2 = iOff1 + nDim, index1 = 1, index2 = 3
        // backward: iOff2 = iOff1 - nDim, index1 = 3, index2 = 1
        unsigned int iOff1 = nDim * nOffset + nCurv + nLocals + 1; // first offset ('i' in u_i)
        unsigned int index1 = 3 - 2 * nJacobian; // index of first offset
        unsigned int iOff2 = iOff1 + nDim * (nJacobian * 2 - 1); // second offset ('i' in u_i)
        unsigned int index2 = 1 + 2 * nJacobian; // index of second offset
        // local offset
        aJacobian(3, index1) = 1.0; // from 1st Offset
        aJacobian(4, index1 + 1) = 1.0;
        for (unsigned int i = 0; i < nDim; ++i) {
            anIndex[index1 + theDimension[i]] = iOff1 + i;
        }

        // local slope and curvature
        if (measDim > 2) {
            SMatrix22 matW, matWJ;
            SVector2 vecWd;
            aPoint.getDerivatives(nJacobian, matW, matWJ, vecWd); // W, W * J, W * d
            double sign = (nJacobian > 0) ? 1. : -1.;
            if (nCurv > 0) {
                aJacobian(0, 0) = 1.0;
                aJacobian.Place_in_col(-sign * vecWd, 1, 0); // from curvature
                anIndex[0] = nLocals + 1;
            }
            aJacobian.Place_at(-sign * matWJ, 1, index1); // from 1st Offset
            aJacobian.Place_at(sign * matW, 1, index2); // from 2nd Offset
            for (unsigned int i = 0; i < nDim; ++i) {
                anIndex[index2 + theDimension[i]] = iOff2 + i;
            }
        }
    }
}

std::pair< std::vector< unsigned int >, TMatrixD > gbl::GblTrajectory::getJacobian ( int  aSignedLabel ) const

private

Get jacobian for transformation from fit to track parameters at point.

Jacobian broken lines (q/p,..,u_i,u_i+1..) to track (q/p,u',u) parameters including additional local parameters.

Parameters

[in]

aSignedLabel

(Signed) label of point for external seed (<0: in front, >0: after point, slope changes at scatterer!)

Returns

List of fit parameters with non zero derivatives and corresponding transformation matrix

Definition at line 343 of file GblTrajectory.cc.

References funct::abs(), getFitToLocalJacobian(), i, j, numCurvature, numLocals, numPoints, numTrajectories, theDimension, and thePoints.

Referenced by getResults(), and prepare().

                                {

    unsigned int nDim = theDimension.size();
    unsigned int nCurv = numCurvature;
    unsigned int nLocals = numLocals;
    unsigned int nBorder = nCurv + nLocals;
    unsigned int nParBRL = nBorder + 2 * nDim;
    unsigned int nParLoc = nLocals + 5;
    std::vector<unsigned int> anIndex;
    anIndex.reserve(nParBRL);
    TMatrixD aJacobian(nParLoc, nParBRL);
    aJacobian.Zero();

    unsigned int aLabel = abs(aSignedLabel);
    unsigned int firstLabel = 1;
    unsigned int lastLabel = 0;
    unsigned int aTrajectory = 0;
    // loop over trajectories
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        aTrajectory = iTraj;
        lastLabel += numPoints[iTraj];
        if (aLabel <= lastLabel)
            break;
        if (iTraj < numTrajectories - 1)
            firstLabel += numPoints[iTraj];
    }
    int nJacobian; // 0: prev, 1: next
    // check consistency of (index, direction)
    if (aSignedLabel > 0) {
        nJacobian = 1;
        if (aLabel >= lastLabel) {
            aLabel = lastLabel;
            nJacobian = 0;
        }
    } else {
        nJacobian = 0;
        if (aLabel <= firstLabel) {
            aLabel = firstLabel;
            nJacobian = 1;
        }
    }
    const GblPoint aPoint = thePoints[aTrajectory][aLabel - firstLabel];
    std::vector<unsigned int> labDer(5);
    SMatrix55 matDer;
    getFitToLocalJacobian(labDer, matDer, aPoint, 5, nJacobian);

    // from local parameters
    for (unsigned int i = 0; i < nLocals; ++i) {
        aJacobian(i + 5, i) = 1.0;
        anIndex.push_back(i + 1);
    }
    // from trajectory parameters
    unsigned int iCol = nLocals;
    for (unsigned int i = 0; i < 5; ++i) {
        if (labDer[i] > 0) {
            anIndex.push_back(labDer[i]);
            for (unsigned int j = 0; j < 5; ++j) {
                aJacobian(j, iCol) = matDer(j, i);
            }
            ++iCol;
        }
    }
    return std::make_pair(anIndex, aJacobian);
}

void gbl::GblTrajectory::getLabels

(

std::vector< unsigned int > &

aLabelList

)

Get (list of) labels of points on (simple) trajectory.

Parameters

[out]

aLabelList

List of labels (aLabelList[i] = i+1)

Definition at line 636 of file GblTrajectory.cc.

References i, and thePoints.

                                                                 {
    unsigned int aLabel = 0;
    unsigned int nPoint = thePoints[0].size();
    aLabelList.resize(nPoint);
    for (unsigned i = 0; i < nPoint; ++i) {
        aLabelList[i] = ++aLabel;
    }
}

void gbl::GblTrajectory::getLabels

(

std::vector< std::vector< unsigned int > > &

aLabelList

)

Get (list of lists of) labels of points on (composed) trajectory.

Parameters

[out]

aLabelList

List of of lists of labels

Definition at line 649 of file GblTrajectory.cc.

References i, numTrajectories, and thePoints.

                                                         {
    unsigned int aLabel = 0;
    aLabelList.resize(numTrajectories);
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        unsigned int nPoint = thePoints[iTraj].size();
        aLabelList[iTraj].resize(nPoint);
        for (unsigned i = 0; i < nPoint; ++i) {
            aLabelList[iTraj][i] = ++aLabel;
        }
    }
}

unsigned int gbl::GblTrajectory::getMeasResults	(	unsigned int	aLabel,
		unsigned int &	numData,
		TVectorD &	aResiduals,
		TVectorD &	aMeasErrors,
		TVectorD &	aResErrors,
		TVectorD &	aDownWeights
	)

Get residuals at point from measurement.

Get (diagonalized) residual, error of measurement and residual and down-weighting factor for measurement at point

Parameters

[in]	aLabel	Label of point on trajectory
[out]	numData	Number of data blocks from measurement at point
[out]	aResiduals	Measurements-Predictions
[out]	aMeasErrors	Errors of Measurements
[out]	aResErrors	Errors of Residuals (including correlations from track fit)
[out]	aDownWeights	Down-Weighting factors

Returns

error code (non-zero if trajectory not fitted successfully)

Definition at line 587 of file GblTrajectory.cc.

References fitOK, getResAndErr(), i, and measDataIndex.

                                                      {
    numData = 0;
    if (not fitOK)
        return 1;

    unsigned int firstData = measDataIndex[aLabel - 1]; // first data block with measurement
    numData = measDataIndex[aLabel] - firstData; // number of data blocks
    for (unsigned int i = 0; i < numData; ++i) {
        getResAndErr(firstData + i, aResiduals[i], aMeasErrors[i],
                aResErrors[i], aDownWeights[i]);
    }
    return 0;
}

unsigned int gbl::GblTrajectory::getNumPoints

(

)

const

Retrieve number of point from trajectory.

Definition at line 225 of file GblTrajectory.cc.

References numAllPoints.

                                               {
    return numAllPoints;
}

void gbl::GblTrajectory::getResAndErr ( unsigned int  aData, double &  aResidual, double &  aMeasError, double &  aResError, double &  aDownWeight  )

private

Get residual and errors from data block.

Get residual, error of measurement and residual and down-weighting factor for (single) data block

Parameters

[in]	aData	Label of data block
[out]	aResidual	Measurement-Prediction
[out]	aMeasError	Error of Measurement
[out]	aResError	Error of Residual (including correlations from track fit)
[out]	aDownWeight	Down-Weighting factor

Definition at line 672 of file GblTrajectory.cc.

References gbl::BorderedBandMatrix::getBlockMatrix(), j, mathSSE::sqrt(), theData, and theMatrix.

Referenced by getMeasResults(), and getScatResults().

                                                                    {

    double aMeasVar;
    std::vector<unsigned int>* indLocal;
    std::vector<double>* derLocal;
    theData[aData].getResidual(aResidual, aMeasVar, aDownWeight, indLocal,
            derLocal);
    unsigned int nParBrl = (*indLocal).size();
    TVectorD aVec(nParBrl); // compressed vector of derivatives
    for (unsigned int j = 0; j < nParBrl; ++j) {
        aVec[j] = (*derLocal)[j];
    }
    TMatrixDSym aMat = theMatrix.getBlockMatrix(*indLocal); // compressed (covariance) matrix
    double aFitVar = aMat.Similarity(aVec); // variance from track fit
    aMeasError = sqrt(aMeasVar); // error of measurement
    aResError = (aFitVar < aMeasVar ? sqrt(aMeasVar - aFitVar) : 0.); // error of residual
}

unsigned int gbl::GblTrajectory::getResults	(	int	aSignedLabel,
		TVectorD &	localPar,
		TMatrixDSym &	localCov
	)		const

Get fit results at point.

Get corrections and covariance matrix for local track and additional parameters in forward or backward direction.

Parameters

[in]	aSignedLabel	(Signed) label of point on trajectory (<0: in front, >0: after point, slope changes at scatterer!)
[out]	localPar	Corrections for local parameters
[out]	localCov	Covariance for local parameters

Returns

error code (non-zero if trajectory not fitted successfully)

Definition at line 557 of file GblTrajectory.cc.

References fitOK, gbl::BorderedBandMatrix::getBlockMatrix(), getJacobian(), i, theMatrix, and theVector.

                                     {
    if (not fitOK)
        return 1;
    std::pair<std::vector<unsigned int>, TMatrixD> indexAndJacobian =
            getJacobian(aSignedLabel);
    unsigned int nParBrl = indexAndJacobian.first.size();
    TVectorD aVec(nParBrl); // compressed vector
    for (unsigned int i = 0; i < nParBrl; ++i) {
        aVec[i] = theVector(indexAndJacobian.first[i] - 1);
    }
    TMatrixDSym aMat = theMatrix.getBlockMatrix(indexAndJacobian.first); // compressed matrix
    localPar = indexAndJacobian.second * aVec;
    localCov = aMat.Similarity(indexAndJacobian.second);
    return 0;
}

unsigned int gbl::GblTrajectory::getScatResults	(	unsigned int	aLabel,
		unsigned int &	numData,
		TVectorD &	aResiduals,
		TVectorD &	aMeasErrors,
		TVectorD &	aResErrors,
		TVectorD &	aDownWeights
	)

Get (kink) residuals at point from scatterer.

Get (diagonalized) residual, error of measurement and residual and down-weighting factor for scatterering kinks at point

Parameters

[in]	aLabel	Label of point on trajectory
[out]	numData	Number of data blocks from scatterer at point
[out]	aResiduals	(kink)Measurements-(kink)Predictions
[out]	aMeasErrors	Errors of (kink)Measurements
[out]	aResErrors	Errors of Residuals (including correlations from track fit)
[out]	aDownWeights	Down-Weighting factors

Returns

error code (non-zero if trajectory not fitted successfully)

Definition at line 616 of file GblTrajectory.cc.

References fitOK, getResAndErr(), i, and scatDataIndex.

                                                      {
    numData = 0;
    if (not fitOK)
        return 1;

    unsigned int firstData = scatDataIndex[aLabel - 1]; // first data block with scatterer
    numData = scatDataIndex[aLabel] - firstData; // number of data blocks
    for (unsigned int i = 0; i < numData; ++i) {
        getResAndErr(firstData + i, aResiduals[i], aMeasErrors[i],
                aResErrors[i], aDownWeights[i]);
    }
    return 0;
}

bool gbl::GblTrajectory::isValid

(

void

)

const

Retrieve validity of trajectory.

Definition at line 220 of file GblTrajectory.cc.

References constructOK.

Referenced by ntuple._Object::_checkIsValid(), and core.AutoHandle.AutoHandle::ReallyLoad().

                                  {
    return constructOK;
}

void gbl::GblTrajectory::milleOut

(

MilleBinary &

aMille

)

Write valid trajectory to Millepede-II binary file.

Definition at line 1004 of file GblTrajectory.cc.

References gbl::MilleBinary::addData(), constructOK, theData, and gbl::MilleBinary::writeRecord().

Referenced by MillePedeAlignmentAlgorithm::addReferenceTrajectory().

                                                {
    double aValue;
    double aErr;
    std::vector<unsigned int>* indLocal;
    std::vector<double>* derLocal;
    std::vector<int>* labGlobal;
    std::vector<double>* derGlobal;

    if (not constructOK)
        return;

//   data: measurements, kinks and external seed
    std::vector<GblData>::iterator itData;
    for (itData = theData.begin(); itData != theData.end(); ++itData) {
        itData->getAllData(aValue, aErr, indLocal, derLocal, labGlobal,
                derGlobal);
        aMille.addData(aValue, aErr, *indLocal, *derLocal, *labGlobal,
                *derGlobal);
    }
    aMille.writeRecord();
}

void gbl::GblTrajectory::predict ( )

private

Calculate predictions for all points.

Definition at line 926 of file GblTrajectory.cc.

References theData, and theVector.

Referenced by fit().

                            {
    std::vector<GblData>::iterator itData;
    for (itData = theData.begin(); itData < theData.end(); ++itData) {
        itData->setPrediction(theVector);
    }
}

void gbl::GblTrajectory::prepare ( )

private

Prepare fit for simple or composed trajectory.

Generate data (blocks) from measurements, kinks, external seed and measurements.

Definition at line 713 of file GblTrajectory.cc.

References gbl::GblData::addDerivatives(), begin, externalDerivatives, externalMeasurements, externalPoint, externalPrecisions, externalSeed, getFitToKinkJacobian(), getFitToLocalJacobian(), getJacobian(), i, diffTreeTool::index, innerTransformations, j, gen::k, measDataIndex, numAllPoints, numCurvature, numInnerTrans, numLocals, numMeasurements, numOffsets, numTrajectories, scatDataIndex, theData, theDimension, thePoints, and mitigatedMETSequence_cff::U.

Referenced by construct().

                            {
    unsigned int nDim = theDimension.size();
    // upper limit
    unsigned int maxData = numMeasurements + nDim * (numOffsets - 2)
            + externalSeed.GetNrows();
    theData.reserve(maxData);
    measDataIndex.resize(numAllPoints + 3); // include external seed and measurements
    scatDataIndex.resize(numAllPoints + 1);
    unsigned int nData = 0;
    std::vector<TMatrixD> innerTransDer;
    std::vector<std::vector<unsigned int> > innerTransLab;
    // composed trajectory ?
    if (numInnerTrans > 0) {
        //std::cout << "composed trajectory" << std::endl;
        for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
            // innermost point
            GblPoint* innerPoint = &thePoints[iTraj].front();
            // transformation fit to local track parameters
            std::vector<unsigned int> firstLabels(5);
            SMatrix55 matFitToLocal, matLocalToFit;
            getFitToLocalJacobian(firstLabels, matFitToLocal, *innerPoint, 5);
            // transformation local track to fit parameters
            int ierr;
            matLocalToFit = matFitToLocal.Inverse(ierr);
            TMatrixD localToFit(5, 5);
            for (unsigned int i = 0; i < 5; ++i) {
                for (unsigned int j = 0; j < 5; ++j) {
                    localToFit(i, j) = matLocalToFit(i, j);
                }
            }
            // transformation external to fit parameters at inner (first) point
            innerTransDer.push_back(localToFit * innerTransformations[iTraj]);
            innerTransLab.push_back(firstLabels);
        }
    }
    // measurements
    SMatrix55 matP;
    // loop over trajectories
    std::vector<GblPoint>::iterator itPoint;
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        for (itPoint = thePoints[iTraj].begin();
                itPoint < thePoints[iTraj].end(); ++itPoint) {
            SVector5 aMeas, aPrec;
            unsigned int nLabel = itPoint->getLabel();
            unsigned int measDim = itPoint->hasMeasurement();
            if (measDim) {
                const TMatrixD localDer = itPoint->getLocalDerivatives();
                const std::vector<int> globalLab = itPoint->getGlobalLabels();
                const TMatrixD globalDer = itPoint->getGlobalDerivatives();
                TMatrixD transDer;
                itPoint->getMeasurement(matP, aMeas, aPrec);
                unsigned int iOff = 5 - measDim; // first active component
                std::vector<unsigned int> labDer(5);
                SMatrix55 matDer, matPDer;
                unsigned int nJacobian =
                        (itPoint < thePoints[iTraj].end() - 1) ? 1 : 0; // last point needs backward propagation
                getFitToLocalJacobian(labDer, matDer, *itPoint, measDim,
                        nJacobian);
                if (measDim > 2) {
                    matPDer = matP * matDer;
                } else { // 'shortcut' for position measurements
                    matPDer.Place_at(
                            matP.Sub<SMatrix22>(3, 3)
                                    * matDer.Sub<SMatrix25>(3, 0), 3, 0);
                }

                if (numInnerTrans > 0) {
                    // transform for external parameters
                    TMatrixD proDer(measDim, 5);
                    // match parameters
                    unsigned int ifirst = 0;
                    unsigned int ilabel = 0;
                    while (ilabel < 5) {
                        if (labDer[ilabel] > 0) {
                            while (innerTransLab[iTraj][ifirst]
                                    != labDer[ilabel] and ifirst < 5) {
                                ++ifirst;
                            }
                            if (ifirst >= 5) {
                                labDer[ilabel] -= 2 * nDim * (iTraj + 1); // adjust label
                            } else {
                                // match
                                labDer[ilabel] = 0; // mark as related to external parameters
                                for (unsigned int k = iOff; k < 5; ++k) {
                                    proDer(k - iOff, ifirst) = matPDer(k,
                                            ilabel);
                                }
                            }
                        }
                        ++ilabel;
                    }
                    transDer.ResizeTo(measDim, numCurvature);
                    transDer = proDer * innerTransDer[iTraj];
                }
                for (unsigned int i = iOff; i < 5; ++i) {
                    if (aPrec(i) > 0.) {
                        GblData aData(nLabel, aMeas(i), aPrec(i));
                        aData.addDerivatives(i, labDer, matPDer, iOff, localDer,
                                globalLab, globalDer, numLocals, transDer);
                        theData.push_back(aData);
                        nData++;
                    }
                }

            }
            measDataIndex[nLabel] = nData;
        }
    }

    // pseudo measurements from kinks
    SMatrix22 matT;
    scatDataIndex[0] = nData;
    scatDataIndex[1] = nData;
    // loop over trajectories
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        for (itPoint = thePoints[iTraj].begin() + 1;
                itPoint < thePoints[iTraj].end() - 1; ++itPoint) {
            SVector2 aMeas, aPrec;
            unsigned int nLabel = itPoint->getLabel();
            if (itPoint->hasScatterer()) {
                itPoint->getScatterer(matT, aMeas, aPrec);
                TMatrixD transDer;
                std::vector<unsigned int> labDer(7);
                SMatrix27 matDer, matTDer;
                getFitToKinkJacobian(labDer, matDer, *itPoint);
                matTDer = matT * matDer;
                if (numInnerTrans > 0) {
                    // transform for external parameters
                    TMatrixD proDer(nDim, 5);
                    // match parameters
                    unsigned int ifirst = 0;
                    unsigned int ilabel = 0;
                    while (ilabel < 7) {
                        if (labDer[ilabel] > 0) {
                            while (innerTransLab[iTraj][ifirst]
                                    != labDer[ilabel] and ifirst < 5) {
                                ++ifirst;
                            }
                            if (ifirst >= 5) {
                                labDer[ilabel] -= 2 * nDim * (iTraj + 1); // adjust label
                            } else {
                                // match
                                labDer[ilabel] = 0; // mark as related to external parameters
                                for (unsigned int k = 0; k < nDim; ++k) {
                                    proDer(k, ifirst) = matTDer(k, ilabel);
                                }
                            }
                        }
                        ++ilabel;
                    }
                    transDer.ResizeTo(nDim, numCurvature);
                    transDer = proDer * innerTransDer[iTraj];
                }
                for (unsigned int i = 0; i < nDim; ++i) {
                    unsigned int iDim = theDimension[i];
                    if (aPrec(iDim) > 0.) {
                        GblData aData(nLabel, aMeas(iDim), aPrec(iDim));
                        aData.addDerivatives(iDim, labDer, matTDer, numLocals,
                                transDer);
                        theData.push_back(aData);
                        nData++;
                    }
                }
            }
            scatDataIndex[nLabel] = nData;
        }
        scatDataIndex[thePoints[iTraj].back().getLabel()] = nData;
    }

    // external seed
    if (externalPoint > 0) {
        std::pair<std::vector<unsigned int>, TMatrixD> indexAndJacobian =
                getJacobian(externalPoint);
        std::vector<unsigned int> externalSeedIndex = indexAndJacobian.first;
        std::vector<double> externalSeedDerivatives(externalSeedIndex.size());
        const TMatrixDSymEigen externalSeedEigen(externalSeed);
        const TVectorD valEigen = externalSeedEigen.GetEigenValues();
        TMatrixD vecEigen = externalSeedEigen.GetEigenVectors();
        vecEigen = vecEigen.T() * indexAndJacobian.second;
        for (int i = 0; i < externalSeed.GetNrows(); ++i) {
            if (valEigen(i) > 0.) {
                for (int j = 0; j < externalSeed.GetNcols(); ++j) {
                    externalSeedDerivatives[j] = vecEigen(i, j);
                }
                GblData aData(externalPoint, 0., valEigen(i));
                aData.addDerivatives(externalSeedIndex, externalSeedDerivatives);
                theData.push_back(aData);
                nData++;
            }
        }
    }
    measDataIndex[numAllPoints + 1] = nData;
    // external measurements
    unsigned int nExt = externalMeasurements.GetNrows();
    if (nExt > 0) {
        std::vector<unsigned int> index(numCurvature);
        std::vector<double> derivatives(numCurvature);
        for (unsigned int iExt = 0; iExt < nExt; ++iExt) {
            for (unsigned int iCol = 0; iCol < numCurvature; ++iCol) {
                index[iCol] = iCol + 1;
                derivatives[iCol] = externalDerivatives(iExt, iCol);
            }
            GblData aData(1U, externalMeasurements(iExt),
                    externalPrecisions(iExt));
            aData.addDerivatives(index, derivatives);
            theData.push_back(aData);
            nData++;
        }
    }
    measDataIndex[numAllPoints + 2] = nData;
}

void gbl::GblTrajectory::printData

(

)

Print GblData blocks for trajectory.

Definition at line 1108 of file GblTrajectory.cc.

References gather_cfg::cout, and theData.

                              {
    std::cout << "GblData blocks " << std::endl;
    std::vector<GblData>::iterator itData;
    for (itData = theData.begin(); itData < theData.end(); ++itData) {
        itData->printData();
    }
}

void gbl::GblTrajectory::printPoints

(

unsigned int

level = 0

)

Print GblPoints on trajectory.

Parameters

[in]

level

print level (0: minimum, >0: more)

Definition at line 1096 of file GblTrajectory.cc.

References begin, gather_cfg::cout, numTrajectories, and thePoints.

                                                  {
    std::cout << "GblPoints " << std::endl;
    for (unsigned int iTraj = 0; iTraj < numTrajectories; ++iTraj) {
        std::vector<GblPoint>::iterator itPoint;
        for (itPoint = thePoints[iTraj].begin();
                itPoint < thePoints[iTraj].end(); ++itPoint) {
            itPoint->printPoint(level);
        }
    }
}

void gbl::GblTrajectory::printTrajectory

(

unsigned int

level = 0

)

Print GblTrajectory.

Parameters

[in]

level

print level (0: minimum, >0: more)

Definition at line 1030 of file GblTrajectory.cc.

References constructOK, gather_cfg::cout, externalDerivatives, externalMeasurements, externalPoint, externalPrecisions, externalSeed, fitOK, i, innerTransformations, numAllPoints, numInnerTrans, numMeasurements, numOffsets, numParameters, gbl::VVector::print(), gbl::BorderedBandMatrix::printMatrix(), theDimension, theMatrix, and theVector.

                                                      {
    if (numInnerTrans) {
        std::cout << "Composed GblTrajectory, " << numInnerTrans
                << " subtrajectories" << std::endl;
    } else {
        std::cout << "Simple GblTrajectory" << std::endl;
    }
    if (theDimension.size() < 2) {
        std::cout << " 2D-trajectory" << std::endl;
    }
    std::cout << " Number of GblPoints          : " << numAllPoints
            << std::endl;
    std::cout << " Number of points with offsets: " << numOffsets << std::endl;
    std::cout << " Number of fit parameters     : " << numParameters
            << std::endl;
    std::cout << " Number of measurements       : " << numMeasurements
            << std::endl;
    if (externalMeasurements.GetNrows()) {
        std::cout << " Number of ext. measurements  : "
                << externalMeasurements.GetNrows() << std::endl;
    }
    if (externalPoint) {
        std::cout << " Label of point with ext. seed: " << externalPoint
                << std::endl;
    }
    if (constructOK) {
        std::cout << " Constructed OK " << std::endl;
    }
    if (fitOK) {
        std::cout << " Fitted OK " << std::endl;
    }
    if (level > 0) {
        if (numInnerTrans) {
            std::cout << " Inner transformations" << std::endl;
            for (unsigned int i = 0; i < numInnerTrans; ++i) {
                innerTransformations[i].Print();
            }
        }
        if (externalMeasurements.GetNrows()) {
            std::cout << " External measurements" << std::endl;
            std::cout << "  Measurements:" << std::endl;
            externalMeasurements.Print();
            std::cout << "  Precisions:" << std::endl;
            externalPrecisions.Print();
            std::cout << "  Derivatives:" << std::endl;
            externalDerivatives.Print();
        }
        if (externalPoint) {
            std::cout << " External seed:" << std::endl;
            externalSeed.Print();
        }
        if (fitOK) {
            std::cout << " Fit results" << std::endl;
            std::cout << "  Parameters:" << std::endl;
            theVector.print();
            std::cout << "  Covariance matrix (bordered band part):"
                    << std::endl;
            theMatrix.printMatrix();
        }
    }
}

Member Data Documentation

bool gbl::GblTrajectory::constructOK

private

Trajectory has been successfully constructed (ready for fit/output)

Definition at line 70 of file GblTrajectory.h.

Referenced by construct(), fit(), isValid(), milleOut(), and printTrajectory().

TMatrixD gbl::GblTrajectory::externalDerivatives

private

Definition at line 80 of file GblTrajectory.h.

Referenced by prepare(), and printTrajectory().

TVectorD gbl::GblTrajectory::externalMeasurements

private

Definition at line 81 of file GblTrajectory.h.

Referenced by prepare(), and printTrajectory().

unsigned int gbl::GblTrajectory::externalPoint

private

Label of external point (or 0)

Definition at line 69 of file GblTrajectory.h.

Referenced by prepare(), and printTrajectory().

TVectorD gbl::GblTrajectory::externalPrecisions

private

Definition at line 82 of file GblTrajectory.h.

Referenced by prepare(), and printTrajectory().

TMatrixDSym gbl::GblTrajectory::externalSeed

private

Precision (inverse covariance matrix) of external seed.

Definition at line 77 of file GblTrajectory.h.

Referenced by prepare(), and printTrajectory().

bool gbl::GblTrajectory::fitOK

private

Trajectory has been successfully fitted (results are valid)

Definition at line 71 of file GblTrajectory.h.

Referenced by construct(), fit(), getMeasResults(), getResults(), getScatResults(), and printTrajectory().

std::vector<TMatrixD> gbl::GblTrajectory::innerTransformations

private

Transformations at innermost points of.

Definition at line 78 of file GblTrajectory.h.

Referenced by GblTrajectory(), prepare(), and printTrajectory().

std::vector<unsigned int> gbl::GblTrajectory::measDataIndex

private

mapping points to data blocks from measurements

Definition at line 75 of file GblTrajectory.h.

Referenced by getMeasResults(), and prepare().

unsigned int gbl::GblTrajectory::numAllPoints

private

Number of all points on trajectory.

Definition at line 60 of file GblTrajectory.h.

Referenced by construct(), GblTrajectory(), getNumPoints(), prepare(), and printTrajectory().

unsigned int gbl::GblTrajectory::numCurvature

private

Number of curvature parameters (0 or 1) or external parameters.

Definition at line 65 of file GblTrajectory.h.

Referenced by buildLinearEquationSystem(), construct(), GblTrajectory(), getFitToKinkJacobian(), getFitToLocalJacobian(), getJacobian(), and prepare().

unsigned int gbl::GblTrajectory::numInnerTrans

private

Number of inner transformations to external parameters.

Definition at line 64 of file GblTrajectory.h.

Referenced by construct(), prepare(), and printTrajectory().

unsigned int gbl::GblTrajectory::numLocals

private

Total number of (additional) local parameters.

Definition at line 67 of file GblTrajectory.h.

Referenced by buildLinearEquationSystem(), construct(), getFitToKinkJacobian(), getFitToLocalJacobian(), getJacobian(), and prepare().

unsigned int gbl::GblTrajectory::numMeasurements

private

Total number of measurements.

Definition at line 68 of file GblTrajectory.h.

Referenced by construct(), prepare(), and printTrajectory().

unsigned int gbl::GblTrajectory::numOffsets

private

Number of (points with) offsets on trajectory.

Definition at line 63 of file GblTrajectory.h.

Referenced by construct(), defineOffsets(), prepare(), and printTrajectory().

unsigned int gbl::GblTrajectory::numParameters

private

Number of fit parameters.

Definition at line 66 of file GblTrajectory.h.

Referenced by buildLinearEquationSystem(), construct(), fit(), and printTrajectory().

std::vector<unsigned int> gbl::GblTrajectory::numPoints

private

Number of points on (sub)trajectory.

Definition at line 61 of file GblTrajectory.h.

Referenced by GblTrajectory(), and getJacobian().

unsigned int gbl::GblTrajectory::numTrajectories

private

Number of trajectories (in composed trajectory)

Definition at line 62 of file GblTrajectory.h.

Referenced by calcJacobians(), construct(), defineOffsets(), getJacobian(), getLabels(), prepare(), and printPoints().

std::vector<unsigned int> gbl::GblTrajectory::scatDataIndex

private

mapping points to data blocks from scatterers

Definition at line 76 of file GblTrajectory.h.

Referenced by getScatResults(), and prepare().

std::vector<GblData> gbl::GblTrajectory::theData

private

List of data blocks.

Definition at line 74 of file GblTrajectory.h.

Referenced by buildLinearEquationSystem(), downWeight(), fit(), getResAndErr(), milleOut(), predict(), prepare(), and printData().

std::vector<unsigned int> gbl::GblTrajectory::theDimension

private

List of active dimensions (0=u1, 1=u2) in fit.

Definition at line 72 of file GblTrajectory.h.

Referenced by construct(), GblTrajectory(), getFitToKinkJacobian(), getFitToLocalJacobian(), getJacobian(), prepare(), and printTrajectory().

BorderedBandMatrix gbl::GblTrajectory::theMatrix

private

(Bordered band) matrix of linear equation system

Definition at line 84 of file GblTrajectory.h.

Referenced by buildLinearEquationSystem(), fit(), getResAndErr(), getResults(), and printTrajectory().

std::vector<std::vector<GblPoint> > gbl::GblTrajectory::thePoints

private

(list of) List of points on trajectory

Definition at line 73 of file GblTrajectory.h.

Referenced by calcJacobians(), construct(), defineOffsets(), GblTrajectory(), getJacobian(), getLabels(), prepare(), and printPoints().

VVector gbl::GblTrajectory::theVector

private

Vector of linear equation system.

Definition at line 83 of file GblTrajectory.h.

Referenced by buildLinearEquationSystem(), fit(), getResults(), predict(), and printTrajectory().