more secure string operations
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibTime.cxx
index e06510d..e25ca39 100644 (file)
@@ -96,12 +96,15 @@ Comments to be written here:
 #include "AliTPCseed.h"
 #include "AliTrackPointArray.h"
 #include "AliTracker.h"
+#include "AliKFVertex.h"
+#include <AliLog.h>
 
 ClassImp(AliTPCcalibTime)
 
 
 AliTPCcalibTime::AliTPCcalibTime() 
-  :AliTPCcalibBase(), 
+  :AliTPCcalibBase(),  
+   fMemoryMode(1), // 0 -do not fill THnSparse with residuals  1- fill only important QA THn 2 - Fill all THnsparse for calibration
    fLaser(0),       // pointer to laser calibration
    fDz(0),          // current delta z
    fCutMaxD(3),        // maximal distance in rfi ditection
@@ -109,10 +112,13 @@ AliTPCcalibTime::AliTPCcalibTime()
    fCutTheta(0.03),    // maximal distan theta
    fCutMinDir(-0.99),  // direction vector products
    fCutTracks(100),
+   fArrayLaserA(0),      //laser  fit parameters C
+   fArrayLaserC(0),      //laser  fit parameters A
    fArrayDz(0),          //NEW! Tmap of V drifts for different triggers
    fAlignITSTPC(0),      //alignemnt array ITS TPC match
    fAlignTRDTPC(0),      //alignemnt array TRD TPC match 
    fAlignTOFTPC(0),      //alignemnt array TOF TPC match
+   fTimeKalmanBin(60*15), //time bin width for kalman - 15 minutes default
    fTimeBins(0),
    fTimeStart(0),
    fTimeEnd(0),
@@ -144,12 +150,27 @@ AliTPCcalibTime::AliTPCcalibTime()
     fResHistoTPCTRD[i]=0;
     fResHistoTPCTOF[i]=0;
     fResHistoTPCvertex[i]=0;
+    fTPCVertex[i]=0;
+  }
+  for (Int_t i=0;i<12;i++) {
+    fTPCVertex[i]=0;
+  }
+  for (Int_t i=0;i<5;i++) {
+    fTPCVertexCorrelation[i]=0;
+  }
+  static Int_t counter=0;
+  if (1) {
+    TTimeStamp s;
+    Int_t time=s;
+    AliInfo(Form("Counter Constructor\t%d\t%d",counter,time));
+    counter++;
   }
 
 }
 
-AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t StartTime, UInt_t EndTime, Int_t deltaIntegrationTimeVdrift)
+AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t StartTime, UInt_t EndTime, Int_t deltaIntegrationTimeVdrift, Int_t memoryMode)
   :AliTPCcalibBase(),
+   fMemoryMode(memoryMode), // 0 -do not fill THnSparse with residuals  1- fill only important QA THn 2 - Fill all THnsparse for calibration
    fLaser(0),            // pointer to laser calibration
    fDz(0),               // current delta z
    fCutMaxD(5*0.5356),   // maximal distance in rfi ditection
@@ -157,10 +178,13 @@ AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t
    fCutTheta(5*0.004644),// maximal distan theta
    fCutMinDir(-0.99),    // direction vector products
    fCutTracks(100),
+   fArrayLaserA(new TObjArray(1000)),      //laser  fit parameters C
+   fArrayLaserC(new TObjArray(1000)),      //laser  fit parameters A
    fArrayDz(0),            //Tmap of V drifts for different triggers
    fAlignITSTPC(0),      //alignemnt array ITS TPC match
    fAlignTRDTPC(0),      //alignemnt array TRD TPC match 
    fAlignTOFTPC(0),      //alignemnt array TOF TPC match
+   fTimeKalmanBin(60*15), //time bin width for kalman - 15 minutes default
    fTimeBins(0),
    fTimeStart(0),
    fTimeEnd(0),
@@ -266,13 +290,27 @@ AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t
   fCosmiMatchingHisto[7]=new TH1F("Cosmics matching","p2-isPair",100,-10*0.01134 ,10*0.01134 );
   fCosmiMatchingHisto[8]=new TH1F("Cosmics matching","p3-isPair",100,-10*0.004644,10*0.004644);
   fCosmiMatchingHisto[9]=new TH1F("Cosmics matching","p4-isPair",100,-10*0.03773 ,10*0.03773 );
+  for (Int_t i=0;i<12;i++) {
+    fTPCVertex[i]=0;
+  }
+  for (Int_t i=0;i<5;i++) {
+    fTPCVertexCorrelation[i]=0;
+  }
   BookDistortionMaps();
+  
 }
 
 AliTPCcalibTime::~AliTPCcalibTime(){
   //
   // Virtual Destructor
   //
+  static Int_t counter=0;
+  if (1) {
+    TTimeStamp s;
+    Int_t time=s;
+    AliInfo(Form("Counter Destructor\t%s\t%d\t%d",GetName(),counter,time));
+    counter++;
+  }
   for(Int_t i=0;i<3;i++){
     if(fHistVdriftLaserA[i]){
       delete fHistVdriftLaserA[i];
@@ -309,7 +347,13 @@ AliTPCcalibTime::~AliTPCcalibTime(){
     fResHistoTPCvertex[i]=0;
   }
 
-
+  if (fTPCVertex) {
+    for (Int_t i=0;i<12;i++)  delete fTPCVertex[i];
+  }
+  if (fTPCVertexCorrelation) {
+    for (Int_t i=0;i<5;i++)  delete fTPCVertexCorrelation[i];
+  }
+  
   fAlignITSTPC->SetOwner(kTRUE);
   fAlignTRDTPC->SetOwner(kTRUE);
   fAlignTOFTPC->SetOwner(kTRUE);
@@ -444,6 +488,20 @@ void AliTPCcalibTime::ProcessLaser(AliESDEvent *event){
   npointsC= TMath::Nint(vdriftC[3]);
   chi2C= vdriftC[4];
 
+  if (npointsA>kMinTracksSide || npointsC>kMinTracksSide){
+    TVectorD *fitA = new TVectorD(6);
+    TVectorD *fitC = new TVectorD(6);
+    for (Int_t ipar=0; ipar<5; ipar++){
+      (*fitA)[ipar]=vdriftA[ipar];
+      (*fitC)[ipar]=vdriftC[ipar];
+    }
+    (*fitA)[5]=fTime;
+    (*fitC)[5]=fTime;
+    fArrayLaserA->AddLast(fitA);
+    fArrayLaserC->AddLast(fitC);
+  }
+  //
+  
   TTimeStamp tstamp(fTime);
   Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
   Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
@@ -473,33 +531,6 @@ void AliTPCcalibTime::ProcessLaser(AliESDEvent *event){
     fHistVdriftLaserA[icalib]->Fill(vecDriftLaserA);
     fHistVdriftLaserC[icalib]->Fill(vecDriftLaserC);
   }
-
-//   THnSparse* curHist=new THnSparseF("","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
-//   TString shortName=curHist->ClassName();
-//   shortName+="_MEAN_DRIFT_LASER_";
-//   delete curHist;
-//   curHist=NULL;
-//   TString name="";
-
-//   name=shortName;
-//   name+=event->GetFiredTriggerClasses();
-//   name.ToUpper();
-//   curHist=(THnSparseF*)fArrayDz->FindObject(name);
-//   if(!curHist){
-//     curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
-//     fArrayDz->AddLast(curHist);
-//   }
-//   curHist->Fill(vecDrift);
-         
-//   name=shortName;
-//   name+="ALL";
-//   name.ToUpper();
-//   curHist=(THnSparseF*)fArrayDz->FindObject(name);
-//   if(!curHist){
-//     curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
-//     fArrayDz->AddLast(curHist);
-//   }
-//   curHist->Fill(vecDrift);
 }
 
 void AliTPCcalibTime::ProcessCosmic(const AliESDEvent *const event){
@@ -745,10 +776,232 @@ void AliTPCcalibTime::ProcessCosmic(const AliESDEvent *const event){
   if (GetDebugLevel()>20) printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
 }
 
-void AliTPCcalibTime::ProcessBeam(const AliESDEvent *const /*event*/){
+void AliTPCcalibTime::ProcessBeam(const AliESDEvent *const event){
   //
-  // Not special treatment yet - the same for cosmic and physic event
   //
+  // 
+  const Int_t kMinClusters  =80;
+  const Int_t kMinTracks    =2;      // minimal number of tracks to define the vertex
+  const Int_t kMinTracksVertex=30;   // minimal number of tracks to define the cumulative vertex
+  const Double_t kMaxTgl    =1.2;    // maximal Tgl (z angle)
+  const Double_t kMinPt     =0.2;    // minimal pt
+  const Double_t kMaxD0     =5.;     // cut on distance to the primary vertex first guess
+  const Double_t kMaxZ0     =20; 
+  const Double_t kMaxD      =2.5;    // cut on distance to the primary vertex 
+  const Double_t kMaxZ      =4;      // maximal z distance between tracks form the same side
+  const Double_t kMaxChi2   =15;     // maximal chi2 of the TPCvertex 
+  const Double_t kCumulCovarXY=0.003; //increase the error of cumul vertex 30 microns profile
+  const Double_t kCumulCovarZ=250.;  //increase the error of cumul vertex
+  const Double_t kMaxDvertex = 1.0;  // cut to accept the vertex; 
+  //
+  Int_t  flags=0;
+  const  Int_t  kBuffSize=100;
+  static Double_t deltaZ[kBuffSize]={0};
+  static Int_t counterZ=0;
+  static AliKFVertex cumulVertexA, cumulVertexC, cumulVertexAC; // cumulative vertex 
+  AliKFVertex vertexA, vertexC;
+
+  Float_t  dca0[2]={0,0};
+  Double_t dcaVertex[2]={0,0};
+  Int_t ntracks=event->GetNumberOfTracks();
+  if (ntracks==0) return;
+  //
+  AliESDfriend *esdFriend=(AliESDfriend*)(((AliESDEvent*)event)->FindListObject("AliESDfriend"));
+  //
+  // Divide tracks to A and C side tracks - using the cluster indexes
+  TObjArray tracksA(ntracks);  
+  TObjArray tracksC(ntracks);  
+  //
+  AliESDVertex *vertexSPD =  (AliESDVertex *)event->GetPrimaryVertexSPD();
+  AliESDVertex *vertex    =  (AliESDVertex *)event->GetPrimaryVertex();
+  AliESDVertex *vertexTracks =  (AliESDVertex *)event->GetPrimaryVertexTracks();
+  Double_t vertexZA[10000], vertexZC[10000];
+  //
+  Int_t ntracksA= 0;
+  Int_t ntracksC= 0;
+  //
+  for (Int_t itrack=0;itrack<ntracks;itrack++) {
+    AliESDtrack *track = event->GetTrack(itrack);
+    AliESDfriendTrack *friendTrack = esdFriend->GetTrack(itrack);
+    if (!friendTrack) continue;
+    if (TMath::Abs(track->GetTgl())>kMaxTgl) continue;
+    if (TMath::Abs(track->Pt())<kMinPt) continue;
+    const AliExternalTrackParam * trackIn  = track->GetInnerParam();
+    TObject *calibObject=0;
+    AliTPCseed *seed = 0;
+    Int_t nA=0, nC=0;
+    for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
+    if (seed) {
+      for (Int_t irow=159;irow>0;irow--) {
+       AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+       if (!cl) continue;
+       if ((cl->GetDetector()%36)<18) nA++;
+       if ((cl->GetDetector()%36)>=18) nC++;
+      }
+      if ((nA>kMinClusters || nC>kMinClusters) && (nA*nC==0) ){
+       track->GetImpactParameters(dca0[0],dca0[1]);
+       if (TMath::Abs(dca0[0])>kMaxD0) continue;
+       if (TMath::Abs(dca0[1])>kMaxZ0) continue;
+       AliExternalTrackParam pTPCvertex(*trackIn);
+       if (!AliTracker::PropagateTrackToBxByBz(&pTPCvertex,4.+4.*TMath::Abs(dca0[0]),0.1,2,kTRUE)) continue;
+       pTPCvertex.PropagateToDCA(vertex,AliTracker::GetBz(), kMaxD, dcaVertex,0);
+       if (TMath::Abs(dcaVertex[0])>kMaxD) continue;
+       if (nA>kMinClusters &&nC==0) { tracksA.AddLast(pTPCvertex.Clone()); vertexZA[ntracksA++] = pTPCvertex.GetZ();}
+       if (nC>kMinClusters &&nA==0) {tracksC.AddLast(pTPCvertex.Clone());  vertexZC[ntracksC++] = pTPCvertex.GetZ();}
+      }
+    }
+  }
+  Double_t medianZA=TMath::Median(ntracksA, vertexZA);  // tracks median
+  Double_t medianZC=TMath::Median(ntracksC, vertexZC);  // tracks median
+  //
+  ntracksA= tracksA.GetEntriesFast();
+  ntracksC= tracksC.GetEntriesFast();
+  if (ntracksA>kMinTracks && ntracksC>kMinTracks){
+    deltaZ[counterZ%kBuffSize]=medianZA-medianZC;
+    counterZ+=1;
+    Double_t medianDelta=(counterZ>=kBuffSize)? TMath::Median(kBuffSize, deltaZ): TMath::Median(counterZ, deltaZ);
+    if (TMath::Abs(medianDelta-(medianZA-medianZC))>kMaxZ) flags+=16;
+    // increse the error of cumulative vertex at the beginning of event
+    cumulVertexA.Covariance(0,0)+=kCumulCovarXY*kCumulCovarXY;
+    cumulVertexA.Covariance(1,1)+=kCumulCovarXY*kCumulCovarXY;
+    cumulVertexA.Covariance(2,2)+=kCumulCovarZ*kCumulCovarZ;
+    cumulVertexC.Covariance(0,0)+=kCumulCovarXY*kCumulCovarXY;
+    cumulVertexC.Covariance(1,1)+=kCumulCovarXY*kCumulCovarXY;
+    cumulVertexC.Covariance(2,2)+=kCumulCovarZ*kCumulCovarZ;
+    cumulVertexAC.Covariance(0,0)+=kCumulCovarXY*kCumulCovarXY;
+    cumulVertexAC.Covariance(1,1)+=kCumulCovarXY*kCumulCovarXY;
+    cumulVertexAC.Covariance(2,2)+=kCumulCovarZ*kCumulCovarZ;
+    //
+    for (Int_t iA=0; iA<ntracksA; iA++){
+      if (flags!=0) continue;
+      AliExternalTrackParam *aliTrack =  (AliExternalTrackParam *)tracksA.At(iA);
+      if (TMath::Abs(aliTrack->GetZ()-medianZA)>kMaxZ) continue;
+      AliKFParticle part(*aliTrack,211);
+      vertexA+=part;
+    }  
+    for (Int_t iC=0; iC<ntracksC; iC++){
+      if (flags!=0) continue;
+      AliExternalTrackParam *aliTrack =  (AliExternalTrackParam *)tracksC.At(iC);
+      if (TMath::Abs(aliTrack->GetZ()-medianZC)>kMaxZ) continue;
+      AliKFParticle part(*aliTrack,211);
+      vertexC+=part;
+    }   
+    //
+    if (vertexA.GetNDF()<kMinTracks) flags+=32;
+    if (vertexC.GetNDF()<kMinTracks) flags+=32;
+    if (TMath::Abs(vertexA.Z()-medianZA)>kMaxZ) flags+=1;   //apply cuts
+    if (TMath::Abs(vertexC.Z()-medianZC)>kMaxZ) flags+=2;
+    if (TMath::Abs(vertexA.GetChi2()/vertexA.GetNDF()+vertexC.GetChi2()/vertexC.GetNDF())> kMaxChi2) flags+=4;
+    //
+    if (flags==0){
+      for (Int_t iA=0; iA<ntracksA; iA++){
+       if (flags!=0) continue;
+       AliExternalTrackParam *aliTrack =  (AliExternalTrackParam *)tracksA.At(iA);
+       if (TMath::Abs(aliTrack->GetZ()-medianZA)>kMaxZ) continue;
+       AliKFParticle part(*aliTrack,211);
+       cumulVertexA+=part;
+       cumulVertexAC+=part;
+      }        
+      for (Int_t iC=0; iC<ntracksC; iC++){
+       if (flags!=0) continue;
+       AliExternalTrackParam *aliTrack =  (AliExternalTrackParam *)tracksC.At(iC);
+       if (TMath::Abs(aliTrack->GetZ()-medianZC)>kMaxZ) continue;
+       AliKFParticle part(*aliTrack,211);
+       cumulVertexC+=part;
+       cumulVertexAC+=part;
+      }             
+      //
+      if (TMath::Abs(cumulVertexA.X()-vertexA.X())>kMaxDvertex) flags+=64;
+      if (TMath::Abs(cumulVertexA.Y()-vertexA.Y())>kMaxDvertex) flags+=64;
+      if (TMath::Abs(cumulVertexA.Z()-vertexA.Z())>kMaxDvertex) flags+=64;
+      //
+      if (TMath::Abs(cumulVertexC.X()-vertexC.X())>kMaxDvertex) flags+=64;
+      if (TMath::Abs(cumulVertexC.Y()-vertexC.Y())>kMaxDvertex) flags+=64;
+      if (TMath::Abs(cumulVertexC.Z()-vertexC.Z())>kMaxDvertex) flags+=64;
+      
+      
+      if ( flags==0 && cumulVertexC.GetNDF()>kMinTracksVertex&&cumulVertexA.GetNDF()>kMinTracksVertex){
+       Double_t cont[2]={0,fTime};
+       //
+       cont[0]= cumulVertexA.X();
+       fTPCVertex[0]->Fill(cont);
+       cont[0]= cumulVertexC.X();
+       fTPCVertex[1]->Fill(cont);
+       cont[0]= 0.5*(cumulVertexA.X()-cumulVertexC.X());
+       fTPCVertex[2]->Fill(cont);
+       cont[0]= 0.5*(cumulVertexA.X()+cumulVertexC.X())-vertexSPD->GetX();
+       fTPCVertex[3]->Fill(cont);
+       //
+       cont[0]= cumulVertexA.Y();
+       fTPCVertex[4]->Fill(cont);
+       cont[0]= cumulVertexC.Y();
+       fTPCVertex[5]->Fill(cont);
+       cont[0]= 0.5*(cumulVertexA.Y()-cumulVertexC.Y());
+       fTPCVertex[6]->Fill(cont);
+       cont[0]= 0.5*(cumulVertexA.Y()+cumulVertexC.Y())-vertexSPD->GetY();
+       fTPCVertex[7]->Fill(cont);
+       //
+       //
+       cont[0]= 0.5*(cumulVertexA.Z()+cumulVertexC.Z());
+       fTPCVertex[8]->Fill(cont);
+       cont[0]= 0.5*(cumulVertexA.Z()-cumulVertexC.Z());
+       fTPCVertex[9]->Fill(cont);
+       cont[0]= 0.5*(cumulVertexA.Z()-cumulVertexC.Z());
+       fTPCVertex[10]->Fill(cont);      
+       cont[0]= 0.5*(cumulVertexA.Z()+cumulVertexC.Z())-vertexSPD->GetZ();
+       fTPCVertex[11]->Fill(cont);
+       //
+       Double_t correl[2]={0,0};
+       //
+       correl[0]=cumulVertexC.Z();
+       correl[1]=cumulVertexA.Z();
+       fTPCVertexCorrelation[0]->Fill(correl);   // fill A side :TPC
+       correl[0]=cumulVertexA.Z();
+       correl[1]=cumulVertexC.Z(); 
+       fTPCVertexCorrelation[1]->Fill(correl);   // fill C side :TPC
+       //
+       correl[0]=vertexSPD->GetZ();
+       correl[1]=cumulVertexA.Z()-correl[0];
+       fTPCVertexCorrelation[2]->Fill(correl);   // fill A side :ITS
+       correl[1]=cumulVertexC.Z()-correl[0]; 
+       fTPCVertexCorrelation[3]->Fill(correl);   // fill C side :ITS
+       correl[1]=0.5*(cumulVertexA.Z()+cumulVertexC.Z())-correl[0]; 
+       fTPCVertexCorrelation[4]->Fill(correl);   // fill C side :ITS
+      }
+    }        
+    TTreeSRedirector *cstream = GetDebugStreamer();
+    if (cstream){
+      /*
+       TCut cutChi2= "sqrt(vA.fChi2/vA.fNDF+vC.fChi2/vC.fNDF)<10";  // chi2 Cut e.g 10         
+       TCut cutXY= "sqrt((vA.fP[0]-vC.fP[0])^2+(vA.fP[0]-vC.fP[1])^2)<5";   // vertex Cut      
+       TCut cutZ= "abs(vA.fP[2]-mZA)<3&&abs(vC.fP[2]-mZC)<5";           // vertex Cut  
+       tree->Draw("sqrt(vA.fChi2/vA.fNDF)","sqrt(vA.fChi2/vA.fNDF)<100","")
+       
+      */
+      //vertexA.Print();
+      //vertexC.Print();      
+      (*cstream)<<"vertexTPC"<<
+       "flags="<<flags<<        // rejection flags
+       "vSPD.="<<vertexSPD<<    // SPD vertex
+       "vT.="<<vertexTracks<<   // track vertex
+       "v.="<<vertex<<          // esd vertex
+       "mZA="<<medianZA<<       // median Z position at vertex A side
+       "mZC="<<medianZC<<       // median Z position at vertex C side
+       "mDelta="<<medianDelta<< // median delta A side -C side
+       "counter="<<counterZ<<    // counter Z
+       //
+       "vA.="<<&vertexA<<       // vertex A side
+       "vC.="<<&vertexC<<       // vertex C side
+       "cvA.="<<&cumulVertexA<<       // cumulative vertex A side
+       "cvC.="<<&cumulVertexC<<       // cumulative vertex C side
+       "cvAC.="<<&cumulVertexAC<<       // cumulative vertex A+C side
+       "nA="<<ntracksA<<        // contributors
+       "nC="<<ntracksC<<        // contributors
+       "\n";
+    }      
+  }
+  tracksA.Delete();
+  tracksC.Delete();
 }
 
 void AliTPCcalibTime::Analyze(){
@@ -850,7 +1103,7 @@ Long64_t AliTPCcalibTime::Merge(TCollection *const li) {
   //
   TIterator* iter = li->MakeIterator();
   AliTPCcalibTime* cal = 0;
-
+  //
   while ((cal = (AliTPCcalibTime*)iter->Next())) {
     if (!cal->InheritsFrom(AliTPCcalibTime::Class())) {
       Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
@@ -863,35 +1116,61 @@ Long64_t AliTPCcalibTime::Merge(TCollection *const li) {
       }
     }
     //
-    for (Int_t imeas=0; imeas<5; imeas++){
-      if (cal->GetResHistoTPCCE(imeas) && cal->GetResHistoTPCCE(imeas)){
-       fResHistoTPCCE[imeas]->Add(cal->fResHistoTPCCE[imeas]);
-      }else{
-        fResHistoTPCCE[imeas]=(THnSparse*)cal->fResHistoTPCCE[imeas]->Clone();
+    if (fTPCVertexCorrelation && cal->fTPCVertexCorrelation){
+      for (Int_t imeas=0; imeas<5; imeas++){
+       if (fTPCVertexCorrelation[imeas] && cal->fTPCVertexCorrelation[imeas]) fTPCVertexCorrelation[imeas]->Add(cal->fTPCVertexCorrelation[imeas]);
+      }
+    }
+      
+    if (fTPCVertex && cal->fTPCVertex) 
+      for (Int_t imeas=0; imeas<12; imeas++){
+       if (fTPCVertex[imeas] && cal->fTPCVertex[imeas]) fTPCVertex[imeas]->Add(cal->fTPCVertex[imeas]);
       }
-      if (cal->GetResHistoTPCITS(imeas) && cal->GetResHistoTPCITS(imeas)){
-       fResHistoTPCITS[imeas]->Add(cal->fResHistoTPCITS[imeas]);
-       fResHistoTPCvertex[imeas]->Add(cal->fResHistoTPCvertex[imeas]);
+    
+    if (fMemoryMode>0) for (Int_t imeas=0; imeas<5; imeas++){
+      if (fMemoryMode>1){
+       if ( cal->GetResHistoTPCCE(imeas) && cal->GetResHistoTPCCE(imeas)){
+         fResHistoTPCCE[imeas]->Add(cal->fResHistoTPCCE[imeas]);
+       }else{
+         fResHistoTPCCE[imeas]=(THnSparse*)cal->fResHistoTPCCE[imeas]->Clone();
+       }
+      }
+      //
+      if ((fMemoryMode>0) &&cal->GetResHistoTPCITS(imeas) && cal->GetResHistoTPCITS(imeas)){
+       if (fMemoryMode>1 || (imeas%2)==1) fResHistoTPCITS[imeas]->Add(cal->fResHistoTPCITS[imeas]);
+       if (fMemoryMode>1) fResHistoTPCvertex[imeas]->Add(cal->fResHistoTPCvertex[imeas]);
       }
-      if (cal->fResHistoTPCTRD[imeas]){
+      //
+      if ((fMemoryMode>1) && cal->fResHistoTPCTRD[imeas]){
        if (fResHistoTPCTRD[imeas])
          fResHistoTPCTRD[imeas]->Add(cal->fResHistoTPCTRD[imeas]);
        else
          fResHistoTPCTRD[imeas]=(THnSparse*)cal->fResHistoTPCTRD[imeas]->Clone();
       }
-      if  (cal->fResHistoTPCTOF[imeas]){
+      //
+      if  ((fMemoryMode>1) && cal->fResHistoTPCTOF[imeas]){
        if (fResHistoTPCTOF[imeas])
          fResHistoTPCTOF[imeas]->Add(cal->fResHistoTPCTOF[imeas]);
        else
          fResHistoTPCTOF[imeas]=(THnSparse*)cal->fResHistoTPCTOF[imeas]->Clone();      
       }
+      //
+      if (cal->fArrayLaserA){
+       fArrayLaserA->Expand(fArrayLaserA->GetEntriesFast()+cal->fArrayLaserA->GetEntriesFast());
+       fArrayLaserC->Expand(fArrayLaserC->GetEntriesFast()+cal->fArrayLaserC->GetEntriesFast());
+       for (Int_t ical=0; ical<cal->fArrayLaserA->GetEntriesFast(); ical++){
+         if (cal->fArrayLaserA->UncheckedAt(ical)) fArrayLaserA->AddLast(cal->fArrayLaserA->UncheckedAt(ical)->Clone());
+         if (cal->fArrayLaserC->UncheckedAt(ical)) fArrayLaserC->AddLast(cal->fArrayLaserC->UncheckedAt(ical)->Clone());
+       }
+      }
+
     }
     TObjArray* addArray=cal->GetHistoDrift();
     if(!addArray) return 0;
     TIterator* iterator = addArray->MakeIterator();
     iterator->Reset();
     THnSparse* addHist=NULL;
-    while((addHist=(THnSparseF*)iterator->Next())){
+    if ((fMemoryMode>1)) while((addHist=(THnSparseF*)iterator->Next())){
       if(!addHist) continue;
       addHist->Print();
       THnSparse* localHist=(THnSparseF*)fArrayDz->FindObject(addHist->GetName());
@@ -1073,8 +1352,8 @@ void  AliTPCcalibTime::ProcessSame(AliESDtrack *const track, AliESDfriendTrack *
   Double_t xyz[3]={0,0.,0.0};  
   Double_t bz   =0;
   Int_t nclIn=0,nclOut=0;
-  trackIn.ResetCovariance(10.);
-  trackOut.ResetCovariance(10.);
+  trackIn.ResetCovariance(1000.);
+  trackOut.ResetCovariance(1000.);
   //
   //2.a Refit inner
   // 
@@ -1209,41 +1488,40 @@ void  AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTra
   const Int_t    kMinITS  = 3;     // minimal number of ITS cluster
   const Double_t kMinZ    = 10;    // maximal dz distance
   const Double_t kMaxDy   = 2.;    // maximal dy distance
-  const Double_t kMaxAngle= 0.015;  // maximal angular distance
+  const Double_t kMaxAngle= 0.07;  // maximal angular distance
   const Double_t kSigmaCut= 5;     // maximal sigma distance to median
   const Double_t kVdErr   = 0.1;  // initial uncertainty of the vd correction 
+  const Double_t kT0Err   = 3.;  // initial uncertainty of the T0 time
   const Double_t kVdYErr  = 0.05;  // initial uncertainty of the vd correction 
   const Double_t kOutCut  = 1.0;   // outlyer cut in AliRelAlgnmentKalman
   const Double_t kMinPt   = 0.3;   // minimal pt
+  const Double_t kMax1Pt=0.5;        //maximal 1/pt distance
   const  Int_t     kN=50;         // deepnes of history
   static Int_t     kglast=0;
   static Double_t* kgdP[4]={new Double_t[kN], new Double_t[kN], new Double_t[kN], new Double_t[kN]};
-  /*
-    0. Standrd cuts:
-    TCut cut="abs(pTPC.fP[2]-pITS.fP[2])<0.01&&abs(pTPC.fP[3]-pITS.fP[3])<0.01&&abs(pTPC.fP[2]-pITS.fP[2])<1";
-  */
   //
   // 0. Apply standard cuts
-  //
+  // 
   Int_t dummycl[1000];
   if (track->GetTPCNcls()<kMinTPC) return;  // minimal amount of clusters cut
-  if (track->GetITSclusters(dummycl)<kMinITS) return;  // minimal amount of clusters
   if (!track->IsOn(AliESDtrack::kTPCrefit)) return;
-  if (!friendTrack->GetITSOut()) return;  
   if (!track->GetInnerParam())   return;
   if (!track->GetOuterParam())   return;
   if (track->GetInnerParam()->Pt()<kMinPt)  return;
   // exclude crossing track
   if (track->GetOuterParam()->GetZ()*track->GetInnerParam()->GetZ()<0)   return;
-  if (TMath::Abs(track->GetInnerParam()->GetZ())<kMinZ)   return;
+  if (TMath::Abs(track->GetInnerParam()->GetZ())<kMinZ/3.)   return;
   if (track->GetInnerParam()->GetX()>90)   return;
   //
   AliExternalTrackParam &pTPC=(AliExternalTrackParam &)(*(track->GetInnerParam()));
-  AliExternalTrackParam pITS(*(friendTrack->GetITSOut()));   // ITS standalone if possible
-  AliExternalTrackParam pITS2(*(friendTrack->GetITSOut()));  //TPC-ITS track
-  pITS2.Rotate(pTPC.GetAlpha());
-  //  pITS2.PropagateTo(pTPC.GetX(),fMagF);
-  AliTracker::PropagateTrackToBxByBz(&pITS2,pTPC.GetX(),0.1,0.1,kFALSE);
+  //  
+  AliExternalTrackParam pITS;   // ITS standalone if possible
+  AliExternalTrackParam pITS2;  //TPC-ITS track
+  if (friendTrack->GetITSOut()){
+    pITS2=(*(friendTrack->GetITSOut()));  //TPC-ITS track - snapshot ITS out
+    pITS2.Rotate(pTPC.GetAlpha());
+    AliTracker::PropagateTrackToBxByBz(&pITS2,pTPC.GetX(),0.1,0.1,kFALSE);
+  }
 
   AliESDfriendTrack *itsfriendTrack=0;
   //
@@ -1252,21 +1530,27 @@ void  AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTra
   Bool_t hasAlone=kFALSE;
   Int_t ntracks=event->GetNumberOfTracks();
   for (Int_t i=0; i<ntracks; i++){
-    AliESDtrack *trackS = event->GetTrack(i);
-    if (trackS->GetTPCNcls()>0) continue;   //continue if has TPC info
+    AliESDtrack * trackITS = event->GetTrack(i); 
+    if (!trackITS) continue;
+    if (trackITS->GetITSclusters(dummycl)<kMinITS) continue;  // minimal amount of clusters
     itsfriendTrack = esdFriend->GetTrack(i);
     if (!itsfriendTrack) continue;
     if (!itsfriendTrack->GetITSOut()) continue;
-    if (TMath::Abs(pITS2.GetTgl()-itsfriendTrack->GetITSOut()->GetTgl())> kMaxAngle) continue;
+     
+    if (TMath::Abs(pTPC.GetTgl()-itsfriendTrack->GetITSOut()->GetTgl())> kMaxAngle) continue;
+    if (TMath::Abs(pTPC.GetSigned1Pt()-itsfriendTrack->GetITSOut()->GetSigned1Pt())> kMax1Pt) continue;
     pITS=(*(itsfriendTrack->GetITSOut()));
     //
     pITS.Rotate(pTPC.GetAlpha());
-    //pITS.PropagateTo(pTPC.GetX(),fMagF);
     AliTracker::PropagateTrackToBxByBz(&pITS,pTPC.GetX(),0.1,0.1,kFALSE);
-    if (TMath::Abs(pITS2.GetY()-pITS.GetY())> kMaxDy) continue;
+    if (TMath::Abs(pTPC.GetY()-pITS.GetY())> kMaxDy) continue;
+    if (TMath::Abs(pTPC.GetSnp()-pITS.GetSnp())> kMaxAngle) continue;
     hasAlone=kTRUE;
   }
-  if (!hasAlone) pITS=pITS2;
+  if (!hasAlone) {
+    if (track->GetITSclusters(dummycl)<kMinITS) return;
+    pITS=pITS2;  // use combined track if it has ITS
+  }
   //
   if (TMath::Abs(pITS.GetY()-pTPC.GetY())    >kMaxDy)    return;
   if (TMath::Abs(pITS.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
@@ -1301,7 +1585,7 @@ void  AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTra
   //
   // 3. Update alignment
   //
-  Int_t htime = fTime/3600; //time in hours
+  Int_t htime = (fTime-fTimeKalmanBin/2)/fTimeKalmanBin; //time bins number
   if (fAlignITSTPC->GetEntriesFast()<htime){
     fAlignITSTPC->Expand(htime*2+20);
   }
@@ -1311,6 +1595,7 @@ void  AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTra
     align=new AliRelAlignerKalman(); 
     align->SetRunNumber(fRun);
     (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
+    (*align->GetStateCov())(7,7)=kT0Err*kT0Err;
     (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
     align->SetOutRejSigma(kOutCut+kOutCut*kN);
     align->SetRejectOutliers(kFALSE);
@@ -1320,7 +1605,12 @@ void  AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTra
     fAlignITSTPC->AddAt(align,htime);
   }
   align->AddTrackParams(&pITS,&pTPC);
-  align->SetTimeStamp(fTime);
+  Double_t averageTime =  fTime;
+  if (align->GetTimeStamp()>0&&align->GetNUpdates()>0){
+    averageTime=((Double_t(align->GetTimeStamp())*Double_t(align->GetNUpdates())+Double_t(fTime)))/(Double_t(align->GetNUpdates())+1.);
+  }
+  align->SetTimeStamp(Int_t(averageTime));
+
   align->SetRunNumber(fRun );
   Float_t dca[2],cov[3];
   track->GetImpactParameters(dca,cov);
@@ -1386,6 +1676,7 @@ void  AliTPCcalibTime::ProcessAlignTRD(AliESDtrack *const track, AliESDfriendTra
   const Double_t kMaxAngle= 0.1;  // maximal angular distance
   const Double_t kSigmaCut= 10;     // maximal sigma distance to median
   const Double_t kVdErr   = 0.1;  // initial uncertainty of the vd correction 
+  const Double_t kT0Err   = 3.;  // initial uncertainty of the T0 time
   const Double_t kVdYErr  = 0.05;  // initial uncertainty of the vd correction 
   const Double_t kOutCut  = 1.0;   // outlyer cut in AliRelAlgnmentKalman
   const Double_t kRefX    = 275;   // reference X
@@ -1451,7 +1742,8 @@ void  AliTPCcalibTime::ProcessAlignTRD(AliESDtrack *const track, AliESDfriendTra
   //
   // 3. Update alignment
   //
-  Int_t htime = fTime/3600; //time in hours
+  //Int_t htime = fTime/3600; //time in hours
+  Int_t htime = (Int_t)(fTime-fTimeKalmanBin/2)/fTimeKalmanBin; //time in half hour
   if (fAlignTRDTPC->GetEntriesFast()<htime){
     fAlignTRDTPC->Expand(htime*2+20);
   }
@@ -1461,6 +1753,7 @@ void  AliTPCcalibTime::ProcessAlignTRD(AliESDtrack *const track, AliESDfriendTra
     align=new AliRelAlignerKalman(); 
     align->SetRunNumber(fRun);
     (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
+    (*align->GetStateCov())(7,7)=kT0Err*kT0Err;
     (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
     align->SetOutRejSigma(kOutCut+kOutCut*kN);
     align->SetRejectOutliers(kFALSE);
@@ -1469,7 +1762,16 @@ void  AliTPCcalibTime::ProcessAlignTRD(AliESDtrack *const track, AliESDfriendTra
     fAlignTRDTPC->AddAt(align,htime);
   }
   align->AddTrackParams(&pTRD,&pTPC);
-  align->SetTimeStamp(fTime);
+  //align->SetTimeStamp(fTime);
+  Double_t averageTime =  fTime;
+  if (align->GetTimeStamp()>0 && align->GetNUpdates()>0) {
+    averageTime = (((Double_t)fTime) + ((Double_t)align->GetTimeStamp())*align->GetNUpdates()) / (align->GetNUpdates() + 1.);
+    //printf("align->GetTimeStamp() %d, align->GetNUpdates() %d \n", align->GetTimeStamp(), align->GetNUpdates());
+  }
+  align->SetTimeStamp((Int_t)averageTime);
+
+  //printf("fTime %d, averageTime %d \n", fTime, (Int_t)averageTime);
+
   align->SetRunNumber(fRun );
   Float_t dca[2],cov[3];
   track->GetImpactParameters(dca,cov);
@@ -1527,9 +1829,10 @@ void  AliTPCcalibTime::ProcessAlignTOF(AliESDtrack *const track, AliESDfriendTra
   const Int_t      kMinTPC  = 80;    // minimal number of TPC cluster
   //  const Double_t   kMinZ    = 10;    // maximal dz distance
   const Double_t   kMaxDy   = 5.;    // maximal dy distance
-  const Double_t   kMaxAngle= 0.015;  // maximal angular distance
+  const Double_t   kMaxAngle= 0.05;  // maximal angular distance
   const Double_t   kSigmaCut= 5;     // maximal sigma distance to median
   const Double_t   kVdErr   = 0.1;  // initial uncertainty of the vd correction 
+  const Double_t   kT0Err   = 3.;  // initial uncertainty of the T0 time
   const Double_t   kVdYErr  = 0.05;  // initial uncertainty of the vd correction 
 
   const Double_t   kOutCut  = 1.0;   // outlyer cut in AliRelAlgnmentKalman
@@ -1616,7 +1919,8 @@ void  AliTPCcalibTime::ProcessAlignTOF(AliESDtrack *const track, AliESDfriendTra
   //
   // 3. Update alignment
   //
-  Int_t htime = fTime/3600; //time in hours
+  //Int_t htime = fTime/3600; //time in hours
+  Int_t htime = (Int_t)(fTime-fTimeKalmanBin)/fTimeKalmanBin; //time bin
   if (fAlignTOFTPC->GetEntriesFast()<htime){
     fAlignTOFTPC->Expand(htime*2+20);
   }
@@ -1626,6 +1930,7 @@ void  AliTPCcalibTime::ProcessAlignTOF(AliESDtrack *const track, AliESDfriendTra
     align=new AliRelAlignerKalman(); 
     align->SetRunNumber(fRun);
     (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
+    (*align->GetStateCov())(7,7)=kT0Err*kT0Err;
     (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
     align->SetOutRejSigma(kOutCut+kOutCut*kN);
     align->SetRejectOutliers(kFALSE);
@@ -1639,7 +1944,16 @@ void  AliTPCcalibTime::ProcessAlignTOF(AliESDtrack *const track, AliESDfriendTra
   if (TMath::Abs(dca[0])<kMaxDy){
     FillResHistoTPCTOF(&pTPC,&pTOF);
   }
-  align->SetTimeStamp(fTime);
+  //align->SetTimeStamp(fTime);
+  Double_t averageTime =  fTime;
+  if (align->GetTimeStamp()>0 && align->GetNUpdates()>0) {
+    averageTime = (((Double_t)fTime) + ((Double_t)align->GetTimeStamp())*align->GetNUpdates()) / (align->GetNUpdates() + 1.);
+    //printf("align->GetTimeStamp() %d, align->GetNUpdates() %d \n", align->GetTimeStamp(), align->GetNUpdates());
+  }
+  align->SetTimeStamp((Int_t)averageTime);
+
+  //printf("fTime %d, averageTime %d \n", fTime, (Int_t)averageTime);
+
   align->SetRunNumber(fRun );
   //
   Int_t nupdates=align->GetNUpdates();
@@ -1693,12 +2007,12 @@ void  AliTPCcalibTime::BookDistortionMaps(){
   axisName[0]   ="#Delta";
   axisTitle[0]  ="#Delta";
   //
-  binsTrack[1] =20;
-  xminTrack[1] =-1.5; xmaxTrack[1]=1.5;
+  binsTrack[1] =44;
+  xminTrack[1] =-1.1; xmaxTrack[1]=1.1;
   axisName[1]  ="tanTheta";
   axisTitle[1]  ="tan(#Theta)";
   //
-  binsTrack[2] =90;
+  binsTrack[2] =180;
   xminTrack[2] =-TMath::Pi(); xmaxTrack[2]=TMath::Pi(); 
   axisName[2]  ="phi";
   axisTitle[2]  ="#phi";
@@ -1770,6 +2084,35 @@ void  AliTPCcalibTime::BookDistortionMaps(){
       }
     }
   }
+  //
+  // Book vertex: time histograms
+  //
+  Int_t    binsVertex[2]={500, fTimeBins};
+  Double_t    aminVertex[2]={-5,fTimeStart};
+  Double_t    amaxVertex[2]={5, fTimeEnd};
+  const char* hnames[12]={"TPCXAside", "TPCXCside","TPCXACdiff","TPCXAPCdiff",
+                         "TPCYAside", "TPCYCside","TPCYACdiff","TPCYAPCdiff",
+                         "TPCZAPCside", "TPCZAMCside","TPCZACdiff","TPCZAPCdiff"}; 
+  const char* anames[12]={"x (cm) - A side ", "x (cm) - C side","#Delta_{x} (cm) - TPC-A-C","#Delta_{x} (cm) - TPC-Common",
+                         "y (cm) - A side ", "y (cm) - C side","#Delta_{x} (cm) - TPC-A-C","#Delta_{y} (cm) - TPC-Common",
+                         "z (cm)", "#Delta_{Z} (cm) A-C side","#Delta_{x} (cm) - TPC-A-C","#Delta_{Z} (cm) TPC-common"}; 
+  for (Int_t ihis=0; ihis<12; ihis++) {
+    if (ihis>=8) aminVertex[0]=-20.;
+    if (ihis>=8) amaxVertex[0]=20.;
+    fTPCVertex[ihis]=new THnSparseF(hnames[ihis],hnames[ihis],2,binsVertex,aminVertex,amaxVertex);
+    fTPCVertex[ihis]->GetAxis(1)->SetTitle("Time");
+    fTPCVertex[ihis]->GetAxis(0)->SetTitle(anames[ihis]);
+  }
+  
+  Int_t    binsVertexC[2]={40, 300};
+  Double_t aminVertexC[2]={-20,-30};
+  Double_t amaxVertexC[2]={20,30};
+  const char* hnamesC[5]={"TPCA_TPC","TPCC_TPC","TPCA_ITS","TPCC_ITS","TPC_ITS"};
+  for (Int_t ihis=0; ihis<5; ihis++) {
+    fTPCVertexCorrelation[ihis]=new THnSparseF(hnamesC[ihis],hnamesC[ihis],2,binsVertexC,aminVertexC,amaxVertexC);
+    fTPCVertexCorrelation[ihis]->GetAxis(1)->SetTitle("z (cm)");
+    fTPCVertexCorrelation[ihis]->GetAxis(0)->SetTitle("z (cm)");
+  }
 }
 
 
@@ -1778,6 +2121,7 @@ void        AliTPCcalibTime::FillResHistoTPCCE(const AliExternalTrackParam * pTP
   // fill residual histograms pTPCOut-pTPCin - trac crossing CE
   // Histogram 
   //
+  if (fMemoryMode<2) return;
   Double_t histoX[5];
   Double_t xyz[3];
   pTPCIn->GetXYZ(xyz);
@@ -1823,6 +2167,7 @@ void        AliTPCcalibTime::FillResHistoTPC(const AliESDtrack * pTrack){
   // fill residual histograms pTPC - vertex
   // Histogram is filled only for primary tracks
   //
+  if (fMemoryMode<2) return;
   Double_t histoX[4];
   const AliExternalTrackParam * pTPCIn = pTrack->GetInnerParam();
   AliExternalTrackParam pTPCvertex(*(pTrack->GetInnerParam()));
@@ -1854,6 +2199,7 @@ void        AliTPCcalibTime::FillResHistoTPCTRD(const AliExternalTrackParam * pT
   //
   // fill resuidual histogram TPCout-TRDin
   //
+  if (fMemoryMode<2) return;
   Double_t histoX[4];
   Double_t xyz[3];
   pTPCOut->GetXYZ(xyz);
@@ -1878,6 +2224,7 @@ void        AliTPCcalibTime::FillResHistoTPCTOF(const AliExternalTrackParam * pT
   //
   // fill resuidual histogram TPCout-TOFin
   // track propagated to the TOF position
+  if (fMemoryMode<2) return;
   Double_t histoX[4];
   Double_t xyz[3];