From 0e9efcbe588a524bd6674b225697db60cac8e892 Mon Sep 17 00:00:00 2001 From: marian Date: Thu, 4 Jun 2009 09:40:29 +0000 Subject: [PATCH] libTPCbase.pkg TPCbaseLinkDef.h - Adding kalman fit AliTPCTransformation.h AliTPCTransformation.cxx - Adding new types of TPC transformation + Possibility to enable/disable transformation AliTPCkalmanFit.h AliTPCkalmanFit.cxx - Add function for merging of results + Filtering of transformation - SetStatus + DumpCorrealtion and DumpCalib - Possibility to specify transforamtion mask (Marian) --- TPC/AliTPCTransformation.cxx | 132 ++++++++++++++++++++++++++++++++++- TPC/AliTPCTransformation.h | 25 ++++++- TPC/AliTPCkalmanFit.cxx | 99 +++++++++++++++++++++++++- TPC/AliTPCkalmanFit.h | 12 +++- TPC/TPCbaseLinkDef.h | 3 + TPC/libTPCbase.pkg | 3 +- 6 files changed, 264 insertions(+), 10 deletions(-) diff --git a/TPC/AliTPCTransformation.cxx b/TPC/AliTPCTransformation.cxx index 1b8bcefc7f2..2a8a9bffabb 100644 --- a/TPC/AliTPCTransformation.cxx +++ b/TPC/AliTPCTransformation.cxx @@ -64,6 +64,19 @@ Int_t AliTPCTransformation::BuildBasicFormulas(){ // RegisterFormula("TPCscalingZDr",(GenFuncG)(AliTPCTransformation::TPCscalingZDr)); RegisterFormula("TPCscalingPhiLocal",(GenFuncG)(AliTPCTransformation::TPCscalingPhiLocal)); + // + // TPC Local X and Y misalignment + // + RegisterFormula("TPClocaldLxdGX",(GenFuncG)(AliTPCTransformation::TPClocaldLxdGX)); + RegisterFormula("TPClocaldLxdGY",(GenFuncG)(AliTPCTransformation::TPClocaldLxdGY)); + RegisterFormula("TPClocaldLydGX",(GenFuncG)(AliTPCTransformation::TPClocaldLydGX)); + RegisterFormula("TPClocaldLydGY",(GenFuncG)(AliTPCTransformation::TPClocaldLydGY)); + // + // Z offset + // + RegisterFormula("TPCDeltaZ",(GenFuncG)(AliTPCTransformation::TPCDeltaZ)); + RegisterFormula("TPCDeltaZMediumLong",(GenFuncG)(AliTPCTransformation::TPCDeltaZMediumLong)); + RegisterFormula("TPCTiltingZ",(GenFuncG)(AliTPCTransformation::TPCTiltingZ)); return 0; } @@ -103,6 +116,7 @@ AliTPCTransformation::AliTPCTransformation(): fSigma(0), // uncertainty of the parameter fSigma2Time(0), // change of the error in time - (For kalman filter) fFixedParam(0), // fixed parameters of tranformation + fIsActive(kTRUE), // swith - On/Off // fInit(kFALSE), // initialization flag - set to kTRUE if corresponding formulas found fFormulaX(0), // x formula - pointer to the function @@ -128,6 +142,7 @@ AliTPCTransformation::AliTPCTransformation(const char *name, TBits *mask, const fSigma(0), fSigma2Time(0), // change of sigma in time fFixedParam(0), // fixed parameters of tranformation + fIsActive(kTRUE), // swith - On/Off // fInit(kFALSE), // initialization flag - set to kTRUE if corresponding formulas found fFormulaX(0), // x formula - pointer to the function @@ -210,6 +225,8 @@ Bool_t AliTPCTransformation::Init(){ return isOK; } + + TBits * AliTPCTransformation::BitsSide(Bool_t aside){ // TBits * bits = new TBits(72); @@ -237,13 +254,14 @@ Double_t AliTPCTransformation::GetDeltaXYZ(Int_t coord, Int_t volID, Double_t pa // // // + if (!fIsActive) return 0; if (fBitMask && (!(*fBitMask)[volID])) return 0; - Double_t xyz[4]={x,y,z, param}; + Double_t xyz[5]={x,y,z, param,volID}; if (fCoordSystem==0){ // cartezian system if (coord==0 && fFormulaX) return fFormulaX(xyz,fFixedParam->GetMatrixArray()); if (coord==1 && fFormulaY) return fFormulaY(xyz,fFixedParam->GetMatrixArray()); - if (coord==2 && fFormulaY) return fFormulaY(xyz,fFixedParam->GetMatrixArray()); + if (coord==2 && fFormulaZ) return fFormulaZ(xyz,fFixedParam->GetMatrixArray()); } if (fCoordSystem==1){ // cylindrical system @@ -339,3 +357,113 @@ Double_t AliTPCTransformation::TPCscalingROFC(Double_t *xyz, Double_t * pa Double_t value = TMath::Power(driftM,param[0])/ndistR; return xyz[3]*value; } + + +// +// TPC sector local misalignment +// +Double_t AliTPCTransformation:: TPClocaldLxdGX(Double_t *xyz, Double_t * /*param*/){ + // + // xyz - [0..2] - position + // [3] - scale parameter + // [4] - volID + // return delta in global coordiante system + // + Int_t sector = TMath::Nint(xyz[4]); + Double_t alpha = TMath::Pi()*(sector+0.5)/9; + Double_t ca = TMath::Cos(alpha); + // Double_t sa = TMath::Sin(alpha); + return ca*xyz[3]; +} + +Double_t AliTPCTransformation::TPClocaldLxdGY(Double_t *xyz, Double_t * /*param*/){ + // + // xyz - [0..2] - position + // [3] - scale parameter + // [4] - volID + // return delta in global coordiante system + // + Int_t sector = TMath::Nint(xyz[4]); + Double_t alpha = TMath::Pi()*(sector+0.5)/9; + //Double_t ca = TMath::Cos(alpha); + Double_t sa = TMath::Sin(alpha); + return sa*xyz[3]; +} + +Double_t AliTPCTransformation:: TPClocaldLydGX(Double_t *xyz, Double_t * /*param*/){ + // + // xyz - [0..2] - position + // [3] - scale parameter + // [4] - volID + // return delta in global coordiante system + // + Int_t sector = TMath::Nint(xyz[4]); + Double_t alpha = TMath::Pi()*(sector+0.5)/9; + //Double_t ca = TMath::Cos(alpha); + Double_t sa = TMath::Sin(alpha); + return -sa*xyz[3]; +} + +Double_t AliTPCTransformation::TPClocaldLydGY(Double_t *xyz, Double_t * /*param*/){ + // + // xyz - [0..2] - position + // [3] - scale parameter + // [4] - volID + // return delta in global coordiante system + // + Int_t sector = TMath::Nint(xyz[4]); + Double_t alpha = TMath::Pi()*(sector+0.5)/9; + Double_t ca = TMath::Cos(alpha); + //Double_t sa = TMath::Sin(alpha); + return ca*xyz[3]; +} + + +Double_t AliTPCTransformation::TPCDeltaZ(Double_t *xyz, Double_t * /*param*/){ + // + // xyz - [0..2] - position + // [3] - scale parameter + // [4] - volID + // return delta in global coordiante system + // + Int_t sector = TMath::Nint(xyz[4]); + Double_t signZ = (sector%36<18) ? 1: -1; // drift direction + return signZ*xyz[3]; // IROC shift +} + +Double_t AliTPCTransformation::TPCDeltaZMediumLong(Double_t *xyz, Double_t * /*param*/){ + // + // xyz - [0..2] - position + // [3] - scale parameter + // [4] - volID + // return delta in global coordinate system + // + Int_t sector = TMath::Nint(xyz[4]); + Double_t signZ = (sector%36<18) ? 1: -1; // drift direction + if (sector<36) return 0; + // + Double_t radius = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])); + const Double_t radiusLong = 198.1; + Double_t sign = (radius0) deltaR *= TMath::Cos(param[0]*alpha); + if (param[1]>0) deltaR *= TMath::Sin(param[1]*alpha); + return deltaR*xyz[3]; +} + diff --git a/TPC/AliTPCTransformation.h b/TPC/AliTPCTransformation.h index ca134031e30..e486418703c 100644 --- a/TPC/AliTPCTransformation.h +++ b/TPC/AliTPCTransformation.h @@ -22,6 +22,8 @@ public: // void SetParams(Double_t param, Double_t sigma, Double_t sigma2Time, TVectorD* fixedParams); Bool_t Init(); + void SetActive(Bool_t flag){ fIsActive = flag;} + Bool_t IsActive(){return fIsActive;} // virtual Double_t GetDeltaXYZ(Int_t coord, Int_t volID, Double_t param, Double_t x, Double_t y, Double_t z); // @@ -42,7 +44,7 @@ public: Double_t fSigma; // error of the parameter Double_t fSigma2Time; // change of the error in time (per hour) - (For kalman filter) TVectorD *fFixedParam; // fixed parameters of tranformation - + Bool_t fIsActive; // switch - is transformation active // // predefined formulas // @@ -50,9 +52,30 @@ public: static Double_t TPCscalingRPol(Double_t *xyz, Double_t * param); static Double_t TPCscalingZDr(Double_t *xyz, Double_t * param); static Double_t TPCscalingPhiLocal(Double_t *xyz, Double_t * param); + // + // TPC Field cage + ROC misalingment induced distortion + // static Double_t TPCscalingRIFC(Double_t *xyz, Double_t * param); // inner field cage r distorion static Double_t TPCscalingROFC(Double_t *xyz, Double_t * param); // outer field cage r distorion // + // TPC local misalignment + // + static Double_t TPClocaldLxdGX(Double_t *xyz, Double_t * param); + static Double_t TPClocaldLxdGY(Double_t *xyz, Double_t * param); + static Double_t TPClocaldLydGX(Double_t *xyz, Double_t * param); + static Double_t TPClocaldLydGY(Double_t *xyz, Double_t * param); + // + // TPC quadrant misalignment + // + // static Double_t TPCQuadrantDr(Double_t *xyz, Double_t * param){return 0;} + //static Double_t TPCQuadrantDrphi(Double_t *xyz, Double_t * param){return 0;} + // + // Z shift - + // + static Double_t TPCDeltaZ(Double_t *xyz, Double_t * param); + static Double_t TPCDeltaZMediumLong(Double_t *xyz, Double_t * param); + static Double_t TPCTiltingZ(Double_t *xyz, Double_t * param); + // Bool_t fInit; // initialization flag GenFuncG fFormulaX; //! x formula GenFuncG fFormulaY; //! y formula diff --git a/TPC/AliTPCkalmanFit.cxx b/TPC/AliTPCkalmanFit.cxx index 541ca9f198f..8cce32f74e1 100644 --- a/TPC/AliTPCkalmanFit.cxx +++ b/TPC/AliTPCkalmanFit.cxx @@ -32,7 +32,7 @@ e.g: kalmanFit0->SetInstance(kalmanFit0); // kalmanFit0->InitTransformation(); // - TF2 fxRZdz("fxRZdz","AliTPCkalmanFit::SGetTPCDeltaXYZ(0,-1,x,0,y)-AliTPCkalmanFit::SGetTPCDeltaXYZ(0,-1,x,0,y-1)",85,245,-250,250); + TF2 fxRZdz("fxRZdz","sign(y)*1000*(AliTPCkalmanFit::SGetTPCDeltaXYZ(0,-1,x,0,y)-AliTPCkalmanFit::SGetTPCDeltaXYZ(0,-1,x,0,y-1))",85,245,-250,250); fxRZdz->Draw(""); @@ -93,6 +93,23 @@ void AliTPCkalmanFit::InitTransformation(){ } } +void AliTPCkalmanFit::Add(const AliTPCkalmanFit * kalman){ + // + // + // + Update(kalman); + for (Int_t i=0;i<8;i++){ + if (fLinearTrackDelta[i] && kalman->fLinearTrackDelta[i]){ + fLinearTrackDelta[i]->Add(kalman->fLinearTrackDelta[i]); + } + if (fLinearTrackPull[i] && kalman->fLinearTrackPull[i]){ + fLinearTrackPull[i]->Add(kalman->fLinearTrackPull[i]); + } + } + +} + + void AliTPCkalmanFit::Init(){ // // Initialize parameter vector and covariance matrix @@ -164,6 +181,69 @@ void AliTPCkalmanFit::Init(){ } +void AliTPCkalmanFit::SetStatus(const char * mask, Bool_t setOn, Bool_t isOr){ + // + // 0. To activate all transforamtion call SetStatus(0,kTRUE) + // 1. To disable everything SetStatus(0,kFALSE) + // 2. To activate/desactivate SetStatus("xxx",kTRUE/kFALSE,kFALSE) + Int_t ncalibs = fCalibration->GetEntries(); + if (mask==0) { + for (Int_t i=0; iAt(i); + transform->SetActive(setOn); + } + } + else{ + for (Int_t i=0; iAt(i); + TString strName(transform->GetName()); + if (strName.Contains(mask)){ + if (!isOr) transform->SetActive( transform->IsActive() && setOn); + if (isOr) transform->SetActive( transform->IsActive() || setOn); + } + } + } +} + + +void AliTPCkalmanFit::Update(const AliTPCkalmanFit * kalman){ + // + // Update Kalman filter + // + Int_t ncalibs = fCalibration->GetEntries(); + TMatrixD vecXk=*fCalibParam; // X vector + TMatrixD covXk=*fCalibCovar; // X covariance + TMatrixD &vecZk = *(kalman->fCalibParam); + TMatrixD &measR = *(kalman->fCalibCovar); + + TMatrixD matHk(ncalibs,ncalibs); // vector to mesurement + TMatrixD vecYk(ncalibs,1); // Innovation or measurement residual + TMatrixD matHkT(ncalibs,ncalibs); // helper matrix Hk transpose + TMatrixD matSk(ncalibs,ncalibs); // Innovation (or residual) covariance + TMatrixD matKk(ncalibs,ncalibs); // Optimal Kalman gain + TMatrixD covXk2(ncalibs,ncalibs); // helper matrix + TMatrixD covXk3(ncalibs,ncalibs); // helper matrix + // + for (Int_t i=0;ithreshold){ AliTPCTransformation * trans0 = GetTransformation(irow); AliTPCTransformation * trans1 = GetTransformation(icol); + TString strName0(trans0->GetName()); + if (mask){ + if (!strName0.Contains(mask)) continue; + } + TString strName1(trans1->GetName()); + if (mask){ + if (!strName1.Contains(mask)) continue; + } printf("%d\t%d\t%s\t%s\t%f\t%f\t%f\n", irow,icol, trans0->GetName(), trans1->GetName(), TMath::Sqrt(mat(irow,irow)), TMath::Sqrt(mat(icol,icol)), corr0); } @@ -711,14 +799,19 @@ Bool_t AliTPCkalmanFit::DumpCorelation(Double_t threshold){ return (nrow>0); } -Bool_t AliTPCkalmanFit::DumpCalib(){ +Bool_t AliTPCkalmanFit::DumpCalib(const char *mask){ // // Print calibration entries - name, value, error // TMatrixD &mat = *fCalibCovar; Int_t nrow= mat.GetNrows(); + TString strMask(mask); for (Int_t irow=0; irowGetName()); + if (mask){ + if (!strName.Contains(mask)) continue; + } printf("%d\t%s\t%f\t%f\n", irow, trans0->GetName(), diff --git a/TPC/AliTPCkalmanFit.h b/TPC/AliTPCkalmanFit.h index dee9463f2a8..986426d1a65 100644 --- a/TPC/AliTPCkalmanFit.h +++ b/TPC/AliTPCkalmanFit.h @@ -20,8 +20,12 @@ public: AliTPCkalmanFit(); void Init(); void InitTransformation(); + void Add(const AliTPCkalmanFit * kalman); + void AddCalibration(AliTPCTransformation * calib); - AliTPCTransformation * GetTransformation(Int_t i){return (fCalibration)? (AliTPCTransformation *)fCalibration->At(i):0;} + AliTPCTransformation * GetTransformation(Int_t i){return (fCalibration)? (AliTPCTransformation *)fCalibration->At(i):0;} + // + void SetStatus(const char * name, Bool_t setOn, Bool_t isOr=kTRUE); // void FitTrackLinear(AliTrackPointArray& points, Int_t step=1, TTreeSRedirector *debug=0); void DumpTrackLinear(AliTrackPointArray& points, TTreeSRedirector *debug); @@ -29,12 +33,14 @@ public: void Propagate(TTreeSRedirector *debug=0); void PropagateTime(Int_t time); + void Update(const AliTPCkalmanFit * kalman); static AliTrackPointArray * MakePointArrayLinear(Double_t alpha, Double_t y0, Double_t z0, Double_t ky, Double_t kz, Double_t err=0.02); void ApplyCalibration(AliTrackPointArray *array, Double_t csign); Bool_t CheckCovariance(TMatrixD &covar, Float_t maxEl); - Bool_t DumpCorelation(Double_t threshold); - Bool_t DumpCalib(); + + Bool_t DumpCorelation(Double_t threshold, const char *mask=0); + Bool_t DumpCalib(const char *mask=0); // Double_t GetTPCDeltaXYZ(Int_t coord, Int_t volID, Double_t x, Double_t y, Double_t z); static Double_t SGetTPCDeltaXYZ(Int_t coord, Int_t volID, Double_t x, Double_t y, Double_t z); diff --git a/TPC/TPCbaseLinkDef.h b/TPC/TPCbaseLinkDef.h index 2d593b23e64..6688b1c3289 100644 --- a/TPC/TPCbaseLinkDef.h +++ b/TPC/TPCbaseLinkDef.h @@ -67,8 +67,11 @@ #pragma link C++ class AliTPCExB+; #pragma link C++ class AliTransform+; #pragma link C++ class AliTPCTransform+; +#pragma link C++ class AliTPCTransformation+; #pragma link C++ class AliTPCCalibTCF+; #pragma link C++ class AliTPCAltroEmulator+; +#pragma link C++ class AliTPCkalmanFit+; + #pragma link C++ class AliTPCAlign; #pragma link C++ class AliTPCdataQA; diff --git a/TPC/libTPCbase.pkg b/TPC/libTPCbase.pkg index 04d13114587..7ac4641d6b9 100644 --- a/TPC/libTPCbase.pkg +++ b/TPC/libTPCbase.pkg @@ -20,7 +20,8 @@ SRCS:= AliSegmentID.cxx AliSegmentArray.cxx AliDigits.cxx AliH2F.cxx \ AliTPCCalibTCF.cxx AliTPCAltroEmulator.cxx AliTPCeventInfo.cxx \ AliTransform.cxx AliTPCTransform.cxx AliTPCAlign.cxx AliTPCGoofieValues.cxx \ AliTPCdataQA.cxx AliTPCQAChecker.cxx AliTPCConfigDA.cxx AliExternalComparison.cxx \ - AliTPCkalmanTime.cxx AliESDcosmic.cxx AliTPCPointCorrection.cxx + AliTPCkalmanTime.cxx AliESDcosmic.cxx AliTPCPointCorrection.cxx AliTPCTransformation.cxx \ + AliTPCkalmanFit.cxx HDRS:= $(SRCS:.cxx=.h) -- 2.43.0