]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCcalibAlign.cxx
Updates in the ITSsa track QA plotting (Leonardo)
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibAlign.cxx
index 27cd6fdaaec739502ca193dff580ad913901a965..8a95ad320e6fbf77883c5bdcb089a1a8d5c24d9f 100644 (file)
 #include "TH1D.h"
 #include "TH2F.h"
 #include "THnSparse.h"
+#include "THn.h"
 #include "TVectorD.h"
 #include "TTreeStream.h"
 #include "TFile.h"
 #include "AliTPCseed.h"
 #include "AliTracker.h"
 #include "TClonesArray.h"
-#include "AliExternalComparison.h"
 #include "AliLog.h"
 #include "TFile.h"
 #include "TProfile.h"
 
 #include "TTreeStream.h"
 #include "Riostream.h"
+#include "TRandom.h"
 #include <sstream>
 using namespace std;
 
 AliTPCcalibAlign* AliTPCcalibAlign::fgInstance = 0;
+Double_t          AliTPCcalibAlign::fgkMergeEntriesCut=10000000.; //10**7 tracks
 ClassImp(AliTPCcalibAlign)
 
 
@@ -449,6 +451,7 @@ void AliTPCcalibAlign::Process(AliESDEvent *event) {
   //
   // Process pairs of cosmic tracks
   //
+  const Double_t kptDownscale=50; // downscale factor for the low pt particels
   if (!fClusterDelta[0])  MakeResidualHistos();
   if (!fTrackletDelta[0])  MakeResidualHistosTracklet();
   //
@@ -458,6 +461,7 @@ void AliTPCcalibAlign::Process(AliESDEvent *event) {
   const Int_t kminCl = 40;
   AliESDfriend *eESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
   if (!eESDfriend) return;
+  if (eESDfriend->TestSkipBit()) return;
   Int_t ntracks=event->GetNumberOfTracks(); 
   Float_t dca0[2];
   Float_t dca1[2];
@@ -481,7 +485,9 @@ void AliTPCcalibAlign::Process(AliESDEvent *event) {
     fCurrentFriendTrack=friendTrack;
     fCurrentSeed=seed0;
     fCurrentEvent=event;
-    ProcessSeed(seed0);
+    Double_t scalept= TMath::Min(1./TMath::Abs(track0->GetParameter()[4]),2.);
+    Bool_t   isSelected = (TMath::Exp(2*scalept)>kptDownscale*gRandom->Rndm());
+    if (isSelected) ProcessSeed(seed0);
   }
   //
   // process cosmic pairs
@@ -859,6 +865,7 @@ void AliTPCcalibAlign::ProcessTracklets(const AliExternalTrackParam &tp1,
   //
   // Process function to fill fitters
   //
+  if (!seed) return;
   Double_t t1[10],t2[10];
   Double_t &x1=t1[0], &y1=t1[1], &z1=t1[3], &dydx1=t1[2], &dzdx1=t1[4];
   Double_t &x2=t2[0], &y2=t2[1], &z2=t2[3], &dydx2=t2[2], &dzdx2=t2[4];
@@ -913,7 +920,7 @@ void AliTPCcalibAlign::ProcessTracklets(const AliExternalTrackParam &tp1,
   Int_t acceptLinear =   AcceptTracklet(parLine1,parLine2);
 
 
-  if (fStreamLevel>1 && seed){
+  if (fStreamLevel>1){
     TTreeSRedirector *cstream = GetDebugStreamer();
     if (cstream){
       static TVectorD vec1(5);
@@ -954,10 +961,12 @@ void AliTPCcalibAlign::ProcessTracklets(const AliExternalTrackParam &tp1,
     // use Linear fit
     //
     if (nl1>10 && nl2>10 &&(acceptLinear==0)){
-      if (seed) ProcessDiff(tp1,tp2, seed,s1,s2);
+      ProcessDiff(tp1,tp2, seed,s1,s2);
       if (TMath::Abs(parLine1[2])<0.8 &&TMath::Abs(parLine1[2])<0.8 ){ //angular cut
        FillHisto(parLine1,parLine2,s1,s2);  
        ProcessAlign(parLine1,parLine2,s1,s2);
+       FillHisto((AliExternalTrackParam*)&tp1,(AliExternalTrackParam*)&tp2,s1,s2);
+       FillHisto((AliExternalTrackParam*)&tp2,(AliExternalTrackParam*)&tp1,s2,s1);
        //UpdateKalman(s1,s2,par1, cov1, par2, cov2); - OBSOLETE to be removed - 50 % of time here
       }
     }
@@ -969,11 +978,9 @@ void AliTPCcalibAlign::ProcessTracklets(const AliExternalTrackParam &tp1,
     //
     // use Kalman if mag field
     //
-    if (seed) {
-      ProcessDiff(tp1,tp2, seed,s1,s2);
-      FillHisto((AliExternalTrackParam*)&tp1,(AliExternalTrackParam*)&tp2,s1,s2);
-      FillHisto((AliExternalTrackParam*)&tp2,(AliExternalTrackParam*)&tp1,s2,s1);
-    }
+    ProcessDiff(tp1,tp2, seed,s1,s2);
+    FillHisto((AliExternalTrackParam*)&tp1,(AliExternalTrackParam*)&tp2,s1,s2);
+    FillHisto((AliExternalTrackParam*)&tp2,(AliExternalTrackParam*)&tp1,s2,s1);
     FillHisto(t1,t2,s1,s2);  
     ProcessAlign(t1,t2,s1,s2);
   }
@@ -991,88 +998,6 @@ void AliTPCcalibAlign::ProcessAlign(Double_t * t1,
   ++fPoints[GetIndex(s1,s2)];
 }
 
-void AliTPCcalibAlign::ProcessTree(TTree * chainTracklet, AliExternalComparison *comp){
-  //
-  // Process the debug streamer tree
-  // Possible to modify selection criteria
-  // Used with entry list
-  //
-  TTreeSRedirector * cstream = new TTreeSRedirector("aligndump.root");
-
-  AliTPCcalibAlign *align = this;
-  //
-  TVectorD * vec1 = 0;
-  TVectorD * vec2 = 0;
-  AliExternalTrackParam * tp1 = 0;
-  AliExternalTrackParam * tp2 = 0;  
-  Int_t      s1 = 0;
-  Int_t      s2 = 0;                           
-  Int_t npoints =0;
-  {
-    Int_t entries=chainTracklet->GetEntries();
-    for (Int_t i=0; i< entries; i++){
-      chainTracklet->GetBranch("tp1.")->SetAddress(&tp1);
-      chainTracklet->GetBranch("tp2.")->SetAddress(&tp2);
-      chainTracklet->GetBranch("v1.")->SetAddress(&vec1);
-      chainTracklet->GetBranch("v2.")->SetAddress(&vec2);
-      chainTracklet->GetBranch("s1")->SetAddress(&s1);
-      chainTracklet->GetBranch("s2")->SetAddress(&s2);      
-      chainTracklet->GetEntry(i);
-      if (!vec1) continue;
-      if (!vec2) continue;
-      if (!tp1) continue;
-      if (!tp2) continue;
-      if (!vec1->GetMatrixArray()) continue;
-      if (!vec2->GetMatrixArray()) continue;
-      // make a local copy
-      AliExternalTrackParam par1(*tp1);
-      AliExternalTrackParam par2(*tp2);
-      TVectorD svec1(*vec1);
-      TVectorD svec2(*vec2);
-      //
-      if (s1==s2) continue;
-      if (i%100==0) printf("%d\t%d\t%d\t%d\t\n",i, npoints,s1,s2);
-      AliExternalTrackParam  cpar1(par1);
-      AliExternalTrackParam  cpar2(par2);      
-      Constrain1Pt(cpar1,par2,fNoField);
-      Constrain1Pt(cpar2,par1,fNoField);
-      Bool_t acceptComp = kFALSE;
-      if (comp) acceptComp=comp->AcceptPair(&par1,&par2);
-      if (comp) acceptComp&=comp->AcceptPair(&cpar1,&cpar2);
-      //
-      Int_t reject =   align->AcceptTracklet(par1,par2);
-      Int_t rejectC =align->AcceptTracklet(cpar1,cpar2); 
-
-      if (1||fStreamLevel>0){
-       (*cstream)<<"Tracklet"<<
-         "s1="<<s1<<
-         "s2="<<s2<<
-         "reject="<<reject<<
-         "rejectC="<<rejectC<<
-         "acceptComp="<<acceptComp<<
-         "tp1.="<<&par1<<
-         "tp2.="<<&par2<<      
-         "ctp1.="<<&cpar1<<
-         "ctp2.="<<&cpar2<<
-         "v1.="<<&svec1<<
-         "v2.="<<&svec2<<
-         "\n";
-      }
-      //
-      if (fNoField){
-       //
-       //
-      }
-      if (acceptComp) comp->Process(&cpar1,&cpar2);
-      //
-      if (reject>0 || rejectC>0) continue;
-      npoints++;
-      align->ProcessTracklets(cpar1,cpar2,0,s1,s2);
-      align->ProcessTracklets(cpar2,cpar1,0,s2,s1); 
-    }
-  }
-  delete cstream;
-}
 
 
 Int_t AliTPCcalibAlign::AcceptTracklet(const AliExternalTrackParam &p1,
@@ -1671,7 +1596,11 @@ void AliTPCcalibAlign::MakeResidualHistos(){
   // 4 - local   kz
   // 
   axisName[0]="delta";   axisTitle[0]="#Delta (cm)"; 
-  binsTrack[0]=60;       xminTrack[0]=-0.6;        xmaxTrack[0]=0.6; 
+  if (TMath::Abs(AliTracker::GetBz())<0.01){
+    binsTrack[0]=60;       xminTrack[0]=-1.2;        xmaxTrack[0]=1.2; 
+  }else{
+    binsTrack[0]=60;       xminTrack[0]=-0.6;        xmaxTrack[0]=0.6; 
+  }
   //
   axisName[1]="sector";   axisTitle[1]="Sector Number"; 
   binsTrack[1]=180;       xminTrack[1]=0;        xmaxTrack[1]=18; 
@@ -1683,8 +1612,8 @@ void AliTPCcalibAlign::MakeResidualHistos(){
   axisName[3]="kZ";      axisTitle[3]="dz/dx"; 
   binsTrack[3]=36;       xminTrack[3]=-1.8;        xmaxTrack[3]=1.8; 
   //
-  fClusterDelta[0] = new THnSparseS("#Delta_{Y} (cm)","#Delta_{Y} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
-  fClusterDelta[1] = new THnSparseS("#Delta_{Z} (cm)","#Delta_{Z} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
+  fClusterDelta[0] = new THnF("#Delta_{Y} (cm)","#Delta_{Y} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
+  fClusterDelta[1] = new THnF("#Delta_{Z} (cm)","#Delta_{Z} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
   //
   //
   //
@@ -1712,13 +1641,14 @@ void AliTPCcalibAlign::MakeResidualHistosTracklet(){
   // 3 - local   ky
   // 4 - local   kz
   // 5 - sector  1
-  // 5 - sector  0
+  // 6 - sector  0
+  // 7 - z position  0
 
   axisName[0]="delta";   axisTitle[0]="#Delta (cm)"; 
-  binsTrack[0]=60;       xminTrack[0]=-0.6;        xmaxTrack[0]=0.6
+  binsTrack[0]=60;       xminTrack[0]=-0.5;        xmaxTrack[0]=0.5
   //
   axisName[1]="phi";   axisTitle[1]="#phi"; 
-  binsTrack[1]=180;       xminTrack[1]=-TMath::Pi();        xmaxTrack[1]=TMath::Pi(); 
+  binsTrack[1]=90;       xminTrack[1]=-TMath::Pi();        xmaxTrack[1]=TMath::Pi(); 
   //
   axisName[2]="localX";   axisTitle[2]="x (cm)"; 
   binsTrack[2]=10;       xminTrack[2]=120.;        xmaxTrack[2]=200.; 
@@ -1727,28 +1657,34 @@ void AliTPCcalibAlign::MakeResidualHistosTracklet(){
   binsTrack[3]=10;       xminTrack[3]=-0.5;        xmaxTrack[3]=0.5; 
   //
   axisName[4]="kZ";      axisTitle[4]="dz/dx"; 
-  binsTrack[4]=22;       xminTrack[4]=-1.1;        xmaxTrack[4]=1.1; 
+  binsTrack[4]=11;       xminTrack[4]=-1.1;        xmaxTrack[4]=1.1; 
   //
   axisName[5]="is1";      axisTitle[5]="is1"; 
   binsTrack[5]=72;       xminTrack[5]=0;        xmaxTrack[5]=72;
   //
   axisName[6]="is0";      axisTitle[6]="is0"; 
   binsTrack[6]=72;       xminTrack[6]=0;        xmaxTrack[6]=72;
+  //
+  axisName[7]="z";        axisTitle[7]="z(cm)"; 
+  binsTrack[7]=12;        xminTrack[7]=-240;        xmaxTrack[7]=240;
+  //
+  axisName[8]="IsPrimary";        axisTitle[8]="Is Primary"; 
+  binsTrack[8]=2;        xminTrack[8]=-0.1;        xmaxTrack[8]=1.1;
  
   //
-  xminTrack[0]=-0.3;        xmaxTrack[0]=0.3
-  fTrackletDelta[0] = new THnSparseF("#Delta_{Y} (cm)","#Delta_{Y} (cm)", 7, binsTrack,xminTrack, xmaxTrack);
+  xminTrack[0]=-0.25;        xmaxTrack[0]=0.25
+  fTrackletDelta[0] = new THnSparseF("#Delta_{Y} (cm)","#Delta_{Y} (cm)", 9, binsTrack,xminTrack, xmaxTrack);
   xminTrack[0]=-0.5;        xmaxTrack[0]=0.5; 
-  fTrackletDelta[1] = new THnSparseF("#Delta_{Z} (cm)","#Delta_{Z} (cm)", 7, binsTrack,xminTrack, xmaxTrack);
-  xminTrack[0]=-0.005;        xmaxTrack[0]=0.005; 
-  fTrackletDelta[2] = new THnSparseF("#Delta_{kY}","#Delta_{kY}", 7, binsTrack,xminTrack, xmaxTrack);
+  fTrackletDelta[1] = new THnSparseF("#Delta_{Z} (cm)","#Delta_{Z} (cm)", 9, binsTrack,xminTrack, xmaxTrack);
   xminTrack[0]=-0.005;        xmaxTrack[0]=0.005; 
-  fTrackletDelta[3] = new THnSparseF("#Delta_{kZ}","#Delta_{kZ}", 7, binsTrack,xminTrack, xmaxTrack);
+  fTrackletDelta[2] = new THnSparseF("#Delta_{kY}","#Delta_{kY}", 9, binsTrack,xminTrack, xmaxTrack);
+  xminTrack[0]=-0.008;        xmaxTrack[0]=0.008; 
+  fTrackletDelta[3] = new THnSparseF("#Delta_{kZ}","#Delta_{kZ}", 9, binsTrack,xminTrack, xmaxTrack);
   //
   //
   //
   for (Int_t ivar=0;ivar<4;ivar++){
-    for (Int_t ivar2=0;ivar2<7;ivar2++){
+    for (Int_t ivar2=0;ivar2<9;ivar2++){
       fTrackletDelta[ivar]->GetAxis(ivar2)->SetName(axisName[ivar2].Data());
       fTrackletDelta[ivar]->GetAxis(ivar2)->SetTitle(axisName[ivar2].Data());
     }
@@ -1796,8 +1732,9 @@ void AliTPCcalibAlign::FillHisto(AliExternalTrackParam *tp1,
   //
   // Fill residual histograms
   // Track2-Track1
+  if (s2<s1) return;//
   const Double_t kEpsilon=0.001;
-  Double_t x[8]={0,0,0,0,0,0,0,0};
+  Double_t x[9]={0,0,0,0,0,0,0,0,0};
   AliExternalTrackParam p1(*tp1);
   AliExternalTrackParam p2(*tp2);
   if (s1%18==s2%18) {
@@ -1815,7 +1752,11 @@ void AliTPCcalibAlign::FillHisto(AliExternalTrackParam *tp1,
   x[4]=0.5*(p1.GetTgl()+p2.GetTgl());  // mean tgl
   x[5]=s2;
   x[6]=s1;
-  
+  x[7]=0.5*(p1.GetZ()+p2.GetZ());
+  // is primary ?
+  Int_t  isPrimary = (TMath::Abs(p1.GetTgl()-p1.GetZ()/p1.GetX())<0.1) ? 1:0;
+  x[8]= isPrimary;
+  //
   x[0]=p2.GetY()-p1.GetY();
   fTrackletDelta[0]->Fill(x);
   x[0]=p2.GetZ()-p1.GetZ();
@@ -1968,7 +1909,7 @@ TGraphErrors * AliTPCcalibAlign::MakeGraph(Int_t sec0, Int_t sec1, Int_t dsec,
   }
   TGraphErrors *gr = new TGraphErrors(npoints,xsec,ysec,0,0);
   Char_t name[1000];
-  sprintf(name,"Mat[%d,%d]  Type=%d",i0,i1,type);
+  snprintf(name,100,"Mat[%d,%d]  Type=%d",i0,i1,type);
   gr->SetName(name);
   return gr;
 }
@@ -2217,15 +2158,22 @@ void AliTPCcalibAlign::Add(AliTPCcalibAlign * align){
   if (!fClusterDelta[0]) MakeResidualHistos();
 
   for (Int_t i=0; i<2; i++){
-    if (align->fClusterDelta[i]) fClusterDelta[i]->Add(align->fClusterDelta[i]);
+    if (align->fClusterDelta[i]){
+      fClusterDelta[i]->Add(align->fClusterDelta[i]);
+    }
   }
 
+  
   for (Int_t i=0; i<4; i++){
     if (!fTrackletDelta[i] && align->fTrackletDelta[i]) {
       fTrackletDelta[i]= (THnSparse*)(align->fTrackletDelta[i]->Clone());
       continue;
     }
-    if (align->fTrackletDelta[i]) fTrackletDelta[i]->Add(align->fTrackletDelta[i]);
+    if (align->fTrackletDelta[i]) {
+      if (fTrackletDelta[i]->GetEntries()<fgkMergeEntriesCut){
+       fTrackletDelta[i]->Add(align->fTrackletDelta[i]);
+      }
+    }
   }
 
 }
@@ -2638,16 +2586,21 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
   //
   // 1. Apply selection
   // 2. Refit the track - in-out
-  //                    - update the cluster delta in upper part
   // 3. Refit the track - out-in
-  //                    - update the cluster delta histo lower part
+  // 4. Combine In and Out track - - fil cluster residuals
   //
+  if (!fCurrentFriendTrack) return;
+  if (!fCurrentFriendTrack->GetTPCOut()) return;
   const Double_t kPtCut=1.0;    // pt
   const Double_t kSnpCut=0.2; // snp cut
   const Double_t kNclCut=120; //
   const Double_t kVertexCut=1;
   const Double_t kMaxDist=0.5; // max distance between tracks and cluster
   const Double_t kEdgeCut = 2.5;
+  const Double_t kDelta2=0.2*0.2;  // initial increase in covar matrix
+  const Double_t kSigma=0.3;       // error increase towards edges of TPC 
+  const Double_t kSkipBoundary=7.5;  // skip track updates in the boundary IFC,OFC, IO
+  //
   if (!fCurrentTrack) return;
   if (!fCurrentFriendTrack) return;
   Float_t vertexXY=0,vertexZ=0;
@@ -2658,12 +2611,31 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
   if (seed->GetNumberOfClusters()<kNclCut) return;
   if (TMath::Abs(seed->GetSnp())>kSnpCut) return; 
   if (!fClusterDelta[0])  MakeResidualHistos();
-
+  //
+  AliExternalTrackParam fitIn[160];
+  AliExternalTrackParam fitOut[160];
+  AliTPCROC * roc = AliTPCROC::Instance();
+  Double_t xmiddle   = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+  Double_t xDiff     = ( -roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+  Double_t xIFC      = ( roc->GetPadRowRadii(0,0));
+  Double_t xOFC      = ( roc->GetPadRowRadii(36,roc->GetNRows(36)-1));
+  //
   Int_t detector=-1;
   //
   //
   AliExternalTrackParam trackIn  = *(fCurrentTrack->GetInnerParam());
   AliExternalTrackParam trackOut = *(fCurrentFriendTrack->GetTPCOut());
+  trackIn.ResetCovariance(10);
+  trackOut.ResetCovariance(10);
+  Double_t *covarIn = (Double_t*)trackIn.GetCovariance();
+  Double_t *covarOut = (Double_t*)trackOut.GetCovariance();
+  covarIn[0]+=kDelta2; covarIn[2]+=kDelta2; 
+  covarIn[5]+=kDelta2/(100.*100.); covarIn[9]=kDelta2/(100.*100.);
+  covarIn[14]+=kDelta2/(5.*5.);
+  covarOut[0]+=kDelta2; covarOut[2]+=kDelta2; 
+  covarOut[5]+=kDelta2/(100.*100.); covarOut[9]=kDelta2/(100.*100.);
+  covarOut[14]+=kDelta2/(5.*5.);
+  //
   static Double_t mass =    TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
   //  
   Int_t ncl=0;
@@ -2677,7 +2649,6 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
     ncl++;
   }
   if (ncl<kNclCut) return;
-
   Int_t nclIn=0,nclOut=0;
   Double_t xyz[3];
   //
@@ -2690,48 +2661,41 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
     if (detector<0) detector=cl->GetDetector()%36;
     Int_t sector = cl->GetDetector();
     Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha();    
+    if (cl->GetDetector()%36!=detector) 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}; 
-    AliTPCseed::GetError(cl, &trackOut,cov[0],cov[2]);
+    Double_t cov[3]={0.1,0.,0.1}; 
     Double_t dedge =  cl->GetX()*TMath::Tan(TMath::Pi()/18.)-TMath::Abs(trackOut.GetY());
-    cov[0]+=1./(irow+1.); // bigger error at boundary
-    cov[0]+=1./(160.-irow); // bigger error at boundary
-    cov[2]+=1./(irow+1.); // bigger error at boundary
-    cov[2]+=1./(160.-irow); // bigger error at boundary
-    cov[0]+=0.5/dedge;       // bigger error close to the boundary
-    cov[2]+=0.5/dedge;       // bigger error close to the boundary
+    Double_t dmiddle = TMath::Abs(cl->GetX()-xmiddle)/xDiff;
+    dmiddle*=dmiddle;
+    //
+    cov[0]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[0]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[2]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[2]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[0]+=kSigma/dedge;      // bigger error close to the boundary
+    cov[2]+=kSigma/dedge;      // bigger error close to the boundary
     cov[0]*=cov[0];
     cov[2]*=cov[2];
-    if (!AliTracker::PropagateTrackToBxByBz(&trackOut, r[0],mass,1.,kFALSE)) continue;
-    
+    if (!AliTracker::PropagateTrackToBxByBz(&trackOut, r[0],mass,1.,kFALSE)) continue;    
     if (TMath::Abs(dedge)<kEdgeCut) continue;
-
+    //
+    Bool_t doUpdate=kTRUE;
+    if (TMath::Abs(cl->GetX()-xIFC)<kSkipBoundary) doUpdate=kFALSE;
+    if (TMath::Abs(cl->GetX()-xOFC)<kSkipBoundary) doUpdate=kFALSE;
+    if (TMath::Abs(cl->GetX()-fXIO)<kSkipBoundary) doUpdate=kFALSE;
+    //
     if (TMath::Abs(cl->GetY()-trackOut.GetY())<kMaxDist){
       nclOut++;
-      trackOut.Update(&r[1],cov);
+      if (doUpdate) trackOut.Update(&r[1],cov);
     }
-    if (nclOut<kNclCut/2) continue;
-    if (cl->GetDetector()%36!=detector) continue;
-    //
-    // fill residual histogram
-    //
-    Double_t resVector[5]; 
-    trackOut.GetXYZ(xyz);
-    resVector[1]= 9.*TMath::ATan2(xyz[1],xyz[0])/TMath::Pi();
-    if (resVector[1]<0) resVector[1]+=18;
-    resVector[2]= TMath::Sqrt(cl->GetX()*cl->GetX()+cl->GetY()*cl->GetY());
-    resVector[3]= cl->GetZ()/cl->GetX();
-    //
-    resVector[0]= cl->GetY()-trackOut.GetY();
-    fClusterDelta[0]->Fill(resVector);
-    resVector[0]= cl->GetZ()-trackOut.GetZ();
-    fClusterDelta[1]->Fill(resVector);
+    fitOut[irow]=trackOut;
   }
+  
   //
-  // Refit in - store residual maps
+  // Refit In - store residual maps
   //
   for (Int_t irow=159; irow>=0; irow--){
     AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
@@ -2740,48 +2704,65 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
     if (detector<0) detector=cl->GetDetector()%36;
     Int_t sector = cl->GetDetector();
     Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();    
+    if (cl->GetDetector()%36!=detector) 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}; 
-    AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]);
+    Double_t cov[3]={0.1,0.,0.1}; 
     Double_t dedge =  cl->GetX()*TMath::Tan(TMath::Pi()/18.)-TMath::Abs(trackIn.GetY());
-    cov[0]+=1./(irow+1.); // bigger error at boundary
-    cov[0]+=1./(160.-irow); // bigger error at boundary
-    cov[2]+=1./(irow+1.); // bigger error at boundary
-    cov[2]+=1./(160.-irow); // bigger error at boundary
-    cov[0]+=0.5/dedge;       // bigger error close to the boundary +-
-    cov[2]+=0.5/dedge;       // bigger error close to the boundary +-
+    Double_t dmiddle = TMath::Abs(cl->GetX()-xmiddle)/xDiff;
+    dmiddle*=dmiddle;
+    //
+    cov[0]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[0]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[2]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[2]+=kSigma*dmiddle;     // bigger error at boundary
+    cov[0]+=kSigma/dedge;      // bigger error close to the boundary
+    cov[2]+=kSigma/dedge;      // bigger error close to the boundary
     cov[0]*=cov[0];
     cov[2]*=cov[2];
-    if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1.,kFALSE)) continue;
+    if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1.,kFALSE)) continue;    
     if (TMath::Abs(dedge)<kEdgeCut) continue;
-
-
+    Bool_t doUpdate=kTRUE;
+    if (TMath::Abs(cl->GetX()-xIFC)<kSkipBoundary) doUpdate=kFALSE;
+    if (TMath::Abs(cl->GetX()-xOFC)<kSkipBoundary) doUpdate=kFALSE;
+    if (TMath::Abs(cl->GetX()-fXIO)<kSkipBoundary) doUpdate=kFALSE;
     if (TMath::Abs(cl->GetY()-trackIn.GetY())<kMaxDist){
       nclIn++;
-      trackIn.Update(&r[1],cov);
+      if (doUpdate) trackIn.Update(&r[1],cov);
     }
-    if (nclIn<kNclCut/2) continue;
-    if (cl->GetDetector()%36!=detector) continue;
+    fitIn[irow]=trackIn;
+  }
+  //
+  //
+  for (Int_t irow=159; irow>=0; irow--){
     //
-    // fill residual histogram
+    // Update kalman - +- direction
+    // Store cluster residuals
+    AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+    if (!cl) continue;
+    if (cl->GetX()<80) continue;
+    if (detector<0) detector=cl->GetDetector()%36;
+    if (cl->GetDetector()%36!=detector) continue;
+    if (fitIn[irow].GetX()<80) continue;
+    if (fitOut[irow].GetX()<80) continue;
+    AliExternalTrackParam trackSmooth = fitIn[irow];
+    AliTrackerBase::UpdateTrack(trackSmooth, fitOut[irow]);
     //
     Double_t resVector[5]; 
-    trackIn.GetXYZ(xyz);
+    trackSmooth.GetXYZ(xyz);
     resVector[1]= 9.*TMath::ATan2(xyz[1],xyz[0])/TMath::Pi();
     if (resVector[1]<0) resVector[1]+=18;
     resVector[2]= TMath::Sqrt(cl->GetX()*cl->GetX()+cl->GetY()*cl->GetY());
-    resVector[3]= cl->GetZ()/cl->GetX();
+    resVector[3]= cl->GetZ()/resVector[2];
     //
-    resVector[0]= cl->GetY()-trackIn.GetY();
+    resVector[0]= cl->GetY()-trackSmooth.GetY();
     fClusterDelta[0]->Fill(resVector);
-    resVector[0]= cl->GetZ()-trackIn.GetZ();
+    resVector[0]= cl->GetZ()-trackSmooth.GetZ();
     fClusterDelta[1]->Fill(resVector);
   }
 
-
 }
 
 
@@ -2792,7 +2773,7 @@ void  AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){
   //
   if (TMath::Abs(AliTracker::GetBz())>0.5) return;
   if (!fClusterDelta[0])  MakeResidualHistos();
-  const Int_t kMinClusterF=40;
+  //  const Int_t kMinClusterF=40;
   const Int_t kMinClusterFit=10;
   const Int_t kMinClusterQ=10;
   //
@@ -2842,12 +2823,12 @@ void  AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){
       }
       if (TMath::Abs(x[0])<10){
        fyf.AddPoint(x,c->GetY(),0.1); //use only middle rows+-10cm
+       fzf.AddPoint(x,c->GetZ(),0.1);            
       }
-      fzf.AddPoint(x,c->GetZ(),0.1);      
     }
     nf = fyf.GetNpoints();
     if (fyf.GetNpoints()<kMinClusterFit) return;   // not enough points - skip 
-    if (fzf.GetNpoints()<kMinClusterF) return;   // not enough points - skip 
+    if (fzf.GetNpoints()<kMinClusterFit) return;   // not enough points - skip 
     fyf.Eval(); 
     fyf.GetParameters(pyf); 
     fyf.GetErrors(peyf);
@@ -2964,7 +2945,7 @@ void  AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){
       Double_t resVector[5];
       resVector[1]= 9.*gphi/TMath::Pi();
       resVector[2]= TMath::Sqrt(c->GetX()*c->GetX()+c->GetY()*c->GetY());
-      resVector[3]= c->GetZ()/c->GetX();
+      resVector[3]= c->GetZ()/resVector[2];
       //
       //
       resVector[0]= c->GetY()-yfitC;
@@ -3817,3 +3798,4 @@ void AliTPCcalibAlign::MakeReportDyPhi(TFile */*output*/){
   }
 }
 
+