X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=TPC%2FAliTPCcalibCalib.cxx;h=30358727a28934dedd4c6c178624863bc31c1acf;hb=02d1b157cc65df0dfa2fba9e7dceb296f216e2fc;hp=706077dff66d2e2f8fb32f854510c31ff1936f5d;hpb=2cfb8d90e47f9aa22ead60da13dca90eb67af48c;p=u%2Fmrichter%2FAliRoot.git diff --git a/TPC/AliTPCcalibCalib.cxx b/TPC/AliTPCcalibCalib.cxx index 706077dff66..30358727a28 100644 --- a/TPC/AliTPCcalibCalib.cxx +++ b/TPC/AliTPCcalibCalib.cxx @@ -31,8 +31,10 @@ gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros"); gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+") AliXRDPROOFtoolkit tool; - TChain * chain = tool.MakeChain("esd.txt","esdTree",0,5000); - chain->Lookup(); + TChain * chainCl = tool.MakeChain("calib.txt","Clusters",0,1); + chainCl->Lookup(); + TChain * chainTr = tool.MakeChain("calib.txt","Tracks",0,1); + chainTr->Lookup(); @@ -52,16 +54,29 @@ #include "AliESDfriend.h" #include "AliESDtrack.h" #include "AliTracker.h" +#include "AliTPCClusterParam.h" +#include "AliTPCParam.h" #include "AliTPCcalibDB.h" #include "AliTPCTransform.h" +#include "AliTPCRecoParam.h" #include "AliTPCclusterMI.h" #include "AliTPCseed.h" - +#include "AliTPCPointCorrection.h" +#include +#include +#include "TDatabasePDG.h" ClassImp(AliTPCcalibCalib) AliTPCcalibCalib::AliTPCcalibCalib(): - AliTPCcalibBase() +AliTPCcalibBase(), + fApplyExBCorrection(1), // apply ExB correction + fApplyTOFCorrection(1), // apply TOF correction + fApplyPositionCorrection(1), // apply position correction + fApplySectorAlignment(1), // apply sector alignment + fApplyRPhiCorrection(1), // apply R-Phi correction + fApplyRCorrection(1) // apply Radial correction + { // // Constructor @@ -70,7 +85,13 @@ AliTPCcalibCalib::AliTPCcalibCalib(): AliTPCcalibCalib::AliTPCcalibCalib(const Text_t *name, const Text_t *title) - :AliTPCcalibBase() + :AliTPCcalibBase(), + fApplyExBCorrection(1), // apply ExB correction + fApplyTOFCorrection(1), // apply TOF correction + fApplyPositionCorrection(1), // apply position correction + fApplySectorAlignment(1), // apply sector alignment + fApplyRPhiCorrection(1), // apply R-Phi correction + fApplyRCorrection(1) // apply Radial correction { SetName(name); SetTitle(title); @@ -78,7 +99,14 @@ AliTPCcalibCalib::AliTPCcalibCalib(const Text_t *name, const Text_t *title) AliTPCcalibCalib::AliTPCcalibCalib(const AliTPCcalibCalib&calib): - AliTPCcalibBase(calib) + AliTPCcalibBase(calib), + fApplyExBCorrection(calib.GetApplyExBCorrection()), + fApplyTOFCorrection(calib.GetApplyTOFCorrection()), + fApplyPositionCorrection(calib.GetApplyPositionCorrection()), + fApplySectorAlignment(calib.GetApplySectorAlignment()), + fApplyRPhiCorrection(calib.GetApplyRPhiCorrection()), + fApplyRCorrection(calib.GetApplyRCorrection()) + { // // copy constructor @@ -115,76 +143,209 @@ void AliTPCcalibCalib::Process(AliESDEvent *event){ if (GetDebugLevel()>20) printf("Hallo world: Im here\n"); Int_t ntracks=event->GetNumberOfTracks(); + //AliTPCcalibDB::Instance()->SetExBField(fMagF); + // // // for (Int_t i=0;iGetTrack(i); - const AliExternalTrackParam * trackIn = track->GetInnerParam(); + AliESDtrack *track = event->GetTrack(i); + AliESDfriendTrack *friendTrack = (AliESDfriendTrack*) ESDfriend->GetTrack(i); + if (!friendTrack) continue; + //track->SetFriendTrack(friendTrack); + fCurrentFriendTrack=friendTrack; + const AliExternalTrackParam * trackIn = track->GetInnerParam(); const AliExternalTrackParam * trackOut = track->GetOuterParam(); + AliExternalTrackParam * tpcOut = (AliExternalTrackParam *)friendTrack->GetTPCOut(); if (!trackIn) continue; if (!trackOut) continue; - - AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i); + if (!tpcOut) continue; TObject *calibObject; AliTPCseed *seed = 0; for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) { if ((seed=dynamic_cast(calibObject))) break; } if (!seed) continue; - RefitTrack(track, seed); + RefitTrack(track, seed, event->GetMagneticField()); + (*tpcOut)=*(track->GetOuterParam()); } return; } -Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ +Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed, Float_t magesd){ // // Refit track + // if magesd==0 forget the curvature + // + // 0 - Setup transform object + // + const Double_t kxIFC = 83.; // position of IFC + const Double_t kxOFC = 250.; // position of OFC + const Double_t kaFC = 1.; // amplitude + const Double_t ktFC = 5.0; // slope of error + //cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); + //cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); + + static Int_t streamCounter=0; + streamCounter++; + AliESDfriendTrack *friendTrack = fCurrentFriendTrack; + + AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ; + AliTPCParam *param = AliTPCcalibDB::Instance()->GetParameters(); + transform->SetCurrentRun(fRun); + transform->SetCurrentTimeStamp((UInt_t)fTime); + if(!fApplyExBCorrection) { // disable ExB correction in transform + if(transform->GetCurrentRecoParam()) + transform->GetCurrentRecoParamNonConst()->SetUseExBCorrection(0); + } + if(!fApplyTOFCorrection) { // disable TOF correction in transform + if(transform->GetCurrentRecoParam()) + transform->GetCurrentRecoParamNonConst()->SetUseTOFCorrection(kFALSE); + } // // 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; AliTPCclusterMI cl0(*cluster); Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()}; Int_t i[1]={cluster->GetDetector()}; + transform->Transform(x,i,0,1); // + // get position correction + // + Int_t ipad=0; + if (cluster->GetDetector()>35) ipad=1; + Float_t dy =0;//AliTPCClusterParam::SPosCorrection(0,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax()); + Float_t dz =0;//AliTPCClusterParam::SPosCorrection(1,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax()); + // cluster->SetX(x[0]); cluster->SetY(x[1]); cluster->SetZ(x[2]); - if (fStreamLevel>2){ + + // + // Apply alignemnt + // + if (transform->GetCurrentRecoParam()->GetUseSectorAlignment()){ + if (!param->IsGeoRead()) param->ReadGeoMatrices(); + TGeoHMatrix *mat = param->GetClusterMatrix(cluster->GetDetector()); + //TGeoHMatrix mat; + Double_t pos[3]= {cluster->GetX(),cluster->GetY(),cluster->GetZ()}; + Double_t posC[3]={cluster->GetX(),cluster->GetY(),cluster->GetZ()}; + if (mat) mat->LocalToMaster(pos,posC); + else{ + // chack Loading of Geo matrices from GeoManager - TEMPORARY FIX + } + cluster->SetX(posC[0]); + cluster->SetY(posC[1]); + cluster->SetZ(posC[2]); + } + + + + if (fStreamLevel>2 && streamCounter<20*fStreamLevel ){ + // dump debug info if required TTreeSRedirector *cstream = GetDebugStreamer(); if (cstream){ (*cstream)<<"Clusters"<< + "run="<GetNumberOfClusters(); + // + // + // + Int_t ncl = seed->GetNumberOfClusters(); + const Double_t kResetCov=4.; + const Double_t kSigma=5.; + Double_t covar[15]; + for (Int_t i=0;i<15;i++) covar[i]=0; + covar[0]=kSigma*kSigma; + covar[2]=kSigma*kSigma; + covar[5]=kSigma*kSigma/Float_t(ncl*ncl); + covar[9]=kSigma*kSigma/Float_t(ncl*ncl); + covar[14]=0.2*0.2; + if (TMath::Abs(magesd)<0.05) { + covar[14]=0.025*0.025; + } // // And now do refit // - AliExternalTrackParam * trackInOld = (AliExternalTrackParam*)track->GetTPCInnerParam(); - AliExternalTrackParam * trackOutOld = (AliExternalTrackParam*)track->GetOuterParam(); + AliExternalTrackParam * trackInOld = (AliExternalTrackParam*)track->GetInnerParam(); + AliExternalTrackParam * trackOuter = (AliExternalTrackParam*)track->GetOuterParam(); + AliExternalTrackParam * trackOutOld = (AliExternalTrackParam *)friendTrack->GetTPCOut(); + Double_t mass = TDatabasePDG::Instance()->GetParticle("pi+")->Mass(); + + AliExternalTrackParam trackIn = *trackOutOld; - AliExternalTrackParam trackOut = *trackInOld; - trackIn.ResetCovariance(10.); - trackOut.ResetCovariance(10.); + trackIn.ResetCovariance(kResetCov); + trackIn.AddCovariance(covar); + if (TMath::Abs(magesd)<0.05) { + ((Double_t&)(trackIn.GetParameter()[4]))=0.000000001; + ((Double_t&)(trackIn.GetCovariance()[14]))=covar[14]; // fix the line + } + 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){ + if (!trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break; + } + 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]; + cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); + cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); + trackIn.GetXYZ(xyz); + // Double_t bz = AliTracker::GetBz(xyz); + + // if (!trackIn.PropagateTo(r[0],bz)) continue; + if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1.,kFALSE)) continue; + + if (RejectCluster(cl,&trackIn)) continue; + nclIn++; + trackIn.Update(&r[1],cov); + } + // + AliExternalTrackParam trackOut = trackIn; + trackOut.ResetCovariance(kResetCov); + trackOut.AddCovariance(covar); + if (TMath::Abs(magesd)<0.05) { + ((Double_t&)(trackOut.GetParameter()[4]))=0.000000001; + ((Double_t&)(trackOut.GetCovariance()[14]))=covar[14]; // fix the line + } + // // Refit out // + //Bool_t lastEdge=kFALSE; for (Int_t irow=0; irow<160; irow++){ AliTPCclusterMI *cl=seed->GetClusterPointer(irow); if (!cl) continue; @@ -192,38 +353,70 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ Int_t sector = cl->GetDetector(); Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha(); - if (TMath::Abs(dalpha)>0.01) - trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.)); - if (RejectCluster(cl,&trackOut)) continue; + if (TMath::Abs(dalpha)>0.01){ + if (!trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break; + } Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()}; - Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation + + Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation + AliTPCseed::GetError(cl, &trackOut,cov[0],cov[2]); + cov[0]*=cov[0]; + cov[2]*=cov[2]; + cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); + cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); trackOut.GetXYZ(xyz); - Double_t bz = AliTracker::GetBz(xyz); - if (trackOut.PropagateTo(r[0],bz)) nclOut++; - trackOut.Update(&r[1],cov); - + //Double_t bz = AliTracker::GetBz(xyz); + // if (!trackOut.PropagateTo(r[0],bz)) continue; + if (!AliTracker::PropagateTrackToBxByBz(&trackOut, r[0],mass,1.,kFALSE)) continue; + + if (RejectCluster(cl,&trackOut)) continue; + 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(kResetCov); + if (TMath::Abs(magesd)<0.05) { + ((Double_t&)(trackIn.GetParameter()[4]))=0.000000001; + ((Double_t&)(trackIn.GetCovariance()[14]))=covar[14]; // fix the line + } + // + // Refit in one more time + // 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.)); - if (RejectCluster(cl,&trackIn)) continue; + if (TMath::Abs(dalpha)>0.01){ + if (!trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break; + } Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()}; Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation - trackOut.GetXYZ(xyz); - Double_t bz = AliTracker::GetBz(xyz); + AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]); + cov[0]*=cov[0]; + cov[2]*=cov[2]; + cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); + cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC)); + + trackIn.GetXYZ(xyz); + //Double_t bz = AliTracker::GetBz(xyz); - if (trackIn.PropagateTo(r[0],bz)) nclIn++; + // if (!trackIn.PropagateTo(r[0],bz)) continue; + if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1,kFALSE)) continue; + + if (RejectCluster(cl,&trackIn)) continue; + nclIn++; trackIn.Update(&r[1],cov); } + + trackIn.Rotate(trackInOld->GetAlpha()); trackOut.Rotate(trackOutOld->GetAlpha()); // @@ -236,10 +429,16 @@ Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){ trackOut.PropagateTo(trackOutOld->GetX(),bz); - if (fStreamLevel>0){ + if (fStreamLevel>0 && streamCounter<100*fStreamLevel){ TTreeSRedirector *cstream = GetDebugStreamer(); if (cstream){ (*cstream)<<"Tracks"<< + "run="<GetY()-param->GetY())/TMath::Sqrt(cov[0]*cov[0]+param->GetSigmaY2()); + Double_t pz = (cl->GetZ()-param->GetZ())/TMath::Sqrt(cov[2]*cov[2]+param->GetSigmaZ2()); + // + if ((py*py+pz*pz)>kSigmaCut*kSigmaCut) isReject=kTRUE; + return isReject; } +