From 15e48021eec97437234f06e3d6a39edce18aa6dc Mon Sep 17 00:00:00 2001 From: marian Date: Mon, 2 Mar 2009 19:03:46 +0000 Subject: [PATCH] AliTPCcalibBase.h AliTPCcalibBase.cxx Adding trigger class as string AliTPCPointCorrection.cxx AliTPCPointCorrection.h Adding R distortion correction RPhi distortion correction Quadrant misalignment Inter sector misalignment represented by 6 parameter to be still defined. AliTPCcalibCosmic.cxx AliTPCcalibCosmic.h Adding MakeCombinedTrack - merge 2 halves of the track Removing setting of Gain map - done automatically in the CookdEdx from the OCDB Adding AliExternalComparison - ndim THnsparse for performance study AliTPCcalibTracks.cxx Changing binning of QA histograms AliTPCcalibCalib.cxx Applying the R,RPhi and alignment correction before refit of the track It is temporary sollution before we move the all correction to the AliTPCTransform Corresponding OCDB entry should be back compatible therefore I didn't want to commit it before correction class (AliTPCPointCorrection) is fixed --- TPC/AliTPCPointCorrection.cxx | 257 ++++++++++++++++- TPC/AliTPCPointCorrection.h | 41 ++- TPC/AliTPCcalibBase.cxx | 4 + TPC/AliTPCcalibBase.h | 2 + TPC/AliTPCcalibCalib.cxx | 118 +++++++- TPC/AliTPCcalibCosmic.cxx | 503 +++++----------------------------- TPC/AliTPCcalibCosmic.h | 9 +- TPC/AliTPCcalibTracks.cxx | 17 +- 8 files changed, 482 insertions(+), 469 deletions(-) diff --git a/TPC/AliTPCPointCorrection.cxx b/TPC/AliTPCPointCorrection.cxx index 1b0fb3ea17a..7541e64a518 100644 --- a/TPC/AliTPCPointCorrection.cxx +++ b/TPC/AliTPCPointCorrection.cxx @@ -39,6 +39,7 @@ #include "TVectorD.h" #include "AliLog.h" #include "AliTPCROC.h" +#include "AliTPCClusterParam.h" #include "AliTPCPointCorrection.h" AliTPCPointCorrection* AliTPCPointCorrection::fgInstance = 0; @@ -52,11 +53,27 @@ AliTPCPointCorrection::AliTPCPointCorrection(): 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): @@ -66,11 +83,30 @@ AliTPCPointCorrection::AliTPCPointCorrection(const Text_t *name, const Text_t *t 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(){ @@ -151,7 +187,8 @@ Double_t AliTPCPointCorrection::CorrectionOutR0(Bool_t isGlobal, Bool_t type, D 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.); // @@ -171,9 +208,65 @@ Double_t AliTPCPointCorrection::CorrectionOutR0(Bool_t isGlobal, Bool_t type, D 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 = AliTPCClusterParam::Instance(); + 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 @@ -195,6 +288,7 @@ Double_t AliTPCPointCorrection::CorrectionOutZ0(Bool_t isGlobal, Bool_t type, D 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.); // @@ -217,5 +311,162 @@ Double_t AliTPCPointCorrection::CorrectionOutZ0(Bool_t isGlobal, Bool_t type, D } +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[2]<0) return 0; + if (dedge>10.*params[1]) return 1; + result = 1.-TMath::Exp(-params[0]*(dedge-params[2])); + 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*)(sideCPar.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 (lxfXIO){ + if (lx0) 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); +} + + + + + diff --git a/TPC/AliTPCPointCorrection.h b/TPC/AliTPCPointCorrection.h index cf7c51059a4..8a54efac910 100644 --- a/TPC/AliTPCPointCorrection.h +++ b/TPC/AliTPCPointCorrection.h @@ -9,7 +9,6 @@ #include "TObjArray.h" #include "TVectorD.h" - class AliTPCPointCorrection:public TNamed { public: @@ -32,6 +31,22 @@ public: Double_t CorrectionOutR0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector); Double_t CorrectionOutZ0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector); + Double_t 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 threhsold); + Double_t 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 threhsold); + // + Double_t GetEdgeQ0(Int_t sector, Int_t padrow, Float_t y); + static Double_t SGetEdgeQ0(Int_t sector, Int_t padrow, Float_t y); + // + // IROC -OROC+Quadrant alignment + // + void AddCorrectionSector(TObjArray & sideAPar, TObjArray &sideCPar, TObjArray & sideACov, TObjArray &sideCCov, Bool_t reset); + Double_t GetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant =-1); + static Double_t SGetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant=-1); + // + // Global alignment + // + Double_t GetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz); + static Double_t SGetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz); public: // // Correction out @@ -43,11 +58,33 @@ public: TObjArray fErrorsOutZ; // Parameters for z distortion - outer field cage Int_t fParamOutZVersion; // version of the parameterization // + // Edge rfi + // + // + // Alignment part + // + Double_t fXIO; // OROC-IROC boundary + Double_t fXmiddle; // center of the TPC sector local X + Double_t fXquadrant; // x quadrant + // + // IROC OROC alignment + // + TObjArray fArraySectorIntParam; // array of sector alignment parameters + TObjArray fArraySectorIntCovar; // array of sector alignment covariances + // + // Kalman filter for global alignment + // + TMatrixD *fSectorParam; // Kalman parameter + TMatrixD *fSectorCovar; // Kalman covariance + // + // + // private: + AliTPCPointCorrection(const AliTPCPointCorrection&); AliTPCPointCorrection& operator=(const AliTPCPointCorrection&); static AliTPCPointCorrection* fgInstance; //! Instance of this class (singleton implementation) - ClassDef(AliTPCPointCorrection, 1); + ClassDef(AliTPCPointCorrection, 3); }; #endif diff --git a/TPC/AliTPCcalibBase.cxx b/TPC/AliTPCcalibBase.cxx index b77b3cadd3b..6c1a930a2c2 100644 --- a/TPC/AliTPCcalibBase.cxx +++ b/TPC/AliTPCcalibBase.cxx @@ -63,6 +63,7 @@ AliTPCcalibBase::AliTPCcalibBase(): fTriggerMaskAccept(-1), //trigger mask - accept trigger fHasLaser(kFALSE), //flag the laser is overlayed with given event fRejectLaser(kTRUE), //flag- reject laser + fTriggerClass(), fDebugLevel(0) { // @@ -83,6 +84,7 @@ AliTPCcalibBase::AliTPCcalibBase(const char * name, const char * title): fTriggerMaskAccept(-1), //trigger mask - accept trigger fHasLaser(kFALSE), //flag the laser is overlayed with given event fRejectLaser(kTRUE), //flag- reject laser + fTriggerClass(), fDebugLevel(0) { // @@ -166,7 +168,9 @@ void AliTPCcalibBase::UpdateEventInfo(AliESDEvent * event){ fTime = event->GetTimeStamp(); fTrigger = event->GetTriggerMask(); fMagF = event->GetMagneticField(); + fTriggerClass = event->GetFiredTriggerClasses().Data(); fHasLaser = HasLaser(event); + } diff --git a/TPC/AliTPCcalibBase.h b/TPC/AliTPCcalibBase.h index 78e2d62a2ae..44cb2e1b727 100644 --- a/TPC/AliTPCcalibBase.h +++ b/TPC/AliTPCcalibBase.h @@ -9,6 +9,7 @@ //// #include "TNamed.h" +#include "TObjString.h" class AliTPCseed; class AliESDEvent; class AliESDtrack; @@ -54,6 +55,7 @@ protected: Int_t fTriggerMaskAccept; //trigger mask - accept Bool_t fHasLaser; //flag the laser is overlayed with given event Bool_t fRejectLaser; //flag- reject laser + TObjString fTriggerClass; // trigger class private: Int_t fDebugLevel; // debug level diff --git a/TPC/AliTPCcalibCalib.cxx b/TPC/AliTPCcalibCalib.cxx index c08b82be711..e308d66c609 100644 --- a/TPC/AliTPCcalibCalib.cxx +++ b/TPC/AliTPCcalibCalib.cxx @@ -60,6 +60,7 @@ #include "AliTPCTransform.h" #include "AliTPCclusterMI.h" #include "AliTPCseed.h" +#include "AliTPCPointCorrection.h" ClassImp(AliTPCcalibCalib) @@ -148,11 +149,16 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ // Refit track // + // + // 0 - Setup transform object + // + AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ; + transform->SetCurrentRun(fRun); + transform->SetCurrentTimeStamp((UInt_t)fTime); // // First apply calibration // - - AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform(); + AliTPCPointCorrection * corr = AliTPCPointCorrection::Instance(); for (Int_t irow=0;irow<159;irow++) { AliTPCclusterMI *cluster=seed->GetClusterPointer(irow); if (!cluster) continue; @@ -169,7 +175,37 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ Float_t dz =AliTPCClusterParam::SPosCorrection(1,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax()); //x[1]-=dy; //x[2]-=dz; + // + // Apply sector alignment + // + Double_t dxq = AliTPCPointCorrection::SGetCorrectionSector(0,cluster->GetDetector()%36,cluster->GetX(), + cluster->GetY(),cluster->GetZ()); + Double_t dyq = AliTPCPointCorrection::SGetCorrectionSector(1,cluster->GetDetector()%36,cluster->GetX(), + cluster->GetY(),cluster->GetZ()); + Double_t dzq = AliTPCPointCorrection::SGetCorrectionSector(2,cluster->GetDetector()%36,cluster->GetX(), + cluster->GetY(),cluster->GetZ()); + if (kTRUE){ + x[0]-=dxq; + x[1]-=dyq; + x[2]-=dzq; + } + // + // Apply r-phi correction - To be done on track level- knowing the track angle !!! + // + Double_t corrclY = + corr->RPhiCOGCorrection(cluster->GetDetector(),cluster->GetRow(), cluster->GetPad(), + cluster->GetY(),cluster->GetY(), cluster->GetZ(), 0., cluster->GetMax(),2.5); + // R correction + Double_t corrR = corr->CorrectionOutR0(kFALSE,kFALSE,cluster->GetX(),cluster->GetY(),cluster->GetZ(),cluster->GetDetector()); + + if (kTRUE){ + if (cluster->GetY()>0) x[1]+=corrclY; // rphi correction + if (cluster->GetY()<0) x[1]-=corrclY; // rphi correction + x[0]+=corrR; // radial correction + } + // + // // cluster->SetX(x[0]); cluster->SetY(x[1]); @@ -182,15 +218,23 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ "event="<GetNumberOfClusters(); Double_t covar[15]; for (Int_t i=0;i<15;i++) covar[i]=0; @@ -207,16 +251,43 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ AliExternalTrackParam * trackOutOld = (AliExternalTrackParam*)track->GetOuterParam(); AliExternalTrackParam trackIn = *trackOutOld; - AliExternalTrackParam trackOut = *trackInOld; trackIn.ResetCovariance(200.); - trackOut.ResetCovariance(200.); trackIn.AddCovariance(covar); - trackOut.AddCovariance(covar); Double_t xyz[3]; Int_t nclIn=0,nclOut=0; // + // Refit in + // + + for (Int_t irow=159; irow>0; irow--){ + AliTPCclusterMI *cl=seed->GetClusterPointer(irow); + if (!cl) continue; + if (cl->GetX()<80) continue; + Int_t sector = cl->GetDetector(); + Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha(); + if (TMath::Abs(dalpha)>0.01) + trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.)); + Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()}; + Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation + AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]); + cov[0]*=cov[0]; + cov[2]*=cov[2]; + trackIn.GetXYZ(xyz); + Double_t bz = AliTracker::GetBz(xyz); + + if (!trackIn.PropagateTo(r[0],bz)) continue; + if (RejectCluster(cl,&trackIn)) continue; + nclIn++; + trackIn.Update(&r[1],cov); + } + // + AliExternalTrackParam trackOut = trackIn; + trackOut.ResetCovariance(200.); + trackOut.AddCovariance(covar); + // // Refit out // + //Bool_t lastEdge=kFALSE; for (Int_t irow=0; irow<160; irow++){ AliTPCclusterMI *cl=seed->GetClusterPointer(irow); if (!cl) continue; @@ -234,14 +305,22 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ cov[2]*=cov[2]; trackOut.GetXYZ(xyz); Double_t bz = AliTracker::GetBz(xyz); - if (trackOut.PropagateTo(r[0],bz)) nclOut++; + if (!trackOut.PropagateTo(r[0],bz)) continue; if (RejectCluster(cl,&trackOut)) continue; - trackOut.Update(&r[1],cov); + nclOut++; + trackOut.Update(&r[1],cov); + //if (cl->GetType()<0) lastEdge=kTRUE; + //if (cl->GetType()>=0) lastEdge=kFALSE; } // - // Refit in // - + // + nclIn=0; + trackIn = trackOut; + trackIn.ResetCovariance(10.); + // + // Refit in one more time + // for (Int_t irow=159; irow>0; irow--){ AliTPCclusterMI *cl=seed->GetClusterPointer(irow); if (!cl) continue; @@ -258,10 +337,13 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ trackIn.GetXYZ(xyz); Double_t bz = AliTracker::GetBz(xyz); - if (trackIn.PropagateTo(r[0],bz)) nclIn++; + if (!trackIn.PropagateTo(r[0],bz)) continue; if (RejectCluster(cl,&trackIn)) continue; + nclIn++; trackIn.Update(&r[1],cov); } + + trackIn.Rotate(trackInOld->GetAlpha()); trackOut.Rotate(trackOutOld->GetAlpha()); // @@ -282,6 +364,7 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ "event="<GetNumberOfTracks(); fHistNTracks->Fill(ntracks); if (ntracks==0) return; - + AliESDcosmic cosmicESD; + TTreeSRedirector * cstream = GetDebugStreamer(); + cosmicESD.SetDebugStreamer(cstream); + cosmicESD.ProcessEvent(event); + if (cstream) cosmicESD.DumpToTree(); + + } @@ -205,11 +214,11 @@ void AliTPCcalibCosmic::FindPairs(AliESDEvent *event) { Double_t meanP = 0.5*(trackIn->GetP() + trackOut->GetP()); if (seed && track->GetTPCNcls() > 80 + 60/(1+TMath::Exp(-meanP+5))) { - fDeDx->Fill(meanP, seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap)); + fDeDx->Fill(meanP, seed->CookdEdxNorm(0.0,0.45,0,0,159)); // - if (meanP > 0.4 && meanP < 0.45) fDeDxMIP->Fill(seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap)); + if (meanP > 0.4 && meanP < 0.45) fDeDxMIP->Fill(seed->CookdEdxNorm(0.0,0.45,0,0,159)); // - if (GetDebugLevel()>0&&meanP>0.2&&seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap)>300) { + if (GetDebugLevel()>0&&meanP>0.2&&seed->CookdEdxNorm(0.0,0.45,0,0,159)>300) { TFile *curfile = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile(); if (curfile) printf(">>> p+ in file: %s \t event: %i \t Number of ESD tracks: %i \n", curfile->GetName(), (int)event->GetEventNumberInFile(), (int)ntracks); if (track->GetOuterParam()->GetAlpha()<0) cout << " Polartiy: " << track->GetSign() << endl; @@ -246,14 +255,14 @@ void AliTPCcalibCosmic::FindPairs(AliESDEvent *event) { AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j); if (! seed0) continue; if (! seed1) continue; - Float_t dedx0 = seed0->CookdEdxNorm(0.05,0.55,0,0,159,fGainMap); - Float_t dedx1 = seed1->CookdEdxNorm(0.05,0.55,0,0,159,fGainMap); + Float_t dedx0 = seed0->CookdEdxNorm(0.05,0.55,0,0,159); + Float_t dedx1 = seed1->CookdEdxNorm(0.05,0.55,0,0,159); // - Float_t dedx0I = seed0->CookdEdxNorm(0.05,0.55,0,0,63,fGainMap); - Float_t dedx1I = seed1->CookdEdxNorm(0.05,0.55,0,0,63,fGainMap); + Float_t dedx0I = seed0->CookdEdxNorm(0.05,0.55,0,0,63); + Float_t dedx1I = seed1->CookdEdxNorm(0.05,0.55,0,0,63); // - Float_t dedx0O = seed0->CookdEdxNorm(0.05,0.55,0,64,159,fGainMap); - Float_t dedx1O = seed1->CookdEdxNorm(0.05,0.55,0,64,159,fGainMap); + Float_t dedx0O = seed0->CookdEdxNorm(0.05,0.55,0,64,159); + Float_t dedx1O = seed1->CookdEdxNorm(0.05,0.55,0,64,159); // Float_t dir = (dir0[0]*dir1[0] + dir0[1]*dir1[1] + dir0[2]*dir1[2]); Float_t d0 = track0->GetLinearD(0,0); @@ -298,6 +307,17 @@ void AliTPCcalibCosmic::FindPairs(AliESDEvent *event) { Bool_t isPair = IsPair(¶m0,¶m1); // if (isPair) FillAcordeHist(track0); + if (fComp){ + Bool_t acceptComp = fComp->AcceptPair(¶m0,¶m1); + if (acceptComp) fComp->Process(¶m0,¶m1); + } + // + // combined track params + // + AliExternalTrackParam *par0U=MakeCombinedTrack(¶m0,¶m1); + AliExternalTrackParam *par1U=MakeCombinedTrack(¶m1,¶m0); + + // if (fStreamLevel>0){ TTreeSRedirector * cstream = GetDebugStreamer(); @@ -316,6 +336,7 @@ void AliTPCcalibCosmic::FindPairs(AliESDEvent *event) { "event="<Add(cal->GetHistMIP()); } - + //if (fComp && cal->fComp) fComp->Add(cal->fComp); return 0; } @@ -519,19 +544,6 @@ void AliTPCcalibCosmic::BinLogX(TH1 *h) { } -AliExternalTrackParam *AliTPCcalibCosmic::Invert(AliExternalTrackParam *input) -{ - // - // Invert paramerameter - not covariance yet - // - AliExternalTrackParam *output = new AliExternalTrackParam(*input); - Double_t * param = (Double_t*)output->GetParameter(); - param[0]*=-1; - param[3]*=-1; - param[4]*=-1; - // - return output; -} AliExternalTrackParam *AliTPCcalibCosmic::MakeTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1){ // @@ -539,10 +551,12 @@ AliExternalTrackParam *AliTPCcalibCosmic::MakeTrack(const AliExternalTrackParam // AliExternalTrackParam *par1R= new AliExternalTrackParam(*track1); par1R->Rotate(track0->GetAlpha()); + par1R->PropagateTo(track0->GetX(),AliTracker::GetBz()); // // Double_t * param = (Double_t*)par1R->GetParameter(); Double_t * covar = (Double_t*)par1R->GetCovariance(); + param[0]*=1; //OK param[1]*=1; //OK param[2]*=1; //? @@ -552,15 +566,23 @@ AliExternalTrackParam *AliTPCcalibCosmic::MakeTrack(const AliExternalTrackParam covar[6] *=-1.; covar[7] *=-1.; covar[8] *=-1.; //covar[10]*=-1.; covar[11]*=-1.; covar[12]*=-1.; covar[13]*=-1.; - par1R->PropagateTo(track0->GetX(),0); // bz shold be set - - //if (1){ - // printf("Print param\n"); - // track1->Print(); - // par1R->Print(); - //} return par1R; } +AliExternalTrackParam *AliTPCcalibCosmic::MakeCombinedTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1){ + // + // Make combined track + // + // + AliExternalTrackParam * par1T = MakeTrack(track0,track1); + AliExternalTrackParam * par0U = new AliExternalTrackParam(*track0); + // + UpdateTrack(*par0U,*par1T); + delete par1T; + return par0U; +} + + void AliTPCcalibCosmic::UpdateTrack(AliExternalTrackParam &track1, const AliExternalTrackParam &track2){ // // Update track 1 with track 2 @@ -588,19 +610,21 @@ void AliTPCcalibCosmic::UpdateTrack(AliExternalTrackParam &track1, const AliExte // // copy data to the matrix for (Int_t ipar=0; ipar<5; ipar++){ - vecXk(ipar,0) = param1[ipar]; - vecZk(ipar,0) = param2[ipar]; for (Int_t jpar=0; jpar<5; jpar++){ covXk(ipar,jpar) = covar1[track1.GetIndex(ipar, jpar)]; measR(ipar,jpar) = covar2[track2.GetIndex(ipar, jpar)]; + matHk(ipar,jpar)=0; + mat1(ipar,jpar)=0; } + vecXk(ipar,0) = param1[ipar]; + vecZk(ipar,0) = param2[ipar]; + matHk(ipar,ipar)=1; + mat1(ipar,ipar)=0; } // // // // - matHk(0,0)=1; matHk(1,1)= 1; matHk(2,2)= 1; - matHk(3,3)= 1; matHk(4,4)= 1; // vector to measurement // vecYk = vecZk-matHk*vecXk; // Innovation or measurement residual matHkT=matHk.T(); matHk.T(); @@ -608,7 +632,6 @@ void AliTPCcalibCosmic::UpdateTrack(AliExternalTrackParam &track1, const AliExte matSk.Invert(); matKk = (covXk*matHkT)*matSk; // Optimal Kalman gain vecXk += matKk*vecYk; // updated vector - mat1(0,0)=1; mat1(1,1)=1; mat1(2,2)=1; mat1(3,3)=1; mat1(4,4)=1; covXk2 = (mat1-(matKk*matHk)); covOut = covXk2*covXk; // @@ -703,405 +726,3 @@ void AliTPCcalibCosmic::ProcessTree(TTree * chainTracklet, AliExternalComparison - - - -/* - -void Init(){ - -.x ~/UliStyle.C -.x ~/rootlogon.C -gSystem->Load("libSTAT.so"); -gSystem->Load("libANALYSIS"); -gSystem->Load("libTPCcalib"); -gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros"); - -gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+") -AliXRDPROOFtoolkit tool; -TChain * chainCosmic = tool.MakeChain("cosmic.txt","Track0",0,1000000); -chainCosmic->Lookup(); - -TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.01"); // OK -TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<2"); // OK -TCut cutP1("cutP1","abs(Tr0.fP[1]-Tr1.fP[1])<20"); // OK -TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10"); -TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>50"); -TCut cutM("cutM","abs(mag)>0.01"); -TCut cutA=cutT+cutD+cutPt+cutN+cutP1+"trigger!=16"; - -TCut cuthpt("abs(Tr0.fP[4])+abs(Tr1.fP[4])<0.2"); -TCut cutS("cutS","Orig0.fIp.fP[1]*Orig1.fIp.fP[1]>0"); - -// -chainCosmic->Draw(">>listELP",cutA,"entryList"); -//TEntryList *elist = (TEntryList*)gDirectory->Get("listEL"); -//TEntryList *elist = (TEntryList*)gProof->GetOutputList()->At(1); -chainCosmic->SetEntryList(elist); -// -chainCosmic->Draw(">>listV40Z100","abs(d0)<40&&abs(v01)<100","entryList"); -TEntryList *elistV40Z100 = (TEntryList*)gDirectory->Get("listV40Z100"); -chainCosmic->SetEntryList(elistV40Z100); - -// -// Aliases -// -chainCosmic->SetAlias("side","(-1+(Tr0.fP[1]>0)*2)"); -chainCosmic->SetAlias("hpt","abs(Tr0.fP[4])<0.2"); -chainCosmic->SetAlias("signy","(-1+(Tr0.fP[0]>0)*2)"); - -chainCosmic->SetAlias("dy","Tr0.fP[0]+Tr1.fP[0]"); -chainCosmic->SetAlias("dz","Tr0.fP[1]-Tr1.fP[1]"); -chainCosmic->SetAlias("d1pt","Tr0.fP[4]+Tr1.fP[4]"); -chainCosmic->SetAlias("dtheta","(Tr0.fP[3]+Tr1.fP[3])"); -chainCosmic->SetAlias("dphi","(Tr0.fAlpha-Tr1.fAlpha-pi)"); - -chainCosmic->SetAlias("mtheta","(Tr0.fP[3]-Tr1.fP[3])*0.5") -chainCosmic->SetAlias("sa","(sin(Tr0.fAlpha+0.))"); -chainCosmic->SetAlias("ca","(cos(Tr0.fAlpha+0.))"); - - - -chainCosmic->Draw("dy:sqrt(abs(Tr0.fP[4]))>>hisdyA(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0"); -hisdyA->FitSlicesY(); -hisdyA_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdyA_2->SetYTitle("#sigma_{y}(cm)"); -hisdyA_2->SetTitle("Cosmic - Y matching"); -hisdyA_2->SetMaximum(0.5); - - -chainCosmic->Draw("dy:sqrt(abs(Tr0.fP[4]))>>hisdyC(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0"); -hisdyC->FitSlicesY(); -hisdyC_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdyC_2->SetYTitle("#sigma_{y}(cm)"); -hisdyC_2->SetTitle("Cosmic - Y matching"); -hisdyC_2->SetMaximum(1); -hisdyC_2->SetMinimum(0); -hisdyC_2->SetMarkerStyle(22); -hisdyA_2->SetMarkerStyle(21); -hisdyC_2->SetMarkerSize(1.5); -hisdyA_2->SetMarkerSize(1.5); -hisdyC_2->Draw(); -hisdyA_2->Draw("same"); -gPad->SaveAs("~/Calibration/Cosmic/pic/ymatching.gif") - -chainCosmic->Draw("dz:sqrt(abs(Tr0.fP[4]))>>hisdzA(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0"); -hisdzA->FitSlicesY(); -hisdzA_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdzA_2->SetYTitle("#sigma_{z}(cm)"); -hisdzA_2->SetTitle("Cosmic - Z matching - A side "); -hisdzA_2->SetMaximum(0.5); - -chainCosmic->Draw("dz:sqrt(abs(Tr0.fP[4]))>>hisdzC(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0"); -hisdzC->FitSlicesY(); -hisdzC_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdzC_2->SetYTitle("#sigma_{z}(cm)"); -hisdzC_2->SetTitle("Cosmic - Z matching"); -hisdzC_2->SetMaximum(0.5); -hisdzC_2->SetMarkerStyle(22); -hisdzA_2->SetMarkerStyle(21); -hisdzC_2->SetMarkerSize(1.5); -hisdzA_2->SetMarkerSize(1.5); - -hisdzC_2->Draw(); -hisdzA_2->Draw("same"); - - -// -// PICTURE 1/pt -// -chainCosmic->Draw("d1pt:sqrt(abs(Tr0.fP[4]))>>hisd1ptA(5,0,1,30,-0.1,0.1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0"+cutM); -hisd1ptA->FitSlicesY(); -hisd1ptA_2->SetXTitle("#sqrt{1/p_{t}}"); -hisd1ptA_2->SetYTitle("#sigma_{z}(cm)"); -hisd1ptA_2->SetTitle("Cosmic - Z matching - A side "); -hisd1ptA_2->SetMaximum(0.5); - -chainCosmic->Draw("d1pt:sqrt(abs(Tr0.fP[4]))>>hisd1ptC(5,0,1,30,-0.1,0.1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0"+cutM); -hisd1ptC->FitSlicesY(); -hisd1ptC_2->SetXTitle("#sqrt{1/p_{t}}"); -hisd1ptC_2->SetYTitle("#sigma_{1/pt}(1/GeV)"); -hisd1ptC_2->SetTitle("Cosmic - 1/pt matching"); -hisd1ptC_2->SetMaximum(0.05); -hisd1ptC_2->SetMarkerStyle(22); -hisd1ptA_2->SetMarkerStyle(21); -hisd1ptC_2->SetMarkerSize(1.5); -hisd1ptA_2->SetMarkerSize(1.5); - -hisd1ptC_2->Draw(); -hisd1ptA_2->Draw("same"); -gPad->SaveAs("~/Calibration/Cosmic/pic/1ptmatching.gif") - -// -// Theta -// -chainCosmic->Draw("dtheta:sqrt(abs(Tr0.fP[4]))>>hisdthetaA(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0"); -hisdthetaA->FitSlicesY(); -hisdthetaA_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdthetaA_2->SetYTitle("#sigma_{#theta}(cm)"); -hisdthetaA_2->SetTitle("Cosmic - Z matching - A side "); -hisdthetaA_2->SetMaximum(0.5); - -chainCosmic->Draw("dtheta:sqrt(abs(Tr0.fP[4]))>>hisdthetaC(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0"); -hisdthetaC->FitSlicesY(); -hisdthetaC_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdthetaC_2->SetYTitle("#sigma_{#theta}(rad)"); -hisdthetaC_2->SetTitle("Cosmic - Theta matching"); -hisdthetaC_2->SetMaximum(0.01); -hisdthetaC_2->SetMinimum(0.0); -hisdthetaC_2->SetMarkerStyle(22); -hisdthetaA_2->SetMarkerStyle(21); -hisdthetaC_2->SetMarkerSize(1.5); -hisdthetaA_2->SetMarkerSize(1.5); - -hisdthetaC_2->Draw(); -hisdthetaA_2->Draw("same"); -gPad->SaveAs("~/Calibration/Cosmic/pic/thetamatching.gif") -// -// Phi -// -chainCosmic->Draw("dphi:sqrt(abs(Tr0.fP[4]))>>hisdphiA(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0"); -hisdphiA->FitSlicesY(); -hisdphiA_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdphiA_2->SetYTitle("#sigma_{#phi}(rad)"); -hisdphiA_2->SetTitle("Cosmic - Z matching - A side "); -hisdphiA_2->SetMaximum(0.5); - -chainCosmic->Draw("dphi:sqrt(abs(Tr0.fP[4]))>>hisdphiC(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0"); -hisdphiC->FitSlicesY(); -hisdphiC_2->SetXTitle("#sqrt{1/p_{t}}"); -hisdphiC_2->SetYTitle("#sigma_{#phi}(rad)"); -hisdphiC_2->SetTitle("Cosmic - Phi matching"); -hisdphiC_2->SetMaximum(0.01); -hisdphiC_2->SetMinimum(0.0); -hisdphiC_2->SetMarkerStyle(22); -hisdphiA_2->SetMarkerStyle(21); -hisdphiC_2->SetMarkerSize(1.5); -hisdphiA_2->SetMarkerSize(1.5); - -hisdphiC_2->Draw(); -hisdphiA_2->Draw("same"); -gPad->SaveAs("~/Calibration/Cosmic/pic/phimatching.gif") - - - -} - - -*/ - - -/* -void MatchTheta(){ - -TStatToolkit toolkit; -Double_t chi2=0; -Int_t npoints=0; -TVectorD fitParamA0; -TVectorD fitParamA1; -TVectorD fitParamC0; -TVectorD fitParamC1; -TMatrixD covMatrix; - - -TString fstring=""; -// -fstring+="mtheta++"; -fstring+="ca++"; -fstring+="sa++"; -fstring+="ca*mtheta++"; -fstring+="sa*mtheta++"; -// -fstring+="side++"; -fstring+="side*mtheta++"; -fstring+="side*ca++"; -fstring+="side*sa++"; -fstring+="side*ca*mtheta++"; -fstring+="side*sa*mtheta++"; - - -TString *strTheta0 = toolkit.FitPlane(chain,"dtheta",fstring->Data(), "hpt&&!crossI&&!crossO", chi2,npoints,fitParamA0,covMatrix,0.8); -chainCosmic->SetAlias("dtheta0",strTheta0.Data()); -strTheta0->Tokenize("+")->Print(); - - -//fstring+="mtheta++"; -//fstring+="mtheta^2++"; -//fstring+="ca*mtheta^2++"; -//fstring+="sa*mtheta^2++"; - - - -} - -*/ - - - - -/* - void PosCorrection() - - - - - TStatToolkit toolkit; - Double_t chi2=0; - Int_t npoints=0; - TVectorD fitParam; - TMatrixD covMatrix; - - //Theta -chainCosmic->SetAlias("dthe","(Tr0.fP[3]+Tr1.fP[3])"); -chainCosmic->SetAlias("sign","(-1+(Tr0.fP[1]>0)*2)"); -chainCosmic->SetAlias("di","(1-abs(Tr0.fP[1])/250)"); -// -// -TString strFit=""; -// -strFit+="sign++"; //1 -strFit+="Tr0.fP[3]++"; //2 -// -strFit+="sin(Tr0.fAlpha)*(Tr0.fP[1]>0)++"; //3 -strFit+="sin(Tr0.fAlpha)*(Tr0.fP[1]<0)++"; //4 -// -strFit+="cos(Tr0.fAlpha)*(Tr0.fP[1]>0)++"; //5 -strFit+="cos(Tr0.fAlpha)*(Tr0.fP[1]<0)++"; //6 -// -strFit+="sin(Tr0.fAlpha)*Tr0.fP[3]++"; //7 -strFit+="cos(Tr0.fAlpha)*Tr0.fP[3]++"; //8 - - - // - TString * thetaParam = toolkit.FitPlane(chain,"dthe", strFit.Data(),"1", chi2,npoints,fitParam,covMatrix,0.8,0,10000) - - chainCosmic->SetAlias("corTheta",thetaParam->Data()); - chainCosmic->Draw("dthe:Tr0.fP[1]","","",50000); - - - -*/ - - - -/* - -void AliTPCcalibCosmic::dEdxCorrection(){ - TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.01"); // OK - TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<2"); // OK - TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10"); - TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>100"); - TCut cutS("cutS","Orig0.fIp.fP[1]*Orig1.fIp.fP[1]>0"); - TCut cutA=cutT+cutD+cutPt+cutN+cutS; - - - .x ~/UliStyle.C - gSystem->Load("libANALYSIS"); - gSystem->Load("libTPCcalib"); - gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros"); - gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+") - AliXRDPROOFtoolkit tool; - TChain * chainCosmic = tool.MakeChain("cosmic.txt","Track0",0,1000000); - chainCosmic->Lookup(); - - chainCosmic->Draw(">>listEL",cutA,"entryList"); - TEntryList *elist = (TEntryList*)gDirectory->Get("listEL"); - chainCosmic->SetEntryList(elist); - - .x ~/rootlogon.C - gSystem->Load("libSTAT.so"); - TStatToolkit toolkit; - Double_t chi2=0; - Int_t npoints=0; - TVectorD fitParam; - TMatrixD covMatrix; - - chainCosmic->Draw("Tr0.fP[4]+Tr1.fP[4]","OK"+cutA); - - TString strFit; - strFit+="(Tr0.fP[1]/250)++"; - strFit+="(Tr0.fP[1]/250)^2++"; - strFit+="(Tr0.fP[3])++"; - strFit+="(Tr0.fP[3])^2++"; - - TString * ptParam = TStatToolkit::FitPlane(chain,"Tr0.fP[4]+Tr1.fP[4]", strFit.Data(),"1", chi2,npoints,fitParam,covMatrix) - - - -*/ - - -/* -gSystem->Load("libANALYSIS"); -gSystem->Load("libSTAT"); -gSystem->Load("libTPCcalib"); - -TStatToolkit toolkit; -Double_t chi2; -TVectorD fitParam; -TMatrixD covMatrix; -Int_t npoints; -// -TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03"); // OK -TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5"); // OK -TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.2&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10"); -TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>110"); -TCut cutA=cutT+cutD+cutPt+cutN; - - - -TTree * chainCosmic = Track0; - - -chainCosmic->SetAlias("norm","signalTot0.fElements[3]/signalTot1.fElements[3]"); -// -chainCosmic->SetAlias("dr1","(signalTot0.fElements[1]/signalTot0.fElements[3])"); -chainCosmic->SetAlias("dr2","(signalTot0.fElements[2]/signalTot0.fElements[3])"); -chainCosmic->SetAlias("dr4","(signalTot0.fElements[4]/signalTot0.fElements[3])"); -chainCosmic->SetAlias("dr5","(signalTot0.fElements[5]/signalTot0.fElements[3])"); - -TString fstring=""; -fstring+="dr1++"; -fstring+="dr2++"; -fstring+="dr4++"; -fstring+="dr5++"; -// -fstring+="dr1*dr2++"; -fstring+="dr1*dr4++"; -fstring+="dr1*dr5++"; -fstring+="dr2*dr4++"; -fstring+="dr2*dr5++"; -fstring+="dr4*dr5++"; - - - -TString *strqdedx = toolkit.FitPlane(chain,"norm",fstring->Data(), cutA, chi2,npoints,fitParam,covMatrix,-1,0,200000); - -chainCosmic->SetAlias("corQT",strqdedx->Data()); - -*/ - - -/* - chainCosmic->SetProof(kTRUE); - chainCosmic->Draw("Seed0.CookdEdxNorm(0,0.6,1,0,159,0,kTRUE,kTRUE):Seed0.CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE)",""+cutA,"",100000); - - -chainCosmic->Draw("Seed0.CookdEdxNorm(0,0.6,1,0,159,0,kTRUE,kTRUE)/Seed1.CookdEdxNorm(0,0.6,1,0,159,0,kTRUE,kTRUE)>>his(100,0.5,1.5)","min(Orig0.fTPCncls,Orig1.fTPCncls)>130"+cutA,"",50000); - -*/ - - -/* -chainCosmic->Draw("Tr0.fP[1]-Tr1.fP[1]:sqrt(abs(Tr0.fP[4]))>>hisdzA(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0&&abs(mag)>0.1"+cutA); - -TGraph *grdzA = (TGraph*)gProof->GetOutputList()->At(1)->Clone(); - - - - -*/ - - - - diff --git a/TPC/AliTPCcalibCosmic.h b/TPC/AliTPCcalibCosmic.h index ffca65ad10f..35c3e976c45 100644 --- a/TPC/AliTPCcalibCosmic.h +++ b/TPC/AliTPCcalibCosmic.h @@ -30,12 +30,15 @@ public: // void FindPairs(AliESDEvent *event); Bool_t IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1); - void SetGainMap(AliTPCCalPad *GainMap){fGainMap = GainMap;}; static void CalculateBetheParams(TH2F *hist, Double_t * initialParam); static Double_t CalculateMIPvalue(TH1F * hist); - AliExternalTrackParam *Invert(AliExternalTrackParam *input); AliExternalTrackParam *MakeTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1); + AliExternalTrackParam *MakeCombinedTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1); + void UpdateTrack(AliExternalTrackParam &track0, const AliExternalTrackParam &track1); + void SetComparison(AliExternalComparison * comp) { fComp=comp;} + AliExternalComparison *GetComparison() { return fComp;} + // TH1F * GetHistNTracks(){return fHistNTracks;}; TH1F * GetHistClusters(){return fClusters;}; @@ -55,7 +58,7 @@ private: void FillAcordeHist(AliESDtrack *upperTrack); - AliTPCCalPad *fGainMap; // gain map from Krypton calibration + AliExternalComparison * fComp; // comparison histogram TH1F *fHistNTracks; // histogram showing number of ESD tracks per event TH1F *fClusters; // histogram showing the number of clusters per track TH2F *fModules; // 2d histogram of tracks which are propagated to the ACORDE scintillator array diff --git a/TPC/AliTPCcalibTracks.cxx b/TPC/AliTPCcalibTracks.cxx index 14c86d4bba8..ee6f472738f 100644 --- a/TPC/AliTPCcalibTracks.cxx +++ b/TPC/AliTPCcalibTracks.cxx @@ -350,22 +350,22 @@ AliTPCcalibTracks::AliTPCcalibTracks(const Text_t *name, const Text_t *title, Al // amplitude sprintf(chname,"Amp_Sector%d",i); - his1 = new TH1F(chname,chname,250,0,500); // valgrind + his1 = new TH1F(chname,chname,100,0,32); // valgrind his1->SetXTitle("Max Amplitude (ADC)"); fArrayAmp->AddAt(his1,i); sprintf(chname,"Amp_Sector%d",i+36); - his1 = new TH1F(chname,chname,200,0,600); // valgrind 3 13,408,208 bytes in 229 blocks are still reachable + his1 = new TH1F(chname,chname,100,0,32); // valgrind 3 13,408,208 bytes in 229 blocks are still reachable his1->SetXTitle("Max Amplitude (ADC)"); fArrayAmp->AddAt(his1,i+36); // driftlength sprintf(chname, "driftlengt vs. charge, ROC %i", i); - prof1 = new TProfile(chname, chname, 500, 0, 250); + prof1 = new TProfile(chname, chname, 25, 0, 250); prof1->SetYTitle("Charge"); prof1->SetXTitle("Driftlength"); fArrayChargeVsDriftlength->AddAt(prof1,i); sprintf(chname, "driftlengt vs. charge, ROC %i", i+36); - prof1 = new TProfile(chname, chname, 500, 0, 250); + prof1 = new TProfile(chname, chname, 25, 0, 250); prof1->SetYTitle("Charge"); prof1->SetXTitle("Driftlength"); fArrayChargeVsDriftlength->AddAt(prof1,i+36); @@ -889,16 +889,17 @@ void AliTPCcalibTracks::FillResolutionHistoLocal(AliTPCseed * track){ TProfile *profAmpRow = (TProfile*)fArrayAmpRow->At(sector); if ( irow >= kFirstLargePad) max /= kLargePadSize; - profAmpRow->Fill( (Double_t)cluster0->GetRow(), max ); + Double_t smax = TMath::Sqrt(max); + profAmpRow->Fill( (Double_t)cluster0->GetRow(), smax ); TH1F *hisAmp = (TH1F*)fArrayAmp->At(sector); - hisAmp->Fill(max); + hisAmp->Fill(smax); // remove the following two lines one day: TProfile *profDriftLength = (TProfile*)fArrayChargeVsDriftlength->At(sector); - profDriftLength->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), max ); + profDriftLength->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), smax ); TProfile *profDriftLengthTmp = (TProfile*)(fcalPadRegionChargeVsDriftlength->GetObject(isegment, padSize)); - profDriftLengthTmp->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), max ); + profDriftLengthTmp->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), smax ); fHclusterPerPadrow->Fill(irow); // fill histogram showing clusters per padrow Int_t ipad = TMath::Nint(cluster0->GetPad()); -- 2.39.3