]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCcalibLaser.cxx
adapted macro to QAManager
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibLaser.cxx
index 3978a3ba59f53764e5e017c094a09810253a9dd8..fcc5a6af59c6b41cdc50175153cf789de9137e7e 100644 (file)
   //
   // To make laser scan the user interaction neccessary
   //
-  .x ~/UliStyle.C
+  .x ~/NimStyle.C
   gSystem->Load("libANALYSIS");
   gSystem->Load("libTPCcalib");
   TFile fcalib("CalibObjects.root");
   TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib");
   AliTPCcalibLaser * laser = ( AliTPCcalibLaser *)array->FindObject("laserTPC");
-  laser->DumpMeanInfo(-0,0,10)
+  laser->DumpMeanInfo(-0,0)
   TFile fmean("laserMean.root")
   //
   //  laser track clasification;
@@ -71,6 +71,9 @@
   gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
   gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
   AliXRDPROOFtoolkit tool;
+  TChain * chainDrift = tool.MakeChain("laser.txt","driftv",0,10200);
+  chainDrift->Lookup();
+
   TChain * chain = tool.MakeChain("laser.txt","Residuals",0,10200);
   chain->Lookup();
   TChain * chainFit = tool.MakeChain("laser.txt","FitModels",0,10200);
 #include <iostream>
 #include <sstream>
 #include "AliTPCLaserTrack.h"
+#include "AliTPCcalibDB.h"
+#include "AliTPCParam.h"
+#include "TTimeStamp.h"
+#include "AliDCSSensorArray.h"
+#include "AliDCSSensor.h"
 
 using namespace std;
 
@@ -121,11 +129,13 @@ ClassImp(AliTPCcalibLaser)
 AliTPCcalibLaser::AliTPCcalibLaser():
   AliTPCcalibBase(),
   fESD(0),
-  fESDfriend(),
+  fESDfriend(0),
+  fNtracks(0),
   fTracksMirror(336),
   fTracksEsd(336),
   fTracksEsdParam(336),
   fTracksTPC(336),
+  fFullCalib(kTRUE),
   fDeltaZ(336),
   fDeltaP3(336),
   fDeltaP4(336),
@@ -176,16 +186,19 @@ AliTPCcalibLaser::AliTPCcalibLaser():
   fHisPz2vP2Out(0),  //-> Curv  P2outer - parabola
   fHisPz3vP2IO(0),   //-> Curv  P2outerinner - common parabola
   //
-  fDeltaYres(336),
-  fDeltaZres(336),  
-  fFitAside(new TVectorD(3)),
-  fFitCside(new TVectorD(3)),      
+  fDeltaYres(336),   //->2D histo of residuals
+  fDeltaZres(336),   //->2D histo fo residuals
+  fDeltaYres2(336),   //->2D histo of residuals
+  fDeltaZres2(336),   //->2D histo fo residuals
+  //fDeltaYres3(336),   //->2D histo of residuals
+  //fDeltaZres3(336),   //->2D histo fo residuals
+  fFitAside(new TVectorD(5)),
+  fFitCside(new TVectorD(5)),      
+  fFitACside(new TVectorD(6)),      
   fEdgeXcuts(3),    
   fEdgeYcuts(3),    
   fNClCuts(5),      
-  fNcuts(0),        
-  fRun(0),
-  fEvent(0)
+  fNcuts(0)        
 {
   //
   // Constructor
@@ -193,14 +206,16 @@ AliTPCcalibLaser::AliTPCcalibLaser():
   fTracksEsdParam.SetOwner(kTRUE);
 }
 
-AliTPCcalibLaser::AliTPCcalibLaser(const Text_t *name, const Text_t *title):
+AliTPCcalibLaser::AliTPCcalibLaser(const Text_t *name, const Text_t *title, Bool_t full):
   AliTPCcalibBase(),
   fESD(0),
   fESDfriend(0),
+  fNtracks(0),
   fTracksMirror(336),
   fTracksEsd(336),
   fTracksEsdParam(336),
   fTracksTPC(336),
+  fFullCalib(full),
   //
   fDeltaZ(336),          // array of histograms of delta z for each track
   fDeltaP3(336),          // array of histograms of delta z for each track
@@ -257,14 +272,17 @@ AliTPCcalibLaser::AliTPCcalibLaser(const Text_t *name, const Text_t *title):
   //
   fDeltaYres(336),
   fDeltaZres(336),  
-  fFitAside(new TVectorD(3)),        // drift fit - A side
-  fFitCside(new TVectorD(3)),        // drift fit - C- side
+  fDeltaYres2(336),
+  fDeltaZres2(336),  
+  //  fDeltaYres3(336),
+  //fDeltaZres3(336),  
+  fFitAside(new TVectorD(5)),        // drift fit - A side
+  fFitCside(new TVectorD(5)),        // drift fit - C- side
+  fFitACside(new TVectorD(6)),        // drift fit - AC- side
   fEdgeXcuts(3),       // cuts in local x direction; used in the refit of the laser tracks
   fEdgeYcuts(3),       // cuts in local y direction; used in the refit of the laser tracks
   fNClCuts(5),         // cuts on the number of clusters per tracklet; used in the refit of the laser tracks
-  fNcuts(0),           // number of cuts
-  fRun(0),
-  fEvent(0)
+  fNcuts(0)           // number of cuts
 {
   SetName(name);
   SetTitle(title);
@@ -278,10 +296,12 @@ AliTPCcalibLaser::AliTPCcalibLaser(const AliTPCcalibLaser& calibLaser):
   AliTPCcalibBase(calibLaser), 
   fESD(0),
   fESDfriend(0),
+  fNtracks(0),
   fTracksMirror(336),
   fTracksEsd(336),
   fTracksEsdParam(336),
   fTracksTPC(336),
+  fFullCalib(calibLaser.fFullCalib),
   //
   fDeltaZ(calibLaser.fDeltaZ),          // array of histograms of delta z for each track
   fDeltaP3(((calibLaser.fDeltaP3))),          // array of histograms of delta z for each track
@@ -337,14 +357,17 @@ AliTPCcalibLaser::AliTPCcalibLaser(const AliTPCcalibLaser& calibLaser):
   //
   fDeltaYres(((calibLaser.fDeltaYres))),
   fDeltaZres(((calibLaser.fDeltaZres))),  
-  fFitAside(new TVectorD(3)),        // drift fit - A side
-  fFitCside(new TVectorD(3)),        // drift fit - C- side
+  fDeltaYres2(((calibLaser.fDeltaYres))),
+  fDeltaZres2(((calibLaser.fDeltaZres))),  
+  //  fDeltaYres3(((calibLaser.fDeltaYres))),
+  //fDeltaZres3(((calibLaser.fDeltaZres))),  
+  fFitAside(new TVectorD(5)),        // drift fit - A side
+  fFitCside(new TVectorD(5)),        // drift fit - C- side
+  fFitACside(new TVectorD(6)),        // drift fit - C- side
   fEdgeXcuts(3),       // cuts in local x direction; used in the refit of the laser tracks
   fEdgeYcuts(3),       // cuts in local y direction; used in the refit of the laser tracks
   fNClCuts(5),         // cuts on the number of clusters per tracklet; used in the refit of the laser tracks
-  fNcuts(0),           // number of cuts
-  fRun(0),             // run number
-  fEvent(0)            // current eventnumber
+  fNcuts(0)           // number of cuts
 {
   //
   // copy constructor
@@ -415,6 +438,33 @@ AliTPCcalibLaser::~AliTPCcalibLaser() {
     delete fHisPz2vP2Out;  //-> Curv  P2outer - parabola
     delete fHisPz3vP2IO;   //-> Curv  P2outerinner - common parabola
   }
+  //
+  //
+  //
+  fDeltaZ.SetOwner();          //-> array of histograms of delta z for each track
+  fDeltaP3.SetOwner();         //-> array of histograms of P3      for each track
+  fDeltaP4.SetOwner();         //-> array of histograms of P4      for each track
+  fDeltaPhi.SetOwner();        //-> array of histograms of delta z for each track
+  fDeltaPhiP.SetOwner();       //-> array of histograms of delta z for each track
+  fSignals.SetOwner();         //->Array of dedx signals
+  
+  fDeltaZ.Delete();          //-> array of histograms of delta z for each track
+  fDeltaP3.Delete();         //-> array of histograms of P3      for each track
+  fDeltaP4.Delete();         //-> array of histograms of P4      for each track
+  fDeltaPhi.Delete();        //-> array of histograms of delta z for each track
+  fDeltaPhiP.Delete();       //-> array of histograms of delta z for each track
+  fSignals.Delete();         //->Array of dedx signals
+
+  fDeltaYres.SetOwner();
+  fDeltaYres.Delete();
+  fDeltaZres.SetOwner();
+  fDeltaZres.Delete();
+  fDeltaYres2.SetOwner();
+  fDeltaYres2.Delete();
+  fDeltaZres2.SetOwner();
+  fDeltaZres2.Delete();
+  
+
 }
 
 
@@ -424,6 +474,7 @@ void AliTPCcalibLaser::Process(AliESDEvent * event) {
   //
   // Loop over tracks and call  Process function
   //
+  Int_t kMinTracks=20;
   fESD = event;
   if (!fESD) {
     return;
@@ -432,6 +483,7 @@ void AliTPCcalibLaser::Process(AliESDEvent * event) {
   if (!fESDfriend) {
     return;
   }
+  if (fESD->GetNumberOfTracks()<kMinTracks) return; //not enough tracks
   AliDebug(4,Form("Event number in current file: %d",event->GetEventNumberInFile()));
   fTracksTPC.Clear();
   fTracksEsd.Clear();
@@ -441,29 +493,35 @@ void AliTPCcalibLaser::Process(AliESDEvent * event) {
     fClusterCounter[id]=0;
     fClusterSatur[id]=0;
   }
-  static Bool_t init=kFALSE;
-  if (!init){
-    init = kTRUE;  // way around for PROOF - to be investigated
-    MakeFitHistos();
-  }
   //
   Int_t n=fESD->GetNumberOfTracks();
-  Int_t run = fESD->GetRunNumber();
-  fRun = run;
+  Int_t counter=0;
   for (Int_t i=0;i<n;++i) {
     AliESDfriendTrack *friendTrack=fESDfriend->GetTrack(i);
+    if (!friendTrack) continue;
     AliESDtrack *track=fESD->GetTrack(i);
     TObject *calibObject=0;
     AliTPCseed *seed=0;
     for (Int_t j=0;(calibObject=friendTrack->GetCalibObject(j));++j)
       if ((seed=dynamic_cast<AliTPCseed*>(calibObject)))
        break;
-    if (track&&seed) {
-      FindMirror(track,seed);
+    if (track&&seed &&TMath::Abs(track->Pt()) >1 ) {
+      //filter CE tracks
+      Int_t id = FindMirror(track,seed);
+      if (id>0) counter++;
     }
     //
-  }  
+  } 
+  fNtracks=counter;
+  if (counter<kMinTracks) return;
+
   FitDriftV();
+  if (!fFullCalib) return;
+  static Bool_t init=kFALSE;
+  if (!init){
+    init = kTRUE;  // way around for PROOF - to be investigated
+    MakeFitHistos();
+  }
   //
   for (Int_t id=0; id<336; id++){    
     //
@@ -476,7 +534,6 @@ void AliTPCcalibLaser::Process(AliESDEvent * event) {
     }
   }
 
- fEvent++;
 }
 
 void AliTPCcalibLaser::MakeDistHisto(Int_t id){
@@ -491,38 +548,14 @@ void AliTPCcalibLaser::MakeDistHisto(Int_t id){
     //
     //
     TH1F * hisdz = (TH1F*)fDeltaZ.At(id);
+    if (!hisdz) MakeFitHistos();
+    hisdz = (TH1F*)fDeltaZ.At(id);
     TH1F * hisP3 = (TH1F*)fDeltaP3.At(id);
     TH1F * hisP4 = (TH1F*)fDeltaP4.At(id);
-
     TH1F * hisdphi = (TH1F*)fDeltaPhi.At(id);
     TH1F * hisdphiP = (TH1F*)fDeltaPhiP.At(id);
     TH1F * hisSignal = (TH1F*)fSignals.At(id);
-
-    if (!hisdz){
-      hisdz = new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10);
-      hisdz->SetDirectory(0);
-      fDeltaZ.AddAt(hisdz,id);
-
-      hisP3 = new TH1F(Form("hisPar3v%d",id),Form("hisPar3v%d",id),400,-0.06,0.06);
-      hisP3->SetDirectory(0);
-      fDeltaP3.AddAt(hisP3,id);
-      //
-      hisP4 = new TH1F(Form("hisPar4v%d",id),Form("hisPar4v%d",id),200,-0.06,0.06);
-      hisP4->SetDirectory(0);
-      fDeltaP4.AddAt(hisP4,id);
-
-      //
-      hisdphi = new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1);
-      hisdphi->SetDirectory(0);
-      fDeltaPhi.AddAt(hisdphi,id);
-      //
-      hisdphiP = new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01);
-      hisdphiP->SetDirectory(0);
-      fDeltaPhiP.AddAt(hisdphiP,id);
-      hisSignal = new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),100,0,300);
-      hisSignal->SetDirectory(0);
-      fSignals.AddAt(hisSignal,id);
-    }
+    //
 
     AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
     AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
@@ -553,38 +586,76 @@ void AliTPCcalibLaser::MakeDistHisto(Int_t id){
 
 void AliTPCcalibLaser::FitDriftV(){
   //
-  // Fit drift velocity - linear approximation in the z and global y
-  //
-  Float_t kDistCut           = 4;      // cut on distance to the fitted value
-  const Float_t kZCut        = 200;    // remove the closest laser beam
+  // Fit corrections to the drift velocity - linear approximation in the z and global y
+  //The transfromatiom from the drift time to the z position done in AliTPCTracnsform class
+  // 
+  /*
+    Formulas:
+    
+    z  = s* (z0 - vd*(t-t0))
+    
+    s  - side -1 and +1 
+    t0 - time 0
+    vd - nominal drift velocity
+    zs - miscalibrated position
+    
+    zs = s*(z0 - vd*(1+vr)*(t-(t0+dt))
+    vr  - relative change of the drift velocity
+    dzt - vd*dt
+    dr  = zz0-s*z
+    ..
+    ==>
+    zs ~ z - s*vr*(z0-s*z)+s*dzt
+    --------------------------------
+    1. Correction function vr constant:
+    
+    
+    dz = zs-z = -s*vr *(z0-s*z)+s*dzt         
+    dzs/dl = dz/dl +s*s*vr*dz/dl 
+    d(dz/dl) = vr*dz/dl     
+  */
+
+  const Float_t kZCut        = 240;    // remove the closest laser beam
   const Float_t kSaturCut    = 0.05;   // remove saturated lasers - cut on fraction of saturated 
+  const Float_t kDistCut     = 3;    // distance sigma cut
+  const Float_t kDistCutAbs  = 0.25;  
   const Float_t kMinClusters = 60;  // minimal amount of the clusters
   const Float_t kMinSignal   = 16;  // minimal mean height of the signal
   const Float_t kChi2Cut     = 0.1; // chi2 cut to accept drift fit
   static TLinearFitter fdriftA(3,"hyp2");
   static TLinearFitter fdriftC(3,"hyp2");
-  static TLinearFitter fdriftA1(2,"hyp1");
-  static TLinearFitter fdriftC1(2,"hyp1");
-  TVectorD fitA(3),fitC(3);
-
+  static TLinearFitter fdriftAC(4,"hyp3");
+  TVectorD fitA(3),fitC(3),fitAC(4);
+  
+  AliTPCcalibDB*  calib=AliTPCcalibDB::Instance();
+  AliTPCParam  * tpcparam    = calib->GetParameters(); 
+   
+  //
   for (Int_t id=0; id<336; id++) fFitZ[id]=0;
 
+  Float_t chi2A = 10;
+  Float_t chi2C = 10;
+  Float_t chi2AC = 10;
+  Int_t npointsA=0;
+  Int_t npointsC=0;
+  Int_t npointsAC=0;
+
 
-  for (Int_t iter=0; iter<2; iter++){
+  for (Int_t iter=0; iter<3; iter++){
     fdriftA.ClearPoints();
     fdriftC.ClearPoints();
-    fdriftA1.ClearPoints();
-    fdriftC1.ClearPoints();
+    fdriftAC.ClearPoints();
     //
     for (Int_t id=0; id<336; id++){
       if (!fTracksEsdParam.At(id)) continue;
       if (!AcceptLaser(id)) continue;
-      if (TMath::Abs(fFitZ[id])>kDistCut-iter) continue;
       if ( fClusterSatur[id]>kSaturCut)  continue;
       if ( fClusterCounter[id]<kMinClusters)  continue;
-
+      AliESDtrack   *track    = (AliESDtrack*)fTracksEsd.At(id);
+      if (track->GetTPCsignal()<kMinSignal) continue;
       AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
       AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
+
       Double_t xyz[3];
       Double_t pxyz[3];
       Double_t lxyz[3];
@@ -594,33 +665,69 @@ void AliTPCcalibLaser::FitDriftV(){
       ltrp->GetXYZ(lxyz);
       ltrp->GetPxPyPz(lpxyz);
       if (TMath::Abs(lxyz[2])>kZCut) continue;
-      Double_t xxx[2] = {lxyz[2],lxyz[1]};
+      Float_t sz = (ltrp->GetSide()==0) ? TMath::Sqrt(chi2A): TMath::Sqrt(chi2C);
+      if (npointsAC>0) sz =TMath::Sqrt(chi2AC); 
+      if (TMath::Abs(fFitZ[id])>sz*kDistCut) continue;
+      if (iter>0 && TMath::Abs(fFitZ[id])>kDistCutAbs) continue;
+
+      // drift distance
+      Double_t zlength =  tpcparam->GetZLength(0);
+      Double_t ldrift = zlength-TMath::Abs(lxyz[2]);
+      Double_t mdrift = zlength-TMath::Abs(xyz[2]);
+      Double_t xxx[2] = {ldrift,lxyz[1]*ldrift/(zlength*250.)};
       if (ltrp->GetSide()==0){
-       fdriftA1.AddPoint(xxx,xyz[2],1);
-       fdriftA.AddPoint(xxx,xyz[2],1);
+       fdriftA.AddPoint(xxx,mdrift,1);
       }else{
-       fdriftC1.AddPoint(xxx,xyz[2],1);
-       fdriftC.AddPoint(xxx,xyz[2],1);
+       fdriftC.AddPoint(xxx,mdrift,1);
       }
+      Double_t xxx2[3] = {ltrp->GetSide(), ldrift,lxyz[1]*ldrift/(zlength*250.)};
+      fdriftAC.AddPoint(xxx2,mdrift,1);
     }
-    Float_t chi2A = 0;
-    Float_t chi2C = 0;
-    Int_t npointsA=0;
-    Int_t npointsC=0;
     //
     if (fdriftA.GetNpoints()>10){
+      //
       fdriftA.Eval();
+      if (iter==0) fdriftA.EvalRobust(0.7);
+      else fdriftA.EvalRobust(0.8);
       fdriftA.GetParameters(fitA);
       npointsA= fdriftA.GetNpoints();
       chi2A = fdriftA.GetChisquare()/fdriftA.GetNpoints();
-      if (chi2A<kChi2Cut ||(*fFitAside)[0]==0 ) (*fFitAside) = fitA;
+      if (chi2A<kChi2Cut ||(*fFitAside)[0]==0 ) {
+       if (fFitAside->GetNoElements()<5) fFitAside->ResizeTo(5);
+       (*fFitAside)[0] = fitA[0];
+       (*fFitAside)[1] = fitA[1];
+       (*fFitAside)[2] = fitA[2];
+       (*fFitAside)[3] = fdriftA.GetNpoints();
+       (*fFitAside)[4] = chi2A;        
+      }
     }
     if (fdriftC.GetNpoints()>10){
       fdriftC.Eval();
+      if (iter==0) fdriftC.EvalRobust(0.7);
+      else fdriftC.EvalRobust(0.8);
+      //
       fdriftC.GetParameters(fitC);
       npointsC= fdriftC.GetNpoints();
       chi2C = fdriftC.GetChisquare()/fdriftC.GetNpoints();
-      if (chi2C<kChi2Cut||(*fFitCside)[0]==0) (*fFitCside) = fitC;
+      if (chi2C<kChi2Cut||(*fFitCside)[0]==0) {        
+       if (fFitCside->GetNoElements()<5) fFitCside->ResizeTo(5);
+       (*fFitCside)[0] = fitC[0];
+       (*fFitCside)[1] = fitC[1];
+       (*fFitCside)[2] = fitC[2];
+       (*fFitCside)[3] = fdriftC.GetNpoints();
+       (*fFitCside)[4] = chi2C;        
+      }
+    }
+
+    if (fdriftAC.GetNpoints()>10&&fdriftC.GetNpoints()>10&&fdriftA.GetNpoints()>10){
+      fdriftAC.Eval();
+      if (iter==0) fdriftAC.EvalRobust(0.7);
+      else fdriftAC.EvalRobust(0.8);
+      //
+      fdriftAC.GetParameters(fitAC);
+      npointsAC= fdriftAC.GetNpoints();
+      chi2AC = fdriftAC.GetChisquare()/fdriftAC.GetNpoints();
+      if (chi2AC<kChi2Cut||(*fFitACside)[0]==0) (*fFitACside) = fitAC;
     }
     
     for (Int_t id=0; id<336; id++){
@@ -635,28 +742,66 @@ void AliTPCcalibLaser::FitDriftV(){
       param->GetXYZ(xyz);
       param->GetPxPyPz(pxyz);
       ltrp->GetXYZ(lxyz);
-      ltrp->GetPxPyPz(lpxyz);
+      ltrp->GetPxPyPz(lpxyz); 
+      Double_t zlength =  tpcparam->GetZLength(0);
+      Double_t ldrift = zlength-TMath::Abs(lxyz[2]);
+      Double_t mdrift = zlength-TMath::Abs(xyz[2]);
+
       Float_t fz =0;
       if (ltrp->GetSide()==0){
-       fz = (fitA)[0]+(fitA)[1]*lxyz[2]+(fitA)[2]*lxyz[1];
+       fz = (fitA)[0]+(fitA)[1]*ldrift+(fitA)[2]*lxyz[1]*ldrift/(zlength*250.);
       }else{
-       fz = (fitC)[0]+(fitC)[1]*lxyz[2]+(fitC)[2]*lxyz[1];     
+       fz = (fitC)[0]+(fitC)[1]*ldrift+(fitC)[2]*lxyz[1]*ldrift/(zlength*250.);        
+      }
+      if (npointsAC>10){
+       fz = (fitAC)[0]+(fitAC)[1]*ltrp->GetSide()+(fitAC)[2]*ldrift+(fitAC)[3]*lxyz[1]*ldrift/(zlength*250.);
       }
-      fFitZ[id]=xyz[2]-fz;
+      fFitZ[id]=mdrift-fz;
     }
     if (fStreamLevel>0){
       TTreeSRedirector *cstream = GetDebugStreamer();
-      Int_t time = fESD->GetTimeStamp();
+      TTimeStamp tstamp(fTime);
+      Float_t valuePressure0  = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
+      Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
+      Double_t ptrelative0   = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
+      Double_t ptrelative1   = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
+      Double_t temp0         = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
+      Double_t temp1         = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
+      TVectorD vecGoofie(20);
+      AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
+      if (goofieArray) 
+       for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
+         AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
+         if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
+       }
+
       if (cstream){
-       (*cstream)<<"driftv"<<      
+       (*cstream)<<"driftv"<<
+         "run="<<fRun<<              //  run number
+         "event="<<fEvent<<          //  event number
+         "time="<<fTime<<            //  time stamp of event
+         "trigger="<<fTrigger<<      //  trigger
+         "mag="<<fMagF<<             //  magnetic field
+         // Environment values
+         "press0="<<valuePressure0<<
+         "press1="<<valuePressure1<<
+         "pt0="<<ptrelative0<<
+         "pt1="<<ptrelative1<<
+         "temp0="<<temp0<<
+         "temp1="<<temp1<<
+         "vecGoofie.="<<&vecGoofie<<
+         //
+         //
          "iter="<<iter<<
          "driftA.="<<fFitAside<<
          "driftC.="<<fFitCside<<
+         "driftAC.="<<fFitACside<<
          "chi2A="<<chi2A<<
          "chi2C="<<chi2C<<
+         "chi2AC="<<chi2AC<<
          "nA="<<npointsA<<
          "nC="<<npointsC<<
-         "time="<<time<<
+         "nAC="<<npointsAC<<
          "\n";
       }
     }
@@ -712,7 +857,7 @@ Bool_t  AliTPCcalibLaser::AcceptLaser(Int_t id){
   */
   AliExternalTrackParam *param =(AliExternalTrackParam*)fTracksEsdParam.At(id);
   AliTPCLaserTrack *ltrp  = ( AliTPCLaserTrack*)fTracksMirror.At(id);
-  AliESDtrack   *track    = (AliESDtrack*)fTracksEsd.At(id);
+  //AliESDtrack   *track    = (AliESDtrack*)fTracksEsd.At(id);
   Double_t xyz[3];
   Double_t lxyz[3];
   param->GetXYZ(xyz);
@@ -732,9 +877,11 @@ Int_t  AliTPCcalibLaser::FindMirror(AliESDtrack *track, AliTPCseed *seed){
   //
   // Find corresponding mirror
   // add the corresponding tracks
-  //
-  Float_t kRadius0 = 252;
-  Float_t kRadius  = 253.4;
+
+  Float_t kRadius0  = 252;
+  Float_t kRadius   = 253.4;
+
   if (!track->GetOuterParam()) return -1;
   AliExternalTrackParam param(*(track->GetOuterParam()));
   AliTracker::PropagateTrackTo(&param,kRadius0,0.10566,3,kTRUE);
@@ -776,6 +923,7 @@ Int_t  AliTPCcalibLaser::FindMirror(AliESDtrack *track, AliTPCseed *seed){
   //
   Float_t radius=TMath::Abs(ltrp->GetX());
   AliTracker::PropagateTrackTo(&param,radius,0.10566,0.01,kTRUE);
+  param.Rotate(ltrp->GetAlpha());
   //
   if (!fTracksMirror.At(id)) fTracksMirror.AddAt(ltrp,id);
   Bool_t accept=kTRUE;  
@@ -829,7 +977,13 @@ void AliTPCcalibLaser::DumpLaser(Int_t id) {
     Bool_t accept = AcceptLaser(id);
     if (cstream){
       (*cstream)<<"Track"<<
-       "run="<<fRun<<
+       //
+       "run="<<fRun<<              //  run number
+       "event="<<fEvent<<          //  event number
+       "time="<<fTime<<            //  time stamp of event
+       "trigger="<<fTrigger<<      //  trigger
+       "mag="<<fMagF<<             //  magnetic field
+
        "id="<<id<<
        "accept="<<accept<<
        "driftA.="<<fFitAside<<
@@ -876,6 +1030,19 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
   // store x, y, z information for all models and the cluster to calculate the residuals
   //
   
+  // The clusters which do not fulfill given  criteria are skipped
+  //
+  // Cluters removed from fit
+  //  1. Extended shape     > kShapeCut
+  //  2. In saturation  Max > kMax
+  //  3. Distance to edge   < cutEdge
+  //
+  // Clusters not used for the calibration:
+  //
+  //  1. Extended shape     > kShapeCut
+  //  2. In saturation  Max > kMax
+
+
   AliTPCseed *track      = (AliTPCseed*)fTracksTPC.At(id);
   AliExternalTrackParam *extparam=(AliExternalTrackParam*)fTracksEsdParam.At(id);
   AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id);
@@ -883,15 +1050,23 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
   AliTPCclusterMI dummyCl;
 
   //two tracklets
-  Int_t kMaxTracklets=2;
+  Int_t kMaxTracklets=2;  
+  Float_t kShapeCut = 1.3;
+  Float_t kRatioCut = 2.;
+
+  Float_t kMax      = 900.;
+
+
   //=============================================//
   // Linear Fitters for the Different Approaches //
   //=============================================//
-  //linear fit model in y and z; inner - outer sector
+  //linear fit model in y and z; inner - outer sector, combined with offset
   static TLinearFitter fy1I(2,"hyp1");
   static TLinearFitter fy1O(2,"hyp1");
   static TLinearFitter fz1I(2,"hyp1");
   static TLinearFitter fz1O(2,"hyp1");
+  static TLinearFitter fy1IO(3,"hyp2");
+  static TLinearFitter fz1IO(3,"hyp2");
   //quadratic fit model in y and z; inner - sector
   static TLinearFitter fy2I(3,"hyp2");
   static TLinearFitter fy2O(3,"hyp2");
@@ -913,8 +1088,12 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
   // Loop over all Tracklet Cuts //
   //=============================//
   for (Int_t icut=0; icut<fNcuts; icut++){
-    Float_t xinMin = 250, xinMax=90;
-    Float_t xoutMin=250, xoutMax=90;
+    Float_t xinMin = 2500, xinMax=-90;
+    Float_t xoutMin=2500, xoutMax=-90;
+    Float_t msigmaYIn=0;
+    Float_t msigmaYOut=0;
+    Float_t mqratioIn=0;
+    Float_t mqratioOut=0;
 
 
       AliDebug(4,Form("Processing cut %d for track with ID %d",icut,id));
@@ -932,12 +1111,19 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
       TVectorD vecy1resOuter(2),vecz1resOuter(2); //pol1 fit parameters outer
       TVectorD vecy2resOuter(3),vecz2resOuter(3); //pol2 fit parameters outer
       TVectorD vecy4res(5),vecz4res(5);
+      TVectorD vecy1resIO(3),vecz1resIO(3);
       // cluster and track positions for each row - used for residuals
+      TVectorD vecgX(159);        // global X
+      TVectorD vecgY(159);        // global Y
+      TVectorD vecgZ(159);        // global Z
+
       TVectorD vecX(159);        // x is the same for all (row center)
       TVectorD vecYkalman(159);  // y from kalman fit
       TVectorD vecZkalman(159);  // z from kalman fit
       TVectorD vecY1(159);       // y from pol1 fit per ROC
       TVectorD vecZ1(159);       // z from pol1 fit per ROC
+      TVectorD vecY1IO(159);     // y from pol1 fit per ROC
+      TVectorD vecZ1IO(159);     // z from pol1 fit per ROC
       TVectorD vecY2(159);       // y from pol2 fit per ROC
       TVectorD vecZ2(159);       // z from pol2 fit per ROC
       TVectorD vecY4(159);       // y from sector fit
@@ -945,11 +1131,14 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
       TVectorD vecClY(159);      // y cluster position
       TVectorD vecClZ(159);      // z cluster position
       TVectorD vecSec(159);      // sector for each row
+      TVectorD  isReject(159);    // flag - cluster to be rejected 
       //chi2 of fits
       Double_t chi2I1z=0;       // chi2 of pol1 fit in z (inner)
       Double_t chi2I1y=0;       // chi2 of pol1 fit in y (inner)
       Double_t chi2O1z=0;       // chi2 of pol1 fit in z (outer)
       Double_t chi2O1y=0;       // chi2 of pol1 fit in y (outer)
+      Double_t chi2IO1z=0;       // chi2 of pol1 fit in z (outer)
+      Double_t chi2IO1y=0;       // chi2 of pol1 fit in y (outer)
       Double_t chi2I2z=0;       // chi2 of pol2 fit in z (inner)
       Double_t chi2I2y=0;       // chi2 of pol2 fit in y (inner)
       Double_t chi2O2z=0;       // chi2 of pol2 fit in z (outer)
@@ -1015,6 +1204,10 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
       // fill fitters with cluster information //
       //=======================================//
       AliDebug(3,"Filing Cluster Information");
+
+      //
+      //
+
       for (Int_t irow=158;irow>-1;--irow) {
          AliTPCclusterMI *c=track->GetClusterPointer(irow);
          AliTPCclusterMI &cl = (AliTPCclusterMI&) (*arrCl[irow]);
@@ -1026,13 +1219,13 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
          for (Int_t drow=-1;drow<=1;drow++) {
            if (irow+drow<0) continue;
            if (irow+drow>158) continue;
-           AliTPCclusterMI *c=track->GetClusterPointer(irow);
-           if (!c) continue;
-           Int_t roc = static_cast<Int_t>(c->GetDetector());
+           AliTPCclusterMI *ccurrent=track->GetClusterPointer(irow);
+           if (!ccurrent) continue;
+           Int_t roc = static_cast<Int_t>(ccurrent->GetDetector());
            if ( roc!=innerSector && roc!=outerSector ) continue;
-           if (c->GetX()<10) continue;
-           if (c->GetY()==0) continue;
-           meanY+=c->GetY();
+           if (ccurrent->GetX()<10) continue;
+           if (ccurrent->GetY()==0) continue;
+           meanY+=ccurrent->GetY();
            sumY++;
          }
          if (sumY>0)  meanY/=sumY; 
@@ -1040,7 +1233,7 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
           //
           vecSec[irow]=-1;
          if (!c) continue;
-         Double_t pedgeY = c->GetX()*TMath::DegToRad()*(10)-TMath::Abs(meanY);
+         Double_t pedgeY = c->GetX()*TMath::Tan(TMath::DegToRad()*(10))-TMath::Abs(meanY);
          Double_t pedgeX = TMath::Min((irow)*0.75, (159.-irow)*1.5);
          
           //
@@ -1056,6 +1249,12 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
          vecX[irow]   = c->GetX();
          vecClY[irow] = c->GetY();
          vecClZ[irow] = c->GetZ();
+         //
+//       Float_t gxyz[3];
+//       c->GetGlobalXYZ(gxyz);
+//       vecgX[irow]   = gxyz[0];
+//       vecgY[irow]   = gxyz[1];
+//       vecgZ[irow]   = gxyz[2];
           //
           Double_t x = vecX[irow]-133.4; //reference is between IROC and OROC
           Double_t y = vecClY[irow];
@@ -1063,9 +1262,11 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
           //
          Double_t x2[2]={x,x*x};   //linear and parabolic parameters
          Double_t x4[4]={0,0,0,0}; //hyp4 parameters
+          Double_t xIO[2]={0,x};    //common linear + offset IROC-OROC
          if ( roc == innerSector ) {
              x4[0]=1; //offset inner - outer sector
              x4[1]=x; //slope parameter inner sector
+              xIO[0]=1;
          } else {
              x4[2]=x; //slope parameter outer sector
          }
@@ -1080,8 +1281,12 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              fy2I.AddPoint(x2,y);
              fz2I.AddPoint(x2,z);
               ++nclI;
-             if (x<xinMin) xinMin=x;
-             if (x>xinMax) xinMax=x;
+             if (c->GetX()<xinMin) xinMin=c->GetX();
+             if (c->GetX()>xinMax) xinMax=c->GetX();
+
+             msigmaYIn +=TMath::Sqrt(c->GetSigmaY2());
+             mqratioIn +=c->GetMax()/c->GetQ();
+
          }
          if ( roc==outerSector ){
              fy1O.AddPoint(x2,y);
@@ -1089,11 +1294,24 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              fy2O.AddPoint(x2,y);
              fz2O.AddPoint(x2,z);
               ++nclO;
-             if (x<xoutMin) xoutMin=x;
-             if (x>xoutMax) xoutMax=x;
+             if (c->GetX()<xoutMin) xoutMin=c->GetX();
+             if (c->GetX()>xoutMax) xoutMax=c->GetX();
+             msigmaYOut +=TMath::Sqrt(c->GetSigmaY2());
+             mqratioOut +=c->GetMax()/c->GetQ();
+
          }
          fy4.AddPoint(x4,y);
          fz4.AddPoint(x4,z);
+          fy1IO.AddPoint(xIO,y);
+          fz1IO.AddPoint(xIO,z);
+      }
+      if (nclI>0)  {
+       msigmaYIn/=nclI;
+       mqratioIn/=nclI;
+      }
+      if (nclO>0)  {
+       msigmaYOut/=nclO;
+       mqratioOut/=nclO;
       }
       //======================================//
       // Evaluate and retrieve fit parameters //
@@ -1141,6 +1359,16 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              fz4.GetParameters(vecz4res);
              chi2IOz=fz4.GetChisquare()/(fz4.GetNpoints()-5);
          }
+          if (fy1IO.GetNpoints()>0) {
+            fy1IO.Eval();
+            fy1IO.GetParameters(vecy1resIO);
+            chi2IO1y=fy1IO.GetChisquare()/(fy1IO.GetNpoints()-3);
+          }
+          if (fz1IO.GetNpoints()>0) {
+            fz1IO.Eval();
+            fz1IO.GetParameters(vecz1resIO);
+            chi2IO1z=fz1IO.GetChisquare()/(fz1IO.GetNpoints()-3);
+          }
       }
       //clear points
       fy4.ClearPoints();  fz4.ClearPoints();
@@ -1148,11 +1376,35 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
       fz1I.ClearPoints(); fz1O.ClearPoints();
       fy2I.ClearPoints(); fy2O.ClearPoints();
       fz2I.ClearPoints(); fz2O.ClearPoints();
+      fy1IO.ClearPoints(); fz1IO.ClearPoints();
       //==============================//
       // calculate tracklet positions //
       //==============================//
       AliDebug(4,"Calculate tracklet positions");
       for (Int_t irow=158;irow>-1;--irow) {
+       isReject[irow]=0;
+        AliTPCclusterMI *c=track->GetClusterPointer(irow);
+        if ( vecSec[irow]!=innerSector && vecSec[irow]!=outerSector ) { // no cluster in given sectors
+          isReject[irow]+=1;
+        }
+
+        if (!c) { //no cluster    
+          isReject[irow]+=2;
+        }else{
+          if (c->GetMax()>kMax){   //saturation
+            isReject[irow]+=4;
+          }
+          if ( vecSec[irow] == outerSector ) {   //extended shape
+            if (c->GetMax()/c->GetQ()> mqratioOut*kRatioCut)                isReject[irow]+=8;
+            if (TMath::Sqrt(c->GetSigmaY2())>msigmaYOut*kShapeCut)    isReject[irow]+=16;
+          }else{
+            if (c->GetMax()/c->GetQ()> mqratioIn*kRatioCut)                isReject[irow]+=8;
+            if (TMath::Sqrt(c->GetSigmaY2())>msigmaYIn*kShapeCut)    isReject[irow]+=16;
+          }       
+        }
+
+
+
          if ( vecSec[irow]==-1 ) continue;  //no cluster info
           if ( vecSec[irow]!=innerSector && vecSec[irow]!=outerSector ) continue;
          tr=&dummy;
@@ -1161,6 +1413,8 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
           //
           Double_t yoffInner=0;
           Double_t zoffInner=0;
+          Double_t yoffInner1=0;
+          Double_t zoffInner1=0;
           Double_t yslopeInner=0;
           Double_t yslopeOuter=0;
           Double_t zslopeInner=0;
@@ -1181,10 +1435,14 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
               vecY2[irow]=vecy2resInner[0]+vecy2resInner[1]*xref+vecy2resInner[2]*xref*xref;
               vecZ2[irow]=vecz2resInner[0]+vecz2resInner[1]*xref+vecz2resInner[2]*xref*xref;
               yoffInner=vecy4res[1];
-             zoffInner=vecz4res[1];
+              zoffInner=vecz4res[1];
+              yoffInner1=vecy1resIO[1];
+              zoffInner1=vecz1resIO[1];
               yslopeInner=vecy4res[2];
              zslopeInner=vecz4res[2];
          }
+          vecY1IO[irow]=vecy1resIO[0]+yoffInner1+vecy1resIO[2]*xref;
+          vecZ1IO[irow]=vecz1resIO[0]+zoffInner1+vecz1resIO[2]*xref;
          vecY4[irow]=vecy4res[0]+yoffInner+yslopeInner*xref+yslopeOuter*xref+vecy4res[4]*xref*xref;
          vecZ4[irow]=vecz4res[0]+zoffInner+zslopeInner*xref+zslopeOuter*xref+vecz4res[4]*xref*xref;
           //positions of kalman fits
@@ -1200,6 +1458,10 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              vecYkalman[irow]=xyz[1];
              vecZkalman[irow]=xyz[2];
          }
+         //
+         //
+         //
+        
       }
       //=====================================================================//
       // write results from the different tracklet fits with debug streamers //
@@ -1209,6 +1471,12 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
          if (cstream){
            Float_t dedx = track->GetdEdx();
            (*cstream)<<"FitModels"<<
+             "run="<<fRun<<              //  run number
+             "event="<<fEvent<<          //  event number
+             "time="<<fTime<<            //  time stamp of event
+             "trigger="<<fTrigger<<      //  trigger
+             "mag="<<fMagF<<             //  magnetic field          
+             //
              "cutNr="      << icut <<
              "edgeCutX="   << edgeCutX <<
              "edgeCutY="   << edgeCutY <<
@@ -1218,9 +1486,11 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              "dEdx="       << dedx <<
              "LTr.="       << ltrp <<
              "Tr.="        << extparam <<
-             "yPol1In.="   << &vecy1resInner <<
-             "zPol1In.="   << &vecz1resInner <<
-             "yPol2In.="   << &vecy2resInner <<
+              "yPol1In.="   << &vecy1resInner <<
+              "zPol1In.="   << &vecz1resInner <<
+              "yPol1InOut.="<< &vecy1resIO <<
+              "zPol1InOut.="<< &vecz1resIO <<
+              "yPol2In.="   << &vecy2resInner <<
              "zPol2In.="   << &vecz2resInner <<
              "yPol1Out.="  << &vecy1resOuter <<
              "zPol1Out.="  << &vecz1resOuter <<
@@ -1228,9 +1498,11 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              "zPol2Out.="  << &vecz2resOuter <<
              "yInOut.="    << &vecy4res <<
              "zInOut.="    << &vecz4res <<
-             "chi2y1In="   << chi2I1y <<
-             "chi2z1In="   << chi2I1z <<
-             "chi2y1Out="  << chi2O1y <<
+              "chi2y1In="   << chi2I1y <<
+              "chi2z1In="   << chi2I1z <<
+              "chi2y1InOut="<< chi2IO1y <<
+              "chi2z1InOut="<< chi2IO1z <<
+              "chi2y1Out="  << chi2O1y <<
              "chi2z1Out="  << chi2O1z <<
              "chi2y2In="   << chi2I2y <<
              "chi2z2In="   << chi2I2z <<
@@ -1246,6 +1518,10 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              "xinMax="     << xinMax<<
              "xoutMin="    << xoutMin<<
              "xoutMax="     << xoutMax<<
+             "msigmaYIn=" <<msigmaYIn<<
+             "msigmaYOut=" <<msigmaYOut<<
+             "mqratioIn=" <<mqratioIn<<
+             "mqratioOut=" <<  mqratioOut      <<
              "\n";
          }
       }
@@ -1256,6 +1532,12 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
          if (cstream){
            Float_t dedx = track->GetdEdx();
            (*cstream)<<"Residuals"<<
+             "run="<<fRun<<              //  run number
+             "event="<<fEvent<<          //  event number
+             "time="<<fTime<<            //  time stamp of event
+             "trigger="<<fTrigger<<      //  trigger
+             "mag="<<fMagF<<             //  magnetic field          
+             //
              "cutNr="      << icut <<
              "edgeCutX="   << edgeCutX <<
              "edgeCutY="   << edgeCutY <<
@@ -1264,15 +1546,21 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              "Tr.="        << extparam<<
              "dEdx="       << dedx <<
              "Cl.="        << &arrCl <<
+             "vX.="        << &vecgX<<   // global x
+             "vY.="        << &vecgY<<   // global y
+             "vZ.="        << &vecgZ<<   // global z
              "TrX.="       << &vecX <<
              "TrYpol1.="   << &vecY1 <<
              "TrZpol1.="   << &vecZ1 <<
              "TrYpol2.="   << &vecY2 <<
              "TrZpol2.="   << &vecZ2 <<
-             "TrYInOut.="  << &vecY4 <<
-             "TrZInOut.="  << &vecZ4 <<
-             "ClY.="       << &vecClY <<
+              "TrYpol1InOut.="<< &vecY1IO <<
+              "TrZpol1InOut.="<< &vecZ1IO <<
+              "TrYInOut.="  << &vecY4 <<
+              "TrZInOut.="  << &vecZ4 <<
+              "ClY.="       << &vecClY <<
              "ClZ.="       << &vecClZ <<
+             "isReject.="  << &isReject<<
              "sec.="       << &vecSec <<
              "nclI="       << nclI <<
              "nclO="       << nclO <<
@@ -1280,18 +1568,37 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
              "xinMax="     << xinMax<<
              "xoutMin="    << xoutMin<<
              "xoutMax="     << xoutMax<<
+             "msigmaYIn=" <<msigmaYIn<<
+             "msigmaYOut=" <<msigmaYOut<<
+             "mqratioIn=" <<mqratioIn<<
+             "mqratioOut=" <<  mqratioOut      <<
              "yInOut.="    << &vecy4res <<
              "zInOut.="    << &vecz4res <<
+             //chi2s
              "chi2y1In="   << chi2I1y <<    //
              "chi2z1In="   << chi2I1z <<
              "chi2y1Out="  << chi2O1y <<
              "chi2z1Out="  << chi2O1z <<
-             "chi2y2In="   << chi2I2y <<
+              "chi2y1InOut="<< chi2IO1y <<
+              "chi2z1InOut="<< chi2IO1z <<
+              "chi2y2In="   << chi2I2y <<
              "chi2z2In="   << chi2I2z <<
              "chi2y2Out="  << chi2O2y <<
              "chi2z2Out="  << chi2O2z <<
              "chi2yInOut=" << chi2IOy <<
              "chi2zInOut=" << chi2IOz <<
+             // fit parameters
+             "yPol1In.="   << &vecy1resInner <<
+             "zPol1In.="   << &vecz1resInner <<
+             "yPol2In.="   << &vecy2resInner <<
+             "zPol2In.="   << &vecz2resInner <<
+             "yPol1Out.="  << &vecy1resOuter <<
+             "zPol1Out.="  << &vecz1resOuter <<
+              "yPol1InOut.="<< &vecy1resIO <<
+              "zPol1InOut.="<< &vecz1resIO <<
+              "yPol2Out.="  << &vecy2resOuter <<
+             "zPol2Out.="  << &vecz2resOuter <<
+
              "\n";
 
          }
@@ -1299,39 +1606,45 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
       //==========================//
       // Fill Residual Histograms //
       //==========================//
-      TProfile *profy = (TProfile*)fDeltaYres.UncheckedAt(id);
-      TProfile *profz = (TProfile*)fDeltaZres.UncheckedAt(id);
-      if (!profy){
-       profy=new TProfile(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d",id),160,0,160);
-       profy->SetDirectory(0);
-       fDeltaYres.AddAt(profy,id);
-      }
-      if (!profz){
-       profz=new TProfile(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d",id),160,0,160);
-       profz->SetDirectory(0);
-       fDeltaZres.AddAt(profz,id);
-      }
+      if (!fHisNclIn) MakeFitHistos(); 
+
+      TH2F *profy = (TH2F*)fDeltaYres.UncheckedAt(id);
+      TH2F *profz = (TH2F*)fDeltaZres.UncheckedAt(id);
+      TH2F *profy2 = (TH2F*)fDeltaYres2.UncheckedAt(id);
+      TH2F *profz2 = (TH2F*)fDeltaZres2.UncheckedAt(id);
+      //      TH2F *profy3 = (TH2F*)fDeltaYres3.UncheckedAt(id);
+      //TH2F *profz3 = (TH2F*)fDeltaZres3.UncheckedAt(id);
+      //
       for (Int_t irow=158;irow>-1;--irow) {
-       if (vecSec[irow]==-1)continue; //no cluster info
-       Double_t x    = vecX[irow];
+       if (vecSec[irow]==-1)continue; // no cluster info
+       if (isReject[irow]>0.5) continue;  // 
+       //Double_t x    = vecX[irow];
        Double_t ycl  = vecClY[irow];
        Double_t yfit = vecY1[irow];
+       Double_t yfit2 = vecY2[irow];
+       //      Double_t yfit3 = vecY2[irow];
        Double_t zcl  = vecClZ[irow];
        Double_t zfit = vecZ1[irow];
+       Double_t zfit2 = vecZ2[irow];
+       //Double_t zfit3 = vecZ2[irow];
+
        if (TMath::Abs(yfit-ycl)<2&&TMath::Abs(zfit-zcl)<2){
-         if (profy) 
-           if (profy->GetEntries()<1000000) 
-             profy->Fill(irow,yfit-ycl);
-         if (profz) 
-           if (profz->GetEntries()<1000000) 
-             profz->Fill(irow,zfit-zcl);
+         if (profy){       
+             profy->Fill(irow,ycl-yfit);
+             profy2->Fill(irow,ycl-yfit2);
+             //              profy3->Fill(irow,ycl-yfit3);
+         }
+         if (profz) {
+             profz->Fill(irow,zcl-zfit);
+             profz2->Fill(irow,zcl-zfit2);
+             //profz3->Fill(irow,zcl-zfit3);
+         }
        }
       }
       //
       //
       // Fill laser fit histograms
       //
-      if (!fHisNclIn) MakeFitHistos(); 
       Float_t dedx = track->GetdEdx();
       if (nclI>20&&nclO>20){
        fHisNclIn->Fill(id,nclI);      //->Number of clusters inner
@@ -1393,7 +1706,7 @@ void AliTPCcalibLaser::RefitLaserJW(Int_t id){
 
 
 
-void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries){
+void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run){
   //
   //  Dump information about laser beams
   //  isOK variable indicates usability of the beam  
@@ -1423,7 +1736,7 @@ void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries)
     TH1F * hisP4    = (TH1F*)laser->fDeltaP4.At(id);
     TH1F * hisS    = (TH1F*)laser->fSignals.At(id);
     //if (!hisphi) continue;
-    Double_t entries = hisphi->GetEntries();
+    Double_t entries = (hisphi==0)? 0: hisphi->GetEntries();
     //if (entries<minEntries) continue;
     //
     AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id);
@@ -1445,22 +1758,27 @@ void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries)
     Float_t rmsphiP = hisphiP->GetRMS();
     Float_t meanZ = hisZ->GetMean();
     Float_t rmsZ = hisZ->GetRMS();
-    hisphi->Fit(&fg,"","",hisphi->GetMean()-4*hisphi->GetRMS(),hisphi->GetMean()+4*hisphi->GetRMS());
+    if (hisphi->GetRMS()>0)
+      hisphi->Fit(&fg,"","",hisphi->GetMean()-4*hisphi->GetRMS(),hisphi->GetMean()+4*hisphi->GetRMS());
     Double_t gphi1 = fg.GetParameter(1);
     Double_t gphi2 = fg.GetParameter(2);
-    hisphiP->Fit(&fg,"","",hisphiP->GetMean()-4*hisphiP->GetRMS(),hisphiP->GetMean()+4*hisphiP->GetRMS());
+    if (hisphiP->GetRMS()>0)
+      hisphiP->Fit(&fg,"","",hisphiP->GetMean()-4*hisphiP->GetRMS(),hisphiP->GetMean()+4*hisphiP->GetRMS());
     Double_t gphiP1 = fg.GetParameter(1);
     Double_t gphiP2 = fg.GetParameter(2);
     //
-    hisZ->Fit(&fg,"","",hisZ->GetMean()-4*hisZ->GetRMS(),hisZ->GetMean()+4*hisZ->GetRMS());
+    if (hisZ->GetRMS()>0)
+      hisZ->Fit(&fg,"","",hisZ->GetMean()-4*hisZ->GetRMS()-0.1,hisZ->GetMean()+4*hisZ->GetRMS()+0.1);
     Double_t gz1 = fg.GetParameter(1);
     Double_t gz2 = fg.GetParameter(2);
     //
-    hisP3->Fit(&fg,"","",hisP3->GetMean()-4*hisP3->GetRMS(),hisP3->GetMean()+4*hisP3->GetRMS());
+    if (hisP3->GetRMS()>0)
+      hisP3->Fit(&fg,"","",hisP3->GetMean()-4*hisP3->GetRMS(),hisP3->GetMean()+4*hisP3->GetRMS());
     Double_t gp31 = fg.GetParameter(1);
     Double_t gp32 = fg.GetParameter(2);
     //
-    hisP4->Fit(&fg,"","",hisP4->GetMean()-4*hisP4->GetRMS(),hisP4->GetMean()+4*hisP4->GetRMS());
+    if (hisP4->GetRMS()>0)
+      hisP4->Fit(&fg,"","",hisP4->GetMean()-4*hisP4->GetRMS(),hisP4->GetMean()+4*hisP4->GetRMS());
     Double_t gp41 = fg.GetParameter(1);
     Double_t gp42 = fg.GetParameter(2);
     //
@@ -1503,9 +1821,6 @@ void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries)
     Float_t mdEdx =  his->GetMean();
     delete his;
     //
-    his = fHisdZfit->ProjectionY("aaa",id+1,id+1);
-    Float_t mdZfit =  his->GetMean();
-    delete his;
     //
     //
     //
@@ -1542,6 +1857,12 @@ void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries)
     //
     // fit res. histos
     // 
+    his = fHisdZfit->ProjectionY("aaa",id+1,id+1);
+    Float_t  edZfit  = his->GetEntries();
+    Float_t  mdZfit =  his->GetMean();
+    Float_t  rdZfit  = his->GetRMS();
+    delete his;
+
     his = fHisPy1vP0->ProjectionY("aaa",id+1,id+1);     //-> delta y   P0outer-P0inner - line
     Float_t  ePy1vP0  = his->GetEntries();
     Float_t  mPy1vP0  = his->GetMean();
@@ -1653,7 +1974,6 @@ void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries)
     delete his;
     
     //
-    if (run<=0) run=fRun;
     (*pcstream)<<"Mean"<<
       "run="<<run<<               //
       "isOK="<<isOK<<             //
@@ -1669,7 +1989,9 @@ void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries)
       "mLclOut="<<mLclOut<<       // mean number of clusters in outer
       "mLclIO="<<mLclIO<<         // mean number of clusters in inner+outer
       "mdEdx="<<mdEdx<<           // mean dEdx
+      "edZfit="<<edZfit<<           // entries z fit
       "mdZfit="<<mdZfit<<           // mean z fit
+      "rdZfit="<<rdZfit<<           // RMS z fit
       //
       //
       "mChi2YIn1="<<mChi2YIn1<<       //->chi2 y inner - line
@@ -1787,7 +2109,7 @@ void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run, Int_t minEntries)
 
 
 
-void AliTPCcalibLaser::DumpScanInfo(TTree * chain){
+void AliTPCcalibLaser::DumpScanInfo(TTree * chain, const char * cutUser){
   //
   //
   //
@@ -1815,10 +2137,11 @@ void AliTPCcalibLaser::DumpScanInfo(TTree * chain){
   Double_t pphi[3];
   Double_t pphiP[3];
   Double_t pmZ[3];
+  
   //
   for (Int_t id=0; id<336; id++){
     // id =205;
-    sprintf(cut,"isOK&&fId==%d",id);
+    sprintf(cut,"fId==%d&&%s",id,cutUser);
     Int_t entries = chain->Draw("bz",cut,"goff");
     if (entries<3) continue;
     AliTPCLaserTrack *ltrp = 0;
@@ -1836,17 +2159,17 @@ void AliTPCcalibLaser::DumpScanInfo(TTree * chain){
     //
     chain->Draw("gphi1",cut,"goff");
     memcpy(mphi, chain->GetV1(), entries*sizeof(Double_t));
-    chain->Draw("0.05*abs(mphi)+gphi2",cut,"goff");
+    chain->Draw("0.05*abs(mphi)+abs(gphi2)*0.5+0.05",cut,"goff");
     memcpy(smphi, chain->GetV1(), entries*sizeof(Double_t));
     //
     chain->Draw("gphiP1",cut,"goff");
     memcpy(mphiP, chain->GetV1(), entries*sizeof(Double_t));
-    chain->Draw("0.05*abs(mphiP)+gphiP2",cut,"goff");
+    chain->Draw("0.05*abs(mphiP)+abs(gphiP2)*0.5+0.001",cut,"goff");
     memcpy(smphiP, chain->GetV1(), entries*sizeof(Double_t));
     //
     chain->Draw("gz1",cut,"goff");
     memcpy(mZ, chain->GetV1(), entries*sizeof(Double_t));
-    chain->Draw("0.01*abs(meanZ)+gz2",cut,"goff");
+    chain->Draw("0.01*abs(meanZ)+abs(gz2)*0.5+0.1",cut,"goff");
     memcpy(smZ, chain->GetV1(), entries*sizeof(Double_t));
     //
     //
@@ -1855,6 +2178,16 @@ void AliTPCcalibLaser::DumpScanInfo(TTree * chain){
     // store data  
     // phi
     f->cd("dirphi");
+    Float_t phiB0 =0;
+    Float_t phiPB0=0;
+    Float_t zB0=0;
+    for (Int_t ientry=0; ientry<entries; ientry++){
+      if (TMath::Abs(bz[ientry])<0.05){
+       phiB0  = mphi[ientry];
+       phiPB0 = mphiP[ientry];
+       zB0    = mZ[ientry];
+      }
+    }
     TGraphErrors *grphi = new TGraphErrors(entries,bz,mphi,sbz,smphi);
     grphi->Draw("a*");
     grphi->Fit(&fp);
@@ -1904,7 +2237,9 @@ void AliTPCcalibLaser::DumpScanInfo(TTree * chain){
 
     gPad->SaveAs(Form("pic/z/z_%s.gif",grnamefull));
     grmZ->Write();
-    
+    //
+    // P4
+    //
 
     for (Int_t ientry=0; ientry<entries; ientry++){
       (*pcstream)<<"Mean"<<
@@ -1919,19 +2254,23 @@ void AliTPCcalibLaser::DumpScanInfo(TTree * chain){
        "lpx1="<<lpxyz[1]<<          // reference y
        "lpx2="<<lpxyz[2]<<          // refernece z            
        //values
+       "phiB0="<<phiB0<<          // position shift at 0 field
+       "phiPB0="<<phiPB0<<        // angular  shift at 0 field
+       "zB0="<<zB0<<              // z shift for 0 field
+       //
        "gphi1="<<mphi[ientry]<< // mean - from gaus fit
        "pphi0="<<pphi[0]<<   // offset
-       "pphi1="<<pphi[1]<<   // mean
+       "pphi1="<<pphi[1]<<   // slope
        "pphi2="<<pphi[2]<<   // norm chi2
        //
        "gphiP1="<<mphiP[ientry]<< // mean - from gaus fit
        "pphiP0="<<pphiP[0]<< // offset
-       "pphiP1="<<pphiP[1]<< // mean
+       "pphiP1="<<pphiP[1]<< // slope
        "pphiP2="<<pphiP[2]<< // norm chi2
        //
        "gz1="<<mZ[ientry]<<
        "pmZ0="<<pmZ[0]<<     // offset
-       "pmZ1="<<pmZ[1]<<     // mean
+       "pmZ1="<<pmZ[1]<<     // slope
        "pmZ2="<<pmZ[2]<<     // norm chi2
        "\n";
     }
@@ -1953,18 +2292,22 @@ Long64_t AliTPCcalibLaser::Merge(TCollection *li) {
 
   TIterator* iter = li->MakeIterator();
   AliTPCcalibLaser* cal = 0;
-
+  static Int_t counter0=0;
   while ((cal = (AliTPCcalibLaser*)iter->Next())) {
     if (!cal->InheritsFrom(AliTPCcalibLaser::Class())) {
-      Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
+      //      Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
       return -1;
     }
+    printf("Marging number %d\n", counter0);
+    counter0++;
     //
     MergeFitHistos(cal);
     TH1F *h=0x0;
     TH1F *hm=0x0;
-    TProfile *hp=0x0;
-    TProfile *hpm=0x0;
+    TH2F *h2=0x0;
+    TH2F *h2m=0x0;
+    //    TProfile *hp=0x0;
+    //TProfile *hpm=0x0;
 
     for (Int_t id=0; id<336; id++){
       // merge fDeltaZ histograms
@@ -2025,24 +2368,30 @@ Long64_t AliTPCcalibLaser::Merge(TCollection *li) {
       if (hm) h->Add(hm);
       //
       //
-      // merge ProfileY histograms
-      hpm = (TProfile*)cal->fDeltaYres.At(id);
-      hp  = (TProfile*)fDeltaYres.At(id);
-      if (!hp) {
-       hp=new TProfile(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d",id),160,0,160);
-       hp->SetDirectory(0);
-       fDeltaYres.AddAt(hp,id);
-      }
-      if (hpm) hp->Add(hpm);
+      // merge ProfileY histograms -0
+      h2m = (TH2F*)cal->fDeltaYres.At(id);
+      h2  = (TH2F*)fDeltaYres.At(id);
+      if (h2m&&h2) h2->Add(h2m);
       //
-      hpm = (TProfile*)cal->fDeltaZres.At(id);
-      hp  = (TProfile*)fDeltaZres.At(id);
-      if (!hp) {
-       hp=new TProfile(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d",id),160,0,160);
-       hp->SetDirectory(0);
-       fDeltaZres.AddAt(hp,id);
-      }
-      if (hpm) hp->Add(hpm);
+      h2m = (TH2F*)cal->fDeltaZres.At(id);
+      h2  = (TH2F*)fDeltaZres.At(id);
+      if (h2m&&h2) h2->Add(h2m);
+      // merge ProfileY histograms - 2
+      h2m = (TH2F*)cal->fDeltaYres2.At(id);
+      h2  = (TH2F*)fDeltaYres2.At(id);
+      if (h2m&&h2) h2->Add(h2m);
+      //
+      h2m = (TH2F*)cal->fDeltaZres2.At(id);
+      h2  = (TH2F*)fDeltaZres2.At(id);
+      if (h2m&&h2) h2->Add(h2m);
+      // merge ProfileY histograms - 3
+      //h2m = (TH2F*)cal->fDeltaYres3.At(id);
+      //h2  = (TH2F*)fDeltaYres3.At(id);
+      //if (h2m) h2->Add(h2m);
+      //
+      //h2m = (TH2F*)cal->fDeltaZres3.At(id);
+      //h2  = (TH2F*)fDeltaZres3.At(id);
+      //if (h2m) h->Add(h2m);
       //
       //
     }
@@ -2069,7 +2418,7 @@ void   AliTPCcalibLaser::MakeFitHistos(){
   fHisLclIO  = new TH2F("HisLclIO","HisLclIO",336,0,336,160,10,160);
   //
   fHisdEdx = new TH2F("HisdEdx","HisdEdx",336,0,336,160,1,50);
-  fHisdZfit = new TH2F("HisdZfit","HisdZfit",336,0,336,300,-1.,1.);
+  fHisdZfit = new TH2F("HisdZfit","HisdZfit",336,0,336,300,-0.3,0.3);
   
   //
   //  Chi2
@@ -2181,7 +2530,76 @@ void   AliTPCcalibLaser::MakeFitHistos(){
   fHisPz3vP2IO->SetDirectory(0);   //-> Curv  P2outerinner - common parabola
 
 
+  //
+  //
+  //
+  for (Int_t id=0; id<336;id++){
+    TH2F *profy = (TH2F*)fDeltaYres.UncheckedAt(id);
+    TH2F *profz = (TH2F*)fDeltaZres.UncheckedAt(id);
+    TH2F *profy2 = (TH2F*)fDeltaYres2.UncheckedAt(id);
+    TH2F *profz2 = (TH2F*)fDeltaZres2.UncheckedAt(id);
+    //    TH2F *profy3 = (TH2F*)fDeltaYres3.UncheckedAt(id);
+    //TH2F *profz3 = (TH2F*)fDeltaZres3.UncheckedAt(id);
+    if (!profy){
+      profy=new TH2F(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d -Linear",id),160,0,160,50,-0.5,0.5);
+      profy->SetDirectory(0);
+      fDeltaYres.AddAt(profy,id);
+      profy2=new TH2F(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d -Parabolic",id),160,0,160,50,-0.5,0.5);
+      profy2->SetDirectory(0);
+      fDeltaYres2.AddAt(profy2,id);
+      //profy3=new TH2F(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d- Parabolic2",id),160,0,160,100,-0.5,0.5);
+      //profy3->SetDirectory(0);
+      //fDeltaYres3.AddAt(profy3,id);
+    }
+    if (!profz){
+      profz=new TH2F(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d Linear",id),160,0,160,50,-0.5,0.5);
+      profz->SetDirectory(0);
+      fDeltaZres.AddAt(profz,id);
+      profz2=new TH2F(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d - Parabolic",id),160,0,160,50,-0.5,0.5);
+      profz2->SetDirectory(0);
+      fDeltaZres2.AddAt(profz2,id);
+      //profz3=new TH2F(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d- Parabolic2",id),160,0,160,100,-0.5,0.5);
+      //profz3->SetDirectory(0);
+      //fDeltaZres3.AddAt(profz3,id);
+    }
+  }
+  //
+  //
+  for (Int_t id=0; id<336;id++){
+    TH1F * hisdz = (TH1F*)fDeltaZ.At(id);
+    TH1F * hisP3 = (TH1F*)fDeltaP3.At(id);
+    TH1F * hisP4 = (TH1F*)fDeltaP4.At(id);
+    
+    TH1F * hisdphi = (TH1F*)fDeltaPhi.At(id);
+    TH1F * hisdphiP = (TH1F*)fDeltaPhiP.At(id);
+    TH1F * hisSignal = (TH1F*)fSignals.At(id);
+
+    if (!hisdz){
+      hisdz = new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10);
+      hisdz->SetDirectory(0);
+      fDeltaZ.AddAt(hisdz,id);
+
+      hisP3 = new TH1F(Form("hisPar3v%d",id),Form("hisPar3v%d",id),400,-0.06,0.06);
+      hisP3->SetDirectory(0);
+      fDeltaP3.AddAt(hisP3,id);
+      //
+      hisP4 = new TH1F(Form("hisPar4v%d",id),Form("hisPar4v%d",id),200,-0.06,0.06);
+      hisP4->SetDirectory(0);
+      fDeltaP4.AddAt(hisP4,id);
 
+      //
+      hisdphi = new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1);
+      hisdphi->SetDirectory(0);
+      fDeltaPhi.AddAt(hisdphi,id);
+      //
+      hisdphiP = new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01);
+      hisdphiP->SetDirectory(0);
+      fDeltaPhiP.AddAt(hisdphiP,id);
+      hisSignal = new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),100,0,300);
+      hisSignal->SetDirectory(0);
+      fSignals.AddAt(hisSignal,id);
+    }
+  }
 }
 
 void AliTPCcalibLaser::MergeFitHistos(AliTPCcalibLaser * laser){