Minor updates in macros to check the SDD calibration files (F. Prino)

[u/mrichter/AliRoot.git] / ITS / AliITStrackV2.cxx

diff --git a/ITS/AliITStrackV2.cxx b/ITS/AliITStrackV2.cxx
index 4cd0a0da069b978f7477e75ced6d6c244ff07050..09d590df68c2e027ebadb0be172d581526b2c5d5 100644 (file)
--- a/ITS/AliITStrackV2.cxx
+++ b/ITS/AliITStrackV2.cxx
@@ -87,13 +87,19 @@ void AliITStrackV2::ResetClusters() {
   SetNumberOfClusters(0);
 } 
 
+//____________________________________________________________________________
 void AliITStrackV2::UpdateESDtrack(ULong_t flags) const {
+  // Update track params
   fESDtrack->UpdateTrackParams(this,flags);
   // copy the module indices
   for(Int_t i=0;i<12;i++) {
     //   printf("     %d\n",GetModuleIndex(i));
     fESDtrack->SetITSModuleIndex(i,GetModuleIndex(i));
   }
+  // copy the 4 dedx samples
+  Double_t sdedx[4]={0.,0.,0.,0.};
+  for(Int_t i=0; i<4; i++) sdedx[i]=fdEdxSample[i];
+  fESDtrack->SetITSdEdxSamples(sdedx);
 }
 
 //____________________________________________________________________________
@@ -238,7 +244,7 @@ Bool_t AliITStrackV2::Update(const AliCluster* c, Double_t chi2, Int_t index)
   //This function updates track parameters
   //------------------------------------------------------------------
   Double_t p[2]={c->GetY(), c->GetZ()};
-  Double_t cov[3]={c->GetSigmaY2(), 0., c->GetSigmaZ2()};
+  Double_t cov[3]={c->GetSigmaY2(), c->GetSigmaYZ(), c->GetSigmaZ2()};
 
   if (!AliExternalTrackParam::Update(p,cov)) return kFALSE;
 
@@ -273,8 +279,8 @@ Bool_t AliITStrackV2::Invariant() const {
   // take into account the misalignment error
   Float_t maxMisalErrY2=0,maxMisalErrZ2=0;
   for (Int_t lay=0; lay<AliITSgeomTGeo::kNLayers; lay++) {
-    maxMisalErrY2 = TMath::Max(maxMisalErrY2,AliITSReconstructor::GetRecoParam()->GetClusterMisalErrorY(lay));
-    maxMisalErrZ2 = TMath::Max(maxMisalErrZ2,AliITSReconstructor::GetRecoParam()->GetClusterMisalErrorZ(lay));
+    maxMisalErrY2 = TMath::Max(maxMisalErrY2,AliITSReconstructor::GetRecoParam()->GetClusterMisalErrorY(lay,GetBz()));
+    maxMisalErrZ2 = TMath::Max(maxMisalErrZ2,AliITSReconstructor::GetRecoParam()->GetClusterMisalErrorZ(lay,GetBz()));
   }
   maxMisalErrY2 *= maxMisalErrY2;
   maxMisalErrZ2 *= maxMisalErrZ2;
@@ -424,36 +430,47 @@ void AliITStrackV2::CookdEdx(Double_t low, Double_t up) {
   //-----------------------------------------------------------------
   // This function calculates dE/dX within the "low" and "up" cuts.
   // Origin: Boris Batyunya, JINR, Boris.Batiounia@cern.ch 
+  // Updated: F. Prino 8-June-2009
   //-----------------------------------------------------------------
-  // The clusters order is: SSD-2, SSD-1, SDD-2, SDD-1, SPD-2, SPD-1
+  // The cluster order is: SDD-1, SDD-2, SSD-1, SSD-2
 
-  Int_t i;
   Int_t nc=0;
-  for (i=0; i<GetNumberOfClusters(); i++) {
-    Int_t idx=GetClusterIndex(i);
-    idx=(idx&0xf0000000)>>28;
-    if (idx>1) nc++; // Take only SSD and SDD
+  Float_t dedx[4];
+  for (Int_t il=0; il<4; il++) { // count good (>0) dE/dx values
+    if(fdEdxSample[il]>0.){
+      dedx[nc]= fdEdxSample[il];
+      nc++;
+    }
+  }
+  if(nc<1){
+    SetdEdx(0.);
+    return;
   }
 
-  Int_t swap;//stupid sorting
+  Int_t swap; // sort in ascending order
   do {
     swap=0;
-    for (i=0; i<nc-1; i++) {
-      if (fdEdxSample[i]<=fdEdxSample[i+1]) continue;
-      Float_t tmp=fdEdxSample[i];
-      fdEdxSample[i]=fdEdxSample[i+1]; fdEdxSample[i+1]=tmp;
+    for (Int_t i=0; i<nc-1; i++) {
+      if (dedx[i]<=dedx[i+1]) continue;
+      Float_t tmp=dedx[i];
+      dedx[i]=dedx[i+1]; 
+      dedx[i+1]=tmp;
       swap++;
     }
   } while (swap);
 
-  Int_t nl=Int_t(low*nc), nu=Int_t(up*nc); //b.b. to take two lowest dEdX
-                                           // values from four ones choose
-                                           // nu=2
-  Float_t dedx=0;
-  for (i=nl; i<nu; i++) dedx += fdEdxSample[i];
-  if (nu-nl>0) dedx /= (nu-nl);
 
-  SetdEdx(dedx);
+  Double_t nl=low*nc, nu =up*nc;
+  Float_t sumamp = 0;
+  Float_t sumweight =0;
+  for (Int_t i=0; i<nc; i++) {
+    Float_t weight =1;
+    if (i<nl+0.1)     weight = TMath::Max(1.-(nl-i),0.);
+    if (i>nu-1)     weight = TMath::Max(nu-i,0.);
+    sumamp+= dedx[i]*weight;
+    sumweight+=weight;
+  }
+  SetdEdx(sumamp/sumweight);
 }
 
 //____________________________________________________________________________
@@ -465,13 +482,25 @@ GetPhiZat(Double_t r, Double_t &phi, Double_t &z) const {
   // The track curvature is neglected.
   //------------------------------------------------------------------
   Double_t d=GetD(0.,0.);
-  if (TMath::Abs(d) > r) return kFALSE; 
+  if (TMath::Abs(d) > r) {
+    if (r>1e-1) return kFALSE;
+    r = TMath::Abs(d);
+  }
 
   Double_t rcurr=TMath::Sqrt(GetX()*GetX() + GetY()*GetY());
-  if (TMath::Abs(d) > rcurr) return kFALSE; 
-  Double_t phicurr=GetAlpha()+TMath::ASin(GetSnp());
+  if (TMath::Abs(d) > rcurr) return kFALSE;
+  Double_t globXYZcurr[3]; GetXYZ(globXYZcurr); 
+  Double_t phicurr=TMath::ATan2(globXYZcurr[1],globXYZcurr[0]);
+
+  if (GetX()>=0.) {
+    phi=phicurr+TMath::ASin(d/r)-TMath::ASin(d/rcurr);
+  } else {
+    phi=phicurr+TMath::ASin(d/r)+TMath::ASin(d/rcurr)-TMath::Pi();
+  }
 
-  phi=phicurr+TMath::ASin(d/r)-TMath::ASin(d/rcurr);
+  // return a phi in [0,2pi[ 
+  if (phi<0.) phi+=2.*TMath::Pi();
+  else if (phi>=2.*TMath::Pi()) phi-=2.*TMath::Pi();
   z=GetZ()+GetTgl()*(TMath::Sqrt((r-d)*(r+d))-TMath::Sqrt((rcurr-d)*(rcurr+d)));
   return kTRUE;
 }
@@ -484,11 +513,20 @@ GetLocalXat(Double_t r,Double_t &xloc) const {
   // The track curvature is neglected.
   //------------------------------------------------------------------
   Double_t d=GetD(0.,0.);
-  if (TMath::Abs(d) > r) return kFALSE; 
+  if (TMath::Abs(d) > r) { 
+    if (r>1e-1) return kFALSE; 
+    r = TMath::Abs(d); 
+  } 
 
   Double_t rcurr=TMath::Sqrt(GetX()*GetX() + GetY()*GetY());
-  Double_t phicurr=GetAlpha()+TMath::ASin(GetSnp());
-  Double_t phi=phicurr+TMath::ASin(d/r)-TMath::ASin(d/rcurr);
+  Double_t globXYZcurr[3]; GetXYZ(globXYZcurr); 
+  Double_t phicurr=TMath::ATan2(globXYZcurr[1],globXYZcurr[0]);
+  Double_t phi;
+  if (GetX()>=0.) {
+    phi=phicurr+TMath::ASin(d/r)-TMath::ASin(d/rcurr);
+  } else {
+    phi=phicurr+TMath::ASin(d/r)+TMath::ASin(d/rcurr)-TMath::Pi();
+  }
 
   xloc=r*(TMath::Cos(phi)*TMath::Cos(GetAlpha())
          +TMath::Sin(phi)*TMath::Sin(GetAlpha()));