Restart the fit from TPC backup parameters in case

[u/mrichter/AliRoot.git] / TPC / AliTPCtrack.cxx

diff --git a/TPC/AliTPCtrack.cxx b/TPC/AliTPCtrack.cxx
index 10741b91b3dd88fb0d8735c468364b2889427d58..d5f388e86dbe860e849e440847e2565947edabd1 100644 (file)
--- a/TPC/AliTPCtrack.cxx
+++ b/TPC/AliTPCtrack.cxx
@@ -28,6 +28,7 @@
 #include "AliTracker.h"
 #include "AliESDtrack.h"
 #include "AliESDVertex.h"
+#include "TTreeStream.h"
 
 ClassImp(AliTPCtrack)
 
@@ -120,7 +121,7 @@ AliTPCtrack::AliTPCtrack(Double_t x, Double_t alpha, const Double_t p[5],
 }
 
 //_____________________________________________________________________________
-AliTPCtrack::AliTPCtrack(const AliESDtrack& t) :
+AliTPCtrack::AliTPCtrack(const AliESDtrack& t, TTreeSRedirector *pcstream) :
   AliKalmanTrack(),
   fdEdx(t.GetTPCsignal()),
   fSdEdx(1e10),
@@ -137,6 +138,7 @@ AliTPCtrack::AliTPCtrack(const AliESDtrack& t) :
   //-----------------------------------------------------------------
   // Conversion AliESDtrack -> AliTPCtrack.
   //-----------------------------------------------------------------
+  const Double_t kmaxC[4]={10,10,0.1,0.1};  // cuts on the rms /fP0,fP1,fP2,fP3
   SetNumberOfClusters(t.GetTPCclusters(fIndex));
   SetLabel(t.GetLabel());
   SetMass(t.GetMass());
@@ -144,8 +146,55 @@ AliTPCtrack::AliTPCtrack(const AliESDtrack& t) :
   for (Int_t i=0; i<12;i++) fKinkPoint[i]=0.;
   for (Int_t i=0; i<3;i++) fKinkIndexes[i]=0;
   for (Int_t i=0; i<3;i++) fV0Indexes[i]=0;
+  //
+  // choose parameters to start
+  //
+  Int_t reject=0;
+  AliExternalTrackParam param(t);
+  const AliExternalTrackParam  *tpcout=(t.GetFriendTrack())? ((AliESDfriendTrack*)(t.GetFriendTrack()))->GetTPCOut():0;
+  const AliExternalTrackParam  *tpcin = t.GetInnerParam();
+  const AliExternalTrackParam  *tpc=(tpcout)?tpcout:tpcin;
+  if (!tpc) tpc=&param;
+  Bool_t isOK=kTRUE;
+  if (param.GetCovariance()[0]>kmaxC[0]*kmaxC[0]) isOK=kFALSE;
+  if (param.GetCovariance()[2]>kmaxC[1]*kmaxC[1]) isOK=kFALSE;
+  if (param.GetCovariance()[5]>kmaxC[2]*kmaxC[2]) isOK=kFALSE;
+  if (param.GetCovariance()[9]>kmaxC[3]*kmaxC[3]) isOK=kFALSE;
+  if (!isOK){
+    param=*tpc;
+    isOK=kTRUE;
+    reject=1;
+  }
+  param.Rotate(tpc->GetAlpha());
+  isOK=AliTracker::PropagateTrackToBxByBz(&param,tpc->GetX(),t.GetMass(),2.,kFALSE);
+  if (param.GetCovariance()[0]>kmaxC[0]*kmaxC[0]) isOK=kFALSE;
+  if (param.GetCovariance()[2]>kmaxC[1]*kmaxC[1]) isOK=kFALSE;
+  if (param.GetCovariance()[5]>kmaxC[2]*kmaxC[2]) isOK=kFALSE;
+  if (param.GetCovariance()[9]>kmaxC[3]*kmaxC[3]) isOK=kFALSE;
+  if (!isOK){
+    param=*tpc;
+    isOK=kTRUE;
+    reject=2;
+  }
+  if (reject>0){
+    param.ResetCovariance(4.);  // reset covariance if start from backup param
+  }
+  //
+  //
+  if (pcstream){
+    AliExternalTrackParam dummy;
+    AliExternalTrackParam *ptpc=(AliExternalTrackParam *)tpc;
+    if (!ptpc) ptpc=&dummy;
+    AliESDtrack *esd= (AliESDtrack *)&t;
+    (*pcstream)<<"trackP"<<
+      "reject="<<reject<<   // flag - rejection of current esd track parameters
+      "esd.="<<esd<<        // original esd track
+      "tr.="<<&param<<      // starting track parameters
+      "out.="<<ptpc<<       // backup tpc parameters
+      "\n";
+  }
 
-  Set(t.GetX(),t.GetAlpha(),t.GetParameter(),t.GetCovariance());
+  Set(param.GetX(),param.GetAlpha(),param.GetParameter(),param.GetCovariance());
 
   if ((t.GetStatus()&AliESDtrack::kTIME) == 0) return;
   StartTimeIntegral();
@@ -239,10 +288,21 @@ Bool_t AliTPCtrack::PropagateTo(Double_t xk,Double_t rho,Double_t x0) {
   //  rho - density of the crossed matrial (g/cm^3)
   //  x0  - radiation length of the crossed material (g/cm^2) 
   //-----------------------------------------------------------------
-  Double_t oldX=GetX(), oldY=GetY(), oldZ=GetZ();
-
+  //
   Double_t bz=GetBz();
-  if (!AliExternalTrackParam::PropagateTo(xk,bz)) return kFALSE;
+  Double_t zat=0;
+  if (!GetZAt(xk, bz,zat)) return kFALSE;
+  if (TMath::Abs(zat)>250.){
+    // Don't propagate track outside of the fiducial volume - material budget not proper one
+    //
+    //AliWarning("Propagate outside of fiducial volume");
+    return kFALSE;
+  }
+
+  Double_t oldX=GetX(), oldY=GetY(), oldZ=GetZ();
+  //if (!AliExternalTrackParam::PropagateTo(xk,bz)) return kFALSE;
+  Double_t b[3]; GetBxByBz(b);
+  if (!AliExternalTrackParam::PropagateToBxByBz(xk,b)) return kFALSE;
 
   Double_t d = TMath::Sqrt((GetX()-oldX)*(GetX()-oldX) + 
                            (GetY()-oldY)*(GetY()-oldY) + 
@@ -250,23 +310,37 @@ Bool_t AliTPCtrack::PropagateTo(Double_t xk,Double_t rho,Double_t x0) {
   if (IsStartedTimeIntegral() && GetX()>oldX) AddTimeStep(d);
 
   if (oldX < xk) d = -d;
-  if (!AliExternalTrackParam::CorrectForMaterial(d*rho/x0,x0,GetMass())) 
-     return kFALSE;
+  if (!AliExternalTrackParam::CorrectForMeanMaterial(d*rho/x0,d*rho,GetMass(),
+      kFALSE,AliExternalTrackParam::BetheBlochGas)) return kFALSE;
 
   return kTRUE;
 }
 
 //_____________________________________________________________________________
 Bool_t 
-AliTPCtrack::PropagateToVertex(const AliESDVertex *v,Double_t d,Double_t x0) 
+AliTPCtrack::PropagateToVertex(const AliESDVertex *v,Double_t rho,Double_t x0) 
 {
   //-----------------------------------------------------------------
-  // This function propagates tracks to the vertex.
+  // This function propagates tracks to the vertex
+  // rho - density of the crossed matrial (g/cm3)
+  // x0  - radiation length of the crossed material (g/cm2) 
   //-----------------------------------------------------------------
-  Double_t bz=GetBz();
-  if (PropagateToDCA(v,bz,kVeryBig))
-   if (AliExternalTrackParam::CorrectForMaterial(d,x0,GetMass())) return kTRUE;
-  return kFALSE;
+  Double_t oldX=GetX(), oldY=GetY(), oldZ=GetZ();
+
+  //Double_t bz=GetBz();
+  //if (!PropagateToDCA(v,bz,kVeryBig)) return kFALSE;
+  Double_t b[3]; GetBxByBz(b);
+  if (!PropagateToDCABxByBz(v,b,kVeryBig)) return kFALSE;
+
+  Double_t d = TMath::Sqrt((GetX()-oldX)*(GetX()-oldX) + 
+                           (GetY()-oldY)*(GetY()-oldY) + 
+                           (GetZ()-oldZ)*(GetZ()-oldZ));
+
+  if (oldX < GetX()) d = -d;
+  if (!AliExternalTrackParam::CorrectForMeanMaterial(d*rho/x0,d*rho,GetMass(),
+      kFALSE,AliExternalTrackParam::BetheBlochGas)) return kFALSE;
+
+  return kTRUE;
 }
 
 //_____________________________________________________________________________
@@ -390,12 +464,3 @@ void  AliTPCtrack::UpdatePoints()
   //
 }
 
-Double_t AliTPCtrack::GetBz() const {
-  //
-  // returns Bz component of the magnetic field (kG)
-  //
-  if (AliTracker::UniformField()) return AliTracker::GetBz();
-  Double_t r[3]; GetXYZ(r);
-  return AliTracker::GetBz(r);
-}
-