AliTPCAnalysisTaskcalib.cxx - Disable creation of debug streamers (Jacek)

AliTPCcalibDB.cxx           - a) Use ITS only if close measurement in time
                            - b) use global y gradient of vd from ITS if available
AliTPCcalibDButil.cxx       - a) Accept one entry graphs
                              b) Filtering onli in case of suficient number of points
AliTPCcalibTime.cxx         - using average time as a reference time
AliTPCTransform.cxx         - a) caching of the global y gradient of drift velocity
                            - b) disabling pad time 0 correction
CalibViewer:

AliTPCcalibSummary.cxx      - Adding attachment to trend graphs
                            - Use constant for eval outside of range

AliTPCcalibAlign.cxx        - increase the error for edge clusters
AliTPCkalmanAlign.h AliTPCkalmanAlign.cxx - update of matrix with one measurement

diff --git a/TPC/AliTPCAnalysisTaskcalib.cxx b/TPC/AliTPCAnalysisTaskcalib.cxx
index aca74e75f710ab577b2246acaa71bc4095601f92..cbbaef4a093801312276cf3ae29671cd329ebe86 100644 (file)
--- a/TPC/AliTPCAnalysisTaskcalib.cxx
+++ b/TPC/AliTPCAnalysisTaskcalib.cxx
@@ -42,7 +42,7 @@ AliTPCAnalysisTaskcalib::AliTPCAnalysisTaskcalib()
    fCalibJobs(0),
    fESD(0),
    fESDfriend(0),
-   fDebugOutputPath()
+   fDebugOutputPath("")
 {
   //
   // default constructor
@@ -56,7 +56,7 @@ AliTPCAnalysisTaskcalib::AliTPCAnalysisTaskcalib(const char *name)
    fCalibJobs(0),
    fESD(0),
    fESDfriend(0),
-   fDebugOutputPath()
+   fDebugOutputPath("")
 {
   //
   // Constructor
@@ -153,7 +153,9 @@ void AliTPCAnalysisTaskcalib::FinishTaskOutput()
   // on the slaves before sending data
   //
   Terminate("slave");
-  RegisterDebugOutput();
+  if(!fDebugOutputPath.IsNull()) { 
+    RegisterDebugOutput();
+  }
   
 }
 

diff --git a/TPC/AliTPCTransform.cxx b/TPC/AliTPCTransform.cxx
index daadb94bd2a59d041133d2270997b469b3974884..6c09db364433df5c5d16f96db04f71199b9fdda8 100755 (executable)
--- a/TPC/AliTPCTransform.cxx
+++ b/TPC/AliTPCTransform.cxx
@@ -164,7 +164,7 @@ void AliTPCTransform::Transform(Double_t *x,Int_t *i,UInt_t /*time*/,
   Double_t xx[3];
   //  Apply Time0 correction - Pad by pad fluctuation
   //
-  x[2]-=time0TPC->GetCalROC(sector)->GetValue(row,pad);
+  //x[2]-=time0TPC->GetCalROC(sector)->GetValue(row,pad);
   //
   // Tranform from pad - time coordinate system to the rotated global (tracking) system
   //
@@ -304,6 +304,7 @@ void AliTPCTransform::Local2RotatedGlobal(Int_t sector, Double_t *x) const {
   //
   // simple caching non thread save
   static Double_t vdcorrectionTime=1;
+  static Double_t vdcorrectionTimeGY=0;
   static Double_t time0corrTime=0;
   static Int_t    lastStampT=-1;
   //
@@ -323,15 +324,13 @@ void AliTPCTransform::Local2RotatedGlobal(Int_t sector, Double_t *x) const {
     }
     //
     if(fCurrentRecoParam&&fCurrentRecoParam->GetUseDriftCorrectionGY()>0) {
-      Float_t xyzPad[3];
-      AliTPCROC::Instance()->GetPositionGlobal(sector, TMath::Nint(x[0]) ,TMath::Nint(x[1]), xyzPad);
       
-      Double_t corrGy= (1+(xyzPad[1])*AliTPCcalibDB::Instance()->
+      Double_t corrGy= AliTPCcalibDB::Instance()->
                        GetVDriftCorrectionGy(fCurrentTimeStamp, 
                                              AliTPCcalibDB::Instance()->GetRun(),
                                              sector%36>=18,
-                                             fCurrentRecoParam->GetUseDriftCorrectionGY()));
-      vdcorrectionTime *=corrGy;
+                                             fCurrentRecoParam->GetUseDriftCorrectionGY());
+      vdcorrectionTimeGY = corrGy;
     }
   }
 
@@ -345,7 +344,9 @@ void AliTPCTransform::Local2RotatedGlobal(Int_t sector, Double_t *x) const {
   const Int_t kNIS=param->GetNInnerSector(), kNOS=param->GetNOuterSector();
   Double_t sign = 1.;
   Double_t zwidth    = param->GetZWidth()*driftCorr;
-  if (AliTPCRecoParam:: GetUseTimeCalibration()) zwidth*=vdcorrectionTime;
+  Float_t xyzPad[3];
+  AliTPCROC::Instance()->GetPositionGlobal(sector, TMath::Nint(x[0]) ,TMath::Nint(x[1]), xyzPad);
+  if (AliTPCRecoParam:: GetUseTimeCalibration()) zwidth*=vdcorrectionTime*(1+xyzPad[1]*vdcorrectionTimeGY);
   Double_t padWidth  = 0;
   Double_t padLength = 0;
   Double_t    maxPad    = 0;

diff --git a/TPC/AliTPCcalibAlign.cxx b/TPC/AliTPCcalibAlign.cxx
index 0486c8e7036351cff77a0c640344297285bf9351..c5238eb622b05edc27dc0a46300efa62b32f1d66 100644 (file)
--- a/TPC/AliTPCcalibAlign.cxx
+++ b/TPC/AliTPCcalibAlign.cxx
@@ -2000,22 +2000,23 @@ void  AliTPCcalibAlign::MakeTree(const char *fname, Int_t minPoints){
 
       //
       //
-      TMatrixD * kpar = fSectorParamA;
-      TMatrixD * kcov = fSectorCovarA;
-      if (s1%36>=18){
-       kpar = fSectorParamC;
-       kcov = fSectorCovarC;
-      }
-      for (Int_t ipar=0;ipar<6;ipar++){
-       Int_t isec1 = s1%18;
-       Int_t isec2 = s2%18;
-       if (s1>35) isec1+=18;
-       if (s2>35) isec2+=18;   
-       param6s1(ipar)=(*kpar)(6*isec1+ipar,0);
-       param6s2(ipar)=(*kpar)(6*isec2+ipar,0);
+      if (fSectorParamA){
+       TMatrixD * kpar = fSectorParamA;
+       TMatrixD * kcov = fSectorCovarA;
+       if (s1%36>=18){
+         kpar = fSectorParamC;
+         kcov = fSectorCovarC;
+       }
+       for (Int_t ipar=0;ipar<6;ipar++){
+         Int_t isec1 = s1%18;
+         Int_t isec2 = s2%18;
+         if (s1>35) isec1+=18;
+         if (s2>35) isec2+=18; 
+         param6s1(ipar)=(*kpar)(6*isec1+ipar,0);
+         param6s2(ipar)=(*kpar)(6*isec2+ipar,0);
+       }
       }
-
-
+       
       Double_t dy=0, dz=0, dphi=0,dtheta=0;
       Double_t sy=0, sz=0, sphi=0,stheta=0;
       Double_t ny=0, nz=0, nphi=0,ntheta=0;
@@ -2641,11 +2642,12 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
   // 3. Refit the track - out-in
   //                    - update the cluster delta histo lower part
   //
-  const Double_t kPtCut=1;    // pt
+  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;
   if (!fCurrentTrack) return;
   if (!fCurrentFriendTrack) return;
   Float_t vertexXY=0,vertexZ=0;
@@ -2662,8 +2664,20 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
   //
   AliExternalTrackParam trackIn  = *(fCurrentTrack->GetInnerParam());
   AliExternalTrackParam trackOut = *(fCurrentFriendTrack->GetTPCOut());
-  Double_t mass =    TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
+  static Double_t mass =    TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
   //  
+  Int_t ncl=0;
+  for (Int_t irow=0; irow<160; irow++){
+    AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+    if (!cl) continue;
+    if (cl->GetX()<80) continue;
+    if (detector<0) detector=cl->GetDetector()%36;
+    if (detector!=cl->GetDetector()%36) return; // cluster from different sectors
+    // skip such tracks
+    ncl++;
+  }
+  if (ncl<kNclCut) return;
+
   Int_t nclIn=0,nclOut=0;
   Double_t xyz[3];
   //
@@ -2682,13 +2696,19 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
     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 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
     cov[0]*=cov[0];
     cov[2]*=cov[2];
     if (!AliTracker::PropagateTrackToBxByBz(&trackOut, r[0],mass,1.,kFALSE)) continue;
+    
+    if (TMath::Abs(dedge)<kEdgeCut) continue;
+
     if (TMath::Abs(cl->GetY()-trackOut.GetY())<kMaxDist){
       nclOut++;
       trackOut.Update(&r[1],cov);
@@ -2726,13 +2746,19 @@ void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
     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 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 +-
     cov[0]*=cov[0];
     cov[2]*=cov[2];
     if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1.,kFALSE)) continue;
+    if (TMath::Abs(dedge)<kEdgeCut) continue;
+
+
     if (TMath::Abs(cl->GetY()-trackIn.GetY())<kMaxDist){
       nclIn++;
       trackIn.Update(&r[1],cov);
@@ -3074,7 +3100,7 @@ void  AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){
            "vPosG.="<<&vPosG<<        // vertex position in global system
            "vPosL.="<<&vPosL<<        // vertex position in local  system
            "vImpact.="<<&vImpact<<   // track impact parameter
-           "tofSignal="<<tofSignal<<   // tof signal
+           //"tofSignal="<<tofSignal<<   // tof signal
            "tpcPosG.="<<&tpcPosG<<   // global position of track at the middle of fXmiddle
            "lphi="<<lphi<<             // expected local phi angle - if vertex at 0
            "gphi="<<gphi<<             // expected global phi if vertex at 0

diff --git a/TPC/AliTPCcalibDB.cxx b/TPC/AliTPCcalibDB.cxx
index aabccafc569ad811bedf1f072215d084081793d0..a308fb417361511e4feac142d4562136d70e76c0 100644 (file)
--- a/TPC/AliTPCcalibDB.cxx
+++ b/TPC/AliTPCcalibDB.cxx
@@ -243,7 +243,7 @@ AliTPCcalibDB::AliTPCcalibDB(const AliTPCcalibDB& ):
   fVoltageArray(0),
   fTemperatureArray(0),   //! array of temperature sensors - per run - Just for calibration studies
   fVdriftArray(0),         //! array of v drift interfaces
-  fDriftCorrectionArray(0),         //! array of v drift interfaces
+  fDriftCorrectionArray(0), //! array of v drift corrections
   fRunList(0),              //! run list - indicates try to get the run param 
   fDButil(0),
   fCTPTimeParams(0)
@@ -323,6 +323,7 @@ void AliTPCcalibDB::Update(){
   AliCDBManager::Instance()->SetCacheFlag(kTRUE); // activate CDB cache
   fDButil = new AliTPCcalibDButil;   
   //
+
   entry          = GetCDBEntry("TPC/Calib/PadGainFactor");
   if (entry){
     //if (fPadGainFactor) delete fPadGainFactor;
@@ -1002,12 +1003,15 @@ void AliTPCcalibDB::UpdateRunInformations( Int_t run, Bool_t force){
     fDButil->FilterSensor(press,kMinP,kMaxP,kMaxdP,kSigmaCut);
     if (press->GetFit()==0) accept=kFALSE;
   }
+
   if (press && temp &&accept){
     AliTPCCalibVdrift * vdrift = new AliTPCCalibVdrift(temp, press,0);
     fVdriftArray.AddAt(vdrift,run);
   }
+
   fDButil->FilterCE(120., 3., 4.,0);
   fDButil->FilterTracks(run, 10.,0);
+
 }
 
 
@@ -1777,6 +1781,9 @@ Double_t AliTPCcalibDB::GetVDriftCorrectionTime(Int_t timeStamp, Int_t run, Int_
   if (mode==1) {
     const Double_t kEpsilon=0.00000000001;
     const Double_t kdeltaT=360.; // 10 minutes
+    if(TMath::Abs(driftITS) < 12*kdeltaT) {
+      result = driftITS;
+    } else {
     wITS  = 64.*kdeltaT/(deltaITS +kdeltaT);
     wLT   = 16.*kdeltaT/(deltaLT  +kdeltaT);
     wP    = 0. *kdeltaT/(deltaP   +kdeltaT);
@@ -1789,6 +1796,9 @@ Double_t AliTPCcalibDB::GetVDriftCorrectionTime(Int_t timeStamp, Int_t run, Int_
     if (TMath::Abs(driftCE)<kEpsilon) wCE=0;  // invalid calibration
     if (wP+wITS+wLT+wCE<kEpsilon) return 0;
     result = (driftP*wP+driftITS*wITS+driftLT*wLT+driftCE*wCE)/(wP+wITS+wLT+wCE);
+   }
+   
+
   }
 
   return result;
@@ -1829,9 +1839,9 @@ Double_t AliTPCcalibDB::GetVDriftCorrectionGy(Int_t timeStamp, Int_t run, Int_t
   //
   // Get global y correction drift velocity correction factor
   // additive factor        vd = vdnom*(1+GetVDriftCorrectionGy *gy)
-  // Value etracted combining the vdrift correction using laser tracks and CE
+  // Value etracted combining the vdrift correction using laser tracks and CE or TPC-ITS
   // Arguments:
-  // mode determines the algorith how to combine the Laser Track, LaserCE
+  // mode determines the algorith how to combine the Laser Track, LaserCE or TPC-ITS
   // timestamp - timestamp
   // run       - run number
   // side      - the drift velocity gy correction per side (CE and Laser tracks)
@@ -1842,10 +1852,25 @@ Double_t AliTPCcalibDB::GetVDriftCorrectionGy(Int_t timeStamp, Int_t run, Int_t
   UpdateRunInformations(run,kFALSE);
   TObjArray *array =AliTPCcalibDB::Instance()->GetTimeVdriftSplineRun(run);
   if (!array) return 0;
+  Double_t result=0;
+
+  // use TPC-ITS if present
+  TGraphErrors *gr= (TGraphErrors*)array->FindObject("ALIGN_ITSB_TPC_VDGY");
+  if(gr) { 
+    result = (gr->Eval(timeStamp));
+
+    // transform from [(cm/mus)/ m] to [1/cm]
+    result /= (fParam->GetDriftV()/1000000.);
+    result /= 100.;
+
+    //printf("result %e \n", result);
+    return result; 
+  }
+
+  // use laser if ITS-TPC not present
   TGraphErrors *laserA= (TGraphErrors*)array->FindObject("GRAPH_MEAN_GLOBALYGRADIENT_LASER_ALL_A");
   TGraphErrors *laserC= (TGraphErrors*)array->FindObject("GRAPH_MEAN_GLOBALYGRADIENT_LASER_ALL_C");
   
-  Double_t result=0;
   if (laserA && laserC){
    result= (laserA->Eval(timeStamp)+laserC->Eval(timeStamp))*0.5;
   }
@@ -1855,6 +1880,8 @@ Double_t AliTPCcalibDB::GetVDriftCorrectionGy(Int_t timeStamp, Int_t run, Int_t
   if (laserC &&side==1){
     result = (laserC->Eval(timeStamp));
   }
+  //printf("laser result %e \n", -result/250.);
+
   return -result/250.; //normalized before
 }
 

diff --git a/TPC/AliTPCcalibDButil.cxx b/TPC/AliTPCcalibDButil.cxx
index e1bc734b198a5b9cf01515a96c4a8db288655a58..f9726f0d8ad8069b2ca323cd0a903b4cdef50426 100644 (file)
--- a/TPC/AliTPCcalibDButil.cxx
+++ b/TPC/AliTPCcalibDButil.cxx
@@ -1607,17 +1607,27 @@ Int_t AliTPCcalibDButil::GetNearest(TGraph *graph, Double_t xref, Double_t &dx,
   //
   // find the closest point to xref  in x  direction
   // return dx and value 
+  dx = 0;
+  y = 0;
+
+  if(!graph) return 0;
+  if(graph->GetN() < 1) return 0;
+
   Int_t index=0;
   index = TMath::BinarySearch(graph->GetN(), graph->GetX(),xref);
   if (index<0) index=0;
-  if (index>=graph->GetN()-1) index=graph->GetN()-2;
-  if (xref-graph->GetX()[index]>graph->GetX()[index]-xref) index++;
-  dx = xref-graph->GetX()[index];
+  if(graph->GetN()==1) {
+    dx = xref-graph->GetX()[index];
+  }
+  else {
+    if (index>=graph->GetN()-1) index=graph->GetN()-2;
+    if (xref-graph->GetX()[index]>graph->GetX()[index]-xref) index++;
+    dx = xref-graph->GetX()[index];
+  }
   y  = graph->GetY()[index];
   return index;
 }
 
-
 Double_t  AliTPCcalibDButil::GetTriggerOffsetTPC(Int_t run, Int_t timeStamp, Double_t deltaT, Double_t deltaTLaser, Int_t valType){
   //
   // Get the correction of the trigger offset
@@ -1721,7 +1731,7 @@ Double_t  AliTPCcalibDButil::GetVDriftTPC(Double_t &dist, Int_t run, Int_t timeS
   AliTPCcalibDButil::GetNearest(cosmicAll,timeStamp,dist,grY);
 
   Double_t t0= AliTPCcalibDButil::GetTriggerOffsetTPC(run,timeStamp, deltaT, deltaTLaser,valType);
-  Double_t vcosmic=  AliTPCcalibDButil::EvalGraphConst(cosmicAll, timeStamp);
+  Double_t vcosmic =  AliTPCcalibDButil::EvalGraphConst(cosmicAll, timeStamp);
   if (timeStamp>cosmicAll->GetX()[cosmicAll->GetN()-1])  vcosmic=cosmicAll->GetY()[cosmicAll->GetN()-1];
   if (timeStamp<cosmicAll->GetX()[0])  vcosmic=cosmicAll->GetY()[0];
   return  vcosmic-t0;
@@ -1831,12 +1841,12 @@ Double_t  AliTPCcalibDButil::GetVDriftTPCLaserTracks(Double_t &dist, Int_t run,
   grlaserA=(TGraphErrors*)array->FindObject("GRAPH_MEAN_DRIFT_LASER_ALL_A");
   grlaserC=(TGraphErrors*)array->FindObject("GRAPH_MEAN_DRIFT_LASER_ALL_C");
   Double_t deltaY;
-  if (grlaserA) {
+  if (grlaserA && grlaserA->GetN()>0) {
     AliTPCcalibDButil::GetNearest(grlaserA,timeStamp,dist,deltaY);
     if (TMath::Abs(dist)>deltaT)  vlaserA= deltaY;
     else  vlaserA = AliTPCcalibDButil::EvalGraphConst(grlaserA,timeStamp);
   }
-  if (grlaserC) {
+  if (grlaserC && grlaserC->GetN()>0) {
     AliTPCcalibDButil::GetNearest(grlaserC,timeStamp,dist,deltaY);
     if (TMath::Abs(dist)>deltaT)  vlaserC= deltaY;
     else  vlaserC = AliTPCcalibDButil::EvalGraphConst(grlaserC,timeStamp);
@@ -1872,6 +1882,7 @@ Double_t  AliTPCcalibDButil::GetVDriftTPCCE(Double_t &dist,Int_t run, Int_t time
   Double_t gry=0;
   Double_t corrA=0, corrC=0;
   Double_t timeA=0, timeC=0;
+  const Double_t kEpsilon = 0.00001;
   TGraph *graphA = (TGraph*)arrT->At(72);
   TGraph *graphC = (TGraph*)arrT->At(73);
   if (!graphA && !graphC) return 0.;
@@ -1880,6 +1891,7 @@ Double_t  AliTPCcalibDButil::GetVDriftTPCCE(Double_t &dist,Int_t run, Int_t time
     timeA   = AliTPCcalibDButil::EvalGraphConst(graphA,timeStamp);
     Int_t mtime   =TMath::Nint((graphA->GetX()[0]+graphA->GetX()[graphA->GetN()-1])*0.5);
     ltime0A       = GetLaserTime0(run,mtime,TMath::Nint(deltaT),0);
+    if(ltime0A < kEpsilon) return 0;
     if (driftCalib) corrPTA =  driftCalib->GetPTRelative(timeStamp,0);
     corrA = (param->GetZLength(36)/(timeA*param->GetTSample()*(1.-ltime0A)-param->GetL1Delay()-0*param->GetZSigma()/param->GetDriftV()))/param->GetDriftV()-1;
     corrA-=corrPTA;
@@ -1889,7 +1901,8 @@ Double_t  AliTPCcalibDButil::GetVDriftTPCCE(Double_t &dist,Int_t run, Int_t time
     timeC=AliTPCcalibDButil::EvalGraphConst(graphC,timeStamp);
     Int_t mtime=TMath::Nint((graphC->GetX()[0]+graphC->GetX()[graphC->GetN()-1])*0.5);
     ltime0C       = GetLaserTime0(run,mtime,TMath::Nint(deltaT),0);
-    if (driftCalib) corrPTC =  driftCalib->GetPTRelative(timeStamp,0);
+    if(ltime0C < kEpsilon) return 0;   
+if (driftCalib) corrPTC =  driftCalib->GetPTRelative(timeStamp,0);
     corrC = (param->GetZLength(54)/(timeC*param->GetTSample()*(1.-ltime0C)-param->GetL1Delay()-0*param->GetZSigma()/param->GetDriftV()))/param->GetDriftV()-1;
     corrC-=corrPTC;
   }
@@ -1911,6 +1924,7 @@ Double_t  AliTPCcalibDButil::GetVDriftTPCITS(Double_t &dist, Int_t run, Int_t ti
   TGraphErrors *graph=0;
   dist=0;
   if (!array) return 0;
+  //array->ls();
   graph = (TGraphErrors*)array->FindObject("ALIGN_ITSB_TPC_DRIFTVD");
   if (!graph) return 0;
   Double_t deltaY;
@@ -2211,13 +2225,23 @@ Double_t AliTPCcalibDButil::EvalGraphConst(TGraph * const graph, Double_t xref){
     printf("AliTPCcalibDButil::EvalGraphConst: 0 pointer\n");
     return 0;
   }
+
   if (graph->GetN()<1){
-    printf("AliTPCcalibDButil::EvalGraphConst: Empty graph");
+    printf("AliTPCcalibDButil::EvalGraphConst: Empty graph \n");
     return 0;
   }
+ 
+
   if (xref<graph->GetX()[0]) return graph->GetY()[0];
   if (xref>graph->GetX()[graph->GetN()-1]) return graph->GetY()[graph->GetN()-1]; 
-  return graph->Eval( xref);
+
+  printf("graph->Eval(graph->GetX()[0]) %f, graph->Eval(xref) %f \n",graph->Eval(graph->GetX()[0]), graph->Eval(xref));
+
+  if(graph->GetN()==1)
+    return graph->Eval(graph->GetX()[0]);
+
+
+  return graph->Eval(xref);
 }
 
 Float_t AliTPCcalibDButil::FilterSensor(AliDCSSensor * sensor, Double_t ymin, Double_t ymax, Double_t maxdy,  Double_t sigmaCut){
@@ -2544,9 +2568,18 @@ void AliTPCcalibDButil::FilterTracks(Int_t run, Double_t cutSigma, TTreeSRedirec
       arrT->AddAt(0,i);
       continue;
     }
-    TGraphErrors *graph2= FilterGraphMedianErr(graph,cutSigma,medianY);
-    if (!graph2) {
-      delete graph; arrT->AddAt(0,i); continue;
+    TGraphErrors *graph2 = NULL;
+    if(graph->GetN()<10) {
+      graph2 = new TGraphErrors(graph->GetN(),graph->GetX(),graph->GetY(),graph->GetEX(),graph->GetEY()); 
+      if (!graph2) {
+        delete graph; arrT->AddAt(0,i); continue;
+      }
+    } 
+    else {
+      graph2= FilterGraphMedianErr(graph,cutSigma,medianY);
+      if (!graph2) {
+        delete graph; arrT->AddAt(0,i); continue;
+      }
     }
     if (graph2->GetN()<1) {
       delete graph; arrT->AddAt(0,i); continue;
@@ -2825,7 +2858,8 @@ TMatrixD* AliTPCcalibDButil::MakeStatRelKalman(TObjArray * const array, Float_t
       (*valArray[ipar])[naccept]=state[ipar];
     naccept++;
   }
-  if (naccept<2) return 0;
+  //if (naccept<2) return 0;
+  if (naccept<1) return 0;
   TMatrixD *pstat=new TMatrixD(9,3);
   TMatrixD &stat=*pstat;
   for (Int_t ipar=0; ipar<9; ipar++){

diff --git a/TPC/AliTPCcalibSummary.cxx b/TPC/AliTPCcalibSummary.cxx
index 5f3b7f81a93a074508177159296d651ce71891e4..dfa0d90525a58fa363ebef5ce837b622e7ddf50a 100644 (file)
--- a/TPC/AliTPCcalibSummary.cxx
+++ b/TPC/AliTPCcalibSummary.cxx
@@ -645,13 +645,13 @@ void AliTPCcalibSummary::ProcessKryptonTime(Int_t run, Int_t timeStamp){
       krErr[isec]=0;
       gr=(TGraphErrors*)krArray->At(isec);
       if (gr) {
-       krMean[isec]=gr->Eval(timeStamp);
+       krMean[isec]=AliTPCcalibDButil::EvalGraphConst(gr,timeStamp);
        AliTPCcalibDButil::GetNearest(gr, timeStamp,deltaT,gry);
        krDist[isec]=deltaT;
       }     
       if (72+isec<krArray->GetEntries()) {
        gr=(TGraphErrors*)krArray->At(72+isec);
-       if (gr) krErr[isec]=gr->Eval(timeStamp);
+       if (gr) krErr[isec]=AliTPCcalibDButil::EvalGraphConst(gr,timeStamp);
       }
     }
     krMean[72]= TMath::Median(36,krMean.GetMatrixArray());
@@ -694,16 +694,24 @@ void AliTPCcalibSummary::ProcessGain(Int_t irun, Int_t timeStamp){
   //
   static Float_t  gainCosmic = 0;
   static Float_t  gainMIP = 0;
+  static Float_t  attachMIP = 0;
+  static Double_t dMIP=0; 
+  Double_t dummy=0;
   TObjArray * gainSplines = fCalibDB->GetTimeGainSplinesRun(irun);
   if (gainSplines) {
     TGraphErrors * graphMIP = (TGraphErrors *) gainSplines->FindObject("TGRAPHERRORS_MEAN_GAIN_BEAM_ALL");
     TGraphErrors * graphCosmic = (TGraphErrors *) gainSplines->FindObject("TGRAPHERRORS_MEAN_GAIN_COSMIC_ALL");
+    TGraphErrors * graphAttach = (TGraphErrors *) gainSplines->FindObject("TGRAPHERRORS_MEAN_ATTACHMENT_BEAM_ALL");
     if (graphMIP) gainMIP = AliTPCcalibDButil::EvalGraphConst(graphMIP,timeStamp);
     if (graphCosmic) gainCosmic = AliTPCcalibDButil::EvalGraphConst(graphCosmic,timeStamp);
+    if (graphAttach) attachMIP = AliTPCcalibDButil::EvalGraphConst(graphAttach,timeStamp);
+    if (graphMIP)  AliTPCcalibDButil::GetNearest(graphMIP, timeStamp, dMIP,dummy);
   }
   // time dependence of gain
   (*fPcstream)<<"dcs"<<
     "gainMIP="<<gainMIP<<
+    "attachMIP="<<attachMIP<<
+    "dMIP="<<dMIP<<
     "gainCosmic="<<gainCosmic;     
 }
 
@@ -733,7 +741,7 @@ void AliTPCcalibSummary::ProcessAlign(Int_t run, Int_t timeStamp){
       }
       values[9*idet+ipar]=0;
       errs[9*idet+ipar]=0;
-      if (gr) values[9*idet+ipar]=gr->Eval(timeStamp);
+      if (gr) values[9*idet+ipar]=AliTPCcalibDButil::EvalGraphConst(gr,timeStamp);
       (*fPcstream)<<"dcs"<<
        Form("%s=",grName.Data())<<values[9*idet+ipar];
       (*fPcstream)<<"align"<<
@@ -764,8 +772,8 @@ void AliTPCcalibSummary::ProcessDriftCERef(){
   static TVectorD vecALx(72);
   static TVectorD vecAChi2(72);
   //
-  static TVectorD vecQ(72);
-  static TVectorD vecQRef(72);
+  static TVectorD vecASide(4);
+  static TVectorD vecCSide(4);
   static Bool_t isCalc=kFALSE;
   
   TFile f("calPads.root");
@@ -789,7 +797,7 @@ void AliTPCcalibSummary::ProcessDriftCERef(){
     fstringG+="ly.fElements++";
     fstringG+="(lx.fElements-134.)++";  
     for (Int_t isec=0; isec<72; isec++){
-      TStatToolkit::FitPlane(tree,"2.64*dt", fstringG.Data(),Form("sector==%d",isec)+cutAll, chi2,npoints,param,covar,-1,0, 10000000, kFALSE);
+      TStatToolkit::FitPlane(tree,"0.264*dt", fstringG.Data(),Form("sector==%d",isec)+cutAll, chi2,npoints,param,covar,-1,0, 10000000, kFALSE);
       if (npoints<3) continue;
       printf("Sector=%d\n",isec);
       vec0[isec]=param[0];
@@ -797,22 +805,31 @@ void AliTPCcalibSummary::ProcessDriftCERef(){
       vecLx[isec]=param[2];
       sec[isec]=isec;
       vecN[isec]=npoints;
+      vecChi2[isec]=TMath::Sqrt(chi2/npoints);
 
-      TStatToolkit::FitPlane(tree,"2.64*CETmean.fElements", fstringG.Data(),Form("sector==%d",isec)+cutAll, chi2,npoints,param,covar,-1,0, 10000000, kFALSE);
+      TStatToolkit::FitPlane(tree,"0.264*CETmean.fElements", fstringG.Data(),Form("sector==%d",isec)+cutAll, chi2,npoints,param,covar,-1,0, 10000000, kFALSE);
       if (npoints<3) continue;
       printf("Sector=%d\n",isec);
       vecA0[isec]=param[0];
       vecALy[isec]=param[1];
       vecALx[isec]=param[2];
       vecAChi2[isec]=TMath::Sqrt(chi2/npoints);
-      tree->Draw("CETmean.fElements",Form("sector==%d",isec)+cutAll);
-      tree->Draw("CETmean.fElements/R.CETmean.fElements",Form("sector==%d",isec)+cutAll);
     }
     isCalc=kTRUE;
+    //
+    TString  fstringRef="";              // global fit
+    fstringRef+="gx.fElements++";
+    fstringRef+="gy.fElements++";  
+    fstringRef+="lx.fElements++";
+    TStatToolkit::FitPlane(tree,"0.264*dt", fstringG.Data(),"(sector%36)<18"+cutAll, chi2,npoints,vecASide,covar,-1,0, 10000000, kFALSE);
+    TStatToolkit::FitPlane(tree,"0.264*dt", fstringG.Data(),"(sector%36)>=18"+cutAll, chi2,npoints,vecCSide,covar,-1,0, 10000000, kFALSE);
+
   }
   (*fPcstream)<<"dcs"<<     // CE information
     "CETSector.="<<&sec<<    // sector numbers
-    //                      // fit in respect to reference
+    "CETRefA.="<<&vecASide<<   // diff to reference A side
+    "CETRefC.="<<&vecCSide<<   // diff to reference C side    
+    //                      // fit in respect to reference data
     "CETRef0.="<<&vec0<<    // offset change
     "CETRefY.="<<&vecLy<<   // slope y change - rad
     "CETRefX.="<<&vecLx<<   // slope x change - rad
@@ -882,6 +899,7 @@ void AliTPCcalibSummary::ProcessPulserRef(){
       vecALx[isec]=param[2];
       vecAChi2[isec]=TMath::Sqrt(chi2/npoints);
     }
+
     isCalc=kTRUE;
   }
   (*fPcstream)<<"dcs"<<     // Pulser information

diff --git a/TPC/AliTPCcalibTime.cxx b/TPC/AliTPCcalibTime.cxx
index e06510db53da3a7c6cc7ccdb970ced90345945fd..a5587617879ca1ab1bcea9741b7ad983a4387d7b 100644 (file)
--- a/TPC/AliTPCcalibTime.cxx
+++ b/TPC/AliTPCcalibTime.cxx
@@ -1301,7 +1301,7 @@ void  AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTra
   //
   // 3. Update alignment
   //
-  Int_t htime = fTime/3600; //time in hours
+  Int_t htime = (fTime-1800/2)/1800; //time bins in half of hours
   if (fAlignITSTPC->GetEntriesFast()<htime){
     fAlignITSTPC->Expand(htime*2+20);
   }
@@ -1320,7 +1320,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=(align->GetTimeStamp()*align->GetNUpdates()+fTime)/(align->GetNUpdates()+1.);
+  }
+  align->SetTimeStamp(Int_t(averageTime));
+
   align->SetRunNumber(fRun );
   Float_t dca[2],cov[3];
   track->GetImpactParameters(dca,cov);
@@ -1392,6 +1397,7 @@ void  AliTPCcalibTime::ProcessAlignTRD(AliESDtrack *const track, AliESDfriendTra
   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]};
+  const Int_t    kTime=900; // time in seconds
   //
   // 0. Apply standard cuts
   //
@@ -1451,7 +1457,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-kTime)/1800; //time in half hour
   if (fAlignTRDTPC->GetEntriesFast()<htime){
     fAlignTRDTPC->Expand(htime*2+20);
   }
@@ -1469,7 +1476,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);
@@ -1536,6 +1552,7 @@ void  AliTPCcalibTime::ProcessAlignTOF(AliESDtrack *const track, AliESDfriendTra
   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]};
+  const Int_t      kTime=900; // time in seconds
   //
   // -1. Make a TOF track-
   //     Clusters are not in friends - use alingment points
@@ -1616,7 +1633,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-kTime)/1800; //time in half hour
   if (fAlignTOFTPC->GetEntriesFast()<htime){
     fAlignTOFTPC->Expand(htime*2+20);
   }
@@ -1639,7 +1657,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();

diff --git a/TPC/AliTPCkalmanAlign.cxx b/TPC/AliTPCkalmanAlign.cxx
index fca08881cb17d023075cfa83be73c28b9d8ab982..513afa63e7a3dccd5ac7732fa107d774f48a0831 100644 (file)
--- a/TPC/AliTPCkalmanAlign.cxx
+++ b/TPC/AliTPCkalmanAlign.cxx
@@ -192,6 +192,60 @@ void AliTPCkalmanAlign::UpdateAlign1D(Double_t delta, Double_t sigma, Int_t s1,
   covXk2= (mat1-(matKk*matHk));
   covXk3 =  covXk2*covXk;          
   covXk = covXk3;  
+  Int_t nrows=covXk3.GetNrows();
+  
+  for (Int_t irow=0; irow<nrows; irow++)
+    for (Int_t icol=0; icol<nrows; icol++){
+      // rounding problems - make matrix again symteric
+      covXk(irow,icol)=(covXk3(irow,icol)+covXk3(icol,irow))*0.5; 
+    }
+}
+
+void AliTPCkalmanAlign::UpdateAlign1D(Double_t delta, Double_t sigma, Int_t s1, TMatrixD &vecXk, TMatrixD &covXk){
+  //
+  // Update 1D kalman filter - with one measurement
+  //
+  const Int_t knMeas=1;
+  const Int_t knElem=72;
+  TMatrixD mat1(72,72);            // update covariance matrix
+  TMatrixD matHk(1,knElem);        // vector to mesurement
+  TMatrixD vecYk(knMeas,1);        // Innovation or measurement residual
+  TMatrixD matHkT(knElem,knMeas);  // helper matrix Hk transpose
+  TMatrixD matSk(knMeas,knMeas);   // Innovation (or residual) covariance
+  TMatrixD matKk(knElem,knMeas);   // Optimal Kalman gain
+  TMatrixD covXk2(knElem,knElem);  // helper matrix
+  TMatrixD covXk3(knElem,knElem);  // helper matrix
+  TMatrixD vecZk(1,1);
+  TMatrixD measR(1,1);
+  vecZk(0,0)=delta;
+  measR(0,0)=sigma*sigma;
+  //
+  // reset matHk
+  for (Int_t iel=0;iel<knElem;iel++) 
+    for (Int_t ip=0;ip<knMeas;ip++) matHk(ip,iel)=0; 
+  //mat1
+  for (Int_t iel=0;iel<knElem;iel++) {
+    for (Int_t jel=0;jel<knElem;jel++) mat1(iel,jel)=0;
+    mat1(iel,iel)=1;
+  }
+  //
+  matHk(0, s1)=1;
+  vecYk = vecZk-matHk*vecXk;               // Innovation or measurement residual
+  matHkT=matHk.T(); matHk.T();
+  matSk = (matHk*(covXk*matHkT))+measR;    // Innovation (or residual) covariance
+  matSk.Invert();
+  matKk = (covXk*matHkT)*matSk;            //  Optimal Kalman gain
+  vecXk += matKk*vecYk;                    //  updated vector 
+  covXk2= (mat1-(matKk*matHk));
+  covXk3 =  covXk2*covXk;          
+  covXk = covXk3;  
+  Int_t nrows=covXk3.GetNrows();
+  
+  for (Int_t irow=0; irow<nrows; irow++)
+    for (Int_t icol=0; icol<nrows; icol++){
+      // rounding problems - make matrix again symteric
+      covXk(irow,icol)=(covXk3(irow,icol)+covXk3(icol,irow))*0.5; 
+    }
 }
 
 

diff --git a/TPC/AliTPCkalmanAlign.h b/TPC/AliTPCkalmanAlign.h
index 16e5d42d2e0c17c8719c862010db8ab4fa7860f6..09e54fcba321fb97228997b9c689a8abfc483905 100644 (file)
--- a/TPC/AliTPCkalmanAlign.h
+++ b/TPC/AliTPCkalmanAlign.h
@@ -23,6 +23,7 @@ public:
   void DrawDeltaAlign();
   void UpdateOCDBTime0( AliTPCCalPad  *pad, Int_t ustartRun, Int_t uendRun,  const char* storagePath );
   static void UpdateAlign1D(Double_t delta, Double_t sigma, Int_t s1, Int_t s2, TMatrixD &param, TMatrixD &covar);
+  static void UpdateAlign1D(Double_t delta, Double_t sigma, Int_t s1, TMatrixD &param, TMatrixD &covar);
   static void BookAlign1D(TMatrixD &param, TMatrixD &covar, Double_t sigma, Double_t mean);
   void DumpOldAlignment(TTreeSRedirector *pcstream);
   void MakeNewAlignment(Bool_t add,TTreeSRedirector *pcstream=0);