[u/mrichter/AliRoot.git] / MUON / AliMillepede.h

#ifndef ALI_MILLEPEDE_H 
#define ALI_MILLEPEDE_H

/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 * See cxx source for full Copyright notice                               */

/* $Id$ */

/// \ingroup rec
/// \class AliMillepede
/// \brief Class for detector alignment - modified C++ implementation of original
/// millepede package
//
//  Author: Javier Castillo

#include <TObject.h>
#include <TArrayI.h>
#include <TArrayD.h>

class AliMillepede : public TObject {

public: 
  /// Standard constructor
  AliMillepede();

  virtual ~AliMillepede( ); /// Destructor

  /// Initialization
  virtual Int_t InitMille(int nglo, int nloc, int nstd,			  
			  double lResCut, double lResCutInit);
  virtual Int_t GlobalFit(double par[], double error[], double pull[]);
  virtual Int_t SetGlobalParameters(double *param);
  virtual Int_t SetGlobalParameter(int iPar, double vPar);
  virtual Int_t SetParSigma(int iPar, double sigma);
  virtual Int_t SetNonLinear(int index);
  virtual Int_t SetGlobalConstraint(double dercs[], double rhs);
  virtual Int_t SetLocalEquation(double dergb[], double derlc[], double rmeas, double sigma);
  virtual Int_t LocalFit(int n, double localParams[], Bool_t bSingleFit);
  virtual Int_t GetNLocalEquations() const {return fNLocalEquations;};
  virtual void  SetNLocalEquations(int value) {fNLocalEquations = value;};
  virtual Int_t PrintGlobalParameters() const;
  virtual Int_t SetIterations (double cutfac);
  virtual Double_t GetParError(int iPar) const;

private:

// Max. dimensions

  static const int fgkMaxGlobalPar  = 1000; // Max. number of global parameters
  static const int fgkMaxLocalPar   = 20;  // Max. number of local parameters
  static const int fgkMaxGloCsts    = 10;  // Max. number of constraint equations
  static const int fgkMaxGloPC	    = 1010; // fgkMaxGlobalPar+fgkMaxGloCsts


// Private methods 

  Double_t GetParCorrelation(int i, int j);

  int   SpmInv(double v[][fgkMaxGloPC], double b[], int n);
  int   SpmInv(double v[][fgkMaxLocalPar], double b[], int n);
  int SpAVAt(double v[][fgkMaxLocalPar], double a[][fgkMaxLocalPar], double w[][fgkMaxGlobalPar], int n, int m);
  int SpAX(double a[][fgkMaxLocalPar], double x[], double y[], int n, int m);
  double Chi2DoFLim(int n, int nd);

// Matrices

  double fMatCGlo[fgkMaxGloPC][fgkMaxGloPC];            // Matrix C global
  double fMatCLoc[fgkMaxLocalPar][fgkMaxLocalPar];      // Matrix C local
  double fMatCGloLoc[fgkMaxGlobalPar][fgkMaxLocalPar];  // Rectangular matrix C g*l
  double fMatCGloCorr[fgkMaxGlobalPar][fgkMaxGlobalPar];// Correction of matrix C global
  double fMatDerConstr[fgkMaxGloCsts][fgkMaxGlobalPar]; // Constrained derivatives

// Vectors and useful variables

  double fDiagCGlo[fgkMaxGloPC];        // Initial diagonal elements of C global matrix
  double fVecBGlo[fgkMaxGloPC];         // Vector B global (parameters) 
  double fVecBGloCorr[fgkMaxGlobalPar]; // Correction of vector B global
  double fVecBLoc[fgkMaxLocalPar];      // Vector B local (parameters) 

  double fInitPar[fgkMaxGlobalPar];    // Initial global parameters
  double fDeltaPar[fgkMaxGlobalPar];   // Variation of global parameters 
  double fSigmaPar[fgkMaxGlobalPar];   // Sigma of allowed variation of global parameter 

  double fLagMult[fgkMaxGloCsts];   // Lagrange multipliers of constrained equations

  Bool_t fIsNonLinear[fgkMaxGlobalPar]; // Flag for non linear parameters
  int   fGlo2CGLRow[fgkMaxGlobalPar];   // Global parameter to row in "used" g*l matrix
  int   fCGLRow2Glo[fgkMaxGlobalPar];   // Row in "used" g*l matrix to global parameter 

  TArrayI fIndexLocEq; // Table of parameter indexes in local equation 
  TArrayD fDerivLocEq; // Table of local equation derivatives wrt. parameter 
  TArrayI fIndexAllEqs; // Index in all loc. eq. storage for iterations
  TArrayD fDerivAllEqs; // derivative in all loc. eq. storage for iterations
  TArrayI fLocEqPlace;  // Loc. Eq. position in AllEqs storage

  Int_t fNIndexLocEq;   // Number of entries in fIndexLocEq
  Int_t fNDerivLocEq;   // Number of entries in fDerivLocEq
  Int_t fNIndexAllEqs;  // Number of entries in fIndexAllEqs
  Int_t fNDerivAllEqs;  // Number of entries in fDerivAllEqs
  Int_t fNLocEqPlace;   // Number of entries in fLocEqPlace
  

  int    fNLocalEquations;       // Number of local equations 
  double fResCutInit;            // Cut in residual for first iterartion
  double fResCut;                // Cut in residual for other iterartiona

  double fChi2CutFactor;         // Cut factor for chi2 cut to accept local fit 
  double fChi2CutRef;            // Reference cut for chi2 cut to accept local fit 

  int fIter;                   // Current iteration
  int fMaxIter;                // Maximum number of iterations
  int fNStdDev;                // Number of standard deviations for chi2 cut 
  int fNGlobalConstraints;     // Number of constraint equations
  int fNLocalFits;             // Number of local fits
  int fNLocalFitsRejected;     // Number of local fits rejected
  int fNGlobalPar;             // Number of global parameters
  int fNLocalPar;              // Number of local parameters

  ClassDef(AliMillepede, 0)  // Millepede Class
};

#endif
Commit	Line	Data
010eb601	1	#ifndef ALI_MILLEPEDE_H
	2	#define ALI_MILLEPEDE_H
	3
	4	/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	5	* See cxx source for full Copyright notice */
	6
	7	/* $Id$ */
	8
	9	/// \ingroup rec
	10	/// \class AliMillepede
	11	/// \brief Class for detector alignment - modified C++ implementation of original
	12	/// millepede package
	13	//
	14	// Author: Javier Castillo
	15
	16	#include <TObject.h>
	17	#include <TArrayI.h>
	18	#include <TArrayD.h>
	19
	20	class AliMillepede : public TObject {
	21
	22	public:
	23	/// Standard constructor
	24	AliMillepede();
	25
	26	virtual ~AliMillepede( ); /// Destructor
	27
	28	/// Initialization
	29	virtual Int_t InitMille(int nglo, int nloc, int nstd,
	30	double lResCut, double lResCutInit);
	31	virtual Int_t GlobalFit(double par[], double error[], double pull[]);
	32	virtual Int_t SetGlobalParameters(double *param);
	33	virtual Int_t SetGlobalParameter(int iPar, double vPar);
	34	virtual Int_t SetParSigma(int iPar, double sigma);
	35	virtual Int_t SetNonLinear(int index);
	36	virtual Int_t SetGlobalConstraint(double dercs[], double rhs);
	37	virtual Int_t SetLocalEquation(double dergb[], double derlc[], double rmeas, double sigma);
	38	virtual Int_t LocalFit(int n, double localParams[], Bool_t bSingleFit);
	39	virtual Int_t GetNLocalEquations() const {return fNLocalEquations;};
	40	virtual void SetNLocalEquations(int value) {fNLocalEquations = value;};
	41	virtual Int_t PrintGlobalParameters() const;
	42	virtual Int_t SetIterations (double cutfac);
	43	virtual Double_t GetParError(int iPar) const;
	44
	45	private:
	46
	47	// Max. dimensions
	48
	49	static const int fgkMaxGlobalPar = 1000; // Max. number of global parameters
	50	static const int fgkMaxLocalPar = 20; // Max. number of local parameters
	51	static const int fgkMaxGloCsts = 10; // Max. number of constraint equations
	52	static const int fgkMaxGloPC = 1010; // fgkMaxGlobalPar+fgkMaxGloCsts
	53
	54
	55	// Private methods
	56
	57	Double_t GetParCorrelation(int i, int j);
	58
	59	int SpmInv(double v[][fgkMaxGloPC], double b[], int n);
	60	int SpmInv(double v[][fgkMaxLocalPar], double b[], int n);
	61	int SpAVAt(double v[][fgkMaxLocalPar], double a[][fgkMaxLocalPar], double w[][fgkMaxGlobalPar], int n, int m);
	62	int SpAX(double a[][fgkMaxLocalPar], double x[], double y[], int n, int m);
	63	double Chi2DoFLim(int n, int nd);
	64
65	// Matrices
66
67	double fMatCGlo[fgkMaxGloPC][fgkMaxGloPC]; // Matrix C global
68	double fMatCLoc[fgkMaxLocalPar][fgkMaxLocalPar]; // Matrix C local
69	double fMatCGloLoc[fgkMaxGlobalPar][fgkMaxLocalPar]; // Rectangular matrix C g*l
70	double fMatCGloCorr[fgkMaxGlobalPar][fgkMaxGlobalPar];// Correction of matrix C global
71	double fMatDerConstr[fgkMaxGloCsts][fgkMaxGlobalPar]; // Constrained derivatives
72
73	// Vectors and useful variables
74
75	double fDiagCGlo[fgkMaxGloPC]; // Initial diagonal elements of C global matrix
76	double fVecBGlo[fgkMaxGloPC]; // Vector B global (parameters)
77	double fVecBGloCorr[fgkMaxGlobalPar]; // Correction of vector B global
78	double fVecBLoc[fgkMaxLocalPar]; // Vector B local (parameters)
79
80	double fInitPar[fgkMaxGlobalPar]; // Initial global parameters
81	double fDeltaPar[fgkMaxGlobalPar]; // Variation of global parameters
82	double fSigmaPar[fgkMaxGlobalPar]; // Sigma of allowed variation of global parameter
83
84	double fLagMult[fgkMaxGloCsts]; // Lagrange multipliers of constrained equations
85
86	Bool_t fIsNonLinear[fgkMaxGlobalPar]; // Flag for non linear parameters
87	int fGlo2CGLRow[fgkMaxGlobalPar]; // Global parameter to row in "used" g*l matrix
88	int fCGLRow2Glo[fgkMaxGlobalPar]; // Row in "used" g*l matrix to global parameter
89
90	TArrayI fIndexLocEq; // Table of parameter indexes in local equation
91	TArrayD fDerivLocEq; // Table of local equation derivatives wrt. parameter
92	TArrayI fIndexAllEqs; // Index in all loc. eq. storage for iterations
93	TArrayD fDerivAllEqs; // derivative in all loc. eq. storage for iterations
94	TArrayI fLocEqPlace; // Loc. Eq. position in AllEqs storage
95
96	Int_t fNIndexLocEq; // Number of entries in fIndexLocEq
97	Int_t fNDerivLocEq; // Number of entries in fDerivLocEq
98	Int_t fNIndexAllEqs; // Number of entries in fIndexAllEqs
99	Int_t fNDerivAllEqs; // Number of entries in fDerivAllEqs
100	Int_t fNLocEqPlace; // Number of entries in fLocEqPlace
101
102
103	int fNLocalEquations; // Number of local equations
104	double fResCutInit; // Cut in residual for first iterartion
105	double fResCut; // Cut in residual for other iterartiona
106
107	double fChi2CutFactor; // Cut factor for chi2 cut to accept local fit
108	double fChi2CutRef; // Reference cut for chi2 cut to accept local fit
109
110	int fIter; // Current iteration
111	int fMaxIter; // Maximum number of iterations
112	int fNStdDev; // Number of standard deviations for chi2 cut
113	int fNGlobalConstraints; // Number of constraint equations
114	int fNLocalFits; // Number of local fits
115	int fNLocalFitsRejected; // Number of local fits rejected
116	int fNGlobalPar; // Number of global parameters
117	int fNLocalPar; // Number of local parameters
118
119	ClassDef(AliMillepede, 0) // Millepede Class
120	};
121
122	#endif
123
124
125