GBL/include/GblTrajectory.h

0001 /*
0002  * GblTrajectory.h
0003  *
0004  *  Created on: Aug 18, 2011
0005  *      Author: kleinwrt
0006  */
0007
0008 /** \file
0009  *  GblTrajectory definition.
0010  *
0011  *  \author Claus Kleinwort, DESY, 2011 (Claus.Kleinwort@desy.de)
0012  *
0013  *  \copyright
0014  *  Copyright (c) 2011 - 2016 Deutsches Elektronen-Synchroton,
0015  *  Member of the Helmholtz Association, (DESY), HAMBURG, GERMANY \n\n
0016  *  This library is free software; you can redistribute it and/or modify
0017  *  it under the terms of the GNU Library General Public License as
0018  *  published by the Free Software Foundation; either version 2 of the
0019  *  License, or (at your option) any later version. \n\n
0020  *  This library is distributed in the hope that it will be useful,
0021  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
0022  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0023  *  GNU Library General Public License for more details. \n\n
0024  *  You should have received a copy of the GNU Library General Public
0025  *  License along with this program (see the file COPYING.LIB for more
0026  *  details); if not, write to the Free Software Foundation, Inc.,
0027  *  675 Mass Ave, Cambridge, MA 02139, USA.
0028  */
0029
0030 #ifndef GBLTRAJECTORY_H_
0031 #define GBLTRAJECTORY_H_
0032
0033 #include "GblPoint.h"
0034 #include "GblData.h"
0035 #include "GblPoint.h"
0036 #include "BorderedBandMatrix.h"
0037 #include "MilleBinary.h"
0038 #include "TMatrixDSymEigen.h"
0039
0040 //! Namespace for the general broken lines package
0041 namespace gbl {
0042
0043 /// GBL trajectory.
0044 /**
0045  * List of GblPoints ordered by arc length.
0046  * Can be fitted and optionally written to MP-II binary file.
0047  */
0048 class GblTrajectory {
0049 public:
0050     explicit GblTrajectory(const std::vector<GblPoint> &aPointList, bool flagCurv = true,
0051             bool flagU1dir = true, bool flagU2dir = true);
0052     GblTrajectory(const std::vector<GblPoint> &aPointList, unsigned int aLabel,
0053             const TMatrixDSym &aSeed, bool flagCurv = true, bool flagU1dir =
0054                     true, bool flagU2dir = true);
0055     explicit GblTrajectory(
0056             const std::vector<std::pair<std::vector<GblPoint>, TMatrixD> > &aPointaAndTransList);
0057     GblTrajectory(
0058             const std::vector<std::pair<std::vector<GblPoint>, TMatrixD> > &aPointaAndTransList,
0059             const TMatrixD &extDerivatives, const TVectorD &extMeasurements,
0060             const TVectorD &extPrecisions);
0061     GblTrajectory(
0062             const std::vector<std::pair<std::vector<GblPoint>, TMatrixD> > &aPointaAndTransList,
0063             const TMatrixD &extDerivatives, const TVectorD &extMeasurements,
0064             const TMatrixDSym &extPrecisions);
0065     virtual ~GblTrajectory();
0066     bool isValid() const;
0067     unsigned int getNumPoints() const;
0068     unsigned int getResults(int aSignedLabel, TVectorD &localPar,
0069             TMatrixDSym &localCov) const;
0070     unsigned int getMeasResults(unsigned int aLabel, unsigned int &numRes,
0071             TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors,
0072             TVectorD &aDownWeights);
0073     unsigned int getScatResults(unsigned int aLabel, unsigned int &numRes,
0074             TVectorD &aResiduals, TVectorD &aMeasErrors, TVectorD &aResErrors,
0075             TVectorD &aDownWeights);
0076     unsigned int getLabels(std::vector<unsigned int> &aLabelList);
0077     unsigned int getLabels(std::vector<std::vector<unsigned int> > &aLabelList);
0078     unsigned int fit(double &Chi2, int &Ndf, double &lostWeight,
0079             std::string optionList = "");
0080     void milleOut(MilleBinary &aMille);
0081     void printTrajectory(unsigned int level = 0);
0082     void printPoints(unsigned int level = 0);
0083     void printData();
0084         std::vector<GblData> getData() {return theData;}
0085
0086 private:
0087     unsigned int numAllPoints; ///< Number of all points on trajectory
0088     std::vector<unsigned int> numPoints; ///< Number of points on (sub)trajectory
0089     unsigned int numTrajectories; ///< Number of trajectories (in composed trajectory)
0090     unsigned int numOffsets; ///< Number of (points with) offsets on trajectory
0091     unsigned int numInnerTrans; ///< Number of inner transformations to external parameters
0092     unsigned int numCurvature; ///< Number of curvature parameters (0 or 1) or external parameters
0093     unsigned int numParameters; ///< Number of fit parameters
0094     unsigned int numLocals; ///< Total number of (additional) local parameters
0095     unsigned int numMeasurements; ///< Total number of measurements
0096     unsigned int externalPoint; ///< Label of external point (or 0)
0097     bool constructOK; ///< Trajectory has been successfully constructed (ready for fit/output)
0098     bool fitOK; ///< Trajectory has been successfully fitted (results are valid)
0099     std::vector<unsigned int> theDimension; ///< List of active dimensions (0=u1, 1=u2) in fit
0100     std::vector<std::vector<GblPoint> > thePoints; ///< (list of) List of points on trajectory
0101     std::vector<GblData> theData; ///< List of data blocks
0102     std::vector<unsigned int> measDataIndex; ///< mapping points to data blocks from measurements
0103     std::vector<unsigned int> scatDataIndex; ///< mapping points to data blocks from scatterers
0104     TMatrixDSym externalSeed; ///< Precision (inverse covariance matrix) of external seed
0105     std::vector<TMatrixD> innerTransformations; ///< Transformations at innermost points of
0106     // composed trajectory (from common external parameters)
0107     TMatrixD externalDerivatives; // Derivatives for external measurements of composed trajectory
0108     TVectorD externalMeasurements; // Residuals for external measurements of composed trajectory
0109     TVectorD externalPrecisions; // Precisions for external measurements of composed trajectory
0110     VVector theVector; ///< Vector of linear equation system
0111     BorderedBandMatrix theMatrix; ///< (Bordered band) matrix of linear equation system
0112
0113     std::pair<std::vector<unsigned int>, TMatrixD> getJacobian(
0114             int aSignedLabel) const;
0115     void getFitToLocalJacobian(std::vector<unsigned int> &anIndex,
0116             SMatrix55 &aJacobian, const GblPoint &aPoint, unsigned int measDim,
0117             unsigned int nJacobian = 1) const;
0118     void getFitToKinkJacobian(std::vector<unsigned int> &anIndex,
0119             SMatrix27 &aJacobian, const GblPoint &aPoint) const;
0120     void construct();
0121     void defineOffsets();
0122     void calcJacobians();
0123     void prepare();
0124     void buildLinearEquationSystem();
0125     void predict();
0126     double downWeight(unsigned int aMethod);
0127     void getResAndErr(unsigned int aData, double &aResidual,
0128             double &aMeadsError, double &aResError, double &aDownWeight);
0129 };
0130 }
0131 #endif /* GBLTRAJECTORY_H_ */