]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG1/TRD/AliTRDresolution.cxx
trigger names changed (Markus)
[u/mrichter/AliRoot.git] / PWG1 / TRD / AliTRDresolution.cxx
index 47463e16efe0888a3798d9eaa1c0a607dda43089..8fb2acbaeb8046738f7ffdba7bcd4dd2579f6f24 100644 (file)
@@ -62,6 +62,7 @@
 #include <TLegend.h>
 #include <TGraphErrors.h>
 #include <TGraphAsymmErrors.h>
+#include <TLinearFitter.h>
 #include <TMath.h>
 #include <TMatrixT.h>
 #include <TVectorT.h>
 #include <TDatabasePDG.h>
 
 #include "AliPID.h"
+#include "AliLog.h"
 #include "AliESDtrack.h"
+#include "AliMathBase.h"
+#include "AliTrackPointArray.h"
 
 #include "AliTRDresolution.h"
 #include "AliTRDgeometry.h"
 #include "AliTRDReconstructor.h"
 #include "AliTRDrecoParam.h"
 #include "AliTRDpidUtil.h"
+#include "AliTRDinfoGen.h"
 
 #include "info/AliTRDclusterInfo.h"
 
 ClassImp(AliTRDresolution)
 
-UChar_t const AliTRDresolution::fgNhistos[kNviews] = {
-  2, 2, 5, 5,
-  2, 5, 12, 2, 14
-};
 UChar_t const AliTRDresolution::fgNproj[kNviews] = {
-  2, 2, 5, 5,
-  4, 7, 12, 2, 14
+  2, 2, 5, 5, 5,
+  2, 5, 11, 11, 11
 };
 Char_t const * AliTRDresolution::fgPerformanceName[kNviews] = {
      "Charge"
     ,"Cluster2Track"
     ,"Tracklet2Track"
-    ,"Tracklet2TPC" 
+    ,"Tracklet2TRDin" 
+    ,"Tracklet2TRDout" 
     ,"Cluster2MC"
     ,"Tracklet2MC"
-    ,"TPC2MC"
-    ,"TOF/HMPID2MC"
+    ,"TRDin2MC"
+    ,"TRDout2MC"
     ,"TRD2MC"
 };
-UChar_t const AliTRDresolution::fgNcomp[kNprojs] = {
-  1, 1, //2, 
-  1, 1, //2, 
-  1, 1, 1, 1, 1, //5, 
-  1, 1, 1, 1, 1, 1, 1, //5,
-// MC
-  1, 1, 1, 1,    //4, 
-  1, 1, 1, 1, 1, //5, 
-  1, 1, 1, 1, 1, 1, 1, 1, 11, 11, 11, 11, //12, 
-  1, 1, //2, 
-  1, 1, 1, 1, 1, 1, 1, 1, 66, 66, 66, 66, 66, 66 //14
+Char_t const * AliTRDresolution::fgParticle[11]={
+  " p bar", " K -", " #pi -", " #mu -", " e -",
+  " No PID",
+  " e +", " #mu +", " #pi +", " K +", " p",
 };
-Char_t const *AliTRDresolution::fgAxTitle[kNprojs][4] = {
-  // Charge
-  {"Impv", "x [cm]", "I_{mpv}", "x/x_{0}"}
- ,{"dI/Impv", "x/x_{0}", "#delta I/I_{mpv}", "x[cm]"}
-  // Clusters to Kalman
- ,{"Pos", "tg(#phi)", "#mu_{y}^{cl} [#mum]", "#sigma_{y}^{cl} [#mum]"}
- ,{"Pulls", "tg(#phi)", "PULL: #mu_{y}^{cl}", "PULL: #sigma_{y}^{cl}"}
-  // TRD tracklet to Kalman fit
- ,{"PosY", "tg(#phi)", "#mu_{y}^{trklt} [#mum]", "#sigma_{y}^{trklt} [#mum]"}
- ,{"PullsY", "tg(#phi)", "PULL: #mu_{y}^{trklt}", "PULL: #sigma_{y}^{trklt}"}
- ,{"PosZ", "tg(#theta)", "#mu_{z}^{trklt} [#mum]", "#sigma_{z}^{trklt} [#mum]"}
- ,{"PullsZ", "tg(#theta)", "PULL: #mu_{z}^{trklt}", "PULL: #sigma_{z}^{trklt}"}
- ,{"Phi", "tg(#phi)", "#mu_{#phi}^{trklt} [mrad]", "#sigma_{#phi}^{trklt} [mrad]"}
-  // TPC track 2 first TRD tracklet
- ,{"PosY", "tg(#phi)", "#mu_{y}^{TPC trklt} [#mum]", "#sigma_{y}^{TPC trklt} [#mum]"}
- ,{"PullsY", "tg(#phi)", "PULL: #mu_{y}^{TPC trklt}", "PULL: #sigma_{y}^{TPC trklt}"}
- ,{"PosZ", "tg(#theta)", "#mu_{z}^{TPC trklt} [#mum]", "#sigma_{z}^{TPC trklt} [#mum]"}
- ,{"PullsZ", "tg(#theta)", "PULL: #mu_{z}^{TPC trklt}", "PULL: #sigma_{z}^{TPC trklt}"}
- ,{"Phi", "tg(#phi)", "#mu_{#phi}^{TPC trklt} [mrad]", "#sigma_{#phi}^{TPC trklt} [mrad]"}
-  // MC cluster
- ,{"Pos1", "tg(#phi)", "MC cluster (p_{t}<1 GeV/c): #mu_{y}[#mum]", "#sigma_{y}[#mum]"}
- ,{"Pos2", "tg(#phi)", "MC cluster (p_{t}<2 GeV/c): #mu_{y} [#mum]", "#sigma_{y} [#mum]"}
- ,{"Pos3", "tg(#phi)", "MC cluster (p_{t}>3 GeV/c): #mu_{y} [#mum]", "#sigma_{y}[#mum]"}
- ,{"Pulls", "tg(#phi)", "MC cluster PULL: #mu_{y}", "#sigma_{y}"}
-  // MC tracklet
- ,{"PosY", "tg(#phi)", "MC tracklet (p_{t}<1 GeV/c): #mu_{y}[#mum]", "#sigma_{y}[#mum]"}
- ,{"PosY1", "tg(#phi)", "MC tracklet (p_{t}<2 GeV/c): #mu_{y}[#mum]", "#sigma_{y}[#mum]"}
- ,{"PosY2", "tg(#phi)", "MC tracklet (p_{t}>3 GeV/c): #mu_{y}[#mum]", "#sigma_{y}[#mum]"}
- ,{"PullsY", "tg(#phi)", "MC tracklet PULL: #mu_{y}", "MC PULL: #sigma_{y}"}
- ,{"PosZ", "tg(#theta)", "MC tracklet: #mu_{z}[#mum]", "#sigma_{z}[#mum]"}
- ,{"PullsZ", "tg(#theta)", "MC tracklet PULL: #mu_{z}", "#sigma_{z}"}
- ,{"Phi", "tg(#phi)", "MC tracklet: #mu_{#phi}[mrad]", "#sigma_{#phi}[mrad]"}
-  // MC track TPC
- ,{"PosY", "tg(#phi)", "MC TPC2TRD0: #mu_{y}[#mum]", "#sigma_{y}[#mum]"}
- ,{"PullsY", "tg(#phi)", "MC TPC2TRD0 PULL: #mu_{y}", "#sigma_{y}"}
- ,{"PosZ", "tg(#theta)", "MC TPC2TRD0: #mu_{z}[#mum]", "#sigma_{z}[#mum]"}
- ,{"PullsZ", "tg(#theta)", "MC TPC2TRD0 PULL: #mu_{z}", "#sigma_{z}"}
- ,{"Phi", "tg(#phi)", "MC: #mu_{#phi}^{TPC} [mrad]", "MC: #sigma_{#phi}^{TPC} [mrad]"}
- ,{"PullsSNP", "tg(#phi)", "MC PULL: #mu_{snp}^{TPC}", "MC PULL: #sigma_{snp}^{TPC}"}
- ,{"Theta", "tg(#theta)", "MC: #mu_{#theta}^{TPC} [mrad]", "MC: #sigma_{#theta}^{TPC} [mrad]"}
- ,{"PullsTGL", "tg(#theta)", "MC PULL: #mu_{tgl}^{TPC}", "MC PULL: #sigma_{tgl}^{TPC}"}
- ,{"Pt", "p_{t}^{MC} [GeV/c]", "MC: #mu^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"}
- ,{"Pulls1Pt", "1/p_{t}^{MC} [c/GeV]", "MC PULL: #mu_{1/p_{t}}^{TPC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"}
- ,{"P", "p^{MC} [GeV/c]", "MC: #mu^{TPC}(#Deltap/p^{MC}) [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"}
- ,{"PullsP", "p^{MC} [GeV/c]", "MC PULL: #mu^{TPC}(#Deltap/#sigma_{p})", "MC PULL: #sigma^{TPC}(#Deltap/#sigma_{p})"}
-  // MC track TOF
- ,{"PosZ", "tg(#theta)", "MC: #mu_{z}^{TOF} [#mum]", "MC: #sigma_{z}^{TOF} [#mum]"}
- ,{"PullsZ", "tg(#theta)", "MC PULL: #mu_{z}^{TOF}", "MC PULL: #sigma_{z}^{TOF}"}
-  // MC track in TRD
- ,{"PosY", "tg(#phi)", "MC: #mu_{y}^{Trk} [#mum]", "MC: #sigma_{y}^{Trk} [#mum]"}
- ,{"PullsY", "tg(#phi)", "MC PULL: #mu_{y}^{Trk}", "MC PULL: #sigma_{y}^{Trk}"}
- ,{"PosZ", "tg(#theta)", "MC: #mu_{z}^{Trk} [#mum]", "MC: #sigma_{z}^{Trk} [#mum]"}
- ,{"PullsZ", "tg(#theta)", "MC PULL: #mu_{z}^{Trk}", "MC PULL: #sigma_{z}^{Trk}"}
- ,{"Phi", "tg(#phi)", "MC: #mu_{#phi}^{Trk} [mrad]", "MC: #sigma_{#phi}^{Trk} [mrad]"}
- ,{"PullsSNP", "tg(#phi)", "MC PULL: #mu_{snp}^{Trk}", "MC PULL: #sigma_{snp}^{Trk}"}
- ,{"Theta", "tg(#theta)", "MC: #mu_{#theta}^{Trk} [mrad]", "MC: #sigma_{#theta}^{Trk} [mrad]"}
- ,{"PullsTGL", "tg(#theta)", "MC PULL: #mu_{tgl}^{Trk}", "MC PULL: #sigma_{tgl}^{Trk}"}
- ,{"MC Pt Layer", "p_{t} [GeV/c]", "#mu(#Deltap_{t}/p_{t}^{MC}) [%]", "#sigma(#Deltap_{t}/p_{t}^{MC}) [%]"}
- ,{"MC Pulls1Pt Layer", "1/p_{t} [c/GeV]", "#mu_{1/p_{t}}", "#sigma_{1/p_{t}}"}
- ,{"MC P Layer", "p [GeV/c]", "#mu(#Deltap/p^{MC}) [%]", "#sigma(#Deltap/p^{MC}) [%]"}
- ,{"Pt - SA", "p_{t}^{MC} [GeV/c]", "MC: #mu^{Trk}(#Deltap_{t}/p_{t}^{MC}) [%]", "MC: #sigma^{Trk}(#Deltap_{t}/p_{t}^{MC}) [%]"}
- ,{"Pulls1Pt - SA", "1/p_{t}^{MC} [c/GeV]", "MC PULL: #mu_{1/p_{t}}^{Trk}", "MC PULL: #sigma_{1/p_{t}}^{Trk}"}
- ,{"P - SA", "p^{MC} [GeV/c]", "MC: #mu^{Trk}(#Deltap/p^{MC}) [%]", "MC: #sigma^{Trk}(#Deltap/p^{MC}) [%]"}
+
+// Configure segmentation for y resolution/residuals
+Int_t const AliTRDresolution::fgkNresYsegm[3] = {
+  AliTRDgeometry::kNsector
+ ,AliTRDgeometry::kNsector*AliTRDgeometry::kNstack
+ ,AliTRDgeometry::kNdet
 };
+Char_t const *AliTRDresolution::fgkResYsegmName[3] = {
+  "Sector", "Stack", "Detector"};
+
 
 //________________________________________________________
 AliTRDresolution::AliTRDresolution()
-  :AliTRDrecoTask("resolution", "Spatial and momentum TRD resolution checker")
-  ,fStatus(0)
+  :AliTRDrecoTask()
+  ,fSegmentLevel(0)
   ,fIdxPlot(0)
   ,fIdxFrame(0)
-  ,fReconstructor(NULL)
-  ,fGeo(NULL)
+  ,fPtThreshold(1.)
+  ,fDyRange(1.5)
   ,fDBPDG(NULL)
   ,fGraphS(NULL)
   ,fGraphM(NULL)
   ,fCl(NULL)
-  ,fTrklt(NULL)
   ,fMCcl(NULL)
-  ,fMCtrklt(NULL)
+/*  ,fTrklt(NULL)
+  ,fMCtrklt(NULL)*/
 {
   //
   // Default constructor
   //
+  SetNameTitle("TRDresolution", "TRD spatial and momentum resolution");
+  SetSegmentationLevel();
+}
 
-  fReconstructor = new AliTRDReconstructor();
-  fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
-  fGeo = new AliTRDgeometry();
+//________________________________________________________
+AliTRDresolution::AliTRDresolution(char* name)
+  :AliTRDrecoTask(name, "TRD spatial and momentum resolution")
+  ,fSegmentLevel(0)
+  ,fIdxPlot(0)
+  ,fIdxFrame(0)
+  ,fPtThreshold(1.)
+  ,fDyRange(1.5)
+  ,fDBPDG(NULL)
+  ,fGraphS(NULL)
+  ,fGraphM(NULL)
+  ,fCl(NULL)
+  ,fMCcl(NULL)
+/*  ,fTrklt(NULL)
+  ,fMCtrklt(NULL)*/
+{
+  //
+  // Default constructor
+  //
 
   InitFunctorList();
+  SetSegmentationLevel();
 
-  DefineOutput(1, TObjArray::Class()); // cluster2track
-  DefineOutput(2, TObjArray::Class()); // tracklet2track
-  DefineOutput(3, TObjArray::Class()); // cluster2mc
-  DefineOutput(4, TObjArray::Class()); // tracklet2mc
+  DefineOutput(kClToTrk, TObjArray::Class()); // cluster2track
+  DefineOutput(kClToMC, TObjArray::Class()); // cluster2mc
+/*  DefineOutput(kTrkltToTrk, TObjArray::Class()); // tracklet2track
+  DefineOutput(kTrkltToMC, TObjArray::Class()); // tracklet2mc*/
 }
 
 //________________________________________________________
@@ -223,52 +184,76 @@ AliTRDresolution::~AliTRDresolution()
 
   if(fGraphS){fGraphS->Delete(); delete fGraphS;}
   if(fGraphM){fGraphM->Delete(); delete fGraphM;}
-  delete fGeo;
-  delete fReconstructor;
-  if(gGeoManager) delete gGeoManager;
   if(fCl){fCl->Delete(); delete fCl;}
-  if(fTrklt){fTrklt->Delete(); delete fTrklt;}
   if(fMCcl){fMCcl->Delete(); delete fMCcl;}
-  if(fMCtrklt){fMCtrklt->Delete(); delete fMCtrklt;}
+/*  if(fTrklt){fTrklt->Delete(); delete fTrklt;}
+  if(fMCtrklt){fMCtrklt->Delete(); delete fMCtrklt;}*/
 }
 
 
 //________________________________________________________
-void AliTRDresolution::CreateOutputObjects()
+void AliTRDresolution::UserCreateOutputObjects()
 {
   // spatial resolution
-  OpenFile(0, "RECREATE");
 
-  fContainer = Histos();
+  AliTRDrecoTask::UserCreateOutputObjects();
+  InitExchangeContainers();
+}
+
+//________________________________________________________
+void AliTRDresolution::InitExchangeContainers()
+{
+// Init containers for subsequent tasks (AliTRDclusterResolution)
 
-  fCl = new TObjArray();
+  fCl = new TObjArray(200);
   fCl->SetOwner(kTRUE);
-  fTrklt = new TObjArray();
-  fTrklt->SetOwner(kTRUE);
   fMCcl = new TObjArray();
   fMCcl->SetOwner(kTRUE);
+/*  fTrklt = new TObjArray();
+  fTrklt->SetOwner(kTRUE);
   fMCtrklt = new TObjArray();
-  fMCtrklt->SetOwner(kTRUE);
+  fMCtrklt->SetOwner(kTRUE);*/
+  PostData(kClToTrk, fCl);
+  PostData(kClToMC, fMCcl);
+/*  PostData(kTrkltToTrk, fTrklt);
+  PostData(kTrkltToMC, fMCtrklt);*/
 }
 
 //________________________________________________________
-void AliTRDresolution::Exec(Option_t *opt)
+void AliTRDresolution::UserExec(Option_t *opt)
 {
   //
   // Execution part
   //
 
   fCl->Delete();
-  fTrklt->Delete();
   fMCcl->Delete();
-  fMCtrklt->Delete();
-
-  AliTRDrecoTask::Exec(opt);
+  AliTRDrecoTask::UserExec(opt);
+}
 
-  PostData(1, fCl);
-  PostData(2, fTrklt);
-  PostData(3, fMCcl);
-  PostData(4, fMCtrklt);
+//________________________________________________________
+Bool_t AliTRDresolution::Pulls(Double_t dyz[2], Double_t cov[3], Double_t tilt) const
+{
+// Helper function to calculate pulls in the yz plane 
+// using proper tilt rotation
+// Uses functionality defined by AliTRDseedV1.
+
+  Double_t t2(tilt*tilt);
+
+  // rotate along pad
+  Double_t cc[3];
+  cc[0] = cov[0] - 2.*tilt*cov[1] + t2*cov[2]; 
+  cc[1] = cov[1]*(1.-t2) + tilt*(cov[0] - cov[2]);
+  cc[2] = t2*cov[0] + 2.*tilt*cov[1] + cov[2];
+  // do sqrt
+  Double_t sqr[3]={0., 0., 0.}; 
+  if(AliTRDseedV1::GetCovSqrt(cc, sqr)) return kFALSE;
+  Double_t invsqr[3]={0., 0., 0.}; 
+  if(AliTRDseedV1::GetCovInv(sqr, invsqr)<1.e-40) return kFALSE;
+  Double_t tmp(dyz[0]);
+  dyz[0] = invsqr[0]*tmp + invsqr[1]*dyz[1];
+  dyz[1] = invsqr[1]*tmp + invsqr[2]*dyz[1];
+  return kTRUE;
 }
 
 //________________________________________________________
@@ -280,11 +265,11 @@ TH1* AliTRDresolution::PlotCharge(const AliTRDtrackV1 *track)
 
   if(track) fkTrack = track;
   if(!fkTrack){
-    AliDebug(2, "No Track defined.");
+    AliDebug(4, "No Track defined.");
     return NULL;
   }
   TObjArray *arr = NULL;
-  if(!(arr = ((TObjArray*)fContainer->At(kCharge)))){
+  if(!fContainer || !(arr = ((TObjArray*)fContainer->At(kCharge)))){
     AliWarning("No output container defined.");
     return NULL;
   }
@@ -296,14 +281,15 @@ TH1* AliTRDresolution::PlotCharge(const AliTRDtrackV1 *track)
     if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
     if(!fTracklet->IsOK()) continue;
     Float_t x0 = fTracklet->GetX0();
-    Float_t dq, dl;
+    Float_t dqdl, dl;
     for(Int_t itb=AliTRDseedV1::kNtb; itb--;){
       if(!(c = fTracklet->GetClusters(itb))){ 
         if(!(c = fTracklet->GetClusters(AliTRDseedV1::kNtb+itb))) continue;
       }
-      dq = fTracklet->GetdQdl(itb, &dl);
+      dqdl = fTracklet->GetdQdl(itb, &dl);
+      if(dqdl<1.e-5) continue;
       dl /= 0.15; // dl/dl0, dl0 = 1.5 mm for nominal vd
-      (h = (TH3S*)arr->At(0))->Fill(dl, x0-c->GetX(), dq);
+      (h = (TH3S*)arr->At(0))->Fill(dl, x0-c->GetX(), dqdl);
     }
 
 //     if(!HasMCdata()) continue;
@@ -325,31 +311,108 @@ TH1* AliTRDresolution::PlotCluster(const AliTRDtrackV1 *track)
 
   if(track) fkTrack = track;
   if(!fkTrack){
-    AliDebug(2, "No Track defined.");
+    AliDebug(4, "No Track defined.");
     return NULL;
   }
   TObjArray *arr = NULL;
-  if(!(arr = ((TObjArray*)fContainer->At(kCluster)))){
+  if(!fContainer || !(arr = ((TObjArray*)fContainer->At(kCluster)))){
     AliWarning("No output container defined.");
     return NULL;
   }
   ULong_t status = fkESD ? fkESD->GetStatus():0;
 
-  Double_t cov[7];
-  Float_t x0, y0, z0, dy, dydx, dzdx;
-  AliTRDseedV1 *fTracklet = NULL;  
+  Int_t sgm[3];
+  Double_t covR[7], cov[3], dy[2], dz[2];
+  Float_t pt, x0, y0, z0, dydx, dzdx, tilt(0.);
+  const AliTRDgeometry *geo(AliTRDinfoGen::Geometry());
+  AliTRDseedV1 *fTracklet(NULL);  TObjArray *clInfoArr(NULL);
+  // do LINEAR track refit if asked by the user
+  // it is the user responsibility to check if B=0T
+  Float_t param[10];  memset(param, 0, 10*sizeof(Float_t));
+  Int_t np(0), nrc(0); AliTrackPoint clusters[300];
+  if(HasTrackRefit()){
+    Bool_t kPrimary(kFALSE);
+    for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
+      if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
+      if(!fTracklet->IsOK()) continue;
+      x0   = fTracklet->GetX0();
+      tilt = fTracklet->GetTilt();
+      AliTRDcluster *c = NULL;
+      fTracklet->ResetClusterIter(kFALSE);
+      while((c = fTracklet->PrevCluster())){
+        Float_t xyz[3] = {c->GetX(), c->GetY(), c->GetZ()};
+        clusters[np].SetCharge(tilt);
+        clusters[np].SetClusterType(0);
+        clusters[np].SetVolumeID(ily);
+        clusters[np].SetXYZ(xyz);
+        np++;
+      }
+      if(fTracklet->IsRowCross()){
+        Float_t xcross(0.), zcross(0.);
+        if(fTracklet->GetEstimatedCrossPoint(xcross, zcross)){
+          clusters[np].SetCharge(tilt);
+          clusters[np].SetClusterType(1);
+          clusters[np].SetVolumeID(ily);
+          clusters[np].SetXYZ(xcross, 0., zcross);
+          np++;
+          nrc++;
+        }
+      }
+      if(fTracklet->IsPrimary()) kPrimary = kTRUE;
+    }
+    if(kPrimary){
+      clusters[np].SetCharge(tilt);
+      clusters[np].SetClusterType(1);
+      clusters[np].SetVolumeID(-1);
+      clusters[np].SetXYZ(0., 0., 0.);
+      np++;
+    }
+    if(!FitTrack(np, clusters, param)) {
+      AliDebug(1, "Linear track Fit failed.");
+      return NULL;
+    }
+    if(HasTrackSelection()){
+      Bool_t kReject(kFALSE);
+      if(fkTrack->GetNumberOfTracklets() != AliTRDgeometry::kNlayer)  kReject = kTRUE; 
+      if(!kReject && !UseTrack(np, clusters, param)) kReject = kTRUE;
+      if(kReject){
+        AliDebug(1, "Reject track for residuals analysis.");
+        return NULL;
+      }
+    }
+  }
   for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
     if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
     if(!fTracklet->IsOK()) continue;
     x0 = fTracklet->GetX0();
-
+    pt = fTracklet->GetPt();
+    sgm[2] = fTracklet->GetDetector();
+    sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
+    sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
+    
     // retrive the track angle with the chamber
-    y0   = fTracklet->GetYref(0);
-    z0   = fTracklet->GetZref(0);
-    dydx = fTracklet->GetYref(1);
-    dzdx = fTracklet->GetZref(1);
-    fTracklet->GetCovRef(cov);
-    Float_t tilt = fTracklet->GetTilt();
+    if(HasTrackRefit()){
+      Float_t par[3];
+      if(!FitTracklet(ily, np, clusters, param, par)) continue;
+      dydx = par[2];//param[3];
+      dzdx = param[4];
+      y0   = par[1] + dydx * (x0 - par[0]);//param[1] + dydx * (x0 - param[0]);
+      z0   = param[2] + dzdx * (x0 - param[0]);
+    } else {
+      y0   = fTracklet->GetYref(0);
+      z0   = fTracklet->GetZref(0);
+      dydx = fTracklet->GetYref(1);
+      dzdx = fTracklet->GetZref(1);
+    }
+    /*printf("RC[%c] Primary[%c]\n"
+           "  Fit : y0[%f] z0[%f] dydx[%f] dzdx[%f]\n"
+           "  Ref:  y0[%f] z0[%f] dydx[%f] dzdx[%f]\n", fTracklet->IsRowCross()?'y':'n', fTracklet->IsPrimary()?'y':'n', y0, z0, dydx, dzdx
+           ,fTracklet->GetYref(0),fTracklet->GetZref(0),fTracklet->GetYref(1),fTracklet->GetZref(1));*/
+    tilt = fTracklet->GetTilt();
+    fTracklet->GetCovRef(covR);
+    Double_t t2(tilt*tilt)
+            ,corr(1./(1. + t2))
+            ,cost(TMath::Sqrt(corr));
     AliTRDcluster *c = NULL;
     fTracklet->ResetClusterIter(kFALSE);
     while((c = fTracklet->PrevCluster())){
@@ -359,40 +422,63 @@ TH1* AliTRDresolution::PlotCluster(const AliTRDtrackV1 *track)
       Float_t dx = x0 - xc; 
       Float_t yt = y0 - dx*dydx;
       Float_t zt = z0 - dx*dzdx; 
-      yc -= tilt*(zc-zt); // tilt correction
-      dy = yt - yc;
-
-      //Float_t sx2 = dydx*c->GetSX(c->GetLocalTimeBin()); sx2*=sx2;
-      Float_t sy2 = c->GetSigmaY2();
-      if(sy2<=0.) continue;
-      ((TH2I*)arr->At(0))->Fill(dydx, dy);
-      ((TH2I*)arr->At(1))->Fill(dydx, dy/TMath::Sqrt(cov[0] /*+ sx2*/ + sy2));
+      dy[0] = yc-yt; dz[0]= zc-zt;
+
+      // rotate along pad
+      dy[1] = cost*(dy[0] - dz[0]*tilt);
+      dz[1] = cost*(dz[0] + dy[0]*tilt);
+      if(pt>fPtThreshold && c->IsInChamber()) ((TH3S*)arr->At(0))->Fill(dydx, dy[1], sgm[fSegmentLevel]);
+
+      // tilt rotation of covariance for clusters
+      Double_t sy2(c->GetSigmaY2()), sz2(c->GetSigmaZ2());
+      cov[0] = (sy2+t2*sz2)*corr;
+      cov[1] = tilt*(sz2 - sy2)*corr;
+      cov[2] = (t2*sy2+sz2)*corr;
+      // sum with track covariance
+      cov[0]+=covR[0]; cov[1]+=covR[1]; cov[2]+=covR[2];
+      Double_t dyz[2]= {dy[1], dz[1]};
+      Pulls(dyz, cov, tilt);
+      ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
   
-      if(DebugLevel()>=2){
-        // Get z-position with respect to anode wire
-        Int_t istk = fGeo->GetStack(c->GetDetector());
-        AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
-        Float_t row0 = pp->GetRow0();
-        Float_t d  =  row0 - zt + pp->GetAnodeWireOffset();
-        d -= ((Int_t)(2 * d)) / 2.0;
-        if (d > 0.25) d  = 0.5 - d;
-
-        AliTRDclusterInfo *clInfo = new AliTRDclusterInfo;
-        fCl->Add(clInfo);
-        clInfo->SetCluster(c);
-        Float_t covR[] = {cov[0], cov[1], cov[2]};
-        clInfo->SetGlobalPosition(yt, zt, dydx, dzdx, covR);
-        clInfo->SetResolution(dy);
-        clInfo->SetAnisochronity(d);
-        clInfo->SetDriftLength(dx);
-        (*DebugStream()) << "ClusterREC"
-          <<"status="  << status
-          <<"clInfo.=" << clInfo
-          << "\n";
+      // Get z-position with respect to anode wire
+      Int_t istk = geo->GetStack(c->GetDetector());
+      AliTRDpadPlane *pp = geo->GetPadPlane(ily, istk);
+      Float_t row0 = pp->GetRow0();
+      Float_t d  =  row0 - zt + pp->GetAnodeWireOffset();
+      d -= ((Int_t)(2 * d)) / 2.0;
+      if (d > 0.25) d  = 0.5 - d;
+
+      AliTRDclusterInfo *clInfo(NULL);
+      clInfo = new AliTRDclusterInfo;
+      clInfo->SetCluster(c);
+      Float_t covF[] = {cov[0], cov[1], cov[2]};
+      clInfo->SetGlobalPosition(yt, zt, dydx, dzdx, covF);
+      clInfo->SetResolution(dy[1]);
+      clInfo->SetAnisochronity(d);
+      clInfo->SetDriftLength(dx);
+      clInfo->SetTilt(tilt);
+      if(fCl) fCl->Add(clInfo);
+      else AliDebug(1, "Cl exchange container missing. Activate by calling \"InitExchangeContainers()\"");
+
+      if(DebugLevel()>=1){
+        if(!clInfoArr){ 
+          clInfoArr=new TObjArray(AliTRDseedV1::kNclusters);
+          clInfoArr->SetOwner(kFALSE);
+        }
+        clInfoArr->Add(clInfo);
       }
     }
+    if(DebugLevel()>=1 && clInfoArr){
+      (*DebugStream()) << "cluster"
+        <<"status="  << status
+        <<"clInfo.=" << clInfoArr
+        << "\n";
+      clInfoArr->Clear();
+    }
   }
-  return (TH2I*)arr->At(0);
+  if(clInfoArr) delete clInfoArr;
+
+  return (TH3S*)arr->At(0);
 }
 
 
@@ -409,54 +495,84 @@ TH1* AliTRDresolution::PlotTracklet(const AliTRDtrackV1 *track)
 // reference position. 
   if(track) fkTrack = track;
   if(!fkTrack){
-    AliDebug(2, "No Track defined.");
+    AliDebug(4, "No Track defined.");
     return NULL;
   }
   TObjArray *arr = NULL;
-  if(!(arr = (TObjArray*)fContainer->At(kTrackTRD ))){
+  if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrack ))){
     AliWarning("No output container defined.");
     return NULL;
   }
 
+  Int_t sgm[3];
   Double_t cov[3], covR[7]/*, sqr[3], inv[3]*/;
-  Float_t x, dx, dy, dz;
-  AliTRDseedV1 *fTracklet = NULL;  
-  for(Int_t il=AliTRDgeometry::kNlayer; il--;){
+  Double_t pt, phi, tht, x, dx, dy[2], dz[2];
+  AliTRDseedV1 *fTracklet(NULL);  
+  for(Int_t il(0); il<AliTRDgeometry::kNlayer; il++){
     if(!(fTracklet = fkTrack->GetTracklet(il))) continue;
     if(!fTracklet->IsOK()) continue;
-    x    = fTracklet->GetX();
-    dx   = fTracklet->GetX0() - x;
-    // compute dy^2 and dz^2
-    dy   = fTracklet->GetYref(0)-dx*fTracklet->GetYref(1) - fTracklet->GetY();
-    dz   = fTracklet->GetZref(0)-dx*fTracklet->GetZref(1) - fTracklet->GetZ();
+    sgm[2] = fTracklet->GetDetector();
+    sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
+    sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
+    x   = fTracklet->GetX();
+    dx  = fTracklet->GetX0() - x;
+    pt  = fTracklet->GetPt();
+    phi = fTracklet->GetYref(1);
+    tht = fTracklet->GetZref(1);
+    // compute dy and dz
+    dy[0]= fTracklet->GetYref(0)-dx*fTracklet->GetYref(1) - fTracklet->GetY();
+    dz[0]= fTracklet->GetZref(0)-dx*fTracklet->GetZref(1) - fTracklet->GetZ();
+    Double_t tilt(fTracklet->GetTilt())
+            ,t2(tilt*tilt)
+            ,corr(1./(1. + t2))
+            ,cost(TMath::Sqrt(corr));
+    Bool_t rc(fTracklet->IsRowCross());
+
+    // calculate residuals using tilt rotation
+    dy[1]= cost*(dy[0] - dz[0]*tilt);
+    dz[1]= cost*(dz[0] + dy[0]*tilt);
+    ((TH3S*)arr->At(0))->Fill(phi, dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]);
+    ((TH3S*)arr->At(2))->Fill(tht, dz[1], rc);
+
     // compute covariance matrix
     fTracklet->GetCovAt(x, cov);
     fTracklet->GetCovRef(covR);
     cov[0] += covR[0]; cov[1] += covR[1]; cov[2] += covR[2]; 
-/*  // Correct PULL calculation by considering off  
-    // diagonal elements in the covariance matrix
-    // compute square root matrix
-    if(AliTRDseedV1::GetCovInv(cov, inv)==0.) continue;
-    if(AliTRDseedV1::GetCovSqrt(inv, sqr)<0.) continue;
-    Double_t y = sqr[0]*dy+sqr[1]*dz;
-    Double_t z = sqr[1]*dy+sqr[2]*dz;
-    ((TH3*)h)->Fill(y, z, fTracklet->GetYref(1));*/
-
-    ((TH2I*)arr->At(0))->Fill(fTracklet->GetYref(1), dy);
-    ((TH2I*)arr->At(1))->Fill(fTracklet->GetYref(1), dy/TMath::Sqrt(cov[0]));
-    ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), TMath::ATan((fTracklet->GetYref(1)-fTracklet->GetYfit(1))/(1-fTracklet->GetYref(1)*fTracklet->GetYfit(1))));
-    if(!fTracklet->IsRowCross()) continue;
-    ((TH2I*)arr->At(2))->Fill(fTracklet->GetZref(1), dz);
-    ((TH2I*)arr->At(3))->Fill(fTracklet->GetZref(1), dz/TMath::Sqrt(cov[2]));
-  }
+    Double_t dyz[2]= {dy[1], dz[1]};
+    Pulls(dyz, cov, tilt);
+    ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
+    ((TH3S*)arr->At(3))->Fill(tht, dyz[1], rc);
 
+    Double_t dphi((phi-fTracklet->GetYfit(1))/(1-phi*fTracklet->GetYfit(1)));
+    Double_t dtht((tht-fTracklet->GetZfit(1))/(1-tht*fTracklet->GetZfit(1)));
+    ((TH2I*)arr->At(4))->Fill(phi, TMath::ATan(dphi));
 
+    if(DebugLevel()>=1){
+      UChar_t err(fTracklet->GetErrorMsg());
+      (*DebugStream()) << "tracklet"
+        <<"pt="  << pt
+        <<"phi=" << phi
+        <<"tht=" << tht
+        <<"det=" << sgm[2]
+        <<"dy0="  << dy[0]
+        <<"dz0="  << dz[0]
+        <<"dy="  << dy[1]
+        <<"dz="  << dz[1]
+        <<"dphi="<< dphi
+        <<"dtht="<< dtht
+        <<"dyp=" << dyz[0]
+        <<"dzp=" << dyz[1]
+        <<"rc="  << rc
+        <<"err=" << err
+        << "\n";
+    }
+  }
   return (TH2I*)arr->At(0);
 }
 
 
 //________________________________________________________
-TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track)
+TH1* AliTRDresolution::PlotTrackIn(const AliTRDtrackV1 *track)
 {
 // Store resolution/pulls of Kalman before updating with the TRD information 
 // at the radial position of the first tracklet. The following points are used 
@@ -469,30 +585,45 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track)
 
   if(track) fkTrack = track;
   if(!fkTrack){
-    AliDebug(2, "No Track defined.");
+    AliDebug(4, "No Track defined.");
+    return NULL;
+  }
+  TObjArray *arr = NULL;
+  if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrackIn))){
+    AliWarning("No output container defined.");
     return NULL;
   }
   AliExternalTrackParam *tin = NULL;
-  if(!(tin = fkTrack->GetTrackLow())){
+  if(!(tin = fkTrack->GetTrackIn())){
     AliWarning("Track did not entered TRD fiducial volume.");
     return NULL;
   }
   TH1 *h = NULL;
   
   Double_t x = tin->GetX();
-  AliTRDseedV1 *tracklet = NULL;  
+  AliTRDseedV1 *fTracklet = NULL;  
   for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
-    if(!(tracklet = fkTrack->GetTracklet(ily))) continue;
+    if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
     break;
   }
-  if(!tracklet || TMath::Abs(x-tracklet->GetX())>1.e-3){
-    AliWarning("Tracklet did not match TRD entrance.");
+  if(!fTracklet || TMath::Abs(x-fTracklet->GetX())>1.e-3){
+    AliWarning("Tracklet did not match Track.");
     return NULL;
   }
+  Int_t sgm[3];
+  sgm[2] = fTracklet->GetDetector();
+  sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
+  sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
+  Double_t tilt(fTracklet->GetTilt())
+          ,t2(tilt*tilt)
+          ,corr(1./(1. + t2))
+          ,cost(TMath::Sqrt(corr));
+  Bool_t rc(fTracklet->IsRowCross());
+
   const Int_t kNPAR(5);
   Double_t parR[kNPAR]; memcpy(parR, tin->GetParameter(), kNPAR*sizeof(Double_t));
   Double_t covR[3*kNPAR]; memcpy(covR, tin->GetCovariance(), 3*kNPAR*sizeof(Double_t));
-  Double_t cov[3]; tracklet->GetCovAt(x, cov);
+  Double_t cov[3]; fTracklet->GetCovAt(x, cov);
 
   // define sum covariances
   TMatrixDSym COV(kNPAR); TVectorD PAR(kNPAR);
@@ -507,30 +638,37 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track)
   //COV.Print(); PAR.Print();
 
   //TODO Double_t dydx =  TMath::Sqrt(1.-parR[2]*parR[2])/parR[2]; 
-  Double_t dy = parR[0] - tracklet->GetY(); 
-  TObjArray *arr = (TObjArray*)fContainer->At(kTrackTPC);
-  ((TH2I*)arr->At(0))->Fill(tracklet->GetYref(1), dy);
-  ((TH2I*)arr->At(1))->Fill(tracklet->GetYref(1), dy/TMath::Sqrt(COV(0,0)+cov[0]));
-  if(tracklet->IsRowCross()){
-    Double_t dz = parR[1] - tracklet->GetZ(); 
-    ((TH2I*)arr->At(2))->Fill(tracklet->GetZref(1), dz);
-    ((TH2I*)arr->At(3))->Fill(tracklet->GetZref(1), dz/TMath::Sqrt(COV(1,1)+cov[2]));
-  }
-  Double_t dphi = TMath::ASin(PAR[2])-TMath::ATan(tracklet->GetYfit(1));  ((TH2I*)arr->At(4))->Fill(tracklet->GetYref(1), dphi);
+  Double_t dy[2]={parR[0] - fTracklet->GetY(), 0.}
+          ,dz[2]={parR[1] - fTracklet->GetZ(), 0.}
+          ,dphi(TMath::ASin(PAR[2])-TMath::ATan(fTracklet->GetYfit(1)));
+  // calculate residuals using tilt rotation
+  dy[1] = cost*(dy[0] - dz[0]*tilt);
+  dz[1] = cost*(dz[0] + dy[0]*tilt);
+
+  if(1./PAR[4]>fPtThreshold) ((TH3S*)arr->At(0))->Fill(fTracklet->GetYref(1), dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]);
+  ((TH3S*)arr->At(2))->Fill(fTracklet->GetZref(1), dz[1], rc);
+  ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), dphi);
+
+  Double_t dyz[2] = {dy[1], dz[1]};
+  Double_t cc[3] = {COV(0,0)+cov[0], COV(0,1)+cov[1], COV(1,1)+cov[2]};
+  Pulls(dyz, cc, tilt);
+  ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
+  ((TH3S*)arr->At(3))->Fill(fTracklet->GetZref(1), dyz[1], rc);
+
 
 
   // register reference histo for mini-task
   h = (TH2I*)arr->At(0);
 
-  if(DebugLevel()>=1){
+  if(DebugLevel()>=2){
     (*DebugStream()) << "trackIn"
       << "x="       << x
       << "P="       << &PAR
       << "C="       << &COV
       << "\n";
 
-    Double_t y = tracklet->GetY(); 
-    Double_t z = tracklet->GetZ(); 
+    Double_t y = fTracklet->GetY(); 
+    Double_t z = fTracklet->GetZ(); 
     (*DebugStream()) << "trackletIn"
       << "y="       << y
       << "z="       << z
@@ -543,7 +681,7 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track)
 
   if(!HasMCdata()) return h;
   UChar_t s;
-  Float_t dx, pt0, x0=tracklet->GetX0(), y0, z0, dydx0, dzdx0;
+  Float_t dx, pt0, x0=fTracklet->GetX0(), y0, z0, dydx0, dzdx0;
   if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) return h;
   Int_t pdg = fkMC->GetPDG(),
         sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index
@@ -569,13 +707,17 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track)
 //   TMatrixD sqrcov(evecs, TMatrixD::kMult, TMatrixD(evalsm, TMatrixD::kMult, evecs.T()));
   
   // fill histos
-  arr = (TObjArray*)fContainer->At(kMCtrackTPC);
+  if(!(arr = (TObjArray*)fContainer->At(kMCtrackIn))) {
+    AliWarning("No MC container defined.");
+    return h;
+  }
+
   // y resolution/pulls
-  ((TH2I*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0]);
-  ((TH2I*)arr->At(1))->Fill(dydx0, (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0)));
+  if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0], sgm[fSegmentLevel]);
+  ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0)), (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)));
   // z resolution/pulls
-  ((TH2I*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1]);
-  ((TH2I*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)));
+  ((TH3S*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1], 0);
+  ((TH3S*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)), 0);
   // phi resolution/snp pulls
   ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ASin(PARMC[2])-TMath::ASin(PAR[2]));
   ((TH2I*)arr->At(5))->Fill(dydx0, (PARMC[2]-PAR[2])/TMath::Sqrt(COV(2,2)));
@@ -588,15 +730,15 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track)
 
   Double_t p0 = TMath::Sqrt(1.+ PARMC[3]*PARMC[3])*pt0, p;
   p = TMath::Sqrt(1.+ PAR[3]*PAR[3])/PAR[4];
-  Float_t sp = 
-    p*p*PAR[4]*PAR[4]*COV(4,4)
-   +2.*PAR[3]*COV(3,4)/PAR[4]
-   +PAR[3]*PAR[3]*COV(3,3)/p/p/PAR[4]/PAR[4]/PAR[4]/PAR[4];
   ((TH3S*)arr->At(10))->Fill(p0, p/p0-1., sign*sIdx);
-  if(sp>0.) ((TH3S*)arr->At(11))->Fill(p0, (p0-p)/TMath::Sqrt(sp), sign*sIdx);
+//   Float_t sp = 
+//     p*p*PAR[4]*PAR[4]*COV(4,4)
+//    +2.*PAR[3]*COV(3,4)/PAR[4]
+//    +PAR[3]*PAR[3]*COV(3,3)/p/p/PAR[4]/PAR[4]/PAR[4]/PAR[4];
+//   if(sp>0.) ((TH3S*)arr->At(11))->Fill(p0, (p0-p)/TMath::Sqrt(sp), sign*sIdx);
 
   // fill debug for MC 
-  if(DebugLevel()>=1){
+  if(DebugLevel()>=3){
     (*DebugStream()) << "trackInMC"
       << "P="   << &PARMC
       << "\n";
@@ -604,6 +746,178 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track)
   return h;
 }
 
+//________________________________________________________
+TH1* AliTRDresolution::PlotTrackOut(const AliTRDtrackV1 *track)
+{
+// Store resolution/pulls of Kalman after last update with the TRD information 
+// at the radial position of the first tracklet. The following points are used 
+// for comparison  
+//  - the (y,z,snp) of the first TRD tracklet
+//  - the (y, z, snp, tgl, pt) of the MC track reference
+// 
+// Additionally the momentum resolution/pulls are calculated for usage in the 
+// PID calculation. 
+
+  if(track) fkTrack = track;
+  if(!fkTrack){
+    AliDebug(4, "No Track defined.");
+    return NULL;
+  }
+  TObjArray *arr = NULL;
+  if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrackOut))){
+    AliWarning("No output container defined.");
+    return NULL;
+  }
+  AliExternalTrackParam *tout = NULL;
+  if(!(tout = fkTrack->GetTrackOut())){
+    AliDebug(2, "Track did not exit TRD.");
+    return NULL;
+  }
+  TH1 *h(NULL);
+  
+  Double_t x = tout->GetX();
+  AliTRDseedV1 *fTracklet(NULL);  
+  for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
+    if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
+    break;
+  }
+  if(!fTracklet || TMath::Abs(x-fTracklet->GetX())>1.e-3){
+    AliWarning("Tracklet did not match Track position.");
+    return NULL;
+  }
+  Int_t sgm[3];
+  sgm[2] = fTracklet->GetDetector();
+  sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
+  sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
+  Double_t tilt(fTracklet->GetTilt())
+          ,t2(tilt*tilt)
+          ,corr(1./(1. + t2))
+          ,cost(TMath::Sqrt(corr));
+  Bool_t rc(fTracklet->IsRowCross());
+
+  const Int_t kNPAR(5);
+  Double_t parR[kNPAR]; memcpy(parR, tout->GetParameter(), kNPAR*sizeof(Double_t));
+  Double_t covR[3*kNPAR]; memcpy(covR, tout->GetCovariance(), 3*kNPAR*sizeof(Double_t));
+  Double_t cov[3]; fTracklet->GetCovAt(x, cov);
+
+  // define sum covariances
+  TMatrixDSym COV(kNPAR); TVectorD PAR(kNPAR);
+  Double_t *pc = &covR[0], *pp = &parR[0];
+  for(Int_t ir=0; ir<kNPAR; ir++, pp++){
+    PAR(ir) = (*pp);
+    for(Int_t ic = 0; ic<=ir; ic++,pc++){ 
+      COV(ir,ic) = (*pc); COV(ic,ir) = (*pc);
+    }
+  }
+  PAR[4] = TMath::Abs(PAR[4]); // remove sign of pt !!
+  //COV.Print(); PAR.Print();
+
+  //TODO Double_t dydx =  TMath::Sqrt(1.-parR[2]*parR[2])/parR[2]; 
+  Double_t dy[3]={parR[0] - fTracklet->GetY(), 0., 0.}
+          ,dz[3]={parR[1] - fTracklet->GetZ(), 0., 0.}
+          ,dphi(TMath::ASin(PAR[2])-TMath::ATan(fTracklet->GetYfit(1)));
+  // calculate residuals using tilt rotation
+  dy[1] = cost*(dy[0] - dz[0]*tilt);
+  dz[1] = cost*(dz[0] + dy[0]*tilt);
+
+  if(1./PAR[4]>fPtThreshold) ((TH3S*)arr->At(0))->Fill(fTracklet->GetYref(1), 1.e2*dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]); // scale to fit general residual range !!!
+  ((TH3S*)arr->At(2))->Fill(fTracklet->GetZref(1), dz[1], rc);
+  ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), dphi);
+
+  Double_t dyz[2] = {dy[1], dz[1]};
+  Double_t cc[3] = {COV(0,0)+cov[0], COV(0,1)+cov[1], COV(1,1)+cov[2]};
+  Pulls(dyz, cc, tilt);
+  ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
+  ((TH3S*)arr->At(3))->Fill(fTracklet->GetZref(1), dyz[1], rc);
+
+  // register reference histo for mini-task
+  h = (TH2I*)arr->At(0);
+
+  if(DebugLevel()>=2){
+    (*DebugStream()) << "trackOut"
+      << "x="       << x
+      << "P="       << &PAR
+      << "C="       << &COV
+      << "\n";
+
+    Double_t y = fTracklet->GetY(); 
+    Double_t z = fTracklet->GetZ(); 
+    (*DebugStream()) << "trackletOut"
+      << "y="       << y
+      << "z="       << z
+      << "Vy="      << cov[0]
+      << "Cyz="     << cov[1]
+      << "Vz="      << cov[2]
+      << "\n";
+  }
+
+
+  if(!HasMCdata()) return h;
+  UChar_t s;
+  Float_t dx, pt0, x0=fTracklet->GetX0(), y0, z0, dydx0, dzdx0;
+  if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) return h;
+  Int_t pdg = fkMC->GetPDG(),
+        sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index
+        sign(0);
+  if(!fDBPDG) fDBPDG=TDatabasePDG::Instance();
+  TParticlePDG *ppdg(fDBPDG->GetParticle(pdg));
+  if(ppdg) sign = ppdg->Charge() > 0. ? 1 : -1;
+
+  // translate to reference radial position
+  dx = x0 - x; y0 -= dx*dydx0; z0 -= dx*dzdx0;
+  Float_t norm = 1./TMath::Sqrt(1.+dydx0*dydx0); // 1/sqrt(1+tg^2(phi))
+  //Fill MC info
+  TVectorD PARMC(kNPAR);
+  PARMC[0]=y0; PARMC[1]=z0;
+  PARMC[2]=dydx0*norm; PARMC[3]=dzdx0*norm;
+  PARMC[4]=1./pt0;
+
+//   TMatrixDSymEigen eigen(COV);
+//   TVectorD evals = eigen.GetEigenValues();
+//   TMatrixDSym evalsm(kNPAR);
+//   for(Int_t ir=0; ir<kNPAR; ir++) for(Int_t ic=0; ic<kNPAR; ic++) evalsm(ir,ic) = (ir==ic ? evals(ir): 0.);
+//   TMatrixD evecs = eigen.GetEigenVectors();
+//   TMatrixD sqrcov(evecs, TMatrixD::kMult, TMatrixD(evalsm, TMatrixD::kMult, evecs.T()));
+  
+  // fill histos
+  if(!(arr = (TObjArray*)fContainer->At(kMCtrackOut))){
+    AliWarning("No MC container defined.");
+    return h;
+  }
+  // y resolution/pulls
+  if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0], sgm[fSegmentLevel]);
+  ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0)), (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)));
+  // z resolution/pulls
+  ((TH3S*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1], 0);
+  ((TH3S*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)), 0);
+  // phi resolution/snp pulls
+  ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ASin(PARMC[2])-TMath::ASin(PAR[2]));
+  ((TH2I*)arr->At(5))->Fill(dydx0, (PARMC[2]-PAR[2])/TMath::Sqrt(COV(2,2)));
+  // theta resolution/tgl pulls
+  ((TH2I*)arr->At(6))->Fill(dzdx0, TMath::ATan((PARMC[3]-PAR[3])/(1-PARMC[3]*PAR[3])));
+  ((TH2I*)arr->At(7))->Fill(dzdx0, (PARMC[3]-PAR[3])/TMath::Sqrt(COV(3,3)));
+  // pt resolution\\1/pt pulls\\p resolution/pull
+  ((TH3S*)arr->At(8))->Fill(pt0, PARMC[4]/PAR[4]-1., sign*sIdx);
+  ((TH3S*)arr->At(9))->Fill(PARMC[4], (PARMC[4]-PAR[4])/TMath::Sqrt(COV(4,4)), sign*sIdx);
+
+  Double_t p0 = TMath::Sqrt(1.+ PARMC[3]*PARMC[3])*pt0, p;
+  p = TMath::Sqrt(1.+ PAR[3]*PAR[3])/PAR[4];
+  ((TH3S*)arr->At(10))->Fill(p0, p/p0-1., sign*sIdx);
+//   Float_t sp = 
+//     p*p*PAR[4]*PAR[4]*COV(4,4)
+//    +2.*PAR[3]*COV(3,4)/PAR[4]
+//    +PAR[3]*PAR[3]*COV(3,3)/p/p/PAR[4]/PAR[4]/PAR[4]/PAR[4];
+//   if(sp>0.) ((TH3S*)arr->At(11))->Fill(p0, (p0-p)/TMath::Sqrt(sp), sign*sIdx);
+
+  // fill debug for MC 
+  if(DebugLevel()>=3){
+    (*DebugStream()) << "trackOutMC"
+      << "P="   << &PARMC
+      << "\n";
+  }
+  return h;
+}
+
 //________________________________________________________
 TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
 {
@@ -612,24 +926,27 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
   //
 
   if(!HasMCdata()){ 
-    AliWarning("No MC defined. Results will not be available.");
+    AliDebug(2, "No MC defined. Results will not be available.");
     return NULL;
   }
   if(track) fkTrack = track;
   if(!fkTrack){
-    AliDebug(2, "No Track defined.");
+    AliDebug(4, "No Track defined.");
+    return NULL;
+  }
+  if(!fContainer){
+    AliWarning("No output container defined.");
     return NULL;
   }
   // retriev track characteristics
   Int_t pdg = fkMC->GetPDG(),
         sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index
         sign(0),
-        det(-1),
+        sgm[3],
         label(fkMC->GetLabel());
   if(!fDBPDG) fDBPDG=TDatabasePDG::Instance();
   TParticlePDG *ppdg(fDBPDG->GetParticle(pdg));
   if(ppdg) sign = ppdg->Charge() > 0. ? 1 : -1;
-  Bool_t kBarrel = Bool_t(fkESD->GetStatus() & AliESDtrack::kTRDin);
 
   TObjArray *arr(NULL);TH1 *h(NULL);
   UChar_t s;
@@ -637,40 +954,46 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
   Float_t pt0, x0, y0, z0, dx, dy, dz, dydx0, dzdx0;
   Double_t covR[7]/*, cov[3]*/;
 
-  if(DebugLevel()>=1){
-    TVectorD dX(12), dY(12), dZ(12), dPt(12), cCOV(12*15);
-    fkMC->PropagateKalman(&dX, &dY, &dZ, &dPt, &cCOV);
+  if(DebugLevel()>=3){
+    TVectorD dX(12), dY(12), dZ(12), vPt(12), dPt(12), cCOV(12*15);
+    fkMC->PropagateKalman(&dX, &dY, &dZ, &vPt, &dPt, &cCOV);
     (*DebugStream()) << "MCkalman"
       << "pdg=" << pdg
       << "dx="  << &dX
       << "dy="  << &dY
       << "dz="  << &dZ
+      << "pt="  << &vPt
       << "dpt=" << &dPt
       << "cov=" << &cCOV
       << "\n";
   }
-
-  AliTRDReconstructor rec;
+  AliTRDgeometry *geo(AliTRDinfoGen::Geometry());
   AliTRDseedV1 *fTracklet(NULL); TObjArray *clInfoArr(NULL);
   for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
     if(!(fTracklet = fkTrack->GetTracklet(ily)))/* ||
        !fTracklet->IsOK())*/ continue;
 
-    det = fTracklet->GetDetector();
+    sgm[2] = fTracklet->GetDetector();
+    sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
+    sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
+    Double_t tilt(fTracklet->GetTilt())
+            ,t2(tilt*tilt)
+            ,corr(1./(1. + t2))
+            ,cost(TMath::Sqrt(corr));
     x0  = fTracklet->GetX0();
     //radial shift with respect to the MC reference (radial position of the pad plane)
     x= fTracklet->GetX();
+    Bool_t rc(fTracklet->IsRowCross());
     if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) continue;
     xAnode  = fTracklet->GetX0();
 
     // MC track position at reference radial position
     dx  = x0 - x;
-    if(DebugLevel()>=1){
+    if(DebugLevel()>=4){
       (*DebugStream()) << "MC"
-        << "det="     << det
+        << "det="     << sgm[2]
         << "pdg="     << pdg
         << "sgn="     << sign
-        << "barrel="  << kBarrel
         << "pt="      << pt0
         << "x="       << x0
         << "y="       << y0
@@ -695,13 +1018,13 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
     pt = TMath::Abs(fTracklet->GetPt());
     fTracklet->GetCovRef(covR);
 
-    arr = (TObjArray*)fContainer->At(kMCtrackTRD);
+    arr = (TObjArray*)((TObjArray*)fContainer->At(kMCtrack))->At(ily);
     // y resolution/pulls
-    ((TH2I*)arr->At(0))->Fill(dydx0, dy);
-    ((TH2I*)arr->At(1))->Fill(dydx0, dy/TMath::Sqrt(covR[0]));
+    if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]);
+    ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(covR[0]), dz/TMath::Sqrt(covR[2]));
     // z resolution/pulls
-    ((TH2I*)arr->At(2))->Fill(dzdx0, dz);
-    ((TH2I*)arr->At(3))->Fill(dzdx0, dz/TMath::Sqrt(covR[2]));
+    ((TH3S*)arr->At(2))->Fill(dzdx0, dz, 0);
+    ((TH3S*)arr->At(3))->Fill(dzdx0, dz/TMath::Sqrt(covR[2]), 0);
     // phi resolution/ snp pulls
     Double_t dtgp = (dydx - dydx0)/(1.- dydx*dydx0);
     ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dtgp));
@@ -716,18 +1039,12 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
     // pt resolution  \\ 1/pt pulls \\ p resolution for PID
     Double_t p0 = TMath::Sqrt(1.+ dzdl0*dzdl0)*pt0,
              p  = TMath::Sqrt(1.+ dzdl*dzdl)*pt;
-    if(kBarrel){
-      ((TH3S*)((TObjArray*)arr->At(8))->At(ily))->Fill(pt0, pt/pt0-1., sign*sIdx);
-      ((TH3S*)((TObjArray*)arr->At(9))->At(ily))->Fill(1./pt0, (1./pt-1./pt0)/TMath::Sqrt(covR[6]), sign*sIdx);
-      ((TH3S*)((TObjArray*)arr->At(10))->At(ily))->Fill(p0, p/p0-1., sign*sIdx);
-    } else {
-      ((TH3S*)((TObjArray*)arr->At(11))->At(ily))->Fill(pt0, pt/pt0-1., sign*sIdx);
-      ((TH3S*)((TObjArray*)arr->At(12))->At(ily))->Fill(1./pt0, (1./pt-1./pt0)/TMath::Sqrt(covR[6]), sign*sIdx);
-      ((TH3S*)((TObjArray*)arr->At(13))->At(ily))->Fill(p0, p/p0-1., sign*sIdx);
-    }
+    ((TH3S*)((TObjArray*)arr->At(8)))->Fill(pt0, pt/pt0-1., sign*sIdx);
+    ((TH3S*)((TObjArray*)arr->At(9)))->Fill(1./pt0, (1./pt-1./pt0)/TMath::Sqrt(covR[6]), sign*sIdx);
+    ((TH3S*)((TObjArray*)arr->At(10)))->Fill(p0, p/p0-1., sign*sIdx);
 
     // Fill Debug stream for Kalman track
-    if(DebugLevel()>=1){
+    if(DebugLevel()>=4){
       (*DebugStream()) << "MCtrack"
         << "pt="      << pt
         << "x="       << x
@@ -744,35 +1061,29 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
     AliTRDseedV1 tt(*fTracklet);
     tt.SetZref(0, z0 - (x0-xAnode)*dzdx0);
     tt.SetZref(1, dzdx0); 
-    tt.SetReconstructor(&rec);
-    tt.Fit(kTRUE, kTRUE);
+    tt.SetReconstructor(AliTRDinfoGen::Reconstructor());
+    tt.Fit(1);
     x= tt.GetX();y= tt.GetY();z= tt.GetZ();
     dydx = tt.GetYfit(1);
     dx = x0 - x;
     ymc = y0 - dx*dydx0;
     zmc = z0 - dx*dzdx0;
-    Bool_t rc = tt.IsRowCross(); 
-    
+    dy = y-ymc;
+    dz = z-zmc;
+    Float_t dphi  = (dydx - dydx0);
+    dphi /= (1.- dydx*dydx0);
+
     // add tracklet residuals for y and dydx
     arr = (TObjArray*)fContainer->At(kMCtracklet);
-    if(!rc){
-      dy    = y-ymc;
 
-      Float_t dphi  = (dydx - dydx0);
-      dphi /= (1.- dydx*dydx0);
-
-      ((TH3S*)arr->At(0))->Fill(dydx0, dy, pt0);
-      if(tt.GetS2Y()>0.) ((TH2I*)arr->At(1))->Fill(dydx0, dy/TMath::Sqrt(tt.GetS2Y()));
-      ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dphi));
-    } else {
-      // add tracklet residuals for z
-      dz = z-zmc;
-      ((TH2I*)arr->At(2))->Fill(dzdl0, dz);
-      if(tt.GetS2Z()>0.) ((TH2I*)arr->At(3))->Fill(dzdl0, dz/TMath::Sqrt(tt.GetS2Z()));
-    }
+    if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]);
+    if(tt.GetS2Y()>0. && tt.GetS2Z()>0.) ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(tt.GetS2Y()), dz/TMath::Sqrt(tt.GetS2Z()));
+    ((TH3S*)arr->At(2))->Fill(dzdl0, dz, rc);
+    if(tt.GetS2Z()>0.) ((TH3S*)arr->At(3))->Fill(dzdl0, dz/TMath::Sqrt(tt.GetS2Z()), rc);
+    ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dphi));
   
     // Fill Debug stream for tracklet
-    if(DebugLevel()>=1){
+    if(DebugLevel()>=4){
       Float_t s2y = tt.GetS2Y();
       Float_t s2z = tt.GetS2Z();
       (*DebugStream()) << "MCtracklet"
@@ -786,26 +1097,26 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
         << "\n";
     }
 
-    Int_t istk = AliTRDgeometry::GetStack(det); 
-    AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
+    AliTRDpadPlane *pp = geo->GetPadPlane(ily, AliTRDgeometry::GetStack(sgm[2]));
     Float_t zr0 = pp->GetRow0() + pp->GetAnodeWireOffset();
-    Float_t tilt = fTracklet->GetTilt();
     //Double_t exb = AliTRDCommonParam::Instance()->GetOmegaTau(1.5);
 
     arr = (TObjArray*)fContainer->At(kMCcluster);
     AliTRDcluster *c = NULL;
-    fTracklet->ResetClusterIter(kFALSE);
-    while((c = fTracklet->PrevCluster())){
+    tt.ResetClusterIter(kFALSE);
+    while((c = tt.PrevCluster())){
       Float_t  q = TMath::Abs(c->GetQ());
       x = c->GetX(); y = c->GetY();z = c->GetZ();
       dx = x0 - x; 
       ymc= y0 - dx*dydx0;
       zmc= z0 - dx*dzdx0;
-      dy = (y - tilt*(z-zmc)) - ymc;
+      dy = cost*(y - ymc - tilt*(z-zmc));
+      dz = cost*(z - zmc + tilt*(y-ymc));
+      
       // Fill Histograms
-      if(q>20. && q<250.){ 
-        ((TH3S*)arr->At(0))->Fill(dydx0, dy, pt0);
-        ((TH2I*)arr->At(1))->Fill(dydx0, dy/TMath::Sqrt(c->GetSigmaY2()));
+      if(q>20. && q<250. && pt0>fPtThreshold && c->IsInChamber()){ 
+        ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]);
+        ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(c->GetSigmaY2()), dz/TMath::Sqrt(c->GetSigmaZ2()));
       }
 
       // Fill calibration container
@@ -818,16 +1129,20 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
       clInfo->SetGlobalPosition(ymc, zmc, dydx0, dzdx0);
       clInfo->SetResolution(dy);
       clInfo->SetAnisochronity(d);
-      clInfo->SetDriftLength(dx-.5*AliTRDgeometry::CamHght());
+      clInfo->SetDriftLength(dx);
       clInfo->SetTilt(tilt);
-      fMCcl->Add(clInfo);
-      if(DebugLevel()>=2){ 
-        if(!clInfoArr) clInfoArr=new TObjArray(AliTRDseedV1::kNclusters);
+      if(fMCcl) fMCcl->Add(clInfo);
+      else AliDebug(1, "MCcl exchange container missing. Activate by calling \"InitExchangeContainers()\"");
+      if(DebugLevel()>=5){ 
+        if(!clInfoArr){ 
+          clInfoArr=new TObjArray(AliTRDseedV1::kNclusters);
+          clInfoArr->SetOwner(kFALSE);
+        }
         clInfoArr->Add(clInfo);
       }
     }
     // Fill Debug Tree
-    if(DebugLevel()>=2 && clInfoArr){
+    if(DebugLevel()>=5 && clInfoArr){
       (*DebugStream()) << "MCcluster"
         <<"clInfo.=" << clInfoArr
         << "\n";
@@ -851,367 +1166,961 @@ Bool_t AliTRDresolution::GetRefFigure(Int_t ifig)
     AliWarning("Please provide a canvas to draw results.");
     return kFALSE;
   }
-  TList *l = NULL; TVirtualPad *pad=NULL;
+  Int_t selection[100], n(0), selStart(0); // 
+  Int_t ly0(0), dly(5);
+  //Int_t ly0(1), dly(2); // used for SA
+  TList *l(NULL); TVirtualPad *pad(NULL); 
+  TGraphErrors *g(NULL);TGraphAsymmErrors *ga(NULL);
   switch(ifig){
-  case kCharge:
+  case 0: // charge resolution
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     ((TVirtualPad*)l->At(0))->cd();
-    ((TGraphAsymmErrors*)((TObjArray*)fGraphM->At(kCharge))->At(0))->Draw("apl");
+    ga=((TGraphAsymmErrors*)((TObjArray*)fGraphM->At(kCharge))->At(0));
+    if(ga->GetN()) ga->Draw("apl");
     ((TVirtualPad*)l->At(1))->cd();
-    ((TGraphErrors*)((TObjArray*)fGraphS->At(kCharge))->At(0))->Draw("apl");
+    g = ((TGraphErrors*)((TObjArray*)fGraphS->At(kCharge))->At(0));
+    if(g->GetN()) g->Draw("apl");
     break;
-  case kCluster:
+  case 1: // cluster2track residuals
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0] = -.3; xy[1] = -100.; xy[2] = .3; xy[3] = 1000.;
+    pad = (TVirtualPad*)l->At(0); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kCluster, 0, 1, n, selection)) break;
+    pad=(TVirtualPad*)l->At(1); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kCluster, 0, 1, n, selection)) break;
+    return kTRUE;
+  case 2: // cluster2track residuals
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0] = -.3; xy[1] = -200.; xy[2] = .3; xy[3] = 1000.;
+    xy[0] = -.3; xy[1] = -100.; xy[2] = .3; xy[3] = 1000.;
+    pad = (TVirtualPad*)l->At(0); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kCluster, 0, 1, n, selection)) break;
+    xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-0.5; xy[3] = 2.5;
+    pad=(TVirtualPad*)l->At(1); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    if(!GetGraphArray(xy, kCluster, 1, 1)) break;
+    return kTRUE;
+  case 3: // kTrack y
+    gPad->Divide(3, 2, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0] = -.3; xy[1] = -20.; xy[2] = .3; xy[3] = 100.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kCluster, 0)) break;
-    xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5;
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrack, 0, 1, n, selection)) break;
+
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kCluster, 1)) break;
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrack, 0, 1, n, selection)) break;
+
+    ((TVirtualPad*)l->At(2))->cd();
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrack, 0, 1, n, selection)) break;
+
+    ((TVirtualPad*)l->At(3))->cd();
+    selStart=fgkNresYsegm[fSegmentLevel]; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrack, 0, 1, n, selection, "[RC]")) break;
+
+    ((TVirtualPad*)l->At(4))->cd();
+    selStart=fgkNresYsegm[fSegmentLevel]/3+fgkNresYsegm[fSegmentLevel]; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrack, 0, 1, n, selection, "[RC]")) break;
+
+    ((TVirtualPad*)l->At(5))->cd();
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3+fgkNresYsegm[fSegmentLevel]; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrack, 0, 1, n, selection, "[RC]")) break;
     return kTRUE;
-  case kTrackTRD :
+  case 4: // kTrack z
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0] = -.3; xy[1] = -500.; xy[2] = .3; xy[3] = 1500.;
+
+    xy[0] = -1.; xy[1] = -150.; xy[2] = 1.; xy[3] = 1000.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kTrackTRD , 0)) break;
-    xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5;
+    selection[0]=1;
+    if(!GetGraphArray(xy, kTrack, 2, 1, 1, selection)) break;
+
+    xy[0] = -1.; xy[1] = -1500.; xy[2] = 1.; xy[3] = 10000.;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kTrackTRD , 1)) break;
+    selection[0]=0;
+    if(!GetGraphArray(xy, kTrack, 2, 1, 1, selection)) break;
+
     return kTRUE;
-  case 3: // kTrackTRD  z
+  case 5: // kTrack  pulls
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.;
+
+    xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kTrackTRD , 2)) break;
+    if(!GetGraphArray(xy, kTrack, 1, 1)) break;
+
     xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kTrackTRD , 3)) break;
+    if(!GetGraphArray(xy, kTrack, 3, 1)) break;
     return kTRUE;
-  case 4: // kTrackTRD  phi
+  case 6: // kTrack  phi
     xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.;
-    if(GetGraphPlot(&xy[0], kTrackTRD , 4)) return kTRUE;
+    if(GetGraph(&xy[0], kTrack , 4)) return kTRUE;
     break;
-  case 5: // kTrackTPC y
+  case 7: // kTrackIn y
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0] = -.3; xy[1] = -1500.; xy[2] = .3; xy[3] = 5000.;
+    pad = ((TVirtualPad*)l->At(0)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrackIn, 0, 1, n, selection)) break;
+    pad=((TVirtualPad*)l->At(1)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrackIn, 0, 1, n, selection)) break;
+    return kTRUE;
+  case 8: // kTrackIn y
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0] = -.3; xy[1] = -500.; xy[2] = .3; xy[3] = 1500.;
+    xy[0] = -.3; xy[1] = -1500.; xy[2] = .3; xy[3] = 5000.;
     pad = ((TVirtualPad*)l->At(0)); pad->cd();
     pad->SetMargin(0.1, 0.1, 0.1, 0.01);
-    if(!GetGraphPlot(&xy[0], kTrackTPC, 0)) break;
-    xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5;
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrackIn, 0, 1, n, selection)) break;
+    xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
     pad=((TVirtualPad*)l->At(1)); pad->cd();
     pad->SetMargin(0.1, 0.1, 0.1, 0.01);
-    if(!GetGraphPlot(&xy[0], kTrackTPC, 1)) break;
+    if(!GetGraphArray(xy, kTrackIn, 1, 1)) break;
     return kTRUE;
-  case 6: // kTrackTPC z
+  case 9: // kTrackIn z
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.;
     pad = ((TVirtualPad*)l->At(0)); pad->cd();
     pad->SetMargin(0.1, 0.1, 0.1, 0.01);
-    if(!GetGraphPlot(&xy[0], kTrackTPC, 2)) break;
+    selection[0]=1;
+    if(!GetGraphArray(xy, kTrackIn, 2, 1, 1, selection)) break;
     xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
     pad = ((TVirtualPad*)l->At(1)); pad->cd();
     pad->SetMargin(0.1, 0.1, 0.1, 0.01);
-    if(!GetGraphPlot(&xy[0], kTrackTPC, 3)) break;
+    if(!GetGraphArray(xy, kTrackIn, 3, 1)) break;
     return kTRUE;
-  case 7: // kTrackTPC phi
+  case 10: // kTrackIn phi
     xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.;
-    if(GetGraphPlot(&xy[0], kTrackTPC, 4)) return kTRUE;
+    if(GetGraph(&xy[0], kTrackIn, 4)) return kTRUE;
     break;
-  case 8: // kMCcluster pt <2 GeV/c
-    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.;
-    ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCcluster, 0)) break;
-    ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCcluster, 1)) break;
-    return kTRUE;
-  case 9: // kMCcluster pt > 3 GeV/c
+  case 11: // kTrackOut y
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.;
-    ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCcluster, 2)) break;
-    xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5;
-    ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCcluster, 3)) break;
+    xy[0] = -.3; xy[1] = -50.; xy[2] = .3; xy[3] = 150.;
+    pad = ((TVirtualPad*)l->At(0)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrackOut, 0, 1, n, selection)) break;
+    pad=((TVirtualPad*)l->At(1)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrackOut, 0, 1, n, selection)) break;
     return kTRUE;
-  case 10: //kMCtracklet [y] pt < 3 GeV/c
+  case 12: // kTrackOut y
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =250.;
-    ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtracklet, 0)) break;
-    ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtracklet, 1)) break;
+    xy[0] = -.3; xy[1] = -50.; xy[2] = .3; xy[3] = 150.;
+    pad = ((TVirtualPad*)l->At(0)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kTrackOut, 0, 1, n, selection)) break;
+    xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
+    pad=((TVirtualPad*)l->At(1)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    if(!GetGraphArray(xy, kTrackOut, 1, 1)) break;
     return kTRUE;
-  case 11: //kMCtracklet [y] pt > 3 GeV/c
+  case 13: // kTrackOut z
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =250.;
-    ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtracklet, 2)) break;
-    xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5;
-    ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtracklet, 3)) break;
+    xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.;
+    pad = ((TVirtualPad*)l->At(0)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    if(!GetGraphArray(xy, kTrackOut, 2, 1)) break;
+    xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
+    pad = ((TVirtualPad*)l->At(1)); pad->cd();
+    pad->SetMargin(0.1, 0.1, 0.1, 0.01);
+    if(!GetGraphArray(xy, kTrackOut, 3, 1)) break;
+    return kTRUE;
+  case 14: // kTrackOut phi
+    xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.;
+    if(GetGraph(&xy[0], kTrackOut, 4)) return kTRUE;
+    break;
+  case 15: // kMCcluster
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCcluster, 0, 1, n, selection)) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCcluster, 0, 1, n, selection)) break;
+    return kTRUE;
+  case 16: // kMCcluster
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCcluster, 0, 1, n, selection)) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
+    if(!GetGraphArray(xy, kMCcluster, 1, 1)) break;
+    return kTRUE;
+  case 17: //kMCtracklet [y]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =500.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtracklet, 0, 1, n, selection)) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtracklet, 0, 1, n, selection)) break;
     return kTRUE;
-  case 12: //kMCtracklet [z]
+  case 18: //kMCtracklet [y]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =500.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtracklet, 0, 1, n, selection)) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
+    if(!GetGraphArray(xy, kMCtracklet, 1, 1)) break;
+    return kTRUE;
+  case 19: //kMCtracklet [z]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     xy[0]=-1.; xy[1]=-100.; xy[2]=1.; xy[3] =2500.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtracklet, 4)) break;
+    if(!GetGraphArray(xy, kMCtracklet, 2)) break;
     xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtracklet, 5)) break;
+    if(!GetGraphArray(xy, kMCtracklet, 3)) break;
     return kTRUE;
-  case 13: //kMCtracklet [phi]
+  case 20: //kMCtracklet [phi]
     xy[0]=-.3; xy[1]=-3.; xy[2]=.3; xy[3] =25.;
-    if(!GetGraphPlot(&xy[0], kMCtracklet, 6)) break;
+    if(!GetGraph(&xy[0], kMCtracklet, 4)) break;
     return kTRUE;
-  case 14: //kMCtrackTRD [y]
+  case 21: //kMCtrack [y] ly [0]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =200.;
+    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 0)) break;
-    xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5;
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*0.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer1")) break;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 1)) break;
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*0.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer1")) break;
     return kTRUE;
-  case 15: //kMCtrackTRD [z]
+  case 22: //kMCtrack [y] ly [1]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-1.; xy[1]=-700.; xy[2]=1.; xy[3] =1500.;
+    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 2)) break;
-    xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*1.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer2")) break;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 3)) break;
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*1.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer2")) break;
     return kTRUE;
-  case 16: //kMCtrackTRD [phi/snp]
+  case 23: //kMCtrack [y] ly [2]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-.2; xy[1]=-0.2; xy[2]=.2; xy[3] =2.;
+    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 4)) break;
-    xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5;
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*2.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer3")) break;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 5)) break;
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*2.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer3")) break;
     return kTRUE;
-  case 17: //kMCtrackTRD [theta/tgl]
+  case 24: //kMCtrack [y] ly [3]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*3.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer4")) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*3.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer4")) break;
+    return kTRUE;
+  case 25: //kMCtrack [y] ly [4]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*4.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer5")) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*4.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer5")) break;
+    return kTRUE;
+  case 26: //kMCtrack [y] ly [5]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*5.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer6")) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    selStart=Int_t(fgkNresYsegm[fSegmentLevel]*5.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer6")) break;
+    return kTRUE;
+  case 27: //kMCtrack [y pulls] 
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 5.5;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=0; for(n=0; n<6; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 1, 1, n, selection)) break;
+    ((TVirtualPad*)l->At(1))->cd();
+    selStart=6; for(n=0; n<6; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrack, 1, 1, n, selection)) break;
+    return kTRUE;
+  case 28: //kMCtrack [z]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-1.; xy[1]=-1500.; xy[2]=1.; xy[3] =6000.;
+    ((TVirtualPad*)l->At(0))->cd();
+    if(!GetGraphArray(xy, kMCtrack, 2)) break;
+    xy[0] = -1.; xy[1] = -1.5; xy[2] = 1.; xy[3] = 5.;
+    ((TVirtualPad*)l->At(1))->cd();
+    if(!GetGraphArray(xy, kMCtrack, 3)) break;
+    return kTRUE;
+  case 29: //kMCtrack [phi/snp]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.2; xy[1]=-0.5; xy[2]=.2; xy[3] =10.;
+    ((TVirtualPad*)l->At(0))->cd();
+    if(!GetGraphArray(xy, kMCtrack, 4)) break;
+    xy[0] = -.2; xy[1] = -1.5; xy[2] = .2; xy[3] = 5.;
+    ((TVirtualPad*)l->At(1))->cd();
+    if(!GetGraphArray(xy, kMCtrack, 5)) break;
+    return kTRUE;
+  case 30: //kMCtrack [theta/tgl]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     xy[0]=-1.; xy[1]=-0.5; xy[2]=1.; xy[3] =5.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 6)) break;
+    if(!GetGraphArray(xy, kMCtrack, 6)) break;
     xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTRD, 7)) break;
+    if(!GetGraphArray(xy, kMCtrack, 7)) break;
     return kTRUE;
-  case 18: //kMCtrackTRD [pt]
-    xy[0] = 0.4; xy[1] = -.7; xy[2] = 12.; xy[3] = 4.;
+  case 31: //kMCtrack [pt]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     pad = (TVirtualPad*)l->At(0); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 8, 2, 7)) break;
-    pad->Modified(); pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // pi selection
+    n=0;
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 2; // pi-
+      selection[n++] = il*11 + 8; // pi+
+    }
+    xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 4.;
+    //xy[0] = 0.2; xy[1] = -1.; xy[2] = 7.; xy[3] = 10.; // SA
+    if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "#pi#pm")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     pad = (TVirtualPad*)l->At(1); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 8, 55+2, 7, kTRUE)) break;
-    pad->Modified(); pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // mu selection
+    n=0;    
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 3; // mu-
+      selection[n++] = il*11 + 7; // mu+
+    }
+    if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "#mu#pm")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     return kTRUE;
-  case 19: //kMCtrackTRD [1/pt] pulls
-    xy[0] = 0.; xy[1] = -1.5; xy[2] = 2.; xy[3] = 2.;
+  case 32: //kMCtrack [pt]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     pad = (TVirtualPad*)l->At(0); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 9, 2, 7)) break;
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // p selection
+    n=0;
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 0; // p bar
+      selection[n++] = il*11 + 10; // p
+    }
+    xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 8.;
+    //xy[0] = 0.2; xy[1] = -1.; xy[2] = 7.; xy[3] = 10.; // SA
+    if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "p&p bar")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     pad = (TVirtualPad*)l->At(1); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 9, 55+2, 7, kTRUE)) break;
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // e selection
+    n=0;    
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 4; // e-
+      selection[n++] = il*11 + 6; // e+
+    }
+    xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 12.;
+    //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 14.; // SA
+    if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "e#pm")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     return kTRUE;
-  case 20: //kMCtrackTRD [p]
-    xy[0] = 0.4; xy[1] = -.7; xy[2] = 12.; xy[3] = 4.;
+  case 33: //kMCtrack [1/pt] pulls
+    xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
+    //xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 4.5; // SA
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     pad = (TVirtualPad*)l->At(0); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 10, 2, 7)) break;
-    pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // pi selection
+    n=0;
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 2; // pi-
+      selection[n++] = il*11 + 8; // pi+
+    }
+    if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "#pi#pm")) break;
     pad = (TVirtualPad*)l->At(1); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 10, 55+2, 7, kTRUE)) break;
-    pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // mu selection
+    n=0;    
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 3; // mu-
+      selection[n++] = il*11 + 7; // mu+
+    }
+    if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "#mu#pm")) break;
     return kTRUE;
-  case 21: //kMCtrackTRD - SA [pt]
-    xy[0] = 0.; xy[1] = -5.; xy[2] = 12.; xy[3] = 7.;
+  case 34: //kMCtrack [1/pt] pulls
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     pad = (TVirtualPad*)l->At(0); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 11, 2, 7)) break;
-    pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // p selection
+    n=0;
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 0; // p bar
+      selection[n++] = il*11 + 10; // p
+    }
+    xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
+    //xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 6.; // SA
+    if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "p & p bar")) break;
     pad = (TVirtualPad*)l->At(1); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 11, 55+2, 7, kTRUE)) break;
-    pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // e selection
+    n=0;    
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 4; // e-
+      selection[n++] = il*11 + 6; // e+
+    }
+    xy[0] = 0.; xy[1] = -2.; xy[2] = 2.; xy[3] = 4.5;
+    if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "e#pm")) break;
     return kTRUE;
-  case 22: //kMCtrackTRD - SA [1/pt] pulls
-    xy[0] = 0.; xy[1] = -1.5; xy[2] = 2.; xy[3] = 2.;
+  case 35: //kMCtrack [p]
+    xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 4.;
+    //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.;
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     pad = (TVirtualPad*)l->At(0); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 12, 2, 7)) break;
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // pi selection
+    n=0;
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 2; // pi-
+      selection[n++] = il*11 + 8; // pi+
+    }
+    if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "#pi#pm")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     pad = (TVirtualPad*)l->At(1); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 12, 55+2, 7, kTRUE)) break;
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // mu selection
+    n=0;    
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 3; // mu-
+      selection[n++] = il*11 + 7; // mu+
+    }
+    if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "#mu#pm")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     return kTRUE;
-  case 23: //kMCtrackTRD - SA [p]
-    xy[0] = 0.; xy[1] = -7.5; xy[2] = 12.; xy[3] = 10.5;
+  case 36: //kMCtrack [p]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     pad = (TVirtualPad*)l->At(0); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 13, 2, 7)) break;
-    pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // p selection
+    n=0;
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 0; // p bar
+      selection[n++] = il*11 + 10; // p
+    }
+    xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 8.;
+    //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 12.; // SA
+    if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "p & p bar")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     pad = (TVirtualPad*)l->At(1); pad->cd();
-    pad->SetMargin(0.07, 0.07, 0.1, 0.1);
-    if(!GetGraphTrack(&xy[0], 13, 55+2, 7, kTRUE)) break;
-    pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    // e selection
+    n=0;    
+    for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
+      selection[n++] = il*11 + 4; // e-
+      selection[n++] = il*11 + 6; // e+
+    }
+    xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 12.;
+    //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 14.; // SA
+    if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "e#pm")) break;
+    pad->Modified(); pad->Update(); pad->SetLogx();
     return kTRUE;
-  case 24: // kMCtrackTPC [y]
+  case 37: // kMCtrackIn [y]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0]=-.25; xy[1]=-50.; xy[2]=.25; xy[3] =800.;
+    xy[0]=-.25; xy[1]=-1000.; xy[2]=.25; xy[3] =3000.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 0)) break;
-    xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 2.5;
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrackIn, 0, 1, n, selection)) break;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 1)) break;
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(&xy[0], kMCtrackIn, 0, 1, n, selection)) break;
     return kTRUE;
-  case 25: // kMCtrackTPC [z]
+  case 38: // kMCtrackIn [y]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.25; xy[1]=-1000.; xy[2]=.25; xy[3] =3000.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrackIn, 0, 1, n, selection)) break;
+    xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
+    ((TVirtualPad*)l->At(1))->cd();
+    if(!GetGraphArray(xy, kMCtrackIn, 1, 1)) break;
+    return kTRUE;
+  case 39: // kMCtrackIn [z]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     xy[0]=-1.; xy[1]=-500.; xy[2]=1.; xy[3] =800.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 2)) break;
+    if(!GetGraphArray(xy, kMCtrackIn, 2, 1)) break;
     xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 3)) break;
+    if(!GetGraphArray(xy, kMCtrackIn, 3, 1)) break;
     return kTRUE;
-  case 26: // kMCtrackTPC [phi|snp]
+  case 40: // kMCtrackIn [phi|snp]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     xy[0]=-.25; xy[1]=-0.5; xy[2]=.25; xy[3] =2.5;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 4)) break;
+    if(!GetGraph(&xy[0], kMCtrackIn, 4)) break;
     xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 1.5;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 5)) break;
+    if(!GetGraph(&xy[0], kMCtrackIn, 5)) break;
     return kTRUE;
-  case 27: // kMCtrackTPC [theta|tgl]
+  case 41: // kMCtrackIn [theta|tgl]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
     xy[0]=-1.; xy[1]=-1.; xy[2]=1.; xy[3] =4.;
     ((TVirtualPad*)l->At(0))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 6)) break;
+    if(!GetGraph(&xy[0], kMCtrackIn, 6)) break;
     xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 1.5;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphPlot(&xy[0], kMCtrackTPC, 7)) break;
+    if(!GetGraph(&xy[0], kMCtrackIn, 7)) break;
+    return kTRUE;
+  case 42: // kMCtrackIn [pt]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
+    //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.; // SA
+    pad=(TVirtualPad*)l->At(0); pad->cd(); pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
+    if(!GetGraphArray(xy, kMCtrackIn, 8, 1, n, selection)) break;
+    pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
+    if(!GetGraphArray(xy, kMCtrackIn, 8, 1, n, selection)) break;
+    return kTRUE;
+  case 43: //kMCtrackIn [1/pt] pulls
+    xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    pad = (TVirtualPad*)l->At(0); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
+    if(!GetGraphArray(xy, kMCtrackIn, 9, 1, n, selection)) break;
+    pad = (TVirtualPad*)l->At(1); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
+    if(!GetGraphArray(xy, kMCtrackIn, 9, 1, n, selection)) break;
+    return kTRUE;
+  case 44: // kMCtrackIn [p]
+    xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
+    //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.;
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    pad = ((TVirtualPad*)l->At(0));pad->cd();pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
+    if(!GetGraphArray(xy, kMCtrackIn, 10, 1, n, selection)) break;
+    pad = ((TVirtualPad*)l->At(1)); pad->cd();pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
+    if(!GetGraphArray(xy, kMCtrackIn, 10, 1,  n, selection)) break;
     return kTRUE;
-  case 28: // kMCtrackTPC [pt]
+  case 45: // kMCtrackOut [y]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0] = 0.4; xy[1] = -.8; xy[2] = 12.; xy[3] = 6.;
-    pad=(TVirtualPad*)l->At(0); pad->cd();
-    if(!GetGraphTrackTPC(xy, 8, 2, 7)) break;
-    pad->SetLogx();
-    xy[0]=0.; xy[1]=-0.5; xy[2]=2.; xy[3] =2.5;
+    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =400.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrackOut, 0, 1, n, selection)) break;
     ((TVirtualPad*)l->At(1))->cd();
-    if(!GetGraphTrackTPC(xy, 9, 2, 7, kTRUE)) break;
+    selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(&xy[0], kMCtrackOut, 0, 1, n, selection)) break;
     return kTRUE;
-  case 29: // kMCtrackTPC [p]
+  case 46: // kMCtrackOut [y]
     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
-    xy[0] = 0.4; xy[1] = -.8; xy[2] = 12.; xy[3] = 6.;
-    pad = ((TVirtualPad*)l->At(0));pad->cd();
-    pad->SetMargin(0.12, 0.12, 0.1, 0.04);
-    if(!GetGraphTrackTPC(xy, 10, 2, 7)) break;
-    pad->SetLogx();
-    xy[0]=0.; xy[1]=-1.5; xy[2]=12.; xy[3] =2.5;
-    pad = ((TVirtualPad*)l->At(1)); pad->cd();
-    pad->SetMargin(0.12, 0.12, 0.1, 0.04);
-    if(!GetGraphTrackTPC(xy, 11, 2, 7, kTRUE)) break;
+    xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =400.;
+    ((TVirtualPad*)l->At(0))->cd();
+    selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
+    if(!GetGraphArray(xy, kMCtrackOut, 0, 1, n, selection)) break;
+    xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
+    ((TVirtualPad*)l->At(1))->cd();
+    if(!GetGraphArray(xy, kMCtrackOut, 1, 1)) break;
     return kTRUE;
-  case 30:  // kMCtrackTOF [z]
+  case 47: // kMCtrackOut [z]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-1.; xy[1]=-500.; xy[2]=1.; xy[3] =1500.;
+    ((TVirtualPad*)l->At(0))->cd();
+    if(!GetGraphArray(xy, kMCtrackOut, 2, 1)) break;
+    xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
+    ((TVirtualPad*)l->At(1))->cd();
+    if(!GetGraphArray(xy, kMCtrackOut, 3, 1)) break;
+    return kTRUE;
+  case 48: // kMCtrackOut [phi|snp]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-.25; xy[1]=-0.5; xy[2]=.25; xy[3] =2.5;
+    ((TVirtualPad*)l->At(0))->cd();
+    if(!GetGraph(&xy[0], kMCtrackOut, 4)) break;
+    xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 1.5;
+    ((TVirtualPad*)l->At(1))->cd();
+    if(!GetGraph(&xy[0], kMCtrackOut, 5)) break;
+    return kTRUE;
+  case 49: // kMCtrackOut [theta|tgl]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0]=-1.; xy[1]=-1.; xy[2]=1.; xy[3] =4.;
+    ((TVirtualPad*)l->At(0))->cd();
+    if(!GetGraph(&xy[0], kMCtrackOut, 6)) break;
+    xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 15.;
+    ((TVirtualPad*)l->At(1))->cd();
+    if(!GetGraph(&xy[0], kMCtrackOut, 7)) break;
+    return kTRUE;
+  case 50: // kMCtrackOut [pt]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
+    pad=(TVirtualPad*)l->At(0); pad->cd(); pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
+    if(!GetGraphArray(xy, kMCtrackOut, 8, 1, n, selection)) break;
+    pad = (TVirtualPad*)l->At(1); pad->cd();pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
+    if(!GetGraphArray(xy, kMCtrackOut, 8, 1, n, selection)) break;
+    return kTRUE;
+  case 51: //kMCtrackOut [1/pt] pulls
+    xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    pad = (TVirtualPad*)l->At(0); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
+    if(!GetGraphArray(xy, kMCtrackOut, 9, 1, n, selection)) break;
+    pad = (TVirtualPad*)l->At(1); pad->cd();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
+    if(!GetGraphArray(xy, kMCtrackOut, 9, 1, n, selection)) break;
+    return kTRUE;
+  case 52: // kMCtrackOut [p]
+    gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); 
+    xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
+    pad = ((TVirtualPad*)l->At(0));pad->cd();pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
+    if(!GetGraphArray(xy, kMCtrackOut, 10, 1, n, selection)) break;
+    pad = ((TVirtualPad*)l->At(1)); pad->cd();pad->SetLogx();
+    pad->SetMargin(0.125, 0.015, 0.1, 0.015);
+    n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
+    if(!GetGraphArray(xy, kMCtrackOut, 10, 1, n, selection)) break;
     return kTRUE;
   }
   AliWarning(Form("Reference plot [%d] missing result", ifig));
   return kFALSE;
 }
 
-Char_t const *fgParticle[11]={
-  " p bar", " K -", " #pi -", " #mu -", " e -",
-  " No PID",
-  " e +", " #mu +", " #pi +", " K +", " p",
-};
+//________________________________________________________
+void AliTRDresolution::MakeSummary()
+{
+// Build summary plots
+
+  if(!fGraphS || !fGraphM){
+    AliError("Missing results");
+    return;
+  }  
+  Float_t xy[4] = {0., 0., 0., 0.};
+  Float_t range[2];
+  TH2 *h2 = new TH2I("h2SF", "", 20, -.2, .2, fgkNresYsegm[fSegmentLevel], -0.5, fgkNresYsegm[fSegmentLevel]-0.5);
+  h2->GetXaxis()->CenterTitle();
+  h2->GetYaxis()->CenterTitle();
+  h2->GetZaxis()->CenterTitle();h2->GetZaxis()->SetTitleOffset(1.4);
+
+  Int_t ih2(0), iSumPlot(0);
+  TCanvas *cOut = new TCanvas(Form("TRDsummary%s_%d", GetName(), iSumPlot++), "Cluster & Tracklet Resolution", 1024, 768);
+  cOut->Divide(3,2, 2.e-3, 2.e-3, kYellow-7);
+  TVirtualPad *p(NULL);
+
+  p=cOut->cd(1); 
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetTitle(Form("Cluster-Track R-Phi Residuals;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphS->At(kCluster))->At(0), h2);
+  GetRange(h2, 1, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(2); 
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetTitle(Form("Cluster-Track R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphM->At(kCluster))->At(0), h2);
+  GetRange(h2, 0, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(3); 
+  p->SetRightMargin(0.06);p->SetTopMargin(0.06);
+  xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
+  GetGraphArray(xy, kCluster, 1, 1);
+
+  p=cOut->cd(4); 
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetTitle(Form("Tracklet-Track R-Phi Residuals;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphS->At(kTrack))->At(0), h2);
+  GetRange(h2, 1, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(5); 
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetTitle(Form("Tracklet-Track R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphM->At(kTrack))->At(0), h2);
+  GetRange(h2, 0, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(6); 
+  p->SetRightMargin(0.06);p->SetTopMargin(0.06);
+  xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
+  GetGraphArray(xy, kTrack, 1, 1);
+
+  cOut->SaveAs(Form("%s.gif", cOut->GetName()));
+
+  if(!HasMCdata() ||
+    (!fGraphS->At(kMCcluster) || !fGraphM->At(kMCcluster) ||
+     !fGraphS->At(kMCtracklet) || !fGraphM->At(kMCtracklet))){
+    delete cOut;
+    return;
+  }
+  cOut->Clear(); cOut->SetName(Form("TRDsummary%s_%d", GetName(), iSumPlot++));
+  cOut->Divide(3, 2, 2.e-3, 2.e-3, kBlue-10);
+
+  p=cOut->cd(1);  
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetTitle(Form("Cluster-MC R-Phi Resolution;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphS->At(kMCcluster))->At(0), h2);
+  GetRange(h2, 1, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(2);  
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetContour(7);
+  h2->SetTitle(Form("Cluster-MC R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphM->At(kMCcluster))->At(0), h2);
+  GetRange(h2, 0, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(3); 
+  p->SetRightMargin(0.06);p->SetTopMargin(0.06);
+  xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
+  GetGraphArray(xy, kMCcluster, 1, 1);
+
+  p=cOut->cd(4); 
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetContour(7);
+  h2->SetTitle(Form("Tracklet-MC R-Phi Resolution;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphS->At(kMCtracklet))->At(0), h2);
+  GetRange(h2, 1, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(5); 
+  p->SetRightMargin(0.16);p->SetTopMargin(0.06);
+  h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
+  h2->SetContour(7);
+  h2->SetTitle(Form("Tracklet-MC R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
+  MakeSummaryPlot((TObjArray*)  ((TObjArray*)fGraphM->At(kMCtracklet))->At(0), h2);
+  GetRange(h2, 0, range);
+  h2->GetZaxis()->SetRangeUser(range[0], range[1]);
+  h2->Draw("colz");
+  h2->SetContour(7);
+
+  p=cOut->cd(6); 
+  p->SetRightMargin(0.06);p->SetTopMargin(0.06);
+  xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
+  GetGraphArray(xy, kMCtracklet, 1, 1);
+
+  cOut->SaveAs(Form("%s.gif", cOut->GetName()));
+  delete cOut;
+}
+
+//________________________________________________________
+void AliTRDresolution::GetRange(TH2 *h2, Char_t mod, Float_t *range)
+{
+// Returns the range of the bulk of data in histogram h2. Removes outliers. 
+// The "range" vector should be initialized with 2 elements
+// Option "mod" can be any of
+//   - 0 : gaussian like distribution 
+//   - 1 : tailed distribution 
+
+  Int_t nx(h2->GetNbinsX())
+      , ny(h2->GetNbinsY())
+      , n(nx*ny);
+  Double_t *data=new Double_t[n];
+  for(Int_t ix(1), in(0); ix<=nx; ix++){
+    for(Int_t iy(1); iy<=ny; iy++)
+      data[in++] = h2->GetBinContent(ix, iy);
+  }
+  Double_t mean, sigm;
+  AliMathBase::EvaluateUni(n, data, mean, sigm, Int_t(n*.8));
+
+  range[0]=mean-3.*sigm; range[1]=mean+3.*sigm;
+  if(mod==1) range[0]=TMath::Max(Float_t(1.e-3), range[0]); 
+  AliDebug(2, Form("h[%s] range0[%f %f]", h2->GetName(), range[0], range[1]));
+  TH1S h1("h1SF0", "", 100, range[0], range[1]);
+  h1.FillN(n,data,0);
+  delete [] data;
+  switch(mod){
+  case 0:// gaussian distribution  
+  {
+    TF1 fg("fg", "gaus", mean-3.*sigm, mean+3.*sigm);
+    h1.Fit(&fg, "QN");
+    mean=fg.GetParameter(1); sigm=fg.GetParameter(2);
+    range[0] = mean-2.5*sigm;range[1] = mean+2.5*sigm;
+    AliDebug(2, Form("     rangeG[%f %f]", range[0], range[1]));
+    break;
+  }
+  case 1:// tailed distribution  
+  {  
+    Int_t bmax(h1.GetMaximumBin());
+    Int_t jBinMin(1), jBinMax(100);
+    for(Int_t ibin(bmax); ibin--;){
+      if(h1.GetBinContent(ibin)<1.){
+        jBinMin=ibin; break;
+      }
+    }
+    for(Int_t ibin(bmax); ibin++;){
+      if(h1.GetBinContent(ibin)<1.){
+        jBinMax=ibin; break;
+      }
+    }
+    range[0]=h1.GetBinCenter(jBinMin); range[1]=h1.GetBinCenter(jBinMax);
+    AliDebug(2, Form("     rangeT[%f %f]", range[0], range[1]));
+    break;
+  }
+  }
+
+  return;
+}
+
+//________________________________________________________
+void AliTRDresolution::MakeSummaryPlot(TObjArray *a, TH2 *h2)
+{
+// Core functionality for MakeSummary function.  
+
+  h2->Reset();  
+  Double_t x,y;
+  TGraphErrors *g(NULL); TAxis *ax(h2->GetXaxis());
+  for(Int_t iseg(0); iseg<fgkNresYsegm[fSegmentLevel]; iseg++){
+    g=(TGraphErrors*)a->At(iseg);
+    for(Int_t in(0); in<g->GetN(); in++){
+      g->GetPoint(in, x, y);
+      h2->SetBinContent(ax->FindBin(x), iseg+1, y);
+    }
+  }
+}
+
+
 //________________________________________________________
 Bool_t AliTRDresolution::PostProcess()
 {
-  //fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
+// Fit, Project, Combine, Extract values from the containers filled during execution
+
+  /*fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));*/
   if (!fContainer) {
     AliError("ERROR: list not available");
     return kFALSE;
   }
-  Int_t nc(0);
+
+  // define general behavior parameters
+  const Color_t fgColorS[11]={
+  kOrange, kOrange-3, kMagenta+1, kViolet, kRed,
+  kGray,
+  kRed, kViolet, kMagenta+1, kOrange-3, kOrange
+  };
+  const Color_t fgColorM[11]={
+  kCyan-5, kAzure-4, kBlue-7, kBlue+2, kViolet+10,
+  kBlack,
+  kViolet+10, kBlue+2, kBlue-7, kAzure-4, kCyan-5
+  };
+  const Marker_t fgMarker[11]={
+  30, 30, 26, 25, 24,
+  28,
+  20, 21, 22, 29, 29
+  };
+
   TGraph *gm= NULL, *gs= NULL;
   if(!fGraphS && !fGraphM){ 
-    TObjArray *aM(NULL), *aS(NULL), *a(NULL);
+    TObjArray *aM(NULL), *aS(NULL);
     Int_t n = fContainer->GetEntriesFast();
     fGraphS = new TObjArray(n); fGraphS->SetOwner();
     fGraphM = new TObjArray(n); fGraphM->SetOwner();
-    for(Int_t ig=0; ig<n; ig++){
+    for(Int_t ig(0), nc(0); ig<n; ig++){
       fGraphM->AddAt(aM = new TObjArray(fgNproj[ig]), ig);
       fGraphS->AddAt(aS = new TObjArray(fgNproj[ig]), ig);
-
-      for(Int_t ic=0; ic<fgNproj[ig]; ic++){
-        //printf("g[%d] c[%d] n[%d]\n", ig, ic, fgNcomp[nc]);
-        if(fgNcomp[nc]>1){
-          aS->AddAt(a = new TObjArray(fgNcomp[nc]), ic); 
-          for(Int_t is=fgNcomp[nc]; is--;){
-            a->AddAt(gs = new TGraphErrors(), is);
-            if(is%11 == 4){
-              gs->SetMarkerStyle(20);
-              gs->SetMarkerColor(kMagenta);
-              gs->SetLineColor(kMagenta);
-           } else if(is%11 == 6){
-              gs->SetMarkerStyle(4);
-              gs->SetMarkerColor(kMagenta);
-              gs->SetLineColor(kMagenta);
+      
+      for(Int_t ic=0; ic<fgNproj[ig]; ic++, nc++){
+        AliDebug(2, Form("building G[%d] P[%d] N[%d]", ig, ic, fNcomp[nc]));
+        if(fNcomp[nc]>1){
+          TObjArray *agS(NULL), *agM(NULL);
+          aS->AddAt(agS = new TObjArray(fNcomp[nc]), ic); 
+          aM->AddAt(agM = new TObjArray(fNcomp[nc]), ic); 
+          for(Int_t is=fNcomp[nc]; is--;){
+            agS->AddAt(gs = new TGraphErrors(), is);
+            Int_t is0(is%11), il0(is/11);
+            gs->SetMarkerStyle(fgMarker[is0]);
+            gs->SetMarkerColor(fgColorS[is0]);
+            gs->SetLineColor(fgColorS[is0]);
+            gs->SetLineStyle(il0);gs->SetLineWidth(2);
+            gs->SetName(Form("s_%d_%02d_%02d", ig, ic, is));
+
+            agM->AddAt(gm = new TGraphErrors(), is);
+            gm->SetMarkerStyle(fgMarker[is0]);
+            gm->SetMarkerColor(fgColorM[is0]);
+            gm->SetLineColor(fgColorM[is0]);
+            gm->SetLineStyle(il0);gm->SetLineWidth(2);
+            gm->SetName(Form("m_%d_%02d_%02d", ig, ic, is));
+            // this is important for labels in the legend 
+            if(ic==0) {
+              gs->SetTitle(Form("%s %02d", fgkResYsegmName[fSegmentLevel], is%fgkNresYsegm[fSegmentLevel]));
+              gm->SetTitle(Form("%s %02d", fgkResYsegmName[fSegmentLevel], is%fgkNresYsegm[fSegmentLevel]));
+            } else if(ic==1) {
+              gs->SetTitle(Form("%s Ly[%d]", is%2 ?"z":"y", is/2));
+              gm->SetTitle(Form("%s Ly[%d]", is%2?"z":"y", is/2));
+            } else if(ic==2||ic==3) {
+              gs->SetTitle(Form("%s Ly[%d]", is%2 ?"RC":"no RC", is/2));
+              gm->SetTitle(Form("%s Ly[%d]", is%2?"RC":"no RC", is/2));
+            } else if(ic<=7) {
+              gs->SetTitle(Form("Layer[%d]", is%AliTRDgeometry::kNlayer));
+              gm->SetTitle(Form("Layer[%d]", is%AliTRDgeometry::kNlayer));
             } else {
-              gs->SetMarkerStyle(7);
-              gs->SetMarkerColor(kRed);
-              gs->SetLineColor(kRed);
-              gs->SetLineStyle(2);
+              gs->SetTitle(Form("%s @ ly[%d]", fgParticle[is0], il0));
+              gm->SetTitle(Form("%s @ ly[%d]", fgParticle[is0], il0));
             }
-            gs->SetNameTitle(Form("s_%d%02d%02d", ig, ic, is), fgParticle[is%11]);
-          }
-          aM->AddAt(a = new TObjArray(fgNcomp[nc]), ic); 
-          for(Int_t is=fgNcomp[nc]; is--;){
-            a->AddAt(gm = new TGraphErrors(), is);
-            if(is%11 == 4){
-              gm->SetMarkerStyle(20);
-              gm->SetMarkerColor(kBlue);
-              gm->SetLineColor(kBlue);
-            } else if(is%11 == 6){
-              gm->SetMarkerStyle(4);
-              gm->SetMarkerColor(kBlue);
-              gm->SetLineColor(kBlue);
-            } else {
-              gm->SetMarkerStyle(7);
-              gm->SetMarkerColor(kBlack);
-              gm->SetLineColor(kBlack);
-              gm->SetLineStyle(2);
-            }
-            gm->SetNameTitle(Form("m_%d%02d%02d", ig, ic, is), "");
           }
         } else {
           aS->AddAt(gs = new TGraphErrors(), ic);
           gs->SetMarkerStyle(23);
           gs->SetMarkerColor(kRed);
           gs->SetLineColor(kRed);
-          gs->SetNameTitle(Form("s_%d%02d", ig, ic), "");
+          gs->SetNameTitle(Form("s_%d_%02d", ig, ic), "sigma");
   
           aM->AddAt(gm = ig ? (TGraph*)new TGraphErrors() : (TGraph*)new TGraphAsymmErrors(), ic);
           gm->SetLineColor(kBlack);
           gm->SetMarkerStyle(7);
           gm->SetMarkerColor(kBlack);
-          gm->SetNameTitle(Form("m_%d%02d", ig, ic), "");
+          gm->SetNameTitle(Form("m_%d_%02d", ig, ic), "mean");
         }
-        nc+=fgNcomp[ic];
       }
     }
   }
@@ -1232,23 +2141,30 @@ Bool_t AliTRDresolution::PostProcess()
   // Charge resolution
   //Process3DL(kCharge, 0, &fl); 
   // Clusters residuals
-  Process2D(kCluster, 0, &fg, 1.e4); 
-  Process2D(kCluster, 1, &fg); 
-  fNRefFigures = 1;
+  Process3D(kCluster, 0, &fg, 1.e4); 
+  Process3Dlinked(kCluster, 1, &fg); 
+  fNRefFigures = 3;
   // Tracklet residual/pulls
-  Process2D(kTrackTRD , 0, &fg, 1.e4); 
-  Process2D(kTrackTRD , 1, &fg); 
-  Process2D(kTrackTRD , 2, &fg, 1.e4); 
-  Process2D(kTrackTRD , 3, &fg); 
-  Process2D(kTrackTRD , 4, &fg, 1.e3); 
-  fNRefFigures = 4;
-  // TPC track residual/pulls
-  Process2D(kTrackTPC, 0, &fg, 1.e4); 
-  Process2D(kTrackTPC, 1, &fg); 
-  Process2D(kTrackTPC, 2, &fg, 1.e4); 
-  Process2D(kTrackTPC, 3, &fg); 
-  Process2D(kTrackTPC, 4, &fg, 1.e3); 
+  Process3D(kTrack , 0, &fg, 1.e4);
+  Process3Dlinked(kTrack , 1, &fg); 
+  Process3D(kTrack , 2, &fg, 1.e4); 
+  Process3D(kTrack , 3, &fg); 
+  Process2D(kTrack , 4, &fg, 1.e3);
   fNRefFigures = 7;
+  // TRDin residual/pulls
+  Process3D(kTrackIn, 0, &fg, 1.e4); 
+  Process3Dlinked(kTrackIn, 1, &fg); 
+  Process3D(kTrackIn, 2, &fg, 1.e4); 
+  Process3D(kTrackIn, 3, &fg); 
+  Process2D(kTrackIn, 4, &fg, 1.e3); 
+  fNRefFigures = 11;
+  // TRDout residual/pulls
+  Process3D(kTrackOut, 0, &fg, 1.e3); // scale to fit - see PlotTrackOut
+  Process3Dlinked(kTrackOut, 1, &fg); 
+  Process3D(kTrackOut, 2, &fg, 1.e4); 
+  Process3D(kTrackOut, 3, &fg); 
+  Process2D(kTrackOut, 4, &fg, 1.e3); 
+  fNRefFigures = 15;
 
   if(!HasMCdata()) return kTRUE;
 
@@ -1256,66 +2172,59 @@ Bool_t AliTRDresolution::PostProcess()
   //PROCESS MC RESIDUAL DISTRIBUTIONS
 
   // CLUSTER Y RESOLUTION/PULLS
-  Process3Drange(kMCcluster, 0, 0, &fg, 1.e4, 1, 1);
-  Process3Drange(kMCcluster, 0, 1, &fg, 1.e4, 2, 2);
-  Process3Drange(kMCcluster, 0, 2, &fg, 1.e4, 3, 12);
-  Process2D(kMCcluster, 1, &fg, 1., 3);
-  fNRefFigures = 10;
+  Process3D(kMCcluster, 0, &fg, 1.e4);
+  Process3Dlinked(kMCcluster, 1, &fg, 1.);
+  fNRefFigures = 17;
 
   // TRACKLET RESOLUTION/PULLS
-  Process3Drange(kMCtracklet, 0, 0, &fg, 1.e4, 1, 1); // y
-  Process3Drange(kMCtracklet, 0, 1, &fg, 1.e4, 2, 2); // y
-  Process3Drange(kMCtracklet, 0, 2, &fg, 1.e4, 3, 12); // y
-  Process2D(kMCtracklet, 1, &fg, 1., 3);       // y pulls
-  Process2D(kMCtracklet, 2, &fg, 1.e4, 4); // z
-  Process2D(kMCtracklet, 3, &fg, 1., 5);       // z pulls
-  Process2D(kMCtracklet, 4, &fg, 1.e3, 6); // phi
-  fNRefFigures = 14;
+  Process3D(kMCtracklet, 0, &fg, 1.e4); // y
+  Process3Dlinked(kMCtracklet, 1, &fg, 1.);   // y pulls
+  Process3D(kMCtracklet, 2, &fg, 1.e4); // z
+  Process3D(kMCtracklet, 3, &fg, 1.);   // z pulls
+  Process2D(kMCtracklet, 4, &fg, 1.e3); // phi
+  fNRefFigures = 21;
 
   // TRACK RESOLUTION/PULLS
-  Process2D(kMCtrackTRD, 0, &fg, 1.e4);   // y
-  Process2D(kMCtrackTRD, 1, &fg);         // y PULL
-  Process2D(kMCtrackTRD, 2, &fg, 1.e4);   // z
-  Process2D(kMCtrackTRD, 3, &fg);         // z PULL
-  Process2D(kMCtrackTRD, 4, &fg, 1.e3);   // phi
-  Process2D(kMCtrackTRD, 5, &fg);         // snp PULL
-  Process2D(kMCtrackTRD, 6, &fg, 1.e3);   // theta
-  Process2D(kMCtrackTRD, 7, &fg);         // tgl PULL
-  Process4D(kMCtrackTRD, 8, &fg, 1.e2);   // pt resolution
-  Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 4);// pt resolution e1- @ L0
-  Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 6);// pt resolution e1+ @ L0
-  Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 55+4);// pt resolution e1- @ L5
-  Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 55+6);// pt resolution e1+ @ L5
-  Process4D(kMCtrackTRD, 9, &fg);         // 1/pt pulls
-  Process4D(kMCtrackTRD, 10, &fg, 1.e2);  // p resolution
-  Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 4);// p resolution e1- @ L0
-  Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 6);// p resolution e1+ @ L0
-  Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 55+4);// p resolution e1- @ L5
-  Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 55+6);// p resolution e1+ @ L5
-  Process4D(kMCtrackTRD, 11, &fg, 1.e2);   // pt resolution stand alone
-  Process4D(kMCtrackTRD, 12, &fg);         // 1/pt pulls  stand alone
-  Process4D(kMCtrackTRD, 13, &fg, 1.e2);  // p resolution  stand alone
-  fNRefFigures = 24;
-
-  // TRACK TPC RESOLUTION/PULLS
-  Process2D(kMCtrackTPC, 0, &fg, 1.e4);// y resolution
-  Process2D(kMCtrackTPC, 1, &fg);      // y pulls
-  Process2D(kMCtrackTPC, 2, &fg, 1.e4);// z resolution
-  Process2D(kMCtrackTPC, 3, &fg);      // z pulls
-  Process2D(kMCtrackTPC, 4, &fg, 1.e3);// phi resolution
-  Process2D(kMCtrackTPC, 5, &fg);      // snp pulls
-  Process2D(kMCtrackTPC, 6, &fg, 1.e3);// theta resolution
-  Process2D(kMCtrackTPC, 7, &fg);      // tgl pulls
-  Process3D(kMCtrackTPC, 8, &fg, 1.e2);// pt resolution
-  Process3D(kMCtrackTPC, 9, &fg);      // 1/pt pulls
-  Process3D(kMCtrackTPC, 10, &fg, 1.e2);// p resolution
-  Process3D(kMCtrackTPC, 11, &fg);      // p pulls
-  fNRefFigures = 30;
-
-  // TRACK HMPID RESOLUTION/PULLS
-  Process2D(kMCtrackTOF, 0, &fg, 1.e4); // z towards TOF
-  Process2D(kMCtrackTOF, 1, &fg);       // z towards TOF
-  fNRefFigures = 31;
+  Process3Darray(kMCtrack, 0, &fg, 1.e4);   // y
+  Process3DlinkedArray(kMCtrack, 1, &fg);   // y PULL
+  Process3Darray(kMCtrack, 2, &fg, 1.e4);   // z
+  Process3Darray(kMCtrack, 3, &fg);         // z PULL
+  Process2Darray(kMCtrack, 4, &fg, 1.e3);   // phi
+  Process2Darray(kMCtrack, 5, &fg);         // snp PULL
+  Process2Darray(kMCtrack, 6, &fg, 1.e3);   // theta
+  Process2Darray(kMCtrack, 7, &fg);         // tgl PULL
+  Process3Darray(kMCtrack, 8, &fg, 1.e2);   // pt resolution
+  Process3Darray(kMCtrack, 9, &fg);         // 1/pt pulls
+  Process3Darray(kMCtrack, 10, &fg, 1.e2);  // p resolution
+  fNRefFigures+=16;
+
+  // TRACK TRDin RESOLUTION/PULLS
+  Process3D(kMCtrackIn, 0, &fg, 1.e4);// y resolution
+  Process3Dlinked(kMCtrackIn, 1, &fg);      // y pulls
+  Process3D(kMCtrackIn, 2, &fg, 1.e4);// z resolution
+  Process3D(kMCtrackIn, 3, &fg);      // z pulls
+  Process2D(kMCtrackIn, 4, &fg, 1.e3);// phi resolution
+  Process2D(kMCtrackIn, 5, &fg);      // snp pulls
+  Process2D(kMCtrackIn, 6, &fg, 1.e3);// theta resolution
+  Process2D(kMCtrackIn, 7, &fg);      // tgl pulls
+  Process3D(kMCtrackIn, 8, &fg, 1.e2);// pt resolution
+  Process3D(kMCtrackIn, 9, &fg);      // 1/pt pulls
+  Process3D(kMCtrackIn, 10, &fg, 1.e2);// p resolution
+  fNRefFigures+=8;
+
+  // TRACK TRDout RESOLUTION/PULLS
+  Process3D(kMCtrackOut, 0, &fg, 1.e4);// y resolution
+  Process3Dlinked(kMCtrackOut, 1, &fg);      // y pulls
+  Process3D(kMCtrackOut, 2, &fg, 1.e4);// z resolution
+  Process3D(kMCtrackOut, 3, &fg);      // z pulls
+  Process2D(kMCtrackOut, 4, &fg, 1.e3);// phi resolution
+  Process2D(kMCtrackOut, 5, &fg);      // snp pulls
+  Process2D(kMCtrackOut, 6, &fg, 1.e3);// theta resolution
+  Process2D(kMCtrackOut, 7, &fg);      // tgl pulls
+  Process3D(kMCtrackOut, 8, &fg, 1.e2);// pt resolution
+  Process3D(kMCtrackOut, 9, &fg);      // 1/pt pulls
+  Process3D(kMCtrackOut, 10, &fg, 1.e2);// p resolution
+  fNRefFigures+=8;
 
   return kTRUE;
 }
@@ -1358,6 +2267,158 @@ void AliTRDresolution::AdjustF1(TH1 *h, TF1 *f)
   f->SetParameter(5, 2.e-1);
 }
 
+//________________________________________________________
+TObjArray* AliTRDresolution::BuildMonitorContainerCluster(const char* name, Bool_t expand)
+{
+// Build performance histograms for AliTRDcluster.vs TRD track or MC
+//  - y reziduals/pulls
+
+  TObjArray *arr = new TObjArray(2);
+  arr->SetName(name); arr->SetOwner();
+  TH1 *h(NULL); char hname[100], htitle[300];
+
+  // tracklet resolution/pull in y direction
+  sprintf(hname, "%s_%s_Y", GetNameId(), name);
+  sprintf(htitle, "Y res for \"%s\" @ %s;tg(#phi);#Delta y [cm];%s", GetNameId(), name, fgkResYsegmName[fSegmentLevel]);
+  if(!(h = (TH3S*)gROOT->FindObject(hname))){
+    Int_t nybins=fgkNresYsegm[fSegmentLevel];
+    if(expand) nybins*=2;
+    h = new TH3S(hname, htitle, 
+                 48, -.48, .48,            // phi
+                 60, -fDyRange, fDyRange,  // dy
+                 nybins, -0.5, nybins-0.5);// segment
+  } else h->Reset();
+  arr->AddAt(h, 0);
+  sprintf(hname, "%s_%s_YZpull", GetNameId(), name);
+  sprintf(htitle, "YZ pull for \"%s\" @ %s;%s;#Delta y  / #sigma_{y};#Delta z  / #sigma_{z}", GetNameId(), name, fgkResYsegmName[fSegmentLevel]);
+  if(!(h = (TH3S*)gROOT->FindObject(hname))){
+    h = new TH3S(hname, htitle, fgkNresYsegm[fSegmentLevel], -0.5, fgkNresYsegm[fSegmentLevel]-0.5, 100, -4.5, 4.5, 100, -4.5, 4.5);
+  } else h->Reset();
+  arr->AddAt(h, 1);
+
+  return arr;
+}
+
+//________________________________________________________
+TObjArray* AliTRDresolution::BuildMonitorContainerTracklet(const char* name, Bool_t expand)
+{
+// Build performance histograms for AliExternalTrackParam.vs TRD tracklet
+//  - y reziduals/pulls
+//  - z reziduals/pulls
+//  - phi reziduals
+  TObjArray *arr = BuildMonitorContainerCluster(name, expand); 
+  arr->Expand(5);
+  TH1 *h(NULL); char hname[100], htitle[300];
+
+  // tracklet resolution/pull in z direction
+  sprintf(hname, "%s_%s_Z", GetNameId(), name);
+  sprintf(htitle, "Z res for \"%s\" @ %s;tg(#theta);#Delta z [cm];row cross", GetNameId(), name);
+  if(!(h = (TH3S*)gROOT->FindObject(hname))){
+    h = new TH3S(hname, htitle, 50, -1., 1., 100, -1.5, 1.5, 2, -0.5, 1.5);
+  } else h->Reset();
+  arr->AddAt(h, 2);
+  sprintf(hname, "%s_%s_Zpull", GetNameId(), name);
+  sprintf(htitle, "Z pull for \"%s\" @ %s;tg(#theta);#Delta z  / #sigma_{z};row cross", GetNameId(), name);
+  if(!(h = (TH3S*)gROOT->FindObject(hname))){
+    h = new TH3S(hname, htitle, 50, -1., 1., 100, -5.5, 5.5, 2, -0.5, 1.5);
+    h->GetZaxis()->SetBinLabel(1, "no RC");
+    h->GetZaxis()->SetBinLabel(2, "RC");
+  } else h->Reset();
+  arr->AddAt(h, 3);
+
+  // tracklet to track phi resolution
+  sprintf(hname, "%s_%s_PHI", GetNameId(), name);
+  sprintf(htitle, "#Phi res for \"%s\" @ %s;tg(#phi);#Delta #phi [rad];entries", GetNameId(), name);
+  if(!(h = (TH2I*)gROOT->FindObject(hname))){
+    h = new TH2I(hname, htitle, 21, -.33, .33, 100, -.5, .5);
+  } else h->Reset();
+  arr->AddAt(h, 4);
+
+  return arr;
+}
+
+//________________________________________________________
+TObjArray* AliTRDresolution::BuildMonitorContainerTrack(const char* name)
+{
+// Build performance histograms for AliExternalTrackParam.vs MC
+//  - y resolution/pulls
+//  - z resolution/pulls
+//  - phi resolution, snp pulls
+//  - theta resolution, tgl pulls
+//  - pt resolution, 1/pt pulls, p resolution
+
+  TObjArray *arr = BuildMonitorContainerTracklet(name); 
+  arr->Expand(11);
+  TH1 *h(NULL); char hname[100], htitle[300];
+  TAxis *ax(NULL);
+
+  // snp pulls
+  sprintf(hname, "%s_%s_SNPpull", GetNameId(), name);
+  sprintf(htitle, "SNP pull for \"%s\" @ %s;tg(#phi);#Delta snp  / #sigma_{snp};entries", GetNameId(), name);
+  if(!(h = (TH2I*)gROOT->FindObject(hname))){
+    h = new TH2I(hname, htitle, 60, -.3, .3, 100, -4.5, 4.5);
+  } else h->Reset();
+  arr->AddAt(h, 5);
+
+  // theta resolution
+  sprintf(hname, "%s_%s_THT", GetNameId(), name);
+  sprintf(htitle, "#Theta res for \"%s\" @ %s;tg(#theta);#Delta #theta [rad];entries", GetNameId(), name);
+  if(!(h = (TH2I*)gROOT->FindObject(hname))){
+    h = new TH2I(hname, htitle, 100, -1., 1., 100, -5e-3, 5e-3);
+  } else h->Reset();
+  arr->AddAt(h, 6);
+  // tgl pulls
+  sprintf(hname, "%s_%s_TGLpull", GetNameId(), name);
+  sprintf(htitle, "TGL pull for \"%s\" @ %s;tg(#theta);#Delta tgl  / #sigma_{tgl};entries", GetNameId(), name);
+  if(!(h = (TH2I*)gROOT->FindObject(hname))){
+    h = new TH2I(hname, htitle, 100, -1., 1., 100, -4.5, 4.5);
+  } else h->Reset();
+  arr->AddAt(h, 7);
+  
+  const Int_t kNpt(14);
+  const Int_t kNdpt(150); 
+  const Int_t kNspc = 2*AliPID::kSPECIES+1;
+  Float_t lPt=0.1, lDPt=-.1, lSpc=-5.5;
+  Float_t binsPt[kNpt+1], binsSpc[kNspc+1], binsDPt[kNdpt+1];
+  for(Int_t i=0;i<kNpt+1; i++,lPt=TMath::Exp(i*.15)-1.) binsPt[i]=lPt;
+  for(Int_t i=0; i<kNspc+1; i++,lSpc+=1.) binsSpc[i]=lSpc;
+  for(Int_t i=0; i<kNdpt+1; i++,lDPt+=2.e-3) binsDPt[i]=lDPt;
+
+  // Pt resolution
+  sprintf(hname, "%s_%s_Pt", GetNameId(), name);
+  sprintf(htitle, "P_{t} res for \"%s\" @ %s;p_{t} [GeV/c];#Delta p_{t}/p_{t}^{MC};SPECIES", GetNameId(), name);
+  if(!(h = (TH3S*)gROOT->FindObject(hname))){
+    h = new TH3S(hname, htitle, 
+                 kNpt, binsPt, kNdpt, binsDPt, kNspc, binsSpc);
+    ax = h->GetZaxis();
+    for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]);
+  } else h->Reset();
+  arr->AddAt(h, 8);
+  // 1/Pt pulls
+  sprintf(hname, "%s_%s_1Pt", GetNameId(), name);
+  sprintf(htitle, "1/P_{t} pull for \"%s\" @ %s;1/p_{t}^{MC} [c/GeV];#Delta(1/p_{t})/#sigma(1/p_{t});SPECIES", GetNameId(), name);
+  if(!(h = (TH3S*)gROOT->FindObject(hname))){
+    h = new TH3S(hname, htitle, 
+                 kNpt, 0., 2., 100, -4., 4., kNspc, -5.5, 5.5);
+    ax = h->GetZaxis();
+    for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]);
+  } else h->Reset();
+  arr->AddAt(h, 9);
+  // P resolution
+  sprintf(hname, "%s_%s_P", GetNameId(), name);
+  sprintf(htitle, "P res for \"%s\" @ %s;p [GeV/c];#Delta p/p^{MC};SPECIES", GetNameId(), name);
+  if(!(h = (TH3S*)gROOT->FindObject(hname))){
+    h = new TH3S(hname, htitle, 
+                 kNpt, binsPt, kNdpt, binsDPt, kNspc, binsSpc);
+    ax = h->GetZaxis();
+    for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]);
+  } else h->Reset();
+  arr->AddAt(h, 10);
+
+  return arr;
+}
+
+
 //________________________________________________________
 TObjArray* AliTRDresolution::Histos()
 {
@@ -1369,13 +2430,19 @@ TObjArray* AliTRDresolution::Histos()
 
   fContainer  = new TObjArray(kNviews);
   //fContainer->SetOwner(kTRUE);
-  TH1 *h = NULL;
-  TObjArray *arr = NULL;
+  TH1 *h(NULL);
+  TObjArray *arr(NULL);
 
-  const Int_t kPbins(12); // binning in momentum range should depend on the statistics analyzed
+  // binnings for plots containing momentum or pt
+  const Int_t kNpt(14), kNphi(48), kNdy(60);
+  Float_t lPhi=-.48, lDy=-.3, lPt=0.1;
+  Float_t binsPhi[kNphi+1], binsDy[kNdy+1], binsPt[kNpt+1];
+  for(Int_t i=0; i<kNphi+1; i++,lPhi+=.02) binsPhi[i]=lPhi;
+  for(Int_t i=0; i<kNdy+1; i++,lDy+=.01) binsDy[i]=lDy;
+  for(Int_t i=0;i<kNpt+1; i++,lPt=TMath::Exp(i*.15)-1.) binsPt[i]=lPt;
 
   // cluster to track residuals/pulls
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kCharge]), kCharge);
+  fContainer->AddAt(arr = new TObjArray(2), kCharge);
   arr->SetName("Charge");
   if(!(h = (TH3S*)gROOT->FindObject("hCharge"))){
     h = new TH3S("hCharge", "Charge Resolution", 20, 1., 2., 24, 0., 3.6, 100, 0., 500.);
@@ -1386,440 +2453,111 @@ TObjArray* AliTRDresolution::Histos()
   arr->AddAt(h, 0);
 
   // cluster to track residuals/pulls
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kCluster]), kCluster);
-  arr->SetName("Cl");
-  if(!(h = (TH2I*)gROOT->FindObject("hCl"))){
-    h = new TH2I("hCl", "Cluster Residuals", 21, -.33, .33, 100, -.5, .5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  if(!(h = (TH2I*)gROOT->FindObject("hClpull"))){
-    h = new TH2I("hClpull", "Cluster Pulls", 21, -.33, .33, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y/#sigma_{y}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 1);
-
-  // tracklet to track residuals/pulls in y direction
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kTrackTRD ]), kTrackTRD );
-  arr->SetName("Trklt");
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkltY"))){
-    h = new TH2I("hTrkltY", "Tracklet Y Residuals", 21, -.33, .33, 100, -.5, .5);
-    h->GetXaxis()->SetTitle("#tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkltYpull"))){
-    h = new TH2I("hTrkltYpull", "Tracklet Y Pulls", 21, -.33, .33, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("#tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y/#sigma_{y}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 1);
-  // tracklet to track residuals/pulls in z direction
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkltZ"))){
-    h = new TH2I("hTrkltZ", "Tracklet Z Residuals", 50, -1., 1., 100, -1.5, 1.5);
-    h->GetXaxis()->SetTitle("#tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 2);
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkltZpull"))){
-    h = new TH2I("hTrkltZpull", "Tracklet Z Pulls", 50, -1., 1., 100, -5.5, 5.5);
-    h->GetXaxis()->SetTitle("#tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z/#sigma_{z}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 3);
-  // tracklet to track phi residuals
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkltPhi"))){
-    h = new TH2I("hTrkltPhi", "Tracklet #phi Residuals", 21, -.33, .33, 100, -.5, .5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta phi [rad]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 4);
-
-
-  // tracklet to TPC track residuals/pulls in y direction
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kTrackTPC]), kTrackTPC);
-  arr->SetName("TrkTPC");
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkTPCY"))){
-    h = new TH2I("hTrkTPCY", "Track[TPC] Y Residuals", 21, -.33, .33, 100, -.5, .5);
-    h->GetXaxis()->SetTitle("#tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkTPCYpull"))){
-    h = new TH2I("hTrkTPCYpull", "Track[TPC] Y Pulls", 21, -.33, .33, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("#tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y/#sigma_{y}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 1);
-  // tracklet to TPC track residuals/pulls in z direction
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkTPCZ"))){
-    h = new TH2I("hTrkTPCZ", "Track[TPC] Z Residuals", 50, -1., 1., 100, -1.5, 1.5);
-    h->GetXaxis()->SetTitle("#tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 2);
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkTPCZpull"))){
-    h = new TH2I("hTrkTPCZpull", "Track[TPC] Z Pulls", 50, -1., 1., 100, -5.5, 5.5);
-    h->GetXaxis()->SetTitle("#tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z/#sigma_{z}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 3);
-  // tracklet to TPC track phi residuals
-  if(!(h = (TH2I*)gROOT->FindObject("hTrkTPCPhi"))){
-    h = new TH2I("hTrkTPCPhi", "Track[TPC] #phi Residuals", 21, -.33, .33, 100, -.5, .5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta phi [rad]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 4);
+  fContainer->AddAt(BuildMonitorContainerCluster("Cl"), kCluster);
+  // tracklet to TRD track
+  fContainer->AddAt(BuildMonitorContainerTracklet("Trk", kTRUE), kTrack);
+  // tracklet to TRDin
+  fContainer->AddAt(BuildMonitorContainerTracklet("TrkIN", kTRUE), kTrackIn);
+  // tracklet to TRDout
+  fContainer->AddAt(BuildMonitorContainerTracklet("TrkOUT"), kTrackOut);
 
 
   // Resolution histos
   if(!HasMCdata()) return fContainer;
 
-  // cluster y resolution [0]
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kMCcluster]), kMCcluster);
-  arr->SetName("McCl");
-  if(!(h = (TH3S*)gROOT->FindObject("hMcCl"))){
-    h = new TH3S("hMcCl", "Cluster Resolution", 48, -.48, .48, 100, -.3, .3, kPbins, 0., 12.);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y [cm]");
-    h->GetZaxis()->SetTitle("p_{t} [GeV/c]");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  if(!(h = (TH2I*)gROOT->FindObject("hMcClPull"))){
-    h = new TH2I("hMcClPull", "Cluster Pulls", 48, -.48, .48, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Deltay/#sigma_{y}");
-    h->GetZaxis()->SetTitle("p_{t} [GeV/c]");
-  } else h->Reset();
-  arr->AddAt(h, 1);
+  // cluster resolution 
+  fContainer->AddAt(BuildMonitorContainerCluster("MCcl"),  kMCcluster);
 
+  // tracklet resolution
+  fContainer->AddAt(BuildMonitorContainerTracklet("MCtracklet"), kMCtracklet);
 
-  // TRACKLET RESOLUTION
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kMCtracklet]), kMCtracklet);
-  arr->SetName("McTrklt");
-  // tracklet y resolution
-  if(!(h = (TH3S*)gROOT->FindObject("hMcTrkltY"))){
-    h = new TH3S("hMcTrkltY", "Tracklet Resolution (Y)", 48, -.48, .48, 100, -.2, .2, kPbins, 0., 12.);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  // tracklet y pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltYPull"))){
-    h = new TH2I("hMcTrkltYPull", "Tracklet Pulls (Y)", 48, -.48, .48, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 1);
-  // tracklet z resolution
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltZ"))){
-    h = new TH2I("hMcTrkltZ", "Tracklet Resolution (Z)", 100, -1., 1., 100, -1., 1.);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 2);
-  // tracklet z pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltZPull"))){
-    h = new TH2I("hMcTrkltZPull", "Tracklet Pulls (Z)", 100, -1., 1., 100, -3.5, 3.5);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 3);
-  // tracklet phi resolution
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltPhi"))){
-    h = new TH2I("hMcTrkltPhi", "Tracklet Resolution (#Phi)", 48, -.48, .48, 100, -.15, .15);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta #phi [rad]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 4);
+  // track resolution
+  fContainer->AddAt(arr = new TObjArray(AliTRDgeometry::kNlayer), kMCtrack);
+  arr->SetName("MCtrk");
+  for(Int_t il(0); il<AliTRDgeometry::kNlayer; il++) arr->AddAt(BuildMonitorContainerTrack(Form("MCtrk_Ly%d", il)), il);
 
+  // TRDin TRACK RESOLUTION
+  fContainer->AddAt(BuildMonitorContainerTrack("MCtrkIN"), kMCtrackIn);
 
-  // KALMAN TRACK RESOLUTION
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kMCtrackTRD]), kMCtrackTRD);
-  arr->SetName("McTrkTRD");
-  // Kalman track y resolution
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkY"))){
-    h = new TH2I("hMcTrkY", "Track Y Resolution", 48, -.48, .48, 100, -.2, .2);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  // Kalman track y pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkYPull"))){
-    h = new TH2I("hMcTrkYPull", "Track Y Pulls", 48, -.48, .48, 100, -4., 4.);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 1);
-  // Kalman track Z
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkZ"))){
-    h = new TH2I("hMcTrkZ", "Track Z Resolution", 100, -1., 1., 100, -1., 1.);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 2);
-  // Kalman track Z pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkZPull"))){
-    h = new TH2I("hMcTrkZPull", "Track Z Pulls", 100, -1., 1., 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 3);
-  // Kalman track SNP
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkSNP"))){
-    h = new TH2I("hMcTrkSNP", "Track Phi Resolution", 60, -.3, .3, 100, -5e-3, 5e-3);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta #phi [rad]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 4);
-  // Kalman track SNP pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkSNPPull"))){
-    h = new TH2I("hMcTrkSNPPull", "Track SNP Pulls", 60, -.3, .3, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta(sin(#phi)) / #sigma_{sin(#phi)}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 5);
-  // Kalman track TGL
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTGL"))){
-    h = new TH2I("hMcTrkTGL", "Track Theta Resolution", 100, -1., 1., 100, -5e-3, 5e-3);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta#theta [rad]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 6);
-  // Kalman track TGL pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTGLPull"))){
-    h = new TH2I("hMcTrkTGLPull", "Track TGL  Pulls", 100, -1., 1., 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta(tg(#theta)) / #sigma_{tg(#theta)}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 7);
-  const Int_t n = 2*AliPID::kSPECIES+1;
-  TObjArray *arr2 = NULL; TH3S* h3=NULL;
-  // Kalman track Pt resolution
-  arr->AddAt(arr2 = new TObjArray(AliTRDgeometry::kNlayer), 8);
-  arr2->SetName("Pt Resolution");
-  for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
-    if(!(h3 = (TH3S*)gROOT->FindObject(Form("hMcTrkPt%d", il)))){
-      h3 = new TH3S(Form("hMcTrkPt%d", il), "Track Pt Resolution", kPbins, 0., 12., 150, -.1, .2, n, -5.5, 5.5);
-      h3->GetXaxis()->SetTitle("p_{t} [GeV/c]");
-      h3->GetYaxis()->SetTitle("#Delta p_{t}/p_{t}^{MC}");
-      h3->GetZaxis()->SetTitle("SPECIES");
-    } else h3->Reset();
-    arr2->AddAt(h3, il);
-  }
-  // Kalman track Pt pulls
-  arr->AddAt(arr2 = new TObjArray(AliTRDgeometry::kNlayer), 9);
-  arr2->SetName("1/Pt Pulls");
-  for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
-    if(!(h3 = (TH3S*)gROOT->FindObject(Form("hMcTrkPtPulls%d", il)))){
-      h3 = new TH3S(Form("hMcTrkPtPulls%d", il), "Track 1/Pt Pulls", kPbins, 0., 2., 100, -4., 4., n, -5.5, 5.5);
-      h3->GetXaxis()->SetTitle("1/p_{t}^{MC} [c/GeV]");
-      h3->GetYaxis()->SetTitle("#Delta(1/p_{t})/#sigma(1/p_{t}) ");
-      h3->GetZaxis()->SetTitle("SPECIES");
-    } else h3->Reset();
-    arr2->AddAt(h3, il);
-  }
-  // Kalman track P resolution
-  arr->AddAt(arr2 = new TObjArray(AliTRDgeometry::kNlayer), 10);
-  arr2->SetName("P Resolution");
-  for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
-    if(!(h3 = (TH3S*)gROOT->FindObject(Form("hMcTrkP%d", il)))){
-      h3 = new TH3S(Form("hMcTrkP%d", il), "Track P Resolution", kPbins, 0., 12., 150, -.15, .35, n, -5.5, 5.5);
-      h3->GetXaxis()->SetTitle("p [GeV/c]");
-      h3->GetYaxis()->SetTitle("#Delta p/p^{MC}");
-      h3->GetZaxis()->SetTitle("SPECIES");
-    } else h3->Reset();
-    arr2->AddAt(h3, il);
-  }
-  // TRD stand-alone track Pt resolution
-  arr->AddAt(arr2 = new TObjArray(AliTRDgeometry::kNlayer), 11);
-  arr2->SetName("Pt Resolution [SA]");
-  for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
-    if(!(h3 = (TH3S*)gROOT->FindObject(Form("hMcSATrkPt%d", il)))){
-      h3 = new TH3S(Form("hMcSATrkPt%d", il), "Track Pt Resolution", kPbins, 0., 12., 150, -10., 1., n, -5.5, 5.5);
-      h3->GetXaxis()->SetTitle("p_{t} [GeV/c]");
-      h3->GetYaxis()->SetTitle("#Delta p_{t}/p_{t}^{MC}");
-      h3->GetZaxis()->SetTitle("SPECIES");
-    } else h3->Reset();
-    arr2->AddAt(h3, il);
-  }
-  // TRD stand-alone track Pt pulls
-  arr->AddAt(arr2 = new TObjArray(AliTRDgeometry::kNlayer), 12);
-  arr2->SetName("1/Pt Pulls [SA]");
-  for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
-    if(!(h3 = (TH3S*)gROOT->FindObject(Form("hMcSATrkPtPulls%d", il)))){
-      h3 = new TH3S(Form("hMcSATrkPtPulls%d", il), "Track 1/Pt Pulls", kPbins, 0., 2., 100, -4., 4., n, -5.5, 5.5);
-      h3->GetXaxis()->SetTitle("1/p_{t}^{MC} [c/GeV]");
-      h3->GetYaxis()->SetTitle("#Delta(1/p_{t})/#sigma(1/p_{t}) ");
-      h3->GetZaxis()->SetTitle("SPECIES");
-    } else h3->Reset();
-    arr2->AddAt(h3, il);
-  }
-  // TRD stand-alone track P resolution
-  arr->AddAt(arr2 = new TObjArray(AliTRDgeometry::kNlayer), 13);
-  arr2->SetName("P Resolution [SA]");
-  for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
-    if(!(h3 = (TH3S*)gROOT->FindObject(Form("hMcSATrkP%d", il)))){
-      h3 = new TH3S(Form("hMcSATrkP%d", il), "Track P Resolution", kPbins, 0., 12., 150, -1.5, 3.5, n, -5.5, 5.5);
-      h3->GetXaxis()->SetTitle("p [GeV/c]");
-      h3->GetYaxis()->SetTitle("#Delta p/p^{MC}");
-      h3->GetZaxis()->SetTitle("SPECIES");
-    } else h3->Reset();
-    arr2->AddAt(h3, il);
-  }
-
-  // TPC TRACK RESOLUTION
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kMCtrackTPC]), kMCtrackTPC);
-  arr->SetName("McTrkTPC");
-  // Kalman track Y
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCY"))){
-    h = new TH2I("hMcTrkTPCY", "Track[TPC] Y Resolution", 60, -.3, .3, 100, -.5, .5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  // Kalman track Y pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCYPull"))){
-    h = new TH2I("hMcTrkTPCYPull", "Track[TPC] Y Pulls", 60, -.3, .3, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 1);
-  // Kalman track Z
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCZ"))){
-    h = new TH2I("hMcTrkTPCZ", "Track[TPC] Z Resolution", 100, -1., 1., 100, -1., 1.);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 2);
-  // Kalman track Z pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCZPull"))){
-    h = new TH2I("hMcTrkTPCZPull", "Track[TPC] Z Pulls", 100, -1., 1., 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 3);
-  // Kalman track SNP
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCSNP"))){
-    h = new TH2I("hMcTrkTPCSNP", "Track[TPC] Phi Resolution", 60, -.3, .3, 100, -5e-3, 5e-3);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta #phi [rad]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 4);
-  // Kalman track SNP pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCSNPPull"))){
-    h = new TH2I("hMcTrkTPCSNPPull", "Track[TPC] SNP Pulls", 60, -.3, .3, 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#phi)");
-    h->GetYaxis()->SetTitle("#Delta(sin(#phi)) / #sigma_{sin(#phi)}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 5);
-  // Kalman track TGL
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCTGL"))){
-    h = new TH2I("hMcTrkTPCTGL", "Track[TPC] Theta Resolution", 100, -1., 1., 100, -5e-3, 5e-3);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta#theta [rad]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 6);
-  // Kalman track TGL pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTPCTGLPull"))){
-    h = new TH2I("hMcTrkTPCTGLPull", "Track[TPC] TGL  Pulls", 100, -1., 1., 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta(tg(#theta)) / #sigma_{tg(#theta)}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 7);
-  // Kalman track Pt resolution
-  if(!(h3 = (TH3S*)gROOT->FindObject("hMcTrkTPCPt"))){
-    h3 = new TH3S("hMcTrkTPCPt", "Track[TPC] Pt Resolution", kPbins, 0., 12., 100, -.1, .2, n, -5.5, 5.5);
-    h3->GetXaxis()->SetTitle("p_{t} [GeV/c]");
-    h3->GetYaxis()->SetTitle("#Delta p_{t}/p_{t}^{MC}");
-    h3->GetZaxis()->SetTitle("SPECIES");
-  } else h3->Reset();
-  arr->AddAt(h3, 8);
-  // Kalman track Pt pulls
-  if(!(h3 = (TH3S*)gROOT->FindObject("hMcTrkTPCPtPulls"))){
-    h3 = new TH3S("hMcTrkTPCPtPulls", "Track[TPC] 1/Pt Pulls", kPbins, 0., 2., 100, -4., 4., n, -5.5, 5.5);
-    h3->GetXaxis()->SetTitle("1/p_{t}^{MC} [c/GeV]");
-    h3->GetYaxis()->SetTitle("#Delta(1/p_{t})/#sigma(1/p_{t}) ");
-    h3->GetZaxis()->SetTitle("SPECIES");
-  } else h3->Reset();
-  arr->AddAt(h3, 9);
-  // Kalman track P resolution
-  if(!(h3 = (TH3S*)gROOT->FindObject("hMcTrkTPCP"))){
-    h3 = new TH3S("hMcTrkTPCP", "Track[TPC] P Resolution", kPbins, 0., 12., 100, -.15, .35, n, -5.5, 5.5);
-    h3->GetXaxis()->SetTitle("p [GeV/c]");
-    h3->GetYaxis()->SetTitle("#Delta p/p^{MC}");
-    h3->GetZaxis()->SetTitle("SPECIES");
-  } else h3->Reset();
-  arr->AddAt(h3, 10);
-  // Kalman track P pulls
-  if(!(h3 = (TH3S*)gROOT->FindObject("hMcTrkTPCPPulls"))){
-    h3 = new TH3S("hMcTrkTPCPPulls", "Track[TPC] P Pulls", kPbins, 0., 12., 100, -5., 5., n, -5.5, 5.5);
-    h3->GetXaxis()->SetTitle("p^{MC} [GeV/c]");
-    h3->GetYaxis()->SetTitle("#Deltap/#sigma_{p}");
-    h3->GetZaxis()->SetTitle("SPECIES");
-  } else h3->Reset();
-  arr->AddAt(h3, 11);
-
-
-
-  // Kalman track Z resolution [TOF]
-  fContainer->AddAt(arr = new TObjArray(fgNhistos[kMCtrackTOF]), kMCtrackTOF);
-  arr->SetName("McTrkTOF");
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTOFZ"))){
-    h = new TH2I("hMcTrkTOFZ", "Track[TOF] Z Resolution", 100, -1., 1., 100, -1., 1.);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z [cm]");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 0);
-  // Kalman track Z pulls
-  if(!(h = (TH2I*)gROOT->FindObject("hMcTrkTOFZPull"))){
-    h = new TH2I("hMcTrkTOFZPull", "Track[TOF] Z Pulls", 100, -1., 1., 100, -4.5, 4.5);
-    h->GetXaxis()->SetTitle("tg(#theta)");
-    h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
-    h->GetZaxis()->SetTitle("entries");
-  } else h->Reset();
-  arr->AddAt(h, 1);
+  // TRDout TRACK RESOLUTION
+  fContainer->AddAt(BuildMonitorContainerTrack("MCtrkOUT"), kMCtrackOut);
 
   return fContainer;
 }
 
+//________________________________________________________
+Bool_t AliTRDresolution::Load(const Char_t *file, const Char_t *dir)
+{
+// Custom load function. Used to guess the segmentation level of the data.
+
+  if(!AliTRDrecoTask::Load(file, dir)) return kFALSE;
+
+  // look for cluster residual plot - always available
+  TH3S* h3((TH3S*)((TObjArray*)fContainer->At(kClToTrk))->At(0));
+  Int_t segmentation(h3->GetNbinsZ()/2);
+  if(segmentation==fgkNresYsegm[0]){ // default segmentation. Nothing to do
+    return kTRUE;
+  } else if(segmentation==fgkNresYsegm[1]){ // stack segmentation.
+    SetSegmentationLevel(1);
+  } else if(segmentation==fgkNresYsegm[2]){ // detector segmentation.
+    SetSegmentationLevel(2);
+  } else {
+    AliError(Form("Unknown segmentation [%d].", h3->GetNbinsZ()));
+    return kFALSE;
+  }
+
+  AliDebug(2, Form("Segmentation set to level \"%s\"", fgkResYsegmName[fSegmentLevel]));
+  return kTRUE;
+}
+
+//________________________________________________________
+Bool_t AliTRDresolution::Process(TH2* const h2, TGraphErrors **g, Int_t stat)
+{
+// Generic function to process sigma/mean for 2D plot dy(x)
+
+  if(!h2) {
+    if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>0) printf("D-AliTRDresolution::Process() : NULL pointer input container.\n");
+    return kFALSE;
+  }
+  if(!Int_t(h2->GetEntries())){
+    if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>0) printf("D-AliTRDresolution::Process() : Empty h[%s - %s].\n", h2->GetName(), h2->GetTitle());
+    return kFALSE;
+  }
+  if(!g || !g[0]|| !g[1]) {
+    if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>0) printf("D-AliTRDresolution::Process() : NULL pointer output container.\n");
+    return kFALSE;
+  }
+  // prepare
+  TAxis *ax(h2->GetXaxis()), *ay(h2->GetYaxis());
+  TF1 f("f", "gaus", ay->GetXmin(), ay->GetXmax());
+  Int_t n(0);
+  if((n=g[0]->GetN())) for(;n--;) g[0]->RemovePoint(n);
+  if((n=g[1]->GetN())) for(;n--;) g[1]->RemovePoint(n);
+  TH1D *h(NULL);
+  if((h=(TH1D*)gROOT->FindObject("py"))) delete h;
+
+  // do actual loop
+  for(Int_t ix = 1, np=0; ix<=ax->GetNbins(); ix++){
+    Double_t x = ax->GetBinCenter(ix);
+    Double_t ex= ax->GetBinWidth(ix)*0.288; // w/sqrt(12)
+    h = h2->ProjectionY("py", ix, ix);
+    if((n=(Int_t)h->GetEntries())<stat){
+      if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>1) printf("I-AliTRDresolution::Process() : Low statistics @ x[%f] stat[%d]=%d [%d].\n", x, ix, n, stat);
+      continue;
+    }
+    f.SetParameter(1, 0.);
+    f.SetParameter(2, 3.e-2);
+    h->Fit(&f, "QN");
+    g[0]->SetPoint(np, x, f.GetParameter(1));
+    g[0]->SetPointError(np, ex, f.GetParError(1));
+    g[1]->SetPoint(np, x, f.GetParameter(2));
+    g[1]->SetPointError(np, ex, f.GetParError(2));
+    np++;
+  }
+  return kTRUE;
+}
+
+
 //________________________________________________________
 Bool_t AliTRDresolution::Process(TH2 * const h2, TF1 *f, Float_t k, TGraphErrors **g)
 {
@@ -1831,21 +2569,30 @@ Bool_t AliTRDresolution::Process(TH2 * const h2, TF1 *f, Float_t k, TGraphErrors
   Int_t n = 0;
   if((n=g[0]->GetN())) for(;n--;) g[0]->RemovePoint(n);
   if((n=g[1]->GetN())) for(;n--;) g[1]->RemovePoint(n);
+  if(Int_t(h2->GetEntries())){ 
+    AliDebug(4, Form("%s: g[%s %s]", pn, g[0]->GetName(), g[0]->GetTitle()));
+  } else {
+    AliDebug(2, Form("%s: g[%s %s]: Missing entries.", pn, g[0]->GetName(), g[0]->GetTitle()));
+    fIdxPlot++;
+    return kTRUE;
+  }
 
-  for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
-    Double_t x = h2->GetXaxis()->GetBinCenter(ibin);
-    TH1D *h = h2->ProjectionY(pn, ibin, ibin);
-    if(h->GetEntries()<100) continue;
-    //AdjustF1(h, f);
-
+  const Int_t kINTEGRAL=1;
+  for(Int_t ibin = 0; ibin < Int_t(h2->GetNbinsX()/kINTEGRAL); ibin++){
+    Int_t abin(ibin*kINTEGRAL+1),bbin(abin+kINTEGRAL-1),mbin(abin+Int_t(kINTEGRAL/2));
+    Double_t x = h2->GetXaxis()->GetBinCenter(mbin);
+    TH1D *h = h2->ProjectionY(pn, abin, bbin);
+    if((n=(Int_t)h->GetEntries())<150){ 
+      AliDebug(4, Form("  x[%f] range[%d %d] stat[%d] low statistics !", x, abin, bbin, n));
+      continue;
+    }
     h->Fit(f, "QN");
-    
     Int_t ip = g[0]->GetN();
+    AliDebug(4, Form("  x_%d[%f] range[%d %d] stat[%d] M[%f] Sgm[%f]", ip, x, abin, bbin, n, f->GetParameter(1), f->GetParameter(2)));
     g[0]->SetPoint(ip, x, k*f->GetParameter(1));
     g[0]->SetPointError(ip, 0., k*f->GetParError(1));
     g[1]->SetPoint(ip, x, k*f->GetParameter(2));
     g[1]->SetPointError(ip, 0., k*f->GetParError(2));
-    
 /*  
     g[0]->SetPoint(ip, x, k*h->GetMean());
     g[0]->SetPointError(ip, 0., k*h->GetMeanError());
@@ -1866,8 +2613,20 @@ Bool_t AliTRDresolution::Process2D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F
   if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
 
   // retrive containers
-  TH2I *h2 = idx<0 ? (TH2I*)(fContainer->At(plot)) : (TH2I*)((TObjArray*)(fContainer->At(plot)))->At(idx);
-  if(!h2) return kFALSE;
+  TH2I *h2(NULL);
+  if(idx<0){
+    if(!(h2= (TH2I*)(fContainer->At(plot)))) return kFALSE; 
+  } else{
+    TObjArray *a0(NULL);
+    if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE;
+    if(!(h2=(TH2I*)a0->At(idx))) return kFALSE;
+  }
+  if(Int_t(h2->GetEntries())){ 
+    AliDebug(2, Form("p[%d] idx[%d] : h[%s] %s", plot, idx, h2->GetName(), h2->GetTitle()));
+  } else {
+    AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx));
+    return kFALSE;
+  }
 
   TGraphErrors *g[2];
   if(gidx<0) gidx=idx;
@@ -1888,8 +2647,20 @@ Bool_t AliTRDresolution::Process3D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F
   if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
 
   // retrive containers
-  TH3S *h3 = idx<0 ? (TH3S*)(fContainer->At(plot)) : (TH3S*)((TObjArray*)(fContainer->At(plot)))->At(idx);
-  if(!h3) return kFALSE;
+  TH3S *h3(NULL);
+  if(idx<0){
+    if(!(h3= (TH3S*)(fContainer->At(plot)))) return kFALSE; 
+  } else{
+    TObjArray *a0(NULL);
+    if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE;
+    if(!(h3=(TH3S*)a0->At(idx))) return kFALSE;
+  }
+  if(Int_t(h3->GetEntries())){ 
+    AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle()));
+  } else {
+    AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx));
+    return kFALSE;
+  }
 
   TObjArray *gm, *gs;
   if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
@@ -1897,18 +2668,19 @@ Bool_t AliTRDresolution::Process3D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F
   TGraphErrors *g[2];
 
   TAxis *az = h3->GetZaxis();
-  for(Int_t iz=1; iz<=az->GetNbins(); iz++){
-    if(!(g[0] = (TGraphErrors*)gm->At(iz-1))) return kFALSE;
-    if(!(g[1] = (TGraphErrors*)gs->At(iz-1))) return kFALSE;
-    az->SetRange(iz, iz);
+  for(Int_t iz(0); iz<gm->GetEntriesFast(); iz++){
+    if(!(g[0] = (TGraphErrors*)gm->At(iz))) return kFALSE;
+    if(!(g[1] = (TGraphErrors*)gs->At(iz))) return kFALSE;
+    az->SetRange(iz+1, iz+1);
     if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
   }
 
   return kTRUE;
 }
 
+
 //________________________________________________________
-Bool_t AliTRDresolution::Process3Drange(ETRDresolutionPlot plot, Int_t hidx, Int_t gidx, TF1 *f, Float_t k, Int_t zbin0, Int_t zbin1)
+Bool_t AliTRDresolution::Process3Dlinked(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
 {
   //
   // Do the processing
@@ -1917,19 +2689,38 @@ Bool_t AliTRDresolution::Process3Drange(ETRDresolutionPlot plot, Int_t hidx, Int
   if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
 
   // retrive containers
-  TH3S *h3 = hidx<0 ? (TH3S*)(fContainer->At(plot)) : (TH3S*)((TObjArray*)(fContainer->At(plot)))->At(hidx);
-  if(!h3) return kFALSE;
+  TH3S *h3(NULL);
+  if(idx<0){
+    if(!(h3= (TH3S*)(fContainer->At(plot)))) return kFALSE; 
+  } else{
+    TObjArray *a0(NULL);
+    if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE;
+    if(!(h3=(TH3S*)a0->At(idx))) return kFALSE;
+  }
+  if(Int_t(h3->GetEntries())){ 
+    AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle()));
+  } else {
+    AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx));
+    return kFALSE;
+  }
 
+  TObjArray *gm, *gs;
+  if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
+  if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
   TGraphErrors *g[2];
-  if(!(g[0] = (TGraphErrors*)((TObjArray*)(fGraphM->At(plot)))->At(gidx))) return kFALSE;
-  if(!(g[1] = (TGraphErrors*)((TObjArray*)(fGraphS->At(plot)))->At(gidx))) return kFALSE;
 
-  TAxis *az = h3->GetZaxis();
-  az->SetRange(zbin0, zbin1);
+  if(!(g[0] = (TGraphErrors*)gm->At(0))) return kFALSE;
+  if(!(g[1] = (TGraphErrors*)gs->At(0))) return kFALSE;
   if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
+
+  if(!(g[0] = (TGraphErrors*)gm->At(1))) return kFALSE;
+  if(!(g[1] = (TGraphErrors*)gs->At(1))) return kFALSE;
+  if(!Process((TH2*)h3->Project3D("zx"), f, k, g)) return kFALSE;
+
   return kTRUE;
 }
 
+
 //________________________________________________________
 Bool_t AliTRDresolution::Process3DL(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
 {
@@ -1942,6 +2733,7 @@ Bool_t AliTRDresolution::Process3DL(ETRDresolutionPlot plot, Int_t idx, TF1 *f,
   // retrive containers
   TH3S *h3 = (TH3S*)((TObjArray*)fContainer->At(plot))->At(idx);
   if(!h3) return kFALSE;
+  AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle()));
 
   TGraphAsymmErrors *gm; 
   TGraphErrors *gs;
@@ -1974,7 +2766,45 @@ Bool_t AliTRDresolution::Process3DL(ETRDresolutionPlot plot, Int_t idx, TF1 *f,
 }
 
 //________________________________________________________
-Bool_t AliTRDresolution::Process4D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k, Int_t n)
+Bool_t AliTRDresolution::Process2Darray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
+{
+  //
+  // Do the processing
+  //
+
+  if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
+
+  // retrive containers
+  TObjArray *arr = (TObjArray*)(fContainer->At(plot));
+  if(!arr) return kFALSE;
+  AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName()));
+
+  TObjArray *gm, *gs;
+  if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
+  if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
+
+  TGraphErrors *g[2]; TH2I *h2(NULL); TObjArray *a0(NULL);
+  for(Int_t ia(0); ia<arr->GetEntriesFast(); ia++){
+    if(!(a0 = (TObjArray*)arr->At(ia))) continue;
+
+    if(!(h2 = (TH2I*)a0->At(idx))) return kFALSE;
+    if(Int_t(h2->GetEntries())){ 
+      AliDebug(4, Form("   idx[%d] h[%s] %s", ia, h2->GetName(), h2->GetTitle()));
+    } else {
+      AliDebug(2, Form("   idx[%d] : Missing entries.", ia));
+      continue;
+    }
+
+    if(!(g[0] = (TGraphErrors*)gm->At(ia))) return kFALSE;
+    if(!(g[1] = (TGraphErrors*)gs->At(ia))) return kFALSE;
+    if(!Process(h2, f, k, g)) return kFALSE;
+  }
+
+  return kTRUE;
+}
+
+//________________________________________________________
+Bool_t AliTRDresolution::Process3Darray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
 {
   //
   // Do the processing
@@ -1984,173 +2814,153 @@ Bool_t AliTRDresolution::Process4D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F
   //printf("Process4D : processing plot[%d] idx[%d]\n", plot, idx);
 
   // retrive containers
-  TObjArray *arr = (TObjArray*)((TObjArray*)(fContainer->At(plot)))->At(idx);
+  TObjArray *arr = (TObjArray*)(fContainer->At(plot));
   if(!arr) return kFALSE;
+  AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName()));
 
   TObjArray *gm, *gs;
   if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
   if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
 
-  TGraphErrors *g[2];
+  TGraphErrors *g[2]; TH3S *h3(NULL); TObjArray *a0(NULL);
+  Int_t in(0);
+  for(Int_t ia(0); ia<arr->GetEntriesFast(); ia++){
+    if(!(a0 = (TObjArray*)arr->At(ia))) continue;
 
-  TH3S *h3(NULL);
-  for(Int_t ix=0, in=0; ix<arr->GetEntriesFast(); ix++){
-    if(!(h3 = (TH3S*)arr->At(ix))) return kFALSE;
+    if(!(h3 = (TH3S*)a0->At(idx))) return kFALSE;
+     if(Int_t(h3->GetEntries())){ 
+       AliDebug(4, Form("   idx[%d] h[%s] %s", ia, h3->GetName(), h3->GetTitle()));
+    } else {
+      AliDebug(2, Form("   idx[%d] : Missing entries.", ia));
+      continue;
+    }
     TAxis *az = h3->GetZaxis();
-    //printf("  process ix[%d] bins[%d] in[%d]\n", ix, az->GetNbins(), in);
     for(Int_t iz=1; iz<=az->GetNbins(); iz++, in++){
-      if(n>=0 && n!=in) continue;
+      if(in >= gm->GetEntriesFast()) break;
       if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE;
       if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE;
-      //printf("    g0[%s] g1[%s]\n", g[0]->GetName(), g[1]->GetName());
       az->SetRange(iz, iz);
       if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
     }
   }
+  AliDebug(2, Form("Projections [%d] from [%d]", in, gs->GetEntriesFast()));
 
   return kTRUE;
 }
 
 //________________________________________________________
-Bool_t AliTRDresolution::GetGraphPlot(Float_t *bb, ETRDresolutionPlot ip, Int_t idx)
+Bool_t AliTRDresolution::Process3DlinkedArray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
 {
   //
-  // Get the graphs
+  // Do the processing
   //
 
-  if(!fGraphS || !fGraphM) return kFALSE;
-
-  //printf("plotting task[%d] gidx[%d]\n", ip, idx);
-  TGraphErrors *gm = idx<0 ? (TGraphErrors*)fGraphM->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphM->At(ip)))->At(idx);
-  if(!gm) return kFALSE;
-  TGraphErrors *gs = idx<0 ? (TGraphErrors*)fGraphS->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphS->At(ip)))->At(idx);
-  if(!gs) return kFALSE;
-  //printf("gs[%s] gm[%s]\n", gs->GetName(), gm->GetName());
-  gs->Draw("apl"); gm->Draw("pl");
-  //return kTRUE;
-  // titles look up
-  Int_t nref = 0;
-  for(Int_t jp=0; jp<(Int_t)ip; jp++) nref+=fgNproj[jp];
-  UChar_t jdx = idx<0?0:idx;
-  for(Int_t jc=0; jc<TMath::Min(jdx,fgNproj[ip]-1); jc++) nref++;
-  const Char_t **at = fgAxTitle[nref];
-
-  Int_t n(0);
-  if((n=gm->GetN())) {
-    PutTrendValue(Form("%s_%s", fgPerformanceName[ip], at[0]), gm->GetMean(2));
-    PutTrendValue(Form("%s_%sRMS", fgPerformanceName[ip], at[0]), gm->GetRMS(2));
-  }
+  if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
+  //printf("Process4D : processing plot[%d] idx[%d]\n", plot, idx);
 
-  if((n=gs->GetN())){
-    gs->Sort(&TGraph::CompareY);
-    PutTrendValue(Form("%s_%sSigMin", fgPerformanceName[ip], at[0]), gs->GetY()[0]);
-    PutTrendValue(Form("%s_%sSigMax", fgPerformanceName[ip], at[0]), gs->GetY()[n-1]);
-    gs->Sort(&TGraph::CompareX); 
-  }
+  // retrive containers
+  TObjArray *arr = (TObjArray*)(fContainer->At(plot));
+  if(!arr) return kFALSE;
+  AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName()));
 
-  // axis range
-  TAxis *ax = NULL;
-  ax = gs->GetHistogram()->GetXaxis();
-  ax->SetRangeUser(bb[0], bb[2]);
-  ax->SetTitle(at[1]);ax->CenterTitle();
+  TObjArray *gm, *gs;
+  if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
+  if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
 
-  ax = gs->GetHistogram()->GetYaxis();
-  ax->SetRangeUser(bb[1], bb[3]);
-  ax->SetTitleOffset(1.1);
-  ax->SetTitle(at[2]);ax->CenterTitle();
+  TGraphErrors *g[2]; TH3S *h3(NULL); TObjArray *a0(NULL);
+  Int_t in(0);
+  for(Int_t ia(0); ia<arr->GetEntriesFast(); ia++){
+    if(!(a0 = (TObjArray*)arr->At(ia))) continue;
+    if(!(h3 = (TH3S*)a0->At(idx))) return kFALSE;
+    if(Int_t(h3->GetEntries())){     
+      AliDebug(4, Form("   idx[%d] h[%s] %s", ia, h3->GetName(), h3->GetTitle()));
+    } else {
+      AliDebug(2, Form("   idx[%d] : Missing entries.", ia));
+      continue;
+    }
+    if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE;
+    if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE;
+    if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
+    in++;
 
-  TGaxis *gax = NULL;
-  gax = new TGaxis(bb[2], bb[1], bb[2], bb[3], bb[1], bb[3], 510, "+U");
-  gax->SetLineColor(kRed);gax->SetLineWidth(2);gax->SetTextColor(kRed);
-  //gax->SetVertical();
-  gax->CenterTitle(); gax->SetTitleOffset(.7);
-  gax->SetTitle(at[3]); gax->Draw();
+    if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE;
+    if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE;
+    if(!Process((TH2*)h3->Project3D("zx"), f, k, g)) return kFALSE;
+    in++;
+  }
+  AliDebug(2, Form("Projections [%d] from [%d]", in, gs->GetEntriesFast()));
 
-  // bounding box
-  TBox *b = new TBox(-.15, bb[1], .15, bb[3]);
-  b->SetFillStyle(3002);b->SetLineColor(0);
-  b->SetFillColor(ip<=Int_t(kMCcluster)?kGreen:kBlue);
-  b->Draw();
   return kTRUE;
 }
 
-Char_t const *fgLegHeader[6]={
-  "Layer 0","Layer 1","Layer 2","Layer 3","Layer 4","Layer 5"
-};
 //________________________________________________________
-Bool_t AliTRDresolution::GetGraphTrack(Float_t *bb, Int_t idx, Int_t il, Int_t n, Bool_t kLEG)
+Bool_t AliTRDresolution::GetGraph(Float_t *bb, ETRDresolutionPlot ip, Int_t idx, Bool_t kLEG, const Char_t *explain)
 {
   //
   // Get the graphs
   //
 
   if(!fGraphS || !fGraphM) return kFALSE;
-
   // axis titles look up
-  Int_t nref(0);
-  for(Int_t jp=0; jp<Int_t(kMCtrackTRD); jp++) nref+=fgNproj[jp];
-  nref+=idx;
-  const Char_t **at = fgAxTitle[nref];
-  printf("nref[%d] ax[%s] x[%f %f] y[%f %f]\n", nref, at[0], bb[0], bb[2], bb[1], bb[3]);
+  Int_t nref = 0;
+  for(Int_t jp=0; jp<(Int_t)ip; jp++) nref+=fgNproj[jp];
+  UChar_t jdx = idx<0?0:idx;
+  for(Int_t jc=0; jc<TMath::Min(jdx,fgNproj[ip]-1); jc++) nref++;
+  Char_t **at = fAxTitle[nref];
 
+  // build legends if requiered
   TLegend *leg(NULL);
   if(kLEG){
-    leg=new TLegend(.65, .7, .95, .99);
-    leg->SetHeader("Species");
-    leg->SetBorderSize(1);
-    leg->SetFillColor(0);
+    leg=new TLegend(.65, .6, .95, .9);
+    leg->SetBorderSize(0);
+    leg->SetFillStyle(0);
   }
+  // build frame
   TH1S *h1(NULL);
-  h1 = new TH1S(Form("h1TF_%d", fIdxFrame++), at[0], 2, bb[0], bb[2]);
+  h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), Form("%s %s;%s;%s", at[0], explain?explain:"", at[1], at[2]), 2, bb[0], bb[2]);
   h1->SetMinimum(bb[1]);h1->SetMaximum(bb[3]);
   h1->SetLineColor(kBlack); h1->SetLineWidth(1);h1->SetLineStyle(2); 
   // axis range
   TAxis *ax = h1->GetXaxis();
-  ax->SetTitle(at[1]);ax->CenterTitle();ax->SetMoreLogLabels();
+  ax->CenterTitle();ax->SetMoreLogLabels();ax->SetTitleOffset(1.2);
   ax = h1->GetYaxis();
   ax->SetRangeUser(bb[1], bb[3]);
-  ax->SetTitleOffset(.5);//ax->SetTitleSize(.06);
-  ax->SetTitle(at[2]);ax->CenterTitle();
+  ax->CenterTitle(); ax->SetTitleOffset(1.4);
   h1->Draw();
-  TGaxis *gax = NULL;
-  gax = new TGaxis(bb[2], bb[1], bb[2], bb[3], bb[1], bb[3], 510, "+U");
-  gax->SetLineColor(kRed);gax->SetLineWidth(2);gax->SetTextColor(kRed);
-  //gax->SetVertical();
-  gax->CenterTitle(); //gax->SetTitleOffset(.5);gax->SetTitleSize(.06);
-  gax->SetTitle(at[3]); gax->Draw();
-
-  TGraphErrors *gm = NULL, *gs = NULL;
-  TObjArray *a0 = NULL, *a1 = NULL;
-  a0 = (TObjArray*)((TObjArray*)fGraphM->At(kMCtrackTRD))->At(idx);
-  a1 = (TObjArray*)((TObjArray*)fGraphS->At(kMCtrackTRD))->At(idx);
-  Int_t nn(0);
-  for(Int_t is=il, is0=0; is<il+n; is++, is0++){
-    if(is0==3) continue;
-    if(!(gs =  (TGraphErrors*)a1->At(is))) return kFALSE;
-    if(!(gm =  (TGraphErrors*)a0->At(is))) return kFALSE;
+  // bounding box
+  TBox *b = new TBox(-.15, bb[1], .15, bb[3]);
+  b->SetFillStyle(3002);b->SetLineColor(0);
+  b->SetFillColor(ip<=Int_t(kMCcluster)?kGreen:kBlue);
+  b->Draw();
 
-    if((nn=gs->GetN())){
-      gs->Draw("pl");
-      if(leg) leg->AddEntry(gs, gs->GetTitle(), "p");
-      //printf("entry[%d]=%s\n", il+is, gs->GetTitle());
-      gs->Sort(&TGraph::CompareY);
-      PutTrendValue(Form("%s_%sSigMin%s", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[0]);
-      PutTrendValue(Form("%s_%sSigMax%s", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[nn-1]);
-      gs->Sort(&TGraph::CompareX); 
-    }
-    if(gm->GetN()){
-      gm->Draw("pl");
-      PutTrendValue(Form("%s_%s_%s", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gm->GetMean(2));
-      PutTrendValue(Form("%s_%s_%sRMS", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gm->GetRMS(2));
-    }
+  TGraphErrors *gm = idx<0 ? (TGraphErrors*)fGraphM->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphM->At(ip)))->At(idx);
+  if(!gm) return kFALSE;
+  TGraphErrors *gs = idx<0 ? (TGraphErrors*)fGraphS->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphS->At(ip)))->At(idx);
+  if(!gs) return kFALSE;
+
+  Int_t n(0), nPlots(0);
+  if((n=gm->GetN())) {
+    nPlots++;
+    gm->Draw("pl"); if(leg) leg->AddEntry(gm, gm->GetTitle(), "pl");
+    PutTrendValue(Form("%s_%s", fgPerformanceName[ip], at[0]), gm->GetMean(2));
+    PutTrendValue(Form("%s_%sRMS", fgPerformanceName[ip], at[0]), gm->GetRMS(2));
   }
+
+  if((n=gs->GetN())){
+    nPlots++;
+    gs->Draw("pl"); if(leg) leg->AddEntry(gs, gs->GetTitle(), "pl");
+    gs->Sort(&TGraph::CompareY);
+    PutTrendValue(Form("%s_%sSigMin", fgPerformanceName[ip], at[0]), gs->GetY()[0]);
+    PutTrendValue(Form("%s_%sSigMax", fgPerformanceName[ip], at[0]), gs->GetY()[n-1]);
+    gs->Sort(&TGraph::CompareX); 
+  }
+  if(!nPlots) return kFALSE;
   if(leg) leg->Draw();
   return kTRUE;
 }
 
-
 //________________________________________________________
-Bool_t AliTRDresolution::GetGraphTrackTPC(Float_t *bb, Int_t idx, Int_t ist, Int_t n, Bool_t kLEG)
+Bool_t AliTRDresolution::GetGraphArray(Float_t *bb, ETRDresolutionPlot ip, Int_t idx, Bool_t kLEG, Int_t n, Int_t *sel, const Char_t *explain)
 {
   //
   // Get the graphs
@@ -2159,60 +2969,221 @@ Bool_t AliTRDresolution::GetGraphTrackTPC(Float_t *bb, Int_t idx, Int_t ist, Int
   if(!fGraphS || !fGraphM) return kFALSE;
 
   // axis titles look up
-  Int_t nref = 0;
-  for(Int_t jp=0; jp<Int_t(kMCtrackTPC); jp++) nref+=fgNproj[jp];
+  Int_t nref(0);
+  for(Int_t jp(0); jp<ip; jp++) nref+=fgNproj[jp];
   nref+=idx;
-  const Char_t **at = fgAxTitle[nref];
+  Char_t **at = fAxTitle[nref];
 
-  TLegend *leg(NULL);
+  // build legends if requiered
+  TLegend *legM(NULL), *legS(NULL);
   if(kLEG){
-    leg=new TLegend(.65, .7, .95, .99);
-    leg->SetHeader("Species");
-    leg->SetBorderSize(1);
-    leg->SetFillColor(0);
+    legM=new TLegend(.35, .6, .65, .9);
+    legM->SetHeader("Mean");
+    legM->SetBorderSize(0);
+    legM->SetFillStyle(0);
+    legS=new TLegend(.65, .6, .95, .9);
+    legS->SetHeader("Sigma");
+    legS->SetBorderSize(0);
+    legS->SetFillStyle(0);
   }
+  // build frame
   TH1S *h1(NULL);
-  h1 = new TH1S(Form("h1TF_%d", fIdxFrame++), at[0], 2, bb[0], bb[2]);
+  h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), Form("%s %s;%s;%s", at[0], explain?explain:"", at[1], at[2]), 2, bb[0], bb[2]);
   h1->SetMinimum(bb[1]);h1->SetMaximum(bb[3]);
   h1->SetLineColor(kBlack); h1->SetLineWidth(1);h1->SetLineStyle(2); 
   // axis range
   TAxis *ax = h1->GetXaxis();
-  ax->SetTitle(at[1]);ax->CenterTitle();ax->SetMoreLogLabels();
+  ax->CenterTitle();ax->SetMoreLogLabels();ax->SetTitleOffset(1.2);
   ax = h1->GetYaxis();
   ax->SetRangeUser(bb[1], bb[3]);
-  ax->SetTitleOffset(.5);//ax->SetTitleSize(.06);
-  ax->SetTitle(at[2]);ax->CenterTitle();
+  ax->CenterTitle(); ax->SetTitleOffset(1.4);
   h1->Draw();
-  TGaxis *gax = NULL;
-  gax = new TGaxis(bb[2], bb[1], bb[2], bb[3], bb[1], bb[3], 510, "+U");
-  gax->SetLineColor(kRed);gax->SetLineWidth(2);gax->SetTextColor(kRed);
-  //gax->SetVertical();
-  gax->CenterTitle(); //gax->SetTitleOffset(.5);gax->SetTitleSize(.06);
-  gax->SetTitle(at[3]); gax->Draw();
-
-  Int_t nn(0);
-  TGraphErrors *gm = NULL, *gs = NULL;
-  TObjArray *a0 = NULL, *a1 = NULL;
-  a0 = (TObjArray*)((TObjArray*)fGraphM->At(kMCtrackTPC))->At(idx);
-  a1 = (TObjArray*)((TObjArray*)fGraphS->At(kMCtrackTPC))->At(idx);
-  for(Int_t is=ist, is0=0; is<ist+n; is++, is0++){
-    if(is0==3) continue;
-    if(!(gs =  (TGraphErrors*)a1->At(is))) return kFALSE;
-    if(!(gm =  (TGraphErrors*)a0->At(is))) return kFALSE;
+
+  TGraphErrors *gm(NULL), *gs(NULL);
+  TObjArray *a0(NULL), *a1(NULL);
+  a0 = (TObjArray*)((TObjArray*)fGraphM->At(ip))->At(idx);
+  a1 = (TObjArray*)((TObjArray*)fGraphS->At(ip))->At(idx);
+  if(!n) n=a0->GetEntriesFast();
+  AliDebug(4, Form("Graph : Ref[%d] Title[%s] Limits{x[%f %f] y[%f %f]} Comp[%d] Selection[%c]", nref, at[0], bb[0], bb[2], bb[1], bb[3], n, sel ? 'y' : 'n'));
+  Int_t nn(0), nPlots(0);
+  for(Int_t is(0), is0(0); is<n; is++){
+    is0 = sel ? sel[is] : is;
+    if(!(gs =  (TGraphErrors*)a1->At(is0))) return kFALSE;
+    if(!(gm =  (TGraphErrors*)a0->At(is0))) return kFALSE;
+
     if((nn=gs->GetN())){
-      gs->Draw("pl");if(leg) leg->AddEntry(gs, gs->GetTitle(), "p");
+      nPlots++;
+      gs->Draw("pc"); 
+      if(legS){ 
+        //printf("LegEntry %s [%s]%s\n", at[0], gs->GetName(), gs->GetTitle());
+        legS->AddEntry(gs, gs->GetTitle(), "pl");
+      }
       gs->Sort(&TGraph::CompareY);
-      PutTrendValue(Form("%s_%sSigMin%s", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[0]);
-      PutTrendValue(Form("%s_%sSigMax%s", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[nn-1]);
+      PutTrendValue(Form("%s_%sSigMin", fgPerformanceName[kMCtrack], at[0]), gs->GetY()[0]);
+      PutTrendValue(Form("%s_%sSigMax", fgPerformanceName[kMCtrack], at[0]), gs->GetY()[nn-1]);
       gs->Sort(&TGraph::CompareX); 
     }
     if(gm->GetN()){
-      gm->Draw("pl");
-      PutTrendValue(Form("%s_%s_%s", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gm->GetMean(2));
-      PutTrendValue(Form("%s_%s_%sRMS", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gm->GetRMS(2));
+      nPlots++;
+      gm->Draw("pc");
+      if(legM){ 
+        //printf("LegEntry %s [%s]%s\n", at[0], gm->GetName(), gm->GetTitle());
+        legM->AddEntry(gm, gm->GetTitle(), "pl");
+      }
+      PutTrendValue(Form("%s_%s", fgPerformanceName[kMCtrack], at[0]), gm->GetMean(2));
+      PutTrendValue(Form("%s_%sRMS", fgPerformanceName[kMCtrack], at[0]), gm->GetRMS(2));
     }
   }
-  if(leg) leg->Draw();
+  if(!nPlots) return kFALSE;
+  if(kLEG){
+    legM->Draw();
+    legS->Draw();
+  }
+  return kTRUE;
+}
+
+//____________________________________________________________________
+Bool_t AliTRDresolution::FitTrack(const Int_t np, AliTrackPoint *points, Float_t param[10])
+{
+//
+// Fit track with a staight line using the "np" clusters stored in the array "points".
+// The following particularities are stored in the clusters from points:
+//   1. pad tilt as cluster charge
+//   2. pad row cross or vertex constrain as fake cluster with cluster type 1
+// The parameters of the straight line fit are stored in the array "param" in the following order :
+//     param[0] - x0 reference radial position
+//     param[1] - y0 reference r-phi position @ x0
+//     param[2] - z0 reference z position @ x0
+//     param[3] - slope dy/dx
+//     param[4] - slope dz/dx
+//
+// Attention :
+// Function should be used to refit tracks for B=0T
+//
+
+  if(np<40){
+    if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>1) printf("D-AliTRDresolution::FitTrack: Not enough clusters to fit a track [%d].", np);
+    return kFALSE;
+  }
+  TLinearFitter yfitter(2, "pol1"), zfitter(2, "pol1");
+
+  Double_t x0(0.);
+  for(Int_t ip(0); ip<np; ip++) x0+=points[ip].GetX();
+  x0/=Float_t(np);
+
+  Double_t x, y, z, dx, tilt(0.);
+  for(Int_t ip(0); ip<np; ip++){
+    x = points[ip].GetX(); z = points[ip].GetZ();
+    dx = x - x0;
+    zfitter.AddPoint(&dx, z, points[ip].GetClusterType()?1.e-3:1.);
+  }
+  if(zfitter.Eval() != 0) return kFALSE;
+
+  Double_t z0    = zfitter.GetParameter(0);
+  Double_t dzdx  = zfitter.GetParameter(1);
+  for(Int_t ip(0); ip<np; ip++){
+    if(points[ip].GetClusterType()) continue;
+    x    = points[ip].GetX();
+    dx   = x - x0;
+    y    = points[ip].GetY();
+    z    = points[ip].GetZ();
+    tilt = points[ip].GetCharge();
+    y -= tilt*(-dzdx*dx + z - z0);
+    Float_t xyz[3] = {x, y, z}; points[ip].SetXYZ(xyz);
+    yfitter.AddPoint(&dx, y, 1.);
+  }
+  if(yfitter.Eval() != 0) return kFALSE;
+  Double_t y0   = yfitter.GetParameter(0);
+  Double_t dydx = yfitter.GetParameter(1);
+
+  param[0] = x0; param[1] = y0; param[2] = z0; param[3] = dydx; param[4] = dzdx;
+  if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>3) printf("D-AliTRDresolution::FitTrack: x0[%f] y0[%f] z0[%f] dydx[%f] dzdx[%f]\n", x0, y0, z0, dydx, dzdx);
+  return kTRUE;
+}
+
+//____________________________________________________________________
+Bool_t AliTRDresolution::FitTracklet(const Int_t ly, const Int_t np, AliTrackPoint *points, const Float_t param[10], Float_t par[3])
+{
+//
+// Fit tracklet with a staight line using the coresponding subset of clusters out of the total "np" clusters stored in the array "points".
+// See function FitTrack for the data stored in the "clusters" array
+
+// The parameters of the straight line fit are stored in the array "param" in the following order :
+//     par[0] - x0 reference radial position
+//     par[1] - y0 reference r-phi position @ x0
+//     par[2] - slope dy/dx
+//
+// Attention :
+// Function should be used to refit tracks for B=0T
+//
+
+  TLinearFitter yfitter(2, "pol1");
+
+  // grep data for tracklet
+  Double_t x0(0.), x[60], y[60], dy[60];
+  Int_t nly(0);
+  for(Int_t ip(0); ip<np; ip++){
+    if(points[ip].GetClusterType()) continue;
+    if(points[ip].GetVolumeID() != ly) continue;
+    Float_t xt(points[ip].GetX())
+           ,yt(param[1] + param[3] * (xt - param[0]));
+    x[nly] = xt;
+    y[nly] = points[ip].GetY();
+    dy[nly]= y[nly]-yt;
+    x0    += xt;
+    nly++;
+  }
+  if(nly<10){
+    if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>1) printf("D-AliTRDresolution::FitTracklet: Not enough clusters to fit a tracklet [%d].", nly);
+    return kFALSE;
+  }
+  // set radial reference for fit
+  x0 /= Float_t(nly);
+
+  // find tracklet core
+  Double_t mean(0.), sig(1.e3);
+  AliMathBase::EvaluateUni(nly, dy, mean, sig, 0);
+
+  // simple cluster error parameterization
+  Float_t kSigCut = TMath::Sqrt(5.e-4 + param[3]*param[3]*0.018);
+
+  // fit tracklet core
+  for(Int_t jly(0); jly<nly; jly++){
+    if(TMath::Abs(dy[jly]-mean)>kSigCut) continue;
+    Double_t dx(x[jly]-x0);
+    yfitter.AddPoint(&dx, y[jly], 1.);
+  }
+  if(yfitter.Eval() != 0) return kFALSE;
+  par[0] = x0;
+  par[1] = yfitter.GetParameter(0);
+  par[2] = yfitter.GetParameter(1);
+  return kTRUE;
+}
+
+//____________________________________________________________________
+Bool_t AliTRDresolution::UseTrack(const Int_t np, AliTrackPoint *points, Float_t param[10])
+{
+//
+// Global selection mechanism of tracksbased on cluster to fit residuals
+// The parameters are the same as used ni function FitTrack().
+
+  const Float_t kS(0.6), kM(0.2);
+  TH1S h("h1", "", 100, -5.*kS, 5.*kS);
+  Float_t dy, dz, s, m;
+  for(Int_t ip(0); ip<np; ip++){
+    if(points[ip].GetClusterType()) continue;
+    Float_t x0(points[ip].GetX())
+          ,y0(param[1] + param[3] * (x0 - param[0]))
+          ,z0(param[2] + param[4] * (x0 - param[0]));
+    dy=points[ip].GetY() - y0; h.Fill(dy);
+    dz=points[ip].GetZ() - z0;
+  }
+  TF1 fg("fg", "gaus", -5.*kS, 5.*kS);
+  fg.SetParameter(1, 0.);
+  fg.SetParameter(2, 2.e-2);
+  h.Fit(&fg, "QN");
+  m=fg.GetParameter(1); s=fg.GetParameter(2);
+  if(s>kS || TMath::Abs(m)>kM) return kFALSE;
   return kTRUE;
 }
 
@@ -2243,8 +3214,96 @@ void AliTRDresolution::GetLandauMpvFwhm(TF1 * const f, Float_t &mpv, Float_t &xm
 
 
 //________________________________________________________
-void AliTRDresolution::SetRecoParam(AliTRDrecoParam *r)
+void AliTRDresolution::SetSegmentationLevel(Int_t l) 
 {
-
-  fReconstructor->SetRecoParam(r);
+// Setting the segmentation level to "l"
+  fSegmentLevel = l;
+
+  UShort_t const lNcomp[kNprojs] = {
+    1,  1, //2, 
+    fgkNresYsegm[fSegmentLevel], 2, //2, 
+    2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5, 
+    2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5,
+    2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5,
+  // MC
+    fgkNresYsegm[fSegmentLevel], 2,          //2, 
+    fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5, 
+    fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, 1, 1, 1, 11, 11, 11, //11
+    fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, 1, 1, 1, 11, 11, 11, //11
+    6*fgkNresYsegm[fSegmentLevel], 6*2, 6*2, 6*2, 6, 6, 6, 6, 6*11, 6*11, 6*11  //11
+  };
+  memcpy(fNcomp, lNcomp, kNprojs*sizeof(UShort_t));
+
+  Char_t const *lAxTitle[kNprojs][4] = {
+    // Charge
+    {"Impv", "x [cm]", "I_{mpv}", "x/x_{0}"}
+  ,{"dI/Impv", "x/x_{0}", "#delta I/I_{mpv}", "x[cm]"}
+    // Clusters to Kalman
+  ,{"Cluster2Track residuals", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
+  ,{"Cluster2Track  YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+    // TRD tracklet to Kalman fit
+  ,{"Tracklet2Track Y residuals", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
+  ,{"Tracklet2Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+  ,{"Tracklet2Track Z residuals", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
+  ,{"Tracklet2Track Z pulls", "tg(#theta)", "z", "#sigma_{z}"}
+  ,{"Tracklet2Track Phi residuals", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
+    // TRDin 2 first TRD tracklet
+  ,{"Tracklet2Track Y residuals @ TRDin", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
+  ,{"Tracklet2Track YZ pulls @ TRDin", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+  ,{"Tracklet2Track Z residuals @ TRDin", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
+  ,{"Tracklet2Track Z pulls @ TRDin", "tg(#theta)", "z", "#sigma_{z}"}
+  ,{"Tracklet2Track Phi residuals @ TRDin", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
+    // TRDout 2 first TRD tracklet
+  ,{"Tracklet2Track Y residuals @ TRDout", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
+  ,{"Tracklet2Track YZ pulls @ TRDout", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+  ,{"Tracklet2Track Z residuals @ TRDout", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
+  ,{"Tracklet2Track Z pulls @ TRDout", "tg(#theta)", "z", "#sigma_{z}"}
+  ,{"Tracklet2Track Phi residuals @ TRDout", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
+    // MC cluster
+  ,{"MC Cluster Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
+  ,{"MC Cluster YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+    // MC tracklet
+  ,{"MC Tracklet Y resolution", "tg(#phi)", "y [#mum]",  "#sigma_{y}[#mum]"}
+  ,{"MC Tracklet YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+  ,{"MC Tracklet Z resolution", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
+  ,{"MC Tracklet Z pulls", "tg(#theta)", "z", "#sigma_{z}"}
+  ,{"MC Tracklet Phi resolution", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
+    // MC track TRDin
+  ,{"MC Y resolution @ TRDin", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"}
+  ,{"MC YZ pulls @ TRDin", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+  ,{"MC Z resolution @ TRDin", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
+  ,{"MC Z pulls @ TRDin", "tg(#theta)", "z", "#sigma_{z}"}
+  ,{"MC #Phi resolution @ TRDin", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
+  ,{"MC SNP pulls @ TRDin", "tg(#phi)", "SNP", "#sigma_{snp}"}
+  ,{"MC #Theta resolution @ TRDin", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"}
+  ,{"MC TGL pulls @ TRDin", "tg(#theta)", "TGL", "#sigma_{tgl}"}
+  ,{"MC P_{t} resolution @ TRDin", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"}
+  ,{"MC 1/P_{t} pulls @ TRDin", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"}
+  ,{"MC P resolution @ TRDin", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"}
+    // MC track TRDout
+  ,{"MC Y resolution @ TRDout", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"}
+  ,{"MC YZ pulls @ TRDout", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+  ,{"MC Z resolution @ TRDout", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
+  ,{"MC Z pulls @ TRDout", "tg(#theta)", "z", "#sigma_{z}"}
+  ,{"MC #Phi resolution @ TRDout", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
+  ,{"MC SNP pulls @ TRDout", "tg(#phi)", "SNP", "#sigma_{snp}"}
+  ,{"MC #Theta resolution @ TRDout", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"}
+  ,{"MC TGL pulls @ TRDout", "tg(#theta)", "TGL", "#sigma_{tgl}"}
+  ,{"MC P_{t} resolution @ TRDout", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"}
+  ,{"MC 1/P_{t} pulls @ TRDout", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"}
+  ,{"MC P resolution @ TRDout", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"}
+    // MC track in TRD
+  ,{"MC Track Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
+  ,{"MC Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
+  ,{"MC Track Z resolution", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
+  ,{"MC Track Z pulls", "tg(#theta)", "z", "#sigma_{z}"}
+  ,{"MC Track #Phi resolution", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
+  ,{"MC Track SNP pulls", "tg(#phi)", "SNP", "#sigma_{snp}"}
+  ,{"MC Track #Theta resolution", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"}
+  ,{"MC Track TGL pulls", "tg(#theta)", "TGL", "#sigma_{tgl}"}
+  ,{"MC P_{t} resolution", "p_{t} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "#sigma(#Deltap_{t}/p_{t}^{MC}) [%]"}
+  ,{"MC 1/P_{t} pulls", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC} - 1/p_{t}^{MC}", "#sigma_{1/p_{t}}"}
+  ,{"MC P resolution", "p [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "#sigma(#Deltap/p^{MC}) [%]"}
+  };
+  memcpy(fAxTitle, lAxTitle, 4*kNprojs*sizeof(Char_t*));
 }