*/
-
+#include "AliTPCcalibDB.h"
#include "TLinearFitter.h"
#include "Riostream.h"
#include "TList.h"
#include "TVectorD.h"
#include "AliLog.h"
#include "AliTPCROC.h"
+#include "AliTPCClusterParam.h"
#include "AliTPCPointCorrection.h"
AliTPCPointCorrection* AliTPCPointCorrection::fgInstance = 0;
fParamOutRVersion(0),
fErrorsOutR(38),
fErrorsOutZ(38),
- fParamOutZVersion(0)
+ fParamOutZVersion(0),
+ //
+ fXIO(0),
+ fXmiddle(0),
+ fXquadrant(0),
+ fArraySectorIntParam(36), // array of sector alignment parameters
+ fArraySectorIntCovar(36), // array of sector alignment covariances
+ //
+ // Kalman filter for global alignment
+ //
+ fSectorParam(0), // Kalman parameter
+ fSectorCovar(0) // Kalman covariance
+
{
//
// Default constructor
//
+ AliTPCROC * roc = AliTPCROC::Instance();
+ fXquadrant = roc->GetPadRowRadii(36,53);
+ fXmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+ fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5;
}
AliTPCPointCorrection::AliTPCPointCorrection(const Text_t *name, const Text_t *title):
fParamOutRVersion(0),
fErrorsOutR(38),
fErrorsOutZ(38),
- fParamOutZVersion(0)
+ fParamOutZVersion(0),
+ //
+ //
+ //
+ fXIO(0),
+ fXmiddle(0),
+ fXquadrant(0),
+ fArraySectorIntParam(36), // array of sector alignment parameters
+ fArraySectorIntCovar(36), // array of sector alignment covariances
+ //
+ // Kalman filter for global alignment
+ //
+ fSectorParam(0), // Kalman parameter for A side
+ fSectorCovar(0) // Kalman covariance for A side
+
{
//
// Non default constructor
//
+ AliTPCROC * roc = AliTPCROC::Instance();
+ fXquadrant = roc->GetPadRowRadii(36,53);
+ fXmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+ fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5;
+
}
AliTPCPointCorrection::~AliTPCPointCorrection(){
Double_t radius = TMath::Sqrt(cx*cx+cy*cy);
AliTPCROC * roc = AliTPCROC::Instance();
Double_t xout = roc->GetPadRowRadiiUp(roc->GetNRows(37)-1);
- Double_t dout = xout-radius;
+ Double_t dout = xout-radius;
+ if (dout<0) return 0;
//drift
Double_t dr = 0.5 - TMath::Abs(cz/250.);
//
result+=(*vec)[5]*eoutp*dr*dr;
result+=(*vec)[6]*eoutm*dr*dr;
return result;
+}
+
+Double_t AliTPCPointCorrection::RPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky,Float_t qMax, Float_t threshold){
+ //
+ // Calculates COG corection in RPHI direction
+ // cluster and track position y is supposed to be corrected before for other effects
+ // (e.g ExB and alignemnt)
+ // Rphi distortion dependeds on the distance to the edge-pad, distance to the wire edge and
+ // relative distance to the center of the pad. Therefore the y position is trnsfromed to the
+ // pad coordiante frame (correction offset (ExB alignemnt) substracted).
+ //
+ // Input parameters:
+ //
+ // sector - sector number - 0-71 - cl.GetDetector()
+ // padrow - padrow number -0-63 -IROC 0-95 OROC - cl->GetRow()
+ // pad - mean pad number - cl->GetPad()
+ // cy - cluster y - cl->GetY()
+ // y - track -or cluster y
+ // qMax - cluster max charge - cl-.GetMax()
+ // threshold - clusterer threshold
+ //
+ AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
+ if (!clparam) {
+ AliFatal("TPC OCDB not initialized");
+ return 0;
+ }
+ Int_t padtype=0;
+ if (sector>=36) padtype = (padrow>64)?2:1;
+ Double_t padwidth=(padtype==0)? 0.4:0.6;
+ Double_t sigma = clparam->GetRMS0(0,padtype,250-TMath::Abs(z),ky)/padwidth;
+ //
+ Int_t npads = AliTPCROC::Instance()->GetNPads(sector,padrow);
+ Float_t yshift = TMath::Abs(cy)-TMath::Abs((-npads*0.5+pad)*padwidth); // the clusters can be corrected before
+ // shift to undo correction
+
+ y -= yshift*((y>0)?1.:-1.); // y in the sector coordinate
+ Double_t y0 = npads*0.5*padwidth;
+ Double_t dy = (TMath::Abs(y0)-TMath::Abs(y))/padwidth-0.5; // distance to the center of pad0
+ //
+ Double_t padcenter = TMath::Nint(dy);
+ Double_t sumw=0;
+ Double_t sumwi=0;
+ for (Double_t ip=padcenter-2.5; ip<=padcenter+2.5;ip++){
+ Double_t weightGaus = qMax*TMath::Exp(-(dy-ip)*(dy-ip)/(2*sigma*sigma));
+ Double_t ypad = (ip-npads*0.5)*padwidth;
+ Double_t weightGain = GetEdgeQ0(sector,padrow,ypad);
+ //
+ Double_t weight = TMath::Nint(weightGaus*weightGain);
+ if (ip<0 &&weight> threshold) weight = threshold; // this is done in cl finder
+ if (weight<0) continue;
+ sumw+=weight;
+ sumwi+=weight*(ip);
+ }
+ Double_t result =0;
+ if (sumw>0) result = sumwi/sumw;
+ result = (result-dy)*padwidth;
+ return result;
}
+
+
+
Double_t AliTPCPointCorrection::CorrectionOutZ0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
//
// return dR correction - for correction version 0
AliTPCROC * roc = AliTPCROC::Instance();
Double_t xout = roc->GetPadRowRadiiUp(roc->GetNRows(37)-1);
Double_t dout = xout-radius;
+ if (dout<0) return 0;
//drift
Double_t dr = 0.5 - TMath::Abs(cz/250.);
//
}
+Double_t AliTPCPointCorrection::GetEdgeQ0(Int_t sector, Int_t padrow, Float_t y){
+ //
+ /* TF1 fexp("fexp","1-exp(-[0]*(x-[1]))",0,20)
+ | param [0] | param [1]
+ IROC | 4.71413 | 1.39558
+ OROC1| 2.11437 | 1.52643
+ OROC2| 2.15082 | 1.53537
+ */
+
+ Double_t params[2]={100,0};
+ if (sector<36){
+ params[0]=4.71413; params[1]=1.39558;
+ }
+ if (sector>=36){
+ if (padrow<64) { params[0]=2.114; params[1]=1.526;}
+ if (padrow>=64){ params[0]=2.15; params[1]=1.535;}
+ }
+ Double_t result= 1;
+ Double_t xrow = AliTPCROC::Instance()->GetPadRowRadii(sector,padrow);
+ Double_t ymax = TMath::Tan(TMath::Pi()/18.)*xrow;
+ Double_t dedge = ymax-TMath::Abs(y);
+ if (dedge-params[1]<0) return 0;
+ if (dedge>10.*params[1]) return 1;
+ result = 1.-TMath::Exp(-params[0]*(dedge-params[1]));
+ return result;
+}
+
+Double_t AliTPCPointCorrection::SRPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky,Float_t qMax, Float_t threshold){
+ return fgInstance->RPhiCOGCorrection(sector, padrow, pad, cy, y,z, ky, qMax, threshold);
+}
+
+Double_t AliTPCPointCorrection::SGetEdgeQ0(Int_t sector, Int_t padrow, Float_t y){
+ //
+ //
+ return fgInstance->GetEdgeQ0(sector, padrow, y);
+}
+
+void AliTPCPointCorrection::AddCorrectionSector(TObjArray & sideAPar, TObjArray &sideCPar, TObjArray & sideACov, TObjArray &sideCCov, Bool_t reset){
+ //
+ //
+ //
+ for (Int_t isec=0;isec<36;isec++){
+ TMatrixD * mat1 = (TMatrixD*)(sideAPar.At(isec));
+ TMatrixD * mat1Covar = (TMatrixD*)(sideACov.At(isec));
+ if (!mat1) continue;
+ TMatrixD * mat0 = (TMatrixD*)(fArraySectorIntParam.At(isec));
+ TMatrixD * mat0Covar = (TMatrixD*)(fArraySectorIntCovar.At(isec));
+ if (!mat0) {
+ fArraySectorIntParam.AddAt(mat1->Clone(),isec);
+ fArraySectorIntCovar.AddAt(mat1Covar->Clone(),isec);
+ continue;
+ }
+ (*mat0Covar)=(*mat1Covar);
+ if (reset) (*mat0)=(*mat1);
+ if (!reset) (*mat0)+=(*mat1);
+ }
+
+ for (Int_t isec=0;isec<36;isec++){
+ TMatrixD * mat1 = (TMatrixD*)(sideCPar.At(isec));
+ TMatrixD * mat1Covar = (TMatrixD*)(sideCCov.At(isec));
+ if (!mat1) continue;
+ TMatrixD * mat0 = (TMatrixD*)(fArraySectorIntParam.At(isec));
+ TMatrixD * mat0Covar = (TMatrixD*)(fArraySectorIntCovar.At(isec));
+ if (!mat0) {
+ fArraySectorIntParam.AddAt(mat1->Clone(),isec);
+ fArraySectorIntCovar.AddAt(mat1Covar->Clone(),isec);
+ continue;
+ }
+ (*mat0Covar)=(*mat1Covar);
+ if (reset) (*mat0)=(*mat1);
+ if (!reset) (*mat0)+=(*mat1);
+ }
+}
+
+
+Double_t AliTPCPointCorrection::GetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t /*lz*/, Int_t quadrant){
+ //
+ // Get position correction for given sector
+ //
+
+ TMatrixD * param = (TMatrixD*)fArraySectorIntParam.At(sector%36);
+ if (!param) return 0;
+ if (quadrant<0){ //recalc quadrant if not specified
+ if (lx<fXIO) quadrant=0;
+ if(lx>fXIO){
+ if (lx<fXquadrant) {
+ if (ly<0) quadrant=1;
+ if (ly>0) quadrant=2;
+ }
+ if (lx>=fXquadrant) {
+ if (ly<0) quadrant=3;
+ if (ly>0) quadrant=4;
+ }
+ }
+ }
+ Double_t a10 = (*param)(quadrant*6+0,0);
+ Double_t a20 = (*param)(quadrant*6+1,0);
+ Double_t a21 = (*param)(quadrant*6+2,0);
+ Double_t dx = (*param)(quadrant*6+3,0);
+ Double_t dy = (*param)(quadrant*6+4,0);
+ Double_t dz = (*param)(quadrant*6+5,0);
+ Double_t deltaX = lx-fXmiddle;
+ if (coord==0) return dx;
+ if (coord==1) return dy+deltaX*a10;
+ if (coord==2) return dz+deltaX*a20+ly*a21;
+ if (coord==3) return a10;
+ if (coord==4) return a20;
+ if (coord==5) return a21;
+ //
+ return 0;
+}
+
+Double_t AliTPCPointCorrection::SGetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant){
+ //
+ //
+ //
+ if (!Instance()) return 0;
+ return Instance()->GetCorrectionSector(coord,sector,lx,ly,lz, quadrant);
+}
+
+
+
+Double_t AliTPCPointCorrection::GetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t /*lz*/){
+ //
+ // Get position correction for given sector
+ //
+ if (!fSectorParam) return 0;
+
+ Double_t a10 = (*fSectorParam)(sector*6+0,0);
+ Double_t a20 = (*fSectorParam)(sector*6+1,0);
+ Double_t a21 = (*fSectorParam)(sector*6+2,0);
+ Double_t dx = (*fSectorParam)(sector*6+3,0);
+ Double_t dy = (*fSectorParam)(sector*6+4,0);
+ Double_t dz = (*fSectorParam)(sector*6+5,0);
+ Double_t deltaX = lx-fXIO;
+ //
+ if (coord==0) return dx;
+ if (coord==1) return dy+deltaX*a10;
+ if (coord==2) return dz+deltaX*a20+ly*a21;
+ if (coord==3) return a10;
+ if (coord==4) return a20;
+ if (coord==5) return a21;
+ return 0;
+}
+
+Double_t AliTPCPointCorrection::SGetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz){
+ //
+ //
+ //
+ if (!Instance()) return 0;
+ return Instance()->GetCorrection(coord,sector,lx,ly,lz);
+}
+
+
+
+
+