bpfit::Lutz Class Reference

Break-point track fitter. More...

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-01-15/BreakPointFitter/func/Lutz.h"

Public Member Functions
	Lutz (unsigned int n, unsigned int N)
	Construct a fitter for n measurements and N scattering planes. More...
void	SetMeasurement (unsigned int i, double z, double x, double sigx)
	Add measurements. More...
void	SetScatteringPlane (unsigned int j, double Z, double sigSJ)
	Set Z location of a scattering plane j. More...
void	Fit (double *a, double *b, double *alpha, double *chi2)
	Do the fit. More...
double	Chi2XIi (unsigned int i) const
	After fit, contribution to chi^2 of the ith measurement. More...
double	Chi2BetaJ (unsigned int j) const
	After fit, contribution to chi^2 of the jth scattering plane. More...
double	Chi2XI () const
double	Chi2Beta () const
double	X (double z) const
	After fit, best-fit x location of track at location z. More...

Private Member Functions
void	CalcZeta ()
void	CalcD ()
bool	CalcAinv ()
void	Calcq ()
void	Calcx ()
void	NullResult ()

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

Detailed Description

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 ...................................................................

Definition at line 20 of file Lutz.h.

Constructor & Destructor Documentation

bpfit::Lutz::Lutz	(	unsigned int	n,
		unsigned int	N
	)

Construct a fitter for n measurements and N scattering planes.

Parameters

n	- number of measurements to be added
N	- number of scattering surfaces to be added

Definition at line 13 of file Lutz.cxx.

References fZeta, MECModelEnuComparisons::i, bpfit::kNn, getGoodRuns4SAM::n, and stan::math::resize().

                                           :
    fn(n),
    fN(N),
    fz(n),
    fXI(n),
    fwx(n),
    fZ(N),
    fwSJ(N),
    fZeta(n),
    fD(N+2,n),
    fAinv(N+2),
    fa(0),
    fb(0),
    fAlpha(fN),
    fx(n)
  {
    if(N >= n)
      throw BPException(__FILE__, __LINE__, kNn);

    for(unsigned int i = 0; i < n; i++)
      fZeta[i].resize(N);
  }

Member Function Documentation

bool bpfit::Lutz::CalcAinv ( )

private

Inverting the matrix Sometimes this operation fails, producing error: matrix is singular Thus we need to put it into try/catch block to avoid the hard exit.

Definition at line 100 of file Lutz.cxx.

References om::cerr, fillBadChanDBTables::det, allTimeWatchdog::endl, fAinv, fn, fN, fwSJ, fwx, fz, fZeta, MECModelEnuComparisons::i, genie::utils::mec::J(), K, and sum.

Referenced by Fit().

  {
    unsigned int i, J, K;
    
    double sum00 = 0;
    double sum01 = 0;
    double sum11 = 0;
    double fwxfz;
    for (i=0; i<fn; ++i) {
      fwxfz  = fwx[i]*fz[i]; 
      sum00 += fwx[i];
      sum01 += fwxfz;
      sum11 += fwxfz*fz[i];
    }
    fAinv[0][0] = sum00;
    fAinv[0][1] = sum01;
    fAinv[1][0] = sum01;
    fAinv[1][1] = sum11;
    
    double sum0J2;
    double sum1J2;
    double fwxZeta;
    for (J=0; J<fN; ++J) {
      sum0J2 = 0;
      sum1J2 = 0;
      for (i=0; i<fn; ++i) {
        fwxZeta = fwx[i]*fZeta[i][J];
        sum0J2 += fwxZeta;
        sum1J2 += fz[i]*fwxZeta;
      }
      fAinv[0][J+2] = sum0J2;
      fAinv[J+2][0] = sum0J2;
      fAinv[1][J+2] = sum1J2;
      fAinv[J+2][1] = sum1J2;
    }
    
    double sum;
    for (J=0; J<fN; ++J) {
      for (K=J; K<fN; ++K) {
        sum = 0;
        for (i=0; i<fn; ++i) sum += fwx[i]*fZeta[i][J]*fZeta[i][K];
        if (J==K) sum += fwSJ[J];
        fAinv[J+2][K+2] = sum;
        fAinv[K+2][J+2] = sum;
      }
    }
    
    double det = 0.0;

    /// Inverting the matrix
    /// Sometimes this operation fails, producing error: matrix is singular
    /// Thus we need to put it into try/catch block to avoid the hard exit.
    try {
      fAinv.InvertFast(&det);
    }
    catch(...) {
      std::cerr << __FILE__ << ":" << __LINE__
                << " Matrix inversion failed" << std::endl;
      return false;
    } // end of inverting the matrix
    
    if (det<=0.0) {
      std::cerr << __FILE__ << ":" << __LINE__ 
                << " Matrix inversion failed" << std::endl;
      return false;
    }
    return true;
  }

void bpfit::Lutz::CalcD ( )

private

Definition at line 84 of file Lutz.cxx.

References fD, fn, fN, fwx, fz, fZeta, MECModelEnuComparisons::i, and calib::j.

Referenced by Fit().

  {
    unsigned int i, j;
    for (i=0; i<fn; ++i) {
      fD[0][i] = fwx[i];
      fD[1][i] = fwx[i]*fz[i];
      for (j=0; j<fN; ++j) {
        fD[j+2][i] = fwx[i]*fZeta[i][j];
      }
    }
  }

void bpfit::Lutz::Calcq ( )

private

Definition at line 173 of file Lutz.cxx.

References fa, fAinv, fAlpha, fb, fD, fN, fXI, and calib::j.

Referenced by Fit().

  {
    unsigned int j;
    //
    // Solve according to equation 35. Don't bother with the second term
    // which is by construction zero.
    //
    TVectorT<double> q(fXI);
    q *= fD;
    q *= fAinv;
    
    //
    // Unpack the resutls
    //
    fa = q[0];
    fb = q[1];
    for (j=0; j<fN; ++j) fAlpha[j] = q[2+j];
  }

void bpfit::Lutz::Calcx ( )

private

Definition at line 194 of file Lutz.cxx.

References fa, fAlpha, fb, fn, fN, fx, fz, fZeta, MECModelEnuComparisons::i, and calib::j.

Referenced by Fit().

  {
    unsigned int i, j;
    for (i=0; i<fn; ++i) {
      fx[i] = fa + fb*fz[i];
      for (j=0; j<fN; ++j) {
        if (fZeta[i][j]>0.0) fx[i] += fAlpha[j]*fZeta[i][j];
      }
    }
  }

void bpfit::Lutz::CalcZeta ( )

private

Definition at line 69 of file Lutz.cxx.

References fn, fN, fz, fZ, fZeta, MECModelEnuComparisons::i, and calib::j.

Referenced by Fit().

  {
    unsigned int i, j;
    for (i=0; i<fn; ++i) {
      for (j=0; j<fN; ++j) {
        fZeta[i][j] = fz[i]-fZ[j];
        if (fZeta[i][j]<0.0) fZeta[i][j] = 0.0;
      }
    }
  }

double bpfit::Lutz::Chi2Beta

(

)

const

Definition at line 287 of file Lutz.cxx.

References Chi2BetaJ(), fN, calib::j, and sum.

Referenced by bpfit::BreakPoint::FillTrackNt(), and Fit().

  {
    unsigned int j;
    double sum = 0.0;
    for (j=0; j<fN; ++j) sum += this->Chi2BetaJ(j);
    return sum;
  }

double bpfit::Lutz::Chi2BetaJ

(

unsigned int

j

)

const

After fit, contribution to chi^2 of the jth scattering plane.

Definition at line 264 of file Lutz.cxx.

References fAlpha, fwSJ, and calib::j.

Referenced by Chi2Beta().

  {
    double chi;
    
    chi = 0.0-fAlpha[j];
    chi *= chi;
    chi *= fwSJ[j];
    
    return chi;
  }

double bpfit::Lutz::Chi2XI

(

)

const

Definition at line 277 of file Lutz.cxx.

References Chi2XIi(), fn, MECModelEnuComparisons::i, and sum.

Referenced by bpfit::BreakPoint::FillTrackNt(), and Fit().

  {
    unsigned int i;
    double sum = 0.0;
    for (i=0; i<fn; ++i) sum += this->Chi2XIi(i);
    return sum;
  }

double bpfit::Lutz::Chi2XIi

(

unsigned int

i

)

const

After fit, contribution to chi^2 of the ith measurement.

Definition at line 251 of file Lutz.cxx.

References fwx, fx, fXI, and MECModelEnuComparisons::i.

Referenced by Chi2XI().

  {
    double chi;
    
    chi  = fXI[i]-fx[i];
    chi *= chi;
    chi *= fwx[i];
    
    return chi;
  }

void bpfit::Lutz::Fit	(	double *	a,
		double *	b,
		double *	alpha,
		double *	chi2
	)

Do the fit.

Parameters

a	- on return, fitted intercept. OK to pass 0.
b	- on return, fitted initial slope. OK to pass 0.
alpha	- on return, fitted scattering slopes. OK to pass 0.
chi2	- on return, the chi^2 result of the fit. OK to pass 0.

Definition at line 221 of file Lutz.cxx.

References CalcAinv(), CalcD(), Calcq(), Calcx(), CalcZeta(), Chi2Beta(), Chi2XI(), fa, fAlpha, fb, fN, calib::j, and NullResult().

  {
    bool ok;
    unsigned int j;
    //
    // Make the fit
    //
    this->CalcZeta();
    this->CalcD();
    ok = this->CalcAinv();
    if (ok) {
      this->Calcq();
      this->Calcx();    
    }
    else this->NullResult();
    //
    // Unpack the results
    //
    if (a!=0) *a = fa;
    if (b!=0) *b = fb;
    if (alpha!=0) {
      for (j=0; j<fN; ++j) alpha[j] = fAlpha[j];
    }
    if (chi2!=0) {
      *chi2 = this->Chi2XI() + this->Chi2Beta();
    }
  }

void bpfit::Lutz::NullResult ( )

private

Definition at line 297 of file Lutz.cxx.

References fa, fAlpha, fb, fn, fN, fx, MECModelEnuComparisons::i, and calib::j.

Referenced by Fit().

  {
    unsigned int i, j;
    fa = 0.0;
    fb = 0.0;
    for (i=0; i<fn; ++i) fx[i]     = 0.0;
    for (j=0; j<fN; ++j) fAlpha[j] = 0.0;
  }

void bpfit::Lutz::SetMeasurement	(	unsigned int	i,
		double	z,
		double	x,
		double	sigx
	)

Add measurements.

Parameters

i	- measurement index [0...n)
z	- z location of measurement (longitudinal location)
x	- x location of measurement (transverse location)
sigx	- measurement error on x position

Definition at line 38 of file Lutz.cxx.

References fn, fwx, fXI, fz, MECModelEnuComparisons::i, bpfit::kINDEX, bpfit::kNAN, bpfit::kSIGX, submit_syst::x, and test::z.

Referenced by bpfit::BreakPoint::FitTracks(), and bpfit::DimuonFitter::FitView().

  {
    if(i >= fn)
      throw BPException(__FILE__, __LINE__, kINDEX);
    if(z != z || x != x || sigx != sigx)
      throw BPException(__FILE__, __LINE__, kNAN);
    if(sigx <= 0.0)
      throw BPException(__FILE__, __LINE__, kSIGX, sigx);

    fz[i]  = z;
    fXI[i] = x;
    fwx[i] = 1.0/(sigx*sigx);
  }

void bpfit::Lutz::SetScatteringPlane	(	unsigned int	j,
		double	Z,
		double	sigSJ
	)

Set Z location of a scattering plane j.

Parameters

j	- index of scattering plane (start from 0)
Z	- longitudinal location of scattering plane
sigSJ	- RMS scattering of track slope at this plane

Definition at line 54 of file Lutz.cxx.

References fN, fwSJ, fZ, calib::j, bpfit::kINDEX, bpfit::kNAN, bpfit::kSIGSJ, and Z.

  {
    if(j >= fN)
      throw BPException(__FILE__, __LINE__, kINDEX);
    if(Z != Z || sigSJ != sigSJ)
      throw BPException(__FILE__, __LINE__, kNAN);
    if(sigSJ <= 0.0)
      throw BPException(__FILE__, __LINE__, kSIGSJ, sigSJ);
    
    fZ[j]   = Z;
    fwSJ[j] = 1.0/(sigSJ*sigSJ);
  }

double bpfit::Lutz::X

(

double

z

)

const

After fit, best-fit x location of track at location z.

Definition at line 207 of file Lutz.cxx.

References fa, fAlpha, fb, fN, fZ, calib::j, submit_syst::x, and test::z.

Referenced by bpfit::BreakPoint::LogResiduals().

  {
    unsigned int j;
    double zeta;
    double x = fa + fb*z;
    for (j=0; j<fN; ++j) {
      zeta = z-fZ[j];
      if (zeta>0.0) x += fAlpha[j]*zeta;
    }
    return x;
  }

Member Data Documentation

double bpfit::Lutz::fa

private

Fitted intercept.

Definition at line 100 of file Lutz.h.

Referenced by Calcq(), Calcx(), Fit(), NullResult(), and X().

TMatrixTSym<double> bpfit::Lutz::fAinv

private

System of equations, Eqn. 29.

Definition at line 99 of file Lutz.h.

Referenced by CalcAinv(), and Calcq().

std::vector<double> bpfit::Lutz::fAlpha

private

Fitted scatters.

Definition at line 102 of file Lutz.h.

Referenced by Calcq(), Calcx(), Chi2BetaJ(), Fit(), NullResult(), and X().

double bpfit::Lutz::fb

private

Fitted slope.

Definition at line 101 of file Lutz.h.

Referenced by Calcq(), Calcx(), Fit(), NullResult(), and X().

TMatrixT<double> bpfit::Lutz::fD

private

Eqn. 31.

Definition at line 98 of file Lutz.h.

Referenced by CalcD(), and Calcq().

unsigned int bpfit::Lutz::fn

private

Number of measurements.

Definition at line 90 of file Lutz.h.

Referenced by CalcAinv(), CalcD(), Calcx(), CalcZeta(), Chi2XI(), NullResult(), and SetMeasurement().

unsigned int bpfit::Lutz::fN

private

Number of scattering planes.

Definition at line 91 of file Lutz.h.

Referenced by CalcAinv(), CalcD(), Calcq(), Calcx(), CalcZeta(), Chi2Beta(), Fit(), NullResult(), SetScatteringPlane(), and X().

std::vector<double> bpfit::Lutz::fwSJ

private

Scattering angle weight (1/sigSJ^2)

Definition at line 96 of file Lutz.h.

Referenced by CalcAinv(), Chi2BetaJ(), and SetScatteringPlane().

std::vector<double> bpfit::Lutz::fwx

private

Measurement weight (1/sigx^2)

Definition at line 94 of file Lutz.h.

Referenced by CalcAinv(), CalcD(), Chi2XIi(), and SetMeasurement().

std::vector<double> bpfit::Lutz::fx

private

Fitted track locations.

Definition at line 103 of file Lutz.h.

Referenced by Calcx(), Chi2XIi(), and NullResult().

TVectorT<double> bpfit::Lutz::fXI

private

Measurements in transverse location.

Definition at line 93 of file Lutz.h.

Referenced by Calcq(), Chi2XIi(), and SetMeasurement().

std::vector<double> bpfit::Lutz::fz

private

Measurement plane longitudinal locations.

Definition at line 92 of file Lutz.h.

Referenced by CalcAinv(), CalcD(), Calcx(), CalcZeta(), and SetMeasurement().

std::vector<double> bpfit::Lutz::fZ

private

Location of scattering planes.

Definition at line 95 of file Lutz.h.

Referenced by CalcZeta(), SetScatteringPlane(), and X().

std::vector< std::vector<double> > bpfit::Lutz::fZeta

private

measurement-scattering plane distances

Definition at line 97 of file Lutz.h.

Referenced by CalcAinv(), CalcD(), Calcx(), CalcZeta(), and Lutz().

---

The documentation for this class was generated from the following files:

• /cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-01-15/BreakPointFitter/func/Lutz.h
• /cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-01-15/BreakPointFitter/func/Lutz.cxx