]>
Commit | Line | Data |
---|---|---|
043badeb | 1 | #ifndef ALIRELALIGNERKALMAN_H |
2 | #define ALIRELALIGNERKALMAN_H | |
3 | ||
4 | /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
5 | * See cxx source for full Copyright notice */ | |
6 | ||
7 | /////////////////////////////////////////////////////////////////////////////// | |
8 | // | |
9 | // Relative alignment of two tracking volumes (default ITS and TPC) | |
10 | // (see AliRelAlignerKalman.cxx for details) | |
11 | // | |
12 | // Origin: Mikolaj Krzewicki, Nikhef, Mikolaj.Krzewicki@cern.ch | |
13 | // | |
14 | ////////////////////////////////////////////////////////////////////////////// | |
15 | ||
782e5230 | 16 | #include <TMath.h> |
55d5da9e | 17 | #include <TVectorD.h> |
782e5230 | 18 | #include <TMatrix.h> |
55d5da9e | 19 | |
9b7d9478 | 20 | class TObject; |
55d5da9e | 21 | class AliExternalTrackParam; |
22 | class AliESDEvent; | |
23 | class AliESDtrack; | |
52d1b7b0 | 24 | class TArrayI; |
fffc83e6 | 25 | class TObjArray; |
043badeb | 26 | |
044eb03e | 27 | class AliRelAlignerKalman : public TObject { |
043badeb | 28 | |
29 | public: | |
30 | AliRelAlignerKalman(); | |
044eb03e | 31 | virtual ~AliRelAlignerKalman(); |
32 | AliRelAlignerKalman& operator= (const AliRelAlignerKalman& a ); | |
33 | AliRelAlignerKalman(const AliRelAlignerKalman& a); | |
043badeb | 34 | |
35 | //User methods: | |
55d5da9e | 36 | Bool_t AddCosmicEvent( const AliESDEvent* pEvent ); |
18f65700 | 37 | Bool_t AddTrackParams( const AliExternalTrackParam* p1, const AliExternalTrackParam* p2 ); |
043badeb | 38 | |
044eb03e | 39 | void Print(Option_t* option="") const; |
043badeb | 40 | |
55d5da9e | 41 | Double_t GetPsi() const {return (*fPX)(0);} |
42 | Double_t GetTheta() const {return (*fPX)(1);} | |
43 | Double_t GetPhi() const {return (*fPX)(2);} | |
44 | Double_t GetX() const {return (*fPX)(3);} | |
45 | Double_t GetY() const {return (*fPX)(4);} | |
46 | Double_t GetZ() const {return (*fPX)(5);} | |
47 | Double_t GetTPCvdCorr() const {return (*fPX)(6);} | |
3ebfbf52 | 48 | Double_t GetTPCt0() const {return (*fPX)(7);} |
49 | Double_t GetTPCvdY() const {if (fgkNSystemParams>8) return (*fPX)(8); else return 0.0;} | |
55d5da9e | 50 | Double_t GetPsiErr() const {return TMath::Sqrt((*fPXcov)(0,0));} |
51 | Double_t GetThetaErr() const {return TMath::Sqrt((*fPXcov)(1,1));} | |
52 | Double_t GetPhiErr() const {return TMath::Sqrt((*fPXcov)(2,2));} | |
53 | Double_t GetXErr() const {return TMath::Sqrt((*fPXcov)(3,3));} | |
54 | Double_t GetYErr() const {return TMath::Sqrt((*fPXcov)(4,4));} | |
55 | Double_t GetZErr() const {return TMath::Sqrt((*fPXcov)(5,5));} | |
56 | Double_t GetTPCvdCorrErr() const {return TMath::Sqrt((*fPXcov)(6,6));} | |
3ebfbf52 | 57 | Double_t GetTPCt0Err() const {return TMath::Sqrt((*fPXcov)(7,7));} |
58 | Double_t GetTPCvdYErr() const {if (fgkNSystemParams>8) return TMath::Sqrt((*fPXcov)(8,8)); else return 0.0;} | |
55d5da9e | 59 | void GetMeasurement( TVectorD& mes ) const { mes = *fPMeasurement; } |
3ebfbf52 | 60 | TVectorD* GetMeasurement() { return fPMeasurement; } |
55d5da9e | 61 | void GetMeasurementCov( TMatrixDSym& cov ) const { cov = *fPMeasurementCov; } |
043badeb | 62 | void SetMeasurement( const TVectorD& mes ) {*fPMeasurement = mes;} |
044eb03e | 63 | void SetMeasurementCov( const TMatrixDSym& cov ) {*fPMeasurementCov = cov;} |
3ebfbf52 | 64 | TMatrixDSym* GetMeasurementCov() const { return fPMeasurementCov; } |
65 | void GetState( TVectorD& state ) const { state = *fPX; } | |
66 | TVectorD* GetState() const { return fPX; } | |
67 | void GetStateCov ( TMatrixDSym& cov ) const { cov = *fPXcov; } | |
043badeb | 68 | void SetState( const TVectorD& param ) {*fPX = param;} |
69 | void SetStateCov (const TMatrixDSym& cov ) {*fPXcov = cov;} | |
3ebfbf52 | 70 | TMatrixDSym* GetStateCov() const { return fPXcov; } |
71 | void GetSeed( TVectorD& seed, TMatrixDSym& seedCov ) const { seed = *fPX; seedCov = *fPXcov; } | |
043badeb | 72 | void SetSeed( const TVectorD& seed, const TMatrixDSym& seedCov ) {*fPX = seed; *fPXcov = seedCov; } |
a80268df | 73 | Bool_t Merge( const AliRelAlignerKalman* al ); |
52d1b7b0 | 74 | Long64_t Merge( TCollection* list ); |
044eb03e | 75 | |
76 | //Expert methods: | |
3ebfbf52 | 77 | Bool_t AddESDevent( const AliESDEvent* pEvent ); |
78 | Bool_t AddESDtrack( const AliESDtrack* pTrack ); | |
79 | void SetMagField( const Double_t f ) { fMagField=f; } | |
80 | Double_t GetMagField() const { return fMagField; } | |
c3b5bfc1 | 81 | Bool_t FindCosmicTrackletNumbersInEvent( TArrayI& outITStracksTArr, TArrayI& outTPCtracksTArr, const AliESDEvent* pEvent ); |
fffc83e6 | 82 | Int_t FindMatchingTracks(TObjArray& arrITS, TObjArray& arrTPC, AliESDEvent* pESD); |
043badeb | 83 | Bool_t Update(); |
044eb03e | 84 | void PrintCorrelationMatrix(); |
3ebfbf52 | 85 | //void PrintCovarianceCorrection(); |
044eb03e | 86 | void PrintSystemMatrix(); |
55d5da9e | 87 | void Reset(); |
c3b5bfc1 | 88 | void ResetCovariance( const Double_t number=0. ); |
55d5da9e | 89 | void ResetTPCparamsCovariance( const Double_t number=0. ); |
90 | Double_t* GetStateArr() const { return fPX->GetMatrixArray(); } | |
91 | Double_t* GetStateCovArr() const { return fPXcov->GetMatrixArray(); } | |
92 | Double_t* GetMeasurementArr() const { return fPMeasurement->GetMatrixArray(); } | |
93 | Double_t* GetMeasurementCovArr() const { return fPMeasurementCov->GetMatrixArray(); } | |
3ebfbf52 | 94 | TMatrixD* GetH() const { return fPH; } |
fffc83e6 | 95 | TVectorD* GetMeasurementPrediction() const {return fPMeasurementPrediction;} |
55d5da9e | 96 | const Double_t* GetDeltaArr() const {return fDelta;} |
3ebfbf52 | 97 | void SetNumericalParanoia( const Bool_t mode=kFALSE ) { fNumericalParanoia=mode; } |
98 | void SetCorrectionMode( const Bool_t mode=kTRUE ) { fCorrectionMode=mode; } | |
044eb03e | 99 | void SetOutRejSigma( const Double_t a=2. ) { fOutRejSigmas = a; } |
100 | void SetRejectOutliers( const Bool_t r=kTRUE ) {fRejectOutliers = r;} | |
fffc83e6 | 101 | void SetRejectOutliersSigma2Median( const Bool_t b=kTRUE ); |
102 | void SetOutRejSigma2Median( const Double_t s ) {fOutRejSigma2Median = s;} | |
3ebfbf52 | 103 | Bool_t SetTrackParams( const AliExternalTrackParam* exparam1, const AliExternalTrackParam* exparam2 ); |
104 | const AliExternalTrackParam* GetTrackParams1() const {return fPTrackParamArr1;} | |
105 | const AliExternalTrackParam* GetTrackParams2() const {return fPTrackParamArr2;} | |
c3b5bfc1 | 106 | void SetMinPointsVol1( const Int_t min ) {fMinPointsVol1=min;} |
107 | void SetMinPointsVol2( const Int_t min ) {fMinPointsVol2=min;} | |
108 | void SetRequireMatchInTPC( const Bool_t s=kTRUE ) {fRequireMatchInTPC = s;} | |
044eb03e | 109 | void SetQ( const Double_t Q = 1e-10 ) { fQ = Q; } |
c3b5bfc1 | 110 | void SetMaxMatchingDistance( const Double_t m ) {fMaxMatchingDistance=m;} |
111 | void SetMaxMatchingAngle( const Double_t m ) {fMaxMatchingAngle=m;} | |
55d5da9e | 112 | void SetTPCvd( const Float_t v ) {fTPCvd=v;} |
113 | void SetTPCZLengthA( const Double_t l ) {fTPCZLengthA=l;} | |
114 | void SetTPCZLengthC( const Double_t l ) {fTPCZLengthC=l;} | |
3ebfbf52 | 115 | Bool_t CorrectTrack( AliExternalTrackParam* tr, const TVectorD& misalignment ) const; |
116 | Bool_t MisalignTrack( AliExternalTrackParam* tr, const TVectorD& misalignment ) const; | |
044eb03e | 117 | static void Angles( TVectorD &angles, const TMatrixD &rotmat ); |
118 | static void RotMat( TMatrixD& R, const TVectorD& angles ); | |
55d5da9e | 119 | static void TMatrixDSymFromTMatrixD( TMatrixDSym& matsym, const TMatrixD& mat ); |
3ebfbf52 | 120 | Bool_t IsPositiveDefinite( const TMatrixD& mat ) const; |
a80268df | 121 | void SetTimeStamp( const UInt_t ts ) { fTimeStamp = ts; } |
122 | UInt_t GetTimeStamp() const {return fTimeStamp;} | |
9b7d9478 | 123 | void SetRunNumber( const Int_t rn ) { fRunNumber = rn; } |
124 | Int_t GetRunNumber() const {return fRunNumber;} | |
329b5744 | 125 | Int_t Compare(const TObject *obj) const; |
126 | Bool_t IsSortable() const { return kTRUE; } | |
16cbbc2c | 127 | Int_t GetNTracks() const {return fNTracks;} |
128 | Int_t GetNUpdates() const {return fNUpdates;} | |
129 | Int_t GetNOutliers() const {return fNOutliers;} | |
fffc83e6 | 130 | Int_t GetNOutliersSigma2Median() const {return fNOutliersSigma2Median;} |
52d1b7b0 | 131 | Int_t GetNMerges() const {return fNMerges;} |
132 | Int_t GetNMergesFailed() const {return fNMergesFailed;} | |
16cbbc2c | 133 | void SetPoint2Track( Bool_t o ); |
043badeb | 134 | |
135 | protected: | |
043badeb | 136 | Bool_t UpdateEstimateKalman(); |
044eb03e | 137 | Bool_t PrepareMeasurement(); |
138 | Bool_t PrepareSystemMatrix(); | |
a80268df | 139 | Bool_t PreparePrediction(); |
043badeb | 140 | Bool_t PredictMeasurement( TVectorD& z, const TVectorD& x ); |
044eb03e | 141 | Bool_t IsOutlier( const TVectorD& update, const TMatrixDSym& covmatrix ); |
fffc83e6 | 142 | Bool_t IsOutlierSigma2Median( const AliExternalTrackParam* pITS, const AliExternalTrackParam* pTPC ); |
043badeb | 143 | |
144 | private: | |
3ebfbf52 | 145 | static const Int_t fgkNTracksPerMeasurement=1; //how many tracks for one update |
3ebfbf52 | 146 | static const Int_t fgkNSystemParams=9; //how many fit parameters |
fffc83e6 | 147 | static const Int_t fgkNtracksSigma2Median=500; //how many sets for median and rms |
043badeb | 148 | |
044eb03e | 149 | //Track parameters |
3ebfbf52 | 150 | AliExternalTrackParam* fPTrackParamArr1; //!local track parameters |
151 | AliExternalTrackParam* fPTrackParamArr2; //!local track parameters | |
52d1b7b0 | 152 | Double_t fMagField; //magnetic field |
044eb03e | 153 | |
043badeb | 154 | //Kalman filter related stuff |
16cbbc2c | 155 | Int_t fNMeasurementParams; //how many measurables |
043badeb | 156 | TVectorD* fPX; //System (fit) parameters (phi, theta, psi, x, y, z, driftcorr, driftoffset ) |
157 | TMatrixDSym* fPXcov; //covariance matrix of system parameters | |
3ebfbf52 | 158 | TMatrixD* fPH; //!System measurement matrix |
159 | Double_t fQ; //!measure for system noise | |
160 | TVectorD* fPMeasurement; //!the measurement vec for Kalman filter (theta,phi,x,z) | |
161 | TMatrixDSym* fPMeasurementCov; //!measurement vec cvariance | |
a80268df | 162 | TVectorD* fPMeasurementPrediction; //!prediction of the measurement |
043badeb | 163 | Double_t fOutRejSigmas; //number of sigmas for outlier rejection |
fffc83e6 | 164 | Double_t fOutRejSigma2Median; //nsigmas to median of input residual distribution |
329b5744 | 165 | Double_t fDelta[fgkNSystemParams]; //array with differentials for calculating derivatives for every parameter(see PrepareSystemMatrix()) |
fffc83e6 | 166 | Double_t* fResArrSigma2Median[4]; //!holds residuals for median based outlier removal |
044eb03e | 167 | |
168 | //Control | |
16cbbc2c | 169 | Bool_t fYZOnly; //whether to consider only yz without directions. |
52d1b7b0 | 170 | Bool_t fNumericalParanoia; //whether to perform additional checks for numerical stability |
171 | Bool_t fRejectOutliers; //whether to do outlier rejection in the Kalman filter | |
fffc83e6 | 172 | Bool_t fRejectOutliersSigma2Median; //whether to reject input based on distance to median |
52d1b7b0 | 173 | Bool_t fRequireMatchInTPC; //when looking for a cosmic in event, require that TPC has 2 matching segments |
174 | Bool_t fCuts; //track cuts? | |
175 | Int_t fMinPointsVol1; //mininum number of points in volume 1 | |
176 | Int_t fMinPointsVol2; //mininum number of points in volume 2 | |
fffc83e6 | 177 | Double_t fMinPt; //min momentum of track for track cuts |
178 | Double_t fMaxPt; //max momentum of track for track cuts | |
52d1b7b0 | 179 | Double_t fMaxMatchingAngle; //cuts |
180 | Double_t fMaxMatchingDistance; //cuts | |
55d5da9e | 181 | Bool_t fCorrectionMode; //calculate corrective transform for TPC (or monitor actual TPC misal params) |
782e5230 | 182 | |
044eb03e | 183 | //Counters |
184 | Int_t fNTracks; //number of processed tracks | |
185 | Int_t fNUpdates; //number of successful Kalman updates | |
186 | Int_t fNOutliers; //number of outliers | |
fffc83e6 | 187 | Int_t fNOutliersSigma2Median; //number of rejected inputs |
55d5da9e | 188 | Int_t fNMatchedCosmics; //number of cosmic events with matching tracklets (good cosmics) |
189 | Int_t fNMatchedTPCtracklets;//number of cosmic events with 2 matching TPC tracklets | |
044eb03e | 190 | Int_t fNProcessedEvents; //number of processed events |
a80268df | 191 | Int_t fTrackInBuffer; //!number of tracks in buffer |
329b5744 | 192 | UInt_t fTimeStamp; //time stamp |
9b7d9478 | 193 | Int_t fRunNumber; //run number |
52d1b7b0 | 194 | Int_t fNMerges; //how many succesful merges |
195 | Int_t fNMergesFailed; //how many merges failed | |
55d5da9e | 196 | |
197 | //TPC stuff | |
198 | Double_t fTPCvd; //TPC drift velocity | |
199 | Double_t fTPCZLengthA; //TPC length side A | |
200 | Double_t fTPCZLengthC; //TPC length side C | |
043badeb | 201 | |
fffc83e6 | 202 | ClassDef(AliRelAlignerKalman,4) //AliRelAlignerKalman class |
043badeb | 203 | }; |
204 | ||
205 | #endif | |
044eb03e | 206 | |
16cbbc2c | 207 |