[u/mrichter/AliRoot.git] / ITS / UPGRADE / AliITSUSeed.h

#ifndef ALIITSUSEED_H
#define ALIITSUSEED_H

#include "AliExternalTrackParam.h"
#include "AliITSUAux.h"
class AliESDtrack;

using namespace AliITSUAux;


class AliITSUSeed: public AliExternalTrackParam
{
 public:
  enum {kKilled=BIT(14),kSave=BIT(15),kMark=BIT(16),kFake=BIT(17)};
  enum {kF02,kF04,kF12,kF13,kF14,kF24, kF44,kNFElem}; // non-trivial elems of propagation matrix
  enum {kK00,kK01,kK10,kK11,kK20,kK21,kK30,kK31,kK40,kK41, kNKElem}; // non-trivial elems of gain matrix
  enum {kS00,kS10,kS11,kS20,kS21,kS22,kS30,kS31,kS32,kS33,kS40,kS41,kS42,kS43,kS44,kNSElem}; // elements of 5x5 sym matrix
  enum {kR00,kR22,kNRElem};                        // non trivial elements of rotation matrix
  //
  AliITSUSeed();
  AliITSUSeed(const AliITSUSeed& src);
  AliITSUSeed &operator=(const AliITSUSeed &src);
  virtual ~AliITSUSeed();
  virtual void    Print(Option_t* option = "") const;
  //
  void            SetLrClusterID(Int_t lr, Int_t cl);
  void            SetLr(Int_t lr)                        {SetLrClusterID(lr,-1);} // lr w/o cluster
  void            SetLrClusterID(UInt_t id)              {fClID = id;}
  void            SetParent(TObject* par)                {fParent = par;}
  void            SetChi2Cl(Double_t v)                  {fChi2Cl= v; v>0 ? fChi2Glo+=v : fChi2Penalty -= v;}
  void            Kill(Bool_t v=kTRUE)                   {SetBit(kKilled, v);}
  void            SetFake(Bool_t v=kTRUE)                {SetBit(kFake, v);}
  void            Save(Bool_t v=kTRUE)                   {SetBit(kSave,v);}
  void            SetMark(Bool_t v=kTRUE)                {SetBit(kMark,v);}
  void            FlagTree(UInt_t bits, Bool_t v=kTRUE);
  void            SetChi2ITSTPC(Float_t v)               {fChi2Match = v;}
  void            SetChi2ITSSA(Float_t v)                {fChi2ITSSA = v;}
  //
  UInt_t          GetLrClusterID()                 const {return fClID;}
  Int_t           GetLrCluster(Int_t &lr)          const {return UnpackCluster(fClID,lr);}
  Int_t           GetLayerID()                     const {return UnpackLayer(fClID);}
  Int_t           GetClusterID()                   const {return UnpackCluster(fClID);}
  Bool_t          HasClusterOnLayer(Int_t lr)      const {return fHitsPattern&(0x1<<lr);}
  Bool_t          HasCluster()                     const {return IsCluster(fClID);}
  Int_t           GetNLayersHit()                  const {return NumberOfBitsSet(fHitsPattern);}
  Int_t           GetNClusters()                   const;
  Int_t           GetClusterIndex(Int_t ind)       const;
  UShort_t        GetHitsPattern()                 const {return fHitsPattern;}
  Float_t         GetChi2Cl()                      const {return fChi2Cl;}
  Float_t         GetChi2Glo()                     const {return fChi2Glo;}
  Float_t         GetChi2Penalty()                 const {return fChi2Penalty;}
  Float_t         GetChi2GloNrm()                  const;
  Float_t         GetChi2ITSTPC()                  const {return fChi2Match;}
  Float_t         GetChi2ITSSA()                   const {return fChi2ITSSA;}
  Float_t         GetQualityVar()                  const;
  Bool_t          IsKilled()                       const {return TestBit(kKilled);}
  Bool_t          IsFake()                         const {return TestBit(kFake);}
  Bool_t          IsSaved()                        const {return TestBit(kSave);}
  Bool_t          IsMarked()                       const {return TestBit(kMark);}
  Bool_t          ContainsFake()                   const;
  Int_t           FetchClusterInfo(Int_t *clIDarr) const;
  //
  Int_t           GetPoolID()                      const {return int(GetUniqueID())-1;}
  void            SetPoolID(Int_t id)                    {SetUniqueID(id+1);}
  //
  TObject*        GetParent()                      const {return fParent;}
  const AliITSUSeed* GetParent(Int_t lr)              const;
  //
  virtual Bool_t  IsSortable()                     const {return kTRUE;}
  virtual Bool_t  IsEqual(const TObject* obj)      const;
  virtual Int_t	  Compare(const TObject* obj)      const;
  //
  // test
  void            InitFromSeed(const AliExternalTrackParam* seed);
  void            ResetFMatrix();
  void            ApplyELoss2FMatrix(Double_t frac, Bool_t beforeProp);
  Bool_t          ApplyMaterialCorrection(Double_t xOverX0, Double_t xTimesRho, Double_t mass, Bool_t beforeProp);
  Bool_t          PropagateToX(Double_t xk, Double_t b);
  Bool_t          RotateToAlpha(Double_t alpha);
  Bool_t          GetTrackingXAtXAlpha(double xOther,double alpOther,double bz, double &x);
  Double_t        GetPredictedChi2(Double_t p[2],Double_t cov[3]);
  Bool_t          Update();
  Bool_t          Smooth(Double_t vecL[5],Double_t matL[15]);
  Double_t*       ProdABA(const double a[15],const double b[15]) const;
  //
  UInt_t          GetNChildren()                                 const {return fNChildren;}
  Int_t           IncChildren()                                        {return ++fNChildren;}
  Int_t           DecChildren()                                        {return --fNChildren;}
  Double_t        GetResidY()                                    const {return fResid[0];}
  Double_t        GetResidZ()                                    const {return fResid[1];}
  //
#ifdef _ITSU_TUNING_MODE_ // this is used only for tuning histo filling
  void                  SetOrdBranch(UShort_t v) {fOrdBranch=v;}
  void                  SetOrdCand(UShort_t v)   {fOrdCand=v;}
  Int_t                 GetOrdBranch()           const {return fOrdBranch;}
  Int_t                 GetOrdCand()             const {return fOrdCand;}
  Double_t              GetPullY()               const {return fResid[0]*Sqrt(fCovIYZ[0]);}
  Double_t              GetPullZ()               const {return fResid[1]*Sqrt(fCovIYZ[2]);}
#endif
  //
 protected:
  //
  UShort_t              fHitsPattern;       // bit pattern of hits
  UShort_t              fNChildren;         // number of children (prolongations)
  UInt_t                fClID;              // packed cluster info (see AliITSUAux::PackCluster)
  Float_t               fChi2Glo;           // current chi2 global (sum of track-cluster chi2's on layers with hit)
  Float_t               fChi2Cl;            // track-cluster chi2 (if >0) or penalty for missing cluster (if < 0)
  Float_t               fChi2Penalty;       // total penalty (e.g. for missing clusters)
  Float_t               fChi2Match;         // ITS/TPC matching chi2 (per NDF)  // RS: to move to separate object of final seed
  Float_t               fChi2ITSSA;         // ITSSA backward fit chi2 (per NDF) // RS: to move to separate object of final seed
  Double_t              fFMatrix[kNFElem];  // matrix of propagation from prev layer (non-trivial elements)
  Double_t              fKMatrix[kNKElem];  // Gain matrix non-trivial elements (note: standard MBF formula uses I-K*H)
  Double_t              fRMatrix[kNRElem];  // rotation matrix non-trivial elements
  Double_t              fCovIYZ[3];         // inverted matrix of propagation + meas errors = [Hi * Pi|i-1 * Hi^T + Ri]^-1
  Double_t              fResid[2];          // residuals vector
  TObject*              fParent;            // parent track (in higher tree hierarchy)
  //
#ifdef _ITSU_TUNING_MODE_ // this is used only for tuning histo filling
  UShort_t              fOrdBranch;         //! order in the seed branching
  UShort_t              fOrdCand;           //! order in candidates tree
#endif
  //
  ClassDef(AliITSUSeed,1)
};

//_________________________________________________________________________
inline void AliITSUSeed::SetLrClusterID(Int_t lr, Int_t cl)
{
  // assign layer, cluster (if -1 - no hit on this layer)
  fClID = PackCluster(lr,cl);
  if (cl>=0) fHitsPattern |= 0x1<<lr;
}

//_________________________________________________________________________
inline void AliITSUSeed::ResetFMatrix()
{
  // reset transport matrix
  fFMatrix[kF02] = fFMatrix[kF04] = fFMatrix[kF12] = fFMatrix[kF13] = fFMatrix[kF14] = fFMatrix[kF24] = 0;
  fFMatrix[kF44] = 1.0;  // this element accumulates eloss 
}

//_________________________________________________________________________
inline Bool_t AliITSUSeed::ApplyMaterialCorrection(Double_t xOverX0, Double_t xTimesRho, Double_t mass, Bool_t beforeProp)
{
  // apply material correction and modify transport matrix
  double pold = Get1P();
  if (!CorrectForMeanMaterial(xOverX0,xTimesRho,mass)) return kFALSE;
  ApplyELoss2FMatrix( Get1P()/pold, beforeProp);
  return kTRUE;
}


//_________________________________________________________________________
inline void AliITSUSeed::ApplyELoss2FMatrix(Double_t frac, Bool_t beforeProp)
{
  // Accounts for the energy loss in the transport matrix
  // equivalent to multiplying Fmatix by E=diag{1,1,1,1,P4new/P4old}, where P4 is the 1/pt param.
  // If beforeProp is true, then it is assumed that the eloss was applied before the transport,
  // i.e. F' = F * E, otherwise, after transport, F' = E * F
  fFMatrix[kF44] *= frac;
  if (beforeProp) {
    fFMatrix[kF04] *= frac;
    fFMatrix[kF14] *= frac;
    fFMatrix[kF24] *= frac;
  }
}

//_________________________________________________________________________
inline const AliITSUSeed* AliITSUSeed::GetParent(Int_t lr) const
{
  // get parent at given layer
  const AliITSUSeed* par=this;
  int lrt;
  while( par && (lrt=par->GetLayerID())>=0 ) {
    if (lrt==lr) break;
    par = dynamic_cast<const AliITSUSeed*>(par->GetParent());
  }
  return par;
}

//__________________________________________________________________
inline Int_t AliITSUSeed::GetNClusters() const
{
  // count number of clusters (some layers may have >1 hit)
  int ncl = 0;
  const AliITSUSeed* seed = this;
  while(seed) {
    if (seed->HasCluster()) ncl++;
    seed = (AliITSUSeed*)seed->GetParent();
  }
  return ncl;
  //
}

//__________________________________________________________________
inline void AliITSUSeed::FlagTree(UInt_t bits, Bool_t v)
{
  // set bits on all tree levels
  AliITSUSeed* seed = this;
  do {seed->SetBit(bits,v);} while ( (seed=(AliITSUSeed*)seed->GetParent()) );
  //
}

//__________________________________________________________________
inline Float_t AliITSUSeed::GetQualityVar() const
{
  // return variable for quality check (the smaller - the better)
  //  return GetChi2ITSSA()/Pt(); // this may change! RS
  return GetChi2GloNrm()/Pt();
}


#endif
Commit	Line	Data
32d38de2	1	#ifndef ALIITSUSEED_H
	2	#define ALIITSUSEED_H
	3
	4	#include "AliExternalTrackParam.h"
f8832015	5	#include "AliITSUAux.h"
3e4e3c23	6	class AliESDtrack;
3e4e3c23	7
f8832015	8	using namespace AliITSUAux;
32d38de2	9
	10
	11	class AliITSUSeed: public AliExternalTrackParam
	12	{
	13	public:
08419930	14	enum {kKilled=BIT(14),kSave=BIT(15),kMark=BIT(16),kFake=BIT(17)};
44785f3e	15	enum {kF02,kF04,kF12,kF13,kF14,kF24, kF44,kNFElem}; // non-trivial elems of propagation matrix
716ccba7	16	enum {kK00,kK01,kK10,kK11,kK20,kK21,kK30,kK31,kK40,kK41, kNKElem}; // non-trivial elems of gain matrix
	17	enum {kS00,kS10,kS11,kS20,kS21,kS22,kS30,kS31,kS32,kS33,kS40,kS41,kS42,kS43,kS44,kNSElem}; // elements of 5x5 sym matrix
	18	enum {kR00,kR22,kNRElem}; // non trivial elements of rotation matrix
3dd9c283	19	//
32d38de2	20	AliITSUSeed();
	21	AliITSUSeed(const AliITSUSeed& src);
	22	AliITSUSeed &operator=(const AliITSUSeed &src);
	23	virtual ~AliITSUSeed();
f8832015	24	virtual void Print(Option_t* option = "") const;
32d38de2	25	//
f8832015	26	void SetLrClusterID(Int_t lr, Int_t cl);
	27	void SetLr(Int_t lr) {SetLrClusterID(lr,-1);} // lr w/o cluster
	28	void SetLrClusterID(UInt_t id) {fClID = id;}
	29	void SetParent(TObject* par) {fParent = par;}
716ccba7	30	void SetChi2Cl(Double_t v) {fChi2Cl= v; v>0 ? fChi2Glo+=v : fChi2Penalty -= v;}
3dd9c283	31	void Kill(Bool_t v=kTRUE) {SetBit(kKilled, v);}
716ccba7	32	void SetFake(Bool_t v=kTRUE) {SetBit(kFake, v);}
6cd80116	33	void Save(Bool_t v=kTRUE) {SetBit(kSave,v);}
08419930	34	void SetMark(Bool_t v=kTRUE) {SetBit(kMark,v);}
6cd80116	35	void FlagTree(UInt_t bits, Bool_t v=kTRUE);
9cdcba2c	36	void SetChi2ITSTPC(Float_t v) {fChi2Match = v;}
9cdcba2c	37	void SetChi2ITSSA(Float_t v) {fChi2ITSSA = v;}
c61e50c3	38	//
f8832015	39	UInt_t GetLrClusterID() const {return fClID;}
	40	Int_t GetLrCluster(Int_t &lr) const {return UnpackCluster(fClID,lr);}
	41	Int_t GetLayerID() const {return UnpackLayer(fClID);}
	42	Int_t GetClusterID() const {return UnpackCluster(fClID);}
	43	Bool_t HasClusterOnLayer(Int_t lr) const {return fHitsPattern&(0x1<<lr);}
e7d83d38	44	Bool_t HasCluster() const {return IsCluster(fClID);}
f8832015	45	Int_t GetNLayersHit() const {return NumberOfBitsSet(fHitsPattern);}
e7d83d38	46	Int_t GetNClusters() const;
e7d83d38	47	Int_t GetClusterIndex(Int_t ind) const;
c8d1f258	48	UShort_t GetHitsPattern() const {return fHitsPattern;}
f8832015	49	Float_t GetChi2Cl() const {return fChi2Cl;}
f8832015	50	Float_t GetChi2Glo() const {return fChi2Glo;}
716ccba7	51	Float_t GetChi2Penalty() const {return fChi2Penalty;}
3dd9c283	52	Float_t GetChi2GloNrm() const;
9cdcba2c	53	Float_t GetChi2ITSTPC() const {return fChi2Match;}
9cdcba2c	54	Float_t GetChi2ITSSA() const {return fChi2ITSSA;}
08419930	55	Float_t GetQualityVar() const;
3dd9c283	56	Bool_t IsKilled() const {return TestBit(kKilled);}
716ccba7	57	Bool_t IsFake() const {return TestBit(kFake);}
6cd80116	58	Bool_t IsSaved() const {return TestBit(kSave);}
08419930	59	Bool_t IsMarked() const {return TestBit(kMark);}
38997c3c	60	Bool_t ContainsFake() const;
9cdcba2c	61	Int_t FetchClusterInfo(Int_t *clIDarr) const;
f8832015	62	//
b8b59e05	63	Int_t GetPoolID() const {return int(GetUniqueID())-1;}
	64	void SetPoolID(Int_t id) {SetUniqueID(id+1);}
	65	//
f8832015	66	TObject* GetParent() const {return fParent;}
716ccba7	67	const AliITSUSeed* GetParent(Int_t lr) const;
c61e50c3	68	//
3dd9c283	69	virtual Bool_t IsSortable() const {return kTRUE;}
	70	virtual Bool_t IsEqual(const TObject* obj) const;
	71	virtual Int_t Compare(const TObject* obj) const;
	72	//
44785f3e	73	// test
9cdcba2c	74	void InitFromSeed(const AliExternalTrackParam* seed);
44785f3e	75	void ResetFMatrix();
	76	void ApplyELoss2FMatrix(Double_t frac, Bool_t beforeProp);
	77	Bool_t ApplyMaterialCorrection(Double_t xOverX0, Double_t xTimesRho, Double_t mass, Bool_t beforeProp);
	78	Bool_t PropagateToX(Double_t xk, Double_t b);
716ccba7	79	Bool_t RotateToAlpha(Double_t alpha);
943e1898	80	Bool_t GetTrackingXAtXAlpha(double xOther,double alpOther,double bz, double &x);
	81	Double_t GetPredictedChi2(Double_t p[2],Double_t cov[3]);
	82	Bool_t Update();
716ccba7	83	Bool_t Smooth(Double_t vecL[5],Double_t matL[15]);
716ccba7	84	Double_t* ProdABA(const double a[15],const double b[15]) const;
44785f3e	85	//
b8b59e05	86	UInt_t GetNChildren() const {return fNChildren;}
	87	Int_t IncChildren() {return ++fNChildren;}
	88	Int_t DecChildren() {return --fNChildren;}
08419930	89	Double_t GetResidY() const {return fResid[0];}
	90	Double_t GetResidZ() const {return fResid[1];}
	91	//
	92	#ifdef _ITSU_TUNING_MODE_ // this is used only for tuning histo filling
	93	void SetOrdBranch(UShort_t v) {fOrdBranch=v;}
	94	void SetOrdCand(UShort_t v) {fOrdCand=v;}
	95	Int_t GetOrdBranch() const {return fOrdBranch;}
	96	Int_t GetOrdCand() const {return fOrdCand;}
	97	Double_t GetPullY() const {return fResid[0]*Sqrt(fCovIYZ[0]);}
	98	Double_t GetPullZ() const {return fResid[1]*Sqrt(fCovIYZ[2]);}
	99	#endif
b8b59e05	100	//
32d38de2	101	protected:
44785f3e	102	//
f8832015	103	UShort_t fHitsPattern; // bit pattern of hits
b8b59e05	104	UShort_t fNChildren; // number of children (prolongations)
c61e50c3	105	UInt_t fClID; // packed cluster info (see AliITSUAux::PackCluster)
716ccba7	106	Float_t fChi2Glo; // current chi2 global (sum of track-cluster chi2's on layers with hit)
	107	Float_t fChi2Cl; // track-cluster chi2 (if >0) or penalty for missing cluster (if < 0)
	108	Float_t fChi2Penalty; // total penalty (e.g. for missing clusters)
9cdcba2c	109	Float_t fChi2Match; // ITS/TPC matching chi2 (per NDF) // RS: to move to separate object of final seed
9cdcba2c	110	Float_t fChi2ITSSA; // ITSSA backward fit chi2 (per NDF) // RS: to move to separate object of final seed
e7d83d38	111	Double_t fFMatrix[kNFElem]; // matrix of propagation from prev layer (non-trivial elements)
	112	Double_t fKMatrix[kNKElem]; // Gain matrix non-trivial elements (note: standard MBF formula uses I-K*H)
	113	Double_t fRMatrix[kNRElem]; // rotation matrix non-trivial elements
	114	Double_t fCovIYZ[3]; // inverted matrix of propagation + meas errors = [Hi * Pi\|i-1 * Hi^T + Ri]^-1
	115	Double_t fResid[2]; // residuals vector
c61e50c3	116	TObject* fParent; // parent track (in higher tree hierarchy)
08419930	117	//
	118	#ifdef _ITSU_TUNING_MODE_ // this is used only for tuning histo filling
	119	UShort_t fOrdBranch; //! order in the seed branching
	120	UShort_t fOrdCand; //! order in candidates tree
	121	#endif
	122	//
32d38de2	123	ClassDef(AliITSUSeed,1)
	124	};
	125
f8832015	126	//_________________________________________________________________________
	127	inline void AliITSUSeed::SetLrClusterID(Int_t lr, Int_t cl)
	128	{
	129	// assign layer, cluster (if -1 - no hit on this layer)
	130	fClID = PackCluster(lr,cl);
134fa001	131	if (cl>=0) fHitsPattern \|= 0x1<<lr;
f8832015	132	}
f8832015	133
44785f3e	134	//_________________________________________________________________________
	135	inline void AliITSUSeed::ResetFMatrix()
	136	{
	137	// reset transport matrix
	138	fFMatrix[kF02] = fFMatrix[kF04] = fFMatrix[kF12] = fFMatrix[kF13] = fFMatrix[kF14] = fFMatrix[kF24] = 0;
	139	fFMatrix[kF44] = 1.0; // this element accumulates eloss
	140	}
	141
	142	//_________________________________________________________________________
	143	inline Bool_t AliITSUSeed::ApplyMaterialCorrection(Double_t xOverX0, Double_t xTimesRho, Double_t mass, Bool_t beforeProp)
	144	{
	145	// apply material correction and modify transport matrix
	146	double pold = Get1P();
	147	if (!CorrectForMeanMaterial(xOverX0,xTimesRho,mass)) return kFALSE;
	148	ApplyELoss2FMatrix( Get1P()/pold, beforeProp);
	149	return kTRUE;
	150	}
	151
	152
	153	//_________________________________________________________________________
	154	inline void AliITSUSeed::ApplyELoss2FMatrix(Double_t frac, Bool_t beforeProp)
	155	{
	156	// Accounts for the energy loss in the transport matrix
	157	// equivalent to multiplying Fmatix by E=diag{1,1,1,1,P4new/P4old}, where P4 is the 1/pt param.
	158	// If beforeProp is true, then it is assumed that the eloss was applied before the transport,
	159	// i.e. F' = F * E, otherwise, after transport, F' = E * F
	160	fFMatrix[kF44] *= frac;
	161	if (beforeProp) {
	162	fFMatrix[kF04] *= frac;
	163	fFMatrix[kF14] *= frac;
	164	fFMatrix[kF24] *= frac;
	165	}
	166	}
	167
716ccba7	168	//_________________________________________________________________________
	169	inline const AliITSUSeed* AliITSUSeed::GetParent(Int_t lr) const
	170	{
	171	// get parent at given layer
	172	const AliITSUSeed* par=this;
	173	int lrt;
	174	while( par && (lrt=par->GetLayerID())>=0 ) {
	175	if (lrt==lr) break;
	176	par = dynamic_cast<const AliITSUSeed*>(par->GetParent());
	177	}
	178	return par;
	179	}
32d38de2	180
e7d83d38	181	//__________________________________________________________________
	182	inline Int_t AliITSUSeed::GetNClusters() const
	183	{
	184	// count number of clusters (some layers may have >1 hit)
	185	int ncl = 0;
	186	const AliITSUSeed* seed = this;
	187	while(seed) {
	188	if (seed->HasCluster()) ncl++;
	189	seed = (AliITSUSeed*)seed->GetParent();
	190	}
	191	return ncl;
	192	//
	193	}
	194
6cd80116	195	//__________________________________________________________________
	196	inline void AliITSUSeed::FlagTree(UInt_t bits, Bool_t v)
	197	{
	198	// set bits on all tree levels
	199	AliITSUSeed* seed = this;
	200	do {seed->SetBit(bits,v);} while ( (seed=(AliITSUSeed*)seed->GetParent()) );
	201	//
	202	}
	203
08419930	204	//__________________________________________________________________
	205	inline Float_t AliITSUSeed::GetQualityVar() const
	206	{
	207	// return variable for quality check (the smaller - the better)
	208	// return GetChi2ITSSA()/Pt(); // this may change! RS
	209	return GetChi2GloNrm()/Pt();
	210	}
	211
6cd80116	212
32d38de2	213	#endif