]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCcalibCalib.cxx
AliTPCcalibDButil and Summary text:
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibCalib.cxx
index e308d66c6091cc1e34bf84c813ad5194c269b03f..cbdc154294bd6c2a0777bab3d70cb441ece579fb 100644 (file)
 #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 <TGeoManager.h>
+#include <TGeoPhysicalNode.h>
+#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
@@ -74,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);
@@ -82,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
@@ -126,91 +150,100 @@ void     AliTPCcalibCalib::Process(AliESDEvent *event){
   //
 
   for (Int_t i=0;i<ntracks;++i) {
-    AliESDtrack *track = event->GetTrack(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<AliTPCseed*>(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
   //
+  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
   //
-  AliTPCPointCorrection * corr =  AliTPCPointCorrection::Instance();
+  //  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 =AliTPCClusterParam::SPosCorrection(0,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
-    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
+    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());
     //
-    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;
-    }
+    cluster->SetX(x[0]);
+    cluster->SetY(x[1]);
+    cluster->SetZ(x[2]);
+    
     //
-    // Apply r-phi correction  - To be done on track level- knowing the track angle !!!
+    // Apply alignemnt
     //
-    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
+    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]);
     }
 
-    //
-    //
-    //
-    cluster->SetX(x[0]);
-    cluster->SetY(x[1]);
-    cluster->SetZ(x[2]);
-    if (fStreamLevel>2){
+
+
+    if (fStreamLevel>2 && streamCounter<20*fStreamLevel ){
+      // dump debug info if required
       TTreeSRedirector *cstream = GetDebugStreamer();
       if (cstream){
        (*cstream)<<"Clusters"<<
@@ -224,10 +257,6 @@ Bool_t  AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
          "cl.="<<cluster<<
          "cy="<<dy<<
          "cz="<<dz<<
-         "cR="<<corrR<<
-         "dxq="<<dxq<<
-         "dyq="<<dyq<<
-         "dzq="<<dzq<<
          "\n";
       }
     }
@@ -236,54 +265,74 @@ Bool_t  AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
   //
   //
   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]=10.*10.;
-  covar[2]=10.*10.;
-  covar[5]=10.*10./(64.*64.);
-  covar[9]=10.*10./(64.*64.);
-  covar[14]=1*1;
-
+  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->GetInnerParam();
-  AliExternalTrackParam * trackOutOld = (AliExternalTrackParam*)track->GetOuterParam();
+  AliExternalTrackParam * trackOuter = (AliExternalTrackParam*)track->GetOuterParam();
+  AliExternalTrackParam * trackOutOld   = (AliExternalTrackParam *)friendTrack->GetTPCOut();
+  Double_t mass =    TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
+
+    
 
   AliExternalTrackParam trackIn  = *trackOutOld;
-  trackIn.ResetCovariance(200.);
+  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)
-      trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.));
+    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];
     trackIn.GetXYZ(xyz);
-    Double_t bz = AliTracker::GetBz(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 (!trackIn.PropagateTo(r[0],bz)) continue;
     if (RejectCluster(cl,&trackIn)) continue;
     nclIn++;
     trackIn.Update(&r[1],cov);    
   }
   //
   AliExternalTrackParam trackOut = trackIn;
-  trackOut.ResetCovariance(200.);
+  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
   //
@@ -295,8 +344,9 @@ 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 (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    
@@ -304,8 +354,10 @@ Bool_t  AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
     cov[0]*=cov[0];
     cov[2]*=cov[2];
     trackOut.GetXYZ(xyz);
-    Double_t bz = AliTracker::GetBz(xyz);
-    if (!trackOut.PropagateTo(r[0],bz)) continue;
+    //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);        
@@ -317,7 +369,11 @@ Bool_t  AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
   //
   nclIn=0;
   trackIn  = trackOut;
-  trackIn.ResetCovariance(10.);
+  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
   //
@@ -327,17 +383,20 @@ Bool_t  AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
     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 (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];
     trackIn.GetXYZ(xyz);
-    Double_t bz = AliTracker::GetBz(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 (!trackIn.PropagateTo(r[0],bz)) continue;
     if (RejectCluster(cl,&trackIn)) continue;
     nclIn++;
     trackIn.Update(&r[1],cov);    
@@ -356,7 +415,7 @@ 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"<<
@@ -382,10 +441,11 @@ Bool_t  AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
  
   (*trackInOld)  = trackIn;
   (*trackOutOld) = trackOut;
+  (*trackOuter) = trackOut;
   AliExternalTrackParam *t = &trackIn;
-  track->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
+  //track->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
   seed->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
-  seed->SetNumberOfClusters((nclIn+nclOut)*0.5);
+  seed->SetNumberOfClusters((nclIn+nclOut)/2);
   return kTRUE;
 }
 
@@ -407,6 +467,10 @@ Bool_t AliTPCcalibCalib::RejectCluster(AliTPCclusterMI* cl, AliExternalTrackPara
 
   Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation    
   AliTPCseed::GetError(cl, param,cov[0],cov[2]);
+  if (param->GetSigmaY2()<0 || param->GetSigmaZ2()<0){
+    AliError("Wrong parameters");
+    return kFALSE;
+  }
   Double_t py = (cl->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());
   //