Fix from Alex for FPE in TRD: ALIROOT-5503 and tracklet calibration: PWGPP-2
authorshahoian <ruben.shahoyan@cern.ch>
Wed, 9 Jul 2014 13:05:27 +0000 (15:05 +0200)
committershahoian <ruben.shahoyan@cern.ch>
Wed, 9 Jul 2014 13:05:27 +0000 (15:05 +0200)
PWGPP/TRD/AliTRDcheckTRK.cxx
TRD/AliTRDrecoParam.cxx
TRD/AliTRDrecoParam.h
TRD/AliTRDseedV1.cxx
TRD/AliTRDseedV1.h
TRD/AliTRDtrackerV1.cxx

index 5074617..11a731d 100644 (file)
@@ -240,11 +240,9 @@ Bool_t AliTRDcheckTRK::PropagateKalman(AliTRDtrackV1 &t, AliExternalTrackParam *
     }
     if(HasTrkltRefit()){
       //      if(!tr->FitRobust(tt.Charge()>0.)) printf("W - AliTRDcheckTRK::PropagateKalman :: FitRobust() failed for Det[%03d]\n", det);
-      if(!tr->FitRobust(AliTRDgeometry::GetPadPlane(det), prod>0., tt.Charge())) printf("W - AliTRDcheckTRK::PropagateKalman :: FitRobust() failed for Det[%03d]\n", det);
-      else {
-       TGeoHMatrix *matrix = AliTRDgeometry::GetClusterMatrix(det);
-       if (matrix) tr->SetXYZ(matrix);
-      }
+      TGeoHMatrix *matrix = AliTRDgeometry::GetClusterMatrix(det);
+      if(!tr->FitRobust(AliTRDgeometry::GetPadPlane(det), matrix, tt.GetBz(), tt.Charge())) printf("W - AliTRDcheckTRK::PropagateKalman :: FitRobust() failed for Det[%03d]\n", det);
+      else tr->SetXYZ(matrix);
     }
     if(!AliTRDtrackerV1::PropagateToX(tt, tr->GetX0(), fgKalmanStep)) continue;
     if(!tt.GetTrackIn()) tt.SetTrackIn();
index f790fea..294d074 100644 (file)
@@ -156,8 +156,8 @@ AliTRDrecoParam::AliTRDrecoParam(const AliTRDrecoParam &ref)
   // tracklet params
   memcpy(fdzdxCorrFactor, ref.fdzdxCorrFactor, 2*sizeof(Double_t));
   memcpy(fdzdxCorrRCbias, ref.fdzdxCorrRCbias, 2*sizeof(Double_t));
-  memcpy(fYcorrTailCancel, ref.fdzdxCorrRCbias, 6*sizeof(Double_t));
-  memcpy(fS2Ycorr, ref.fS2Ycorr, 2*sizeof(Double_t));
+  memcpy(fYcorrTailCancel, ref.fdzdxCorrRCbias, 12*sizeof(Double_t));
+  memcpy(fS2Ycorr, ref.fS2Ycorr, 4*sizeof(Double_t));
 }
 
 //______________________________________________________________
@@ -215,8 +215,8 @@ AliTRDrecoParam& AliTRDrecoParam::operator=(const AliTRDrecoParam &ref)
   // tracklet params
   memcpy(fdzdxCorrFactor, ref.fdzdxCorrFactor, 2*sizeof(Double_t));
   memcpy(fdzdxCorrRCbias, ref.fdzdxCorrRCbias, 2*sizeof(Double_t));
-  memcpy(fYcorrTailCancel, ref.fdzdxCorrRCbias, 6*sizeof(Double_t));
-  memcpy(fS2Ycorr, ref.fS2Ycorr, 2*sizeof(Double_t));
+  memcpy(fYcorrTailCancel, ref.fdzdxCorrRCbias, 12*sizeof(Double_t));
+  memcpy(fS2Ycorr, ref.fS2Ycorr, 4*sizeof(Double_t));
   return *this;
 }
 
@@ -351,13 +351,16 @@ void  AliTRDrecoParam::SetTrackletParams(Double_t *par)
     fdzdxCorrRCbias[1] = par[3];   // dz/dx < 0
     /// correct x_cross for the bias in dzdx
     fdzdxXcrossFactor  = par[4];
-    // y linear q/pt correction due to wrong tail cancellation. 
-    fYcorrTailCancel[0][0] = par[5]; fYcorrTailCancel[0][1] = par[6];  // opposite sign !RC
-    fYcorrTailCancel[1][0] = par[7]; fYcorrTailCancel[1][1] = par[8];  // same sign !RC
-    fYcorrTailCancel[2][0] = par[9]; fYcorrTailCancel[2][1] = par[10]; // RC
+    // y correction due to wrong tail cancellation. 
+    fYcorrTailCancel[0][0] = par[5];fYcorrTailCancel[0][1] = par[6];fYcorrTailCancel[0][2] = par[7]; 
+    fYcorrTailCancel[1][0] = par[8];fYcorrTailCancel[1][1] = par[9];fYcorrTailCancel[1][2] = par[10]; 
+    fYcorrTailCancel[2][0] = par[11];fYcorrTailCancel[2][1] = par[12];fYcorrTailCancel[2][3] = par[13];
+    fYcorrTailCancel[3][0] = par[14];fYcorrTailCancel[3][1] = par[15];fYcorrTailCancel[3][3] = par[16];
     // inflation factor of error parameterization in r-phi due to wrong estimation of residuals. 
-    fS2Ycorr[0] = par[11];  // opposite sign, 
-    fS2Ycorr[1] = par[12];  // same sign
+    fS2Ycorr[0] = par[17];
+    fS2Ycorr[1] = par[18];
+    fS2Ycorr[2] = par[19];
+    fS2Ycorr[3] = par[20];
     
   } else {
     // correct dzdx for the bias in z
@@ -368,12 +371,23 @@ void  AliTRDrecoParam::SetTrackletParams(Double_t *par)
     fdzdxCorrRCbias[1] = -0.012; // dz/dx < 0
     /// correct x_cross for the bias in dzdx
     fdzdxXcrossFactor  = 0.14;
-    // y linear q/pt correction due to wrong tail cancellation. 
-    fYcorrTailCancel[0][0] = 0.;   fYcorrTailCancel[0][1] = 0.027;  // opposite sign !RC
-    fYcorrTailCancel[1][0] = 0.04; fYcorrTailCancel[1][1] = 0.027;  // same sign !RC
-    fYcorrTailCancel[2][0] = 0.013;fYcorrTailCancel[2][1] = 0.018;  // RC
+    // y correction due to wrong tail cancellation. 
+        // bz<0 && !RC
+    fYcorrTailCancel[0][0] = 0.04; fYcorrTailCancel[0][1] = 2.151; fYcorrTailCancel[0][2] = 0.013;
+        // bz>0 && !RC
+    fYcorrTailCancel[1][0] = 0.034; fYcorrTailCancel[1][1] = 1.817; fYcorrTailCancel[1][2] = -0.01;
+        // bz<0 && RC
+    fYcorrTailCancel[2][0] = 0.04; fYcorrTailCancel[2][1] = 2.513; fYcorrTailCancel[2][2] = 0.015;
+        // bz>0 && RC
+    fYcorrTailCancel[3][0] = 0.034; fYcorrTailCancel[3][1] = 2.476; fYcorrTailCancel[3][2] = -0.01;
     // inflation factor of error parameterization in r-phi due to wrong estimation of residuals. 
-    fS2Ycorr[0] = 5.;  // opposite sign, 
-    fS2Ycorr[1] = 10;  // same sign
+        // chg<0 && !RC
+    fS2Ycorr[0] = 5.52; 
+        // chg>0 && !RC
+    fS2Ycorr[1] = 3.61; 
+        // chg<0 && RC
+    fS2Ycorr[2] = 4.84; 
+        // chg>0 && RC
+    fS2Ycorr[3] = 3.24; 
   }
 }
index d4c9cee..9edca2b 100644 (file)
@@ -94,8 +94,8 @@ public:
   Double_t GetCorrDZDXbiasRC(Bool_t dzdx) const { return fdzdxCorrRCbias[dzdx];}
   Double_t GetCorrDZDX(Bool_t rc) const     { return fdzdxCorrFactor[rc];}
   Double_t GetCorrDZDXxcross() const        { return fdzdxXcrossFactor;}
-  void GetYcorrTailCancel(Int_t it, Double_t par[2]) const;
-  Double_t GetS2Ycorr(Bool_t sgn) const     { return fS2Ycorr[sgn];}
+  void     GetYcorrTailCancel(Int_t it, Double_t par[3]) const;
+  Double_t GetS2Ycorr(Bool_t rc, Bool_t chg) const     { return fS2Ycorr[2*rc+chg];}
 
   Bool_t   IsArgon() const                  { return TESTBIT(fFlags, kDriftGas); }
   Bool_t   IsCheckTimeConsistency() const   { return kCheckTimeConsistency;}
@@ -204,8 +204,8 @@ private:
   Double_t  fdzdxCorrFactor[2];      // correction of dzdx estimation due to z bias; [0] for !RC, [1] for RC
   Double_t  fdzdxCorrRCbias[2];      // correction of dzdx estimation bias for RC; [0] for dz/dx>0, [1] for dz/dx<0
   Double_t  fdzdxXcrossFactor;       // bias in dzdx of estimated xcross [RC]
-  Double_t  fYcorrTailCancel[3][2];  // y linear q/pt correction due to wrong tail cancellation. [0] opposite sign !RC, [1] same sign !RC, [2] RC
-  Double_t  fS2Ycorr[2];             // inflation factor of error parameterization in r-phi due to wrong estimation of residuals. [0] opposite sign, [1] same sign
+  Double_t  fYcorrTailCancel[4][3];  // y correction due to wrong tail cancellation. [0] bz<0 && !RC, [1] bz>0 && !RC, [2] bz<0 && RC [3] bz>0 && RC
+  Double_t  fS2Ycorr[4];             // inflation factor of error parameterization in r-phi due to wrong estimation of residuals.
   
   // Clusterization parameter
   Double_t  fMinMaxCutSigma;         // Threshold sigma noise pad middle
@@ -274,10 +274,10 @@ inline void AliTRDrecoParam::SetTCParams(Double_t *par)
 
 
 //___________________________________________________
-inline void AliTRDrecoParam::GetYcorrTailCancel(Int_t it, Double_t par[2]) const
+inline void AliTRDrecoParam::GetYcorrTailCancel(Int_t it, Double_t par[3]) const
 {
-  if(it<0||it>2) return;
-  par[0] = fYcorrTailCancel[it][0]; par[1] = fYcorrTailCancel[it][1]; 
+  if(it<0||it>3) return;
+  memcpy(par, fYcorrTailCancel[it], 3*sizeof(Double_t));
 }
 
 //___________________________________________________
index b718f5e..9b0e506 100644 (file)
@@ -539,7 +539,7 @@ Int_t AliTRDseedV1::GetChargeGaps(Float_t sz[kNtb], Float_t pos[kNtb], Int_t isz
 
 
 //____________________________________________________________________
-Double_t AliTRDseedV1::EstimatedCrossPoint(AliTRDpadPlane *pp)
+Double_t AliTRDseedV1::EstimatedCrossPoint(AliTRDpadPlane *pp, Float_t bz)
 {
 // Algorithm to estimate cross point in the x-z plane for pad row cross tracklets or the z coordinate of pad row without pad row cross in the local chamber coordinates.
 // Returns variance of the radial offset from anode wire in case of raw cross or 0 otherwise.
@@ -633,7 +633,7 @@ Double_t AliTRDseedV1::EstimatedCrossPoint(AliTRDpadPlane *pp)
   fZfit[1] = (fZfit[0]+zoff)/dx; 
 
   // correct dzdx for the bias
-  UnbiasDZDX(IsRowCross());
+  UnbiasDZDX(IsRowCross(), bz);
   if(IsRowCross()){
     // correct x_cross/sigma(x_cross) for the bias in dzdx
     const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam();
@@ -652,32 +652,28 @@ Double_t AliTRDseedV1::EstimatedCrossPoint(AliTRDpadPlane *pp)
 }
 
 //____________________________________________________________________
-void AliTRDseedV1::UnbiasDZDX(Bool_t rc)
+void AliTRDseedV1::UnbiasDZDX(Bool_t rc, Float_t bz)
 {
   // correct dzdx for the bias in z according to MC
   const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam();
   if(!recoParam) return;
-  fZfit[1] *= recoParam->GetCorrDZDX(rc);
+  fZfit[1] *= recoParam->GetCorrDZDX(rc)-(bz>0?0.01:0.);
   if(rc) fZfit[1] += recoParam->GetCorrDZDXbiasRC(fZfit[1]<0);
 }
 
 //____________________________________________________________________
-Double_t AliTRDseedV1::UnbiasY(Bool_t rc, Bool_t sgn, Int_t chg)
+Double_t AliTRDseedV1::UnbiasY(Bool_t rc, Float_t bz)
 {
 // correct y coordinate for tail cancellation. This should be fixed by considering TC as a function of q/pt. 
 //  rc : TRUE if tracklet crosses rows
-// sgn : TRUE if track has same sign with magnetic field
-// chg : -1 for negative particles, +1 for the rest
+// bz : magnetic field z component
   
   const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam();
   if(!recoParam) return 0.;
-  Double_t par[2]={0.};
-  if(rc) recoParam->GetYcorrTailCancel(2, par);
-  else{
-    if(sgn && 1./fPt > 1.5)  recoParam->GetYcorrTailCancel(1, par);
-    else if(!sgn)  recoParam->GetYcorrTailCancel(0, par);
-  }
-  return par[0]+par[1]*chg/fPt;
+  Double_t par[3]={0.};
+  Int_t idx(2*(rc?1:0)+Int_t(bz>0));
+  recoParam->GetYcorrTailCancel(idx, par);
+  return par[0]*TMath::Sin(par[1]*fYref[1])+par[2];
 }
 
 
@@ -2161,7 +2157,7 @@ Bool_t AliTRDseedV1::Fit(UChar_t opt)
 
 
 //____________________________________________________________________
-Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t opt)
+Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, TGeoHMatrix *mDet, Float_t bz, Int_t chg, Int_t opt)
 {
 //
 // Linear fit of the clusters attached to the tracklet
@@ -2199,7 +2195,7 @@ Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t
   Double_t freq(cp?cp->GetSamplingFrequency():10.);
       
   // evaluate locally z and dzdx from TRD only information
-  if(EstimatedCrossPoint(pp)<0.) return kFALSE;
+  if(EstimatedCrossPoint(pp, bz)<0.) return kFALSE;
   
   //printf("D%03d RC[%c] dzdx[%f %f] opt[%d]\n", fDet, IsRowCross()?'y':'n', fZref[1], fZfit[1], opt);
   Double_t //xchmb = 0.5 * AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick(),
@@ -2210,13 +2206,11 @@ Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t
   Double_t xc[kNclusters], yc[kNclusters], dz(0.), dzdx(0.), 
            s2dz(0.), s2dzdx(0.), sy[kNclusters],
            s2x((8.33e-2/freq/freq+1.56e-2)*fVD*fVD),  // error of 1tb + error of mean time (TRF)
-           t2(fPad[2]*fPad[2]),
-           cs(0.);
+           t2(fPad[2]*fPad[2]), loc[3]={0.};
   Int_t n(0),          // clusters used in fit 
-        row[]={-1, -1},// pad row spanned by the tracklet
-        col(-1);       // pad column of current cluster
-  Double_t ycorr(UnbiasY(IsRowCross(), sgn, chg)),
-           kS2Ycorr(recoParam->GetS2Ycorr(sgn));
+        row[]={-1, -1};// pad row spanned by the tracklet
+  Double_t ycorr(UnbiasY(IsRowCross(), bz)),
+           kS2Ycorr(recoParam->GetS2Ycorr(IsRowCross(), chg>0));
         
   AliTRDcluster *c(NULL), **jc = &fClusters[0];
   for(Int_t ic=0; ic<kNtb; ic++, ++jc) {
@@ -2245,14 +2239,12 @@ Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t
           break;
       }    
     }
-    if(col != c->GetPadCol()){
-      col = c->GetPadCol();
-      cs  = pp->GetColSize(col);
-    }
-    //Use local cluster coordinates - the code should be identical with AliTRDtransform::Transform() !!!
-    //A.Bercuci 27.11.13
-    xc[n] = c->GetXloc(fT0, fVD); // c->GetX();
-    yc[n] = c->GetYloc(pp->GetColPos(col) + .5*cs, fS2PRF, cs) - xc[n]*fExB; //c->GetY();
+    //Use local cluster coordinates 
+    //A.Bercuci 27.11.13/30.06.14
+    Double_t trk[] = {c->GetX(), c->GetY(), c->GetZ()};
+    mDet->MasterToLocal(trk, loc);
+    xc[n] = AliTRDgeometry::AnodePos()-loc[0]; //c->GetXloc(fT0, fVD); // c->GetX();
+    yc[n] = loc[1]; //c->GetYloc(pp->GetColPos(col) + .5*cs, fS2PRF, cs) - xc[n]*fExB; //c->GetY();
     yc[n]-= fPad[2]*(dz+xc[n]*dzdx);
     yc[n]-= ycorr;
     if(IsRowCross()){ // estimate closest distance to anode wire
@@ -2293,14 +2285,13 @@ Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t
           break;
       }    
     }  
-    if(col != c->GetPadCol()){
-      col = c->GetPadCol();
-      cs  = pp->GetColSize(col);
-    }
+
     //Use local cluster coordinates - the code should be identical with AliTRDtransform::Transform() !!!
     //A.Bercuci 27.11.13
-    xc[n] = c->GetXloc(fT0, fVD); // c->GetX();
-    yc[n] = c->GetYloc(pp->GetColPos(col) + .5*cs, fS2PRF, cs) - xc[n]*fExB ;
+    Double_t trk[] = {c->GetX(), c->GetY(), c->GetZ()};
+    mDet->MasterToLocal(trk, loc);
+    xc[n] = AliTRDgeometry::AnodePos()-loc[0]; //c->GetXloc(fT0, fVD); // c->GetX();
+    yc[n] = loc[1]; //c->GetYloc(pp->GetColPos(col) + .5*cs, fS2PRF, cs) - xc[n]*fExB; //c->GetY();
     yc[n]-= fPad[2]*(dz+xc[n]*dzdx);
     yc[n]-= ycorr;
 
index 80bc782..a8dd3b4 100644 (file)
@@ -91,8 +91,8 @@ public:
   void      CookLabels();
   Bool_t    CookPID();
   Bool_t    Fit(UChar_t opt=0); // OBSOLETE
-  Bool_t    FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t opt=0);
-  Double_t  EstimatedCrossPoint(AliTRDpadPlane *pp);
+  Bool_t    FitRobust(AliTRDpadPlane *pp, TGeoHMatrix *mdet, Float_t bz, Int_t chg, Int_t opt=0);
+  Double_t  EstimatedCrossPoint(AliTRDpadPlane *pp, Float_t bz);
   Bool_t    Init(const AliTRDtrackV1 *track);
   void      Init(const AliRieman *fit);
   Bool_t    IsEqual(const TObject *inTracklet) const;
@@ -216,8 +216,8 @@ public:
 
 protected:
   void      Copy(TObject &ref) const;
-  void      UnbiasDZDX(Bool_t rc);
-  Double_t  UnbiasY(Bool_t rc, Bool_t sgn, Int_t chg);
+  void      UnbiasDZDX(Bool_t rc, Float_t bz);
+  Double_t  UnbiasY(Bool_t rc, Float_t bz);
 
 private:
   inline void SetN(Int_t n);
@@ -312,11 +312,12 @@ Double_t AliTRDseedV1::GetS2XcrossDZDX(Double_t absdzdx) const
 //____________________________________________________________
 Double_t AliTRDseedV1::GetS2DZDX(Float_t dzdx) const
 {
-//   Double_t p0[] = {0.03925, 0.03178},
-//            p1[] = {0.06316, 0.06669};
+  // Error parametrization for dzdx.
+  // TODO Should be layer dependent 
+  
   Double_t p0[] = {0.02835, 0.03925},
            p1[] = {0.04746, 0.06316};
-  
+           
   Double_t s2(p0[IsRowCross()]+p1[IsRowCross()]*dzdx*dzdx);
   s2*=s2;
   return s2;
index b481899..cca2e27 100644 (file)
@@ -956,7 +956,7 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
     // tilt correction options
     // 0 : no correction
     // 2 : pseudo tilt correction
-    if(!ptrTracklet->FitRobust(fGeom->GetPadPlane(ily, stk), prod>0., t.Charge())){
+    if(!ptrTracklet->FitRobust(fGeom->GetPadPlane(ily, stk), matrix, t.GetBz(), t.Charge())){
       t.SetErrStat(AliTRDtrackV1::kNoFit, ily);
       AliDebug(4, "Failed Tracklet Fit");
       continue;