Lutz.h

Go to the documentation of this file.

///
/// \file   Lutz.h
/// \brief  "Break-point" track fitter
/// \author messier@indiana.edu
///
#include "TVector.h"
#include "TMatrix.h"

namespace bpfit {
  ///
  /// \brief Break-point track fitter
  ///
  /// Implementation of
  ///...................................................................
  ///    OPTIMUM TRACK FITTING IN THE PRESENCE OF MULTIPLE SCATTERING
  ///                              G. LUTZ
  ///         Nuclear Instruments and Methods in Physics Research
  ///                      A273 (1998) pp. 349-361
  ///...................................................................
  class Lutz {
  public:
    ///
    /// \brief Construct a fitter for n measurements and N scattering planes
    /// \param n - number of measurements to be added
    /// \param N - number of scattering surfaces to be added
    ///
    Lutz(unsigned int n, unsigned int N);

    ///
    /// \brief Add measurements
    /// \param i    - measurement index [0...n)
    /// \param z    - z location of measurement (longitudinal location)
    /// \param x    - x location of measurement (transverse location)
    /// \param sigx - measurement error on x position
    ///
    void SetMeasurement(unsigned int i, double z, double x, double sigx);
    
    ///
    /// \brief Set Z location of a scattering plane j
    /// \param j - index of scattering plane (start from 0)
    /// \param Z - longitudinal location of scattering plane
    /// \param sigSJ - RMS scattering of track slope at this plane
    ///
    void SetScatteringPlane(unsigned int j, double Z, double sigSJ);
    
    ///
    /// \brief Do the fit
    /// 
    /// \param a     - on return, fitted intercept. OK to pass 0.
    /// \param b     - on return, fitted initial slope. OK to pass 0.
    /// \param alpha - on return, fitted scattering slopes. OK to pass 0.
    /// \param chi2  - on return, the chi^2 result of the fit. OK to pass 0.
    ///
    void Fit(double *a, double* b, double* alpha, double* chi2);

    ///
    /// \brief After fit, contribution to chi^2 of the ith measurement
    ///
    double Chi2XIi(unsigned int i) const;
    
    ///
    /// \brief After fit, contribution to chi^2 of the jth scattering plane
    ///
    double Chi2BetaJ(unsigned int j) const;
    
    //
    // \brief After fit, contribution from measurements to chi^2
    //
    double Chi2XI() const;

    //
    // \brief After fit, contribution from scattering to chi^2
    //
    double Chi2Beta() const;
    
    ///
    /// \brief After fit, best-fit x location of track at location z
    ///
    double X(double z) const;
    
  private:
    void CalcZeta();
    void CalcD();
    bool CalcAinv();
    void Calcq();
    void Calcx();
    void NullResult();

  private:
    unsigned int        fn;    ///< Number of measurements
    unsigned int        fN;    ///< Number of scattering planes
    std::vector<double> fz;    ///< Measurement plane longitudinal locations
    TVectorT<double>    fXI;   ///< Measurements in transverse location
    std::vector<double> fwx;   ///< Measurement weight (1/sigx^2)
    std::vector<double> fZ;    ///< Location of scattering planes
    std::vector<double> fwSJ;  ///< Scattering angle weight (1/sigSJ^2)
    std::vector< std::vector<double> > fZeta; ///< measurement-scattering plane distances
    TMatrixT<double>    fD;    ///< Eqn. 31
    TMatrixTSym<double> fAinv; ///< System of equations, Eqn. 29
    double              fa;    ///< Fitted intercept
    double              fb;    ///< Fitted slope
    std::vector<double> fAlpha;///< Fitted scatters
    std::vector<double> fx;    ///< Fitted track locations
  };
}

////////////////////////////////////////////////////////////////////////