TRD/qaRec/AliTRDtrackingResolution.cxx

   1 /**************************************************************************
   2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3 *                                                                        *
   4 * Author: The ALICE Off-line Project.                                    *
   5 * Contributors are mentioned in the code where appropriate.              *
   6 *                                                                        *
   7 * Permission to use, copy, modify and distribute this software and its   *
   8 * documentation strictly for non-commercialf purposes is hereby granted   *
   9 * without fee, provided that the above copyright notice appears in all   *
  10 * copies and that both the copyright notice and this permission notice   *
  11 * appear in the supporting documentation. The authors make no claims     *
  12 * about the suitability of this software for any purpose. It is          *
  13 * provided "as is" without express or implied warranty.                  *
  14 **************************************************************************/
  15
  16 /* $Id: AliTRDtrackingResolution.cxx 27496 2008-07-22 08:35:45Z cblume $ */
  17
  18 ////////////////////////////////////////////////////////////////////////////
  19 //                                                                        //
  20 //  TRD tracking resolution                                               //
  21 //
  22 // The class performs resolution and residual studies
  23 // of the TRD tracks for the following quantities :
  24 //   - spatial position (y, [z])
  25 //   - angular (phi) tracklet
  26 //   - momentum at the track level
  27 //
  28 // The class has to be used for regular detector performance checks using the official macros:
  29 //   - $ALICE_ROOT/TRD/qaRec/run.C
  30 //   - $ALICE_ROOT/TRD/qaRec/makeResults.C
  31 //
  32 // For stand alone usage please refer to the following example:
  33 // {
  34 //   gSystem->Load("libANALYSIS.so");
  35 //   gSystem->Load("libTRDqaRec.so");
  36 //   AliTRDtrackingResolution *res = new AliTRDtrackingResolution();
  37 //   //res->SetMCdata();
  38 //   //res->SetVerbose();
  39 //   //res->SetVisual();
  40 //   res->Load("TRD.TaskResolution.root");
  41 //   if(!res->PostProcess()) return;
  42 //   res->GetRefFigure(0);
  43 // }
  44 //
  45 //  Authors:                                                              //
  46 //    Alexandru Bercuci <A.Bercuci@gsi.de>                                //
  47 //    Markus Fasel <M.Fasel@gsi.de>                                       //
  48 //                                                                        //
  49 ////////////////////////////////////////////////////////////////////////////
  50
  51 #include <cstring>
  52
  53 #include <TROOT.h>
  54 #include <TSystem.h>
  55 #include <TPDGCode.h>
  56 #include <TObjArray.h>
  57 #include <TH2.h>
  58 #include <TH1.h>
  59 #include <TF1.h>
  60 #include <TCanvas.h>
  61 #include <TBox.h>
  62 #include <TProfile.h>
  63 #include <TGraphErrors.h>
  64 #include <TMath.h>
  65 #include "TTreeStream.h"
  66 #include "TGeoManager.h"
  67
  68 #include "AliAnalysisManager.h"
  69 #include "AliTrackReference.h"
  70 #include "AliTrackPointArray.h"
  71 #include "AliCDBManager.h"
  72
  73 #include "AliTRDSimParam.h"
  74 #include "AliTRDgeometry.h"
  75 #include "AliTRDpadPlane.h"
  76 #include "AliTRDcluster.h"
  77 #include "AliTRDseedV1.h"
  78 #include "AliTRDtrackV1.h"
  79 #include "AliTRDtrackerV1.h"
  80 #include "AliTRDReconstructor.h"
  81 #include "AliTRDrecoParam.h"
  82
  83 #include "AliTRDtrackInfo/AliTRDclusterInfo.h"
  84 #include "AliTRDtrackInfo/AliTRDtrackInfo.h"
  85 #include "AliTRDtrackingResolution.h"
  86
  87 ClassImp(AliTRDtrackingResolution)
  88
  89 //________________________________________________________
  90 AliTRDtrackingResolution::AliTRDtrackingResolution()
  91   :AliTRDrecoTask("Resolution", "Tracking Resolution")
  92   ,fStatus(0)
  93   ,fReconstructor(0x0)
  94   ,fGeo(0x0)
  95   ,fGraphS(0x0)
  96   ,fGraphM(0x0)
  97   ,fClResiduals(0x0)
  98   ,fTrkltResiduals(0x0)
  99   ,fTrkltPhiResiduals(0x0)
 100   ,fClResolution(0x0)
 101   ,fTrkltResolution(0x0)
 102 {
 103   fReconstructor = new AliTRDReconstructor();
 104   fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
 105   fGeo = new AliTRDgeometry();
 106
 107   InitFunctorList();
 108
 109   DefineOutput(1+kCluster, TObjArray::Class());
 110   DefineOutput(1+kTrackletY, TObjArray::Class());
 111   DefineOutput(1+kTrackletPhi, TObjArray::Class());
 112   DefineOutput(1+kMCcluster, TObjArray::Class());
 113   DefineOutput(1+kMCtrackletY, TObjArray::Class());
 114 }
 115
 116 //________________________________________________________
 117 AliTRDtrackingResolution::~AliTRDtrackingResolution()
 118 {
 119   if(fGraphS){fGraphS->Delete(); delete fGraphS;}
 120   if(fGraphM){fGraphM->Delete(); delete fGraphM;}
 121   delete fGeo;
 122   delete fReconstructor;
 123   if(gGeoManager) delete gGeoManager;
 124   if(fClResiduals){fClResiduals->Delete(); delete fClResiduals;}
 125   if(fTrkltResiduals){fTrkltResiduals->Delete(); delete fTrkltResiduals;}
 126   if(fTrkltPhiResiduals){fTrkltPhiResiduals->Delete(); delete fTrkltPhiResiduals;}
 127   if(fClResolution){
 128     fClResolution->Delete();
 129     delete fClResolution;
 130   }
 131   if(fTrkltResolution){fTrkltResolution->Delete(); delete fTrkltResolution;}
 132 }
 133
 134
 135 //________________________________________________________
 136 void AliTRDtrackingResolution::CreateOutputObjects()
 137 {
 138   // spatial resolution
 139   OpenFile(0, "RECREATE");
 140
 141   fContainer = Histos();
 142
 143   fClResiduals = new TObjArray();
 144   fClResiduals->SetOwner(kTRUE);
 145   fTrkltResiduals = new TObjArray();
 146   fTrkltResiduals->SetOwner(kTRUE);
 147   fTrkltPhiResiduals = new TObjArray();
 148   fTrkltPhiResiduals->SetOwner(kTRUE);
 149   fClResolution = new TObjArray();
 150   fClResolution->SetOwner(kTRUE);
 151   fTrkltResolution = new TObjArray();
 152   fTrkltResolution->SetOwner(kTRUE);
 153 }
 154
 155 //________________________________________________________
 156 void AliTRDtrackingResolution::Exec(Option_t *opt)
 157 {
 158   fClResiduals->Delete();
 159   fTrkltResiduals->Delete();
 160   fTrkltPhiResiduals->Delete();
 161   fClResolution->Delete();
 162   fTrkltResolution->Delete();
 163
 164   AliTRDrecoTask::Exec(opt);
 165
 166   PostData(1+kCluster, fClResiduals);
 167   PostData(1+kTrackletY, fTrkltResiduals);
 168   PostData(1+kTrackletPhi, fTrkltPhiResiduals);
 169   PostData(1+kMCcluster, fClResolution);
 170   PostData(1+kMCtrackletY, fTrkltResolution);
 171 }
 172
 173 //________________________________________________________
 174 TH1* AliTRDtrackingResolution::PlotCluster(const AliTRDtrackV1 *track)
 175 {
 176   if(track) fTrack = track;
 177   if(!fTrack){
 178     AliWarning("No Track defined.");
 179     return 0x0;
 180   }
 181   TH1 *h = 0x0;
 182   if(!(h = ((TH2I*)fContainer->At(kCluster)))){
 183     AliWarning("No output histogram defined.");
 184     return 0x0;
 185   }
 186
 187   Double_t cov[3];
 188   Float_t x0, y0, z0, dy, dydx, dzdx;
 189   AliTRDseedV1 *fTracklet = 0x0;
 190   for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
 191     if(!(fTracklet = fTrack->GetTracklet(ily))) continue;
 192     if(!fTracklet->IsOK()) continue;
 193     x0 = fTracklet->GetX0();
 194
 195     // retrive the track angle with the chamber
 196     y0   = fTracklet->GetYref(0);
 197     z0   = fTracklet->GetZref(0);
 198     dydx = fTracklet->GetYref(1);
 199     dzdx = fTracklet->GetZref(1);
 200     fTracklet->GetCovRef(cov);
 201     Float_t tilt = fTracklet->GetTilt();
 202     AliTRDcluster *c = 0x0;
 203     fTracklet->ResetClusterIter(kFALSE);
 204     while((c = fTracklet->PrevCluster())){
 205       Float_t xc = c->GetX();
 206       Float_t yc = c->GetY();
 207       Float_t zc = c->GetZ();
 208       Float_t dx = x0 - xc;
 209       Float_t yt = y0 - dx*dydx;
 210       Float_t zt = z0 - dx*dzdx;
 211       yc -= tilt*(zc-zt); // tilt correction
 212       dy = yt - yc;
 213
 214       h->Fill(dydx, dy/TMath::Sqrt(cov[0] + c->GetSigmaY2()));
 215
 216       if(fDebugLevel>=1){
 217         // Get z-position with respect to anode wire
 218         AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
 219         Int_t istk = fGeo->GetStack(c->GetDetector());
 220         AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
 221         Float_t row0 = pp->GetRow0();
 222         Float_t d  =  row0 - zt + simParam->GetAnodeWireOffset();
 223         d -= ((Int_t)(2 * d)) / 2.0;
 224         if (d > 0.25) d  = 0.5 - d;
 225
 226 /*        AliTRDclusterInfo *clInfo = new AliTRDclusterInfo;
 227         fClResiduals->Add(clInfo);
 228         clInfo->SetCluster(c);
 229         clInfo->SetGlobalPosition(yt, zt, dydx, dzdx);
 230         clInfo->SetResolution(dy);
 231         clInfo->SetAnisochronity(d);
 232         clInfo->SetDriftLength(dx);
 233         (*fDebugStream) << "ClusterResiduals"
 234           <<"clInfo.=" << clInfo
 235           << "\n";*/
 236       }
 237     }
 238   }
 239   return h;
 240 }
 241
 242
 243 //________________________________________________________
 244 TH1* AliTRDtrackingResolution::PlotTracklet(const AliTRDtrackV1 *track)
 245 {
 246   if(track) fTrack = track;
 247   if(!fTrack){
 248     AliWarning("No Track defined.");
 249     return 0x0;
 250   }
 251   TH1 *h = 0x0;
 252   if(!(h = ((TH2I*)fContainer->At(kTrackletY)))){
 253     AliWarning("No output histogram defined.");
 254     return 0x0;
 255   }
 256
 257   Double_t cov[3], covR[3];
 258   Float_t xref, y0, dx, dy, dydx;
 259   AliTRDseedV1 *fTracklet = 0x0;
 260   for(Int_t il=AliTRDgeometry::kNlayer; il--;){
 261     if(!(fTracklet = fTrack->GetTracklet(il))) continue;
 262     if(!fTracklet->IsOK()) continue;
 263     y0   = fTracklet->GetYref(0);
 264     dydx = fTracklet->GetYref(1);
 265     xref = fTracklet->GetXref();
 266     dx   = fTracklet->GetX0() - xref;
 267     dy   = y0-dx*dydx - fTracklet->GetYat(xref);
 268     fTracklet->GetCovAt(xref, cov);
 269     fTracklet->GetCovRef(covR);
 270     h->Fill(dydx, dy/TMath::Sqrt(cov[0] + covR[0]));
 271   }
 272   return h;
 273 }
 274
 275 //________________________________________________________
 276 TH1* AliTRDtrackingResolution::PlotTrackletPhi(const AliTRDtrackV1 *track)
 277 {
 278   if(track) fTrack = track;
 279   if(!fTrack){
 280     AliWarning("No Track defined.");
 281     return 0x0;
 282   }
 283   TH1 *h = 0x0;
 284   if(!(h = ((TH2I*)fContainer->At(kTrackletPhi)))){
 285     AliWarning("No output histogram defined.");
 286     return 0x0;
 287   }
 288
 289   AliTRDseedV1 *tracklet = 0x0;
 290   for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
 291     if(!(tracklet = fTrack->GetTracklet(il))) continue;
 292     if(!tracklet->IsOK()) continue;
 293     h->Fill(tracklet->GetYref(1), TMath::ATan(tracklet->GetYref(1))-TMath::ATan(tracklet->GetYfit(1)));
 294   }
 295   return h;
 296 }
 297
 298
 299 //________________________________________________________
 300 TH1* AliTRDtrackingResolution::PlotResolution(const AliTRDtrackV1 *track)
 301 {
 302   if(!HasMCdata()){
 303     AliWarning("No MC defined. Results will not be available.");
 304     return 0x0;
 305   }
 306   if(track) fTrack = track;
 307   if(!fTrack){
 308     AliWarning("No Track defined.");
 309     return 0x0;
 310   }
 311   TH1 *h = 0x0;
 312   UChar_t s;
 313   Int_t pdg = fMC->GetPDG(), det=-1;
 314   Int_t label = fMC->GetLabel();
 315   Float_t p, pt, x0, y0, z0, dx, dy, dz, dydx, dzdx;
 316
 317   if(fDebugLevel>=1){
 318     Double_t DX[12], DY[12], DZ[12], DPt[12], COV[12][15];
 319     fMC->PropagateKalman(DX, DY, DZ, DPt, COV);
 320     (*fDebugStream) << "MCkalman"
 321       << "pdg="  << pdg
 322       << "dx0="  << DX[0]
 323       << "dx1="  << DX[1]
 324       << "dx2="  << DX[2]
 325       << "dy0="  << DY[0]
 326       << "dy1="  << DY[1]
 327       << "dy2="  << DY[2]
 328       << "dz0="  << DZ[0]
 329       << "dz1="  << DZ[1]
 330       << "dz2="  << DZ[2]
 331       << "dpt0=" << DPt[0]
 332       << "dpt1=" << DPt[1]
 333       << "dpt2=" << DPt[2]
 334       << "\n";
 335   }
 336
 337   AliTRDseedV1 *fTracklet = 0x0;
 338   for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
 339     if(!(fTracklet = fTrack->GetTracklet(ily)) ||
 340        !fTracklet->IsOK()) continue;
 341
 342     det = fTracklet->GetDetector();
 343     x0  = fTracklet->GetX0();
 344     //radial shift with respect to the MC reference (radial position of the pad plane)
 345     dx  = x0 - fTracklet->GetXref();
 346     if(!fMC->GetDirections(x0, y0, z0, dydx, dzdx, pt, s)) continue;
 347     // MC track position at reference radial position
 348     Float_t yt = y0 - dx*dydx;
 349     Float_t zt = z0 - dx*dzdx;
 350     p = pt*(1.+dzdx*dzdx); // pt -> p
 351
 352     // add Kalman residuals for y, z and pt
 353     Float_t dxr= fTracklet->GetX0() - x0 + dx;
 354     Float_t yr = fTracklet->GetYref(0) - dxr*fTracklet->GetYref(1);
 355     dy = yt - yr;
 356     Float_t zr = fTracklet->GetZref(0) - dxr*fTracklet->GetZref(1);
 357     dz = zt - zr;
 358     Float_t tgl = fTracklet->GetTgl();
 359     Float_t ptr = fTracklet->GetMomentum()/(1.+tgl*tgl);
 360
 361     ((TH2I*)fContainer->At(kMCtrackY))->Fill(dydx, dy);
 362     ((TH2I*)fContainer->At(kMCtrackZ))->Fill(dzdx, dz);
 363     if(pdg!=kElectron && pdg!=kPositron) ((TH2I*)fContainer->At(kMCtrackPt))->Fill(1./pt, ptr-pt);
 364     // Fill Debug stream for Kalman track
 365     if(fDebugLevel>=1){
 366       Float_t dydxr = fTracklet->GetYref(1);
 367       (*fDebugStream) << "MCtrack"
 368         << "det="     << det
 369         << "pdg="     << pdg
 370         << "pt="      << pt
 371         << "yt="      << yt
 372         << "zt="      << zt
 373         << "dydx="    << dydx
 374         << "dzdx="    << dzdx
 375         << "ptr="     << ptr
 376         << "dy="      << dy
 377         << "dz="      << dz
 378         << "dydxr="   << dydxr
 379         << "dzdxr="   << tgl
 380         << "\n";
 381     }
 382
 383     // recalculate tracklet based on the MC info
 384     AliTRDseedV1 tt(*fTracklet);
 385     tt.SetZref(0, z0);
 386     tt.SetZref(1, dzdx);
 387     if(!tt.Fit(kTRUE)) continue;
 388
 389     // add tracklet residuals for y and dydx
 390     Float_t yf = tt.GetYfit(0) - dxr*tt.GetYfit(1);
 391     dy = yt - yf;
 392     Float_t dphi   = (tt.GetYfit(1) - dydx);
 393     dphi /= 1.- tt.GetYfit(1)*dydx;
 394     ((TH2I*)fContainer->At(kMCtrackletY))->Fill(dydx, dy);
 395     ((TH2I*)fContainer->At(kMCtrackletPhi))->Fill(dydx, dphi*TMath::RadToDeg());
 396
 397     Float_t dz = 100.;
 398     Bool_t rc = fTracklet->IsRowCross();
 399     if(rc){
 400       // add tracklet residuals for z
 401       Double_t *xyz = tt.GetCrossXYZ();
 402       dz = xyz[2] - (z0 - (x0 - xyz[0])*dzdx) ;
 403       ((TH2I*)fContainer->At(kMCtrackletZ))->Fill(dzdx, dz);
 404     }
 405
 406     // Fill Debug stream for tracklet
 407     if(fDebugLevel>=1){
 408       (*fDebugStream) << "MCtracklet"
 409         << "det="     << det
 410         << "pdg="     << pdg
 411         << "p="       << p
 412         << "yt="      << yt
 413         << "zt="      << zt
 414         << "dydx="    << dydx
 415         << "dzdx="    << dzdx
 416         << "rowc="    << rc
 417         << "dy="      << dy
 418         << "dz="      << dz
 419         << "dphi="    << dphi
 420         << "\n";
 421     }
 422
 423     Int_t istk = AliTRDgeometry::GetStack(det);
 424     AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
 425     Float_t zr0 = pp->GetRow0() + AliTRDSimParam::Instance()->GetAnodeWireOffset();
 426     Float_t tilt = fTracklet->GetTilt();
 427
 428     Double_t x,y;
 429     AliTRDcluster *c = 0x0;
 430     fTracklet->ResetClusterIter(kFALSE);
 431     while((c = fTracklet->PrevCluster())){
 432       Float_t  q = TMath::Abs(c->GetQ());
 433       //AliTRDseedV1::GetClusterXY(c,x,y);
 434       x = c->GetX(); y = c->GetY();
 435       Float_t xc = x;
 436       Float_t yc = y;
 437       Float_t zc = c->GetZ();
 438       dx = x0 - xc;
 439       Float_t yt = y0 - dx*dydx;
 440       Float_t zt = z0 - dx*dzdx;
 441       dy = yt - (yc - tilt*(zc-zt));
 442
 443       // Fill Histograms
 444       if(q>20. && q<250.) ((TH2I*)fContainer->At(kMCcluster))->Fill(dydx, dy);
 445
 446       // Fill calibration container
 447       Float_t d = zr0 - zt;
 448       d -= ((Int_t)(2 * d)) / 2.0;
 449       if (d > 0.25) d  = 0.5 - d;
 450       AliTRDclusterInfo *clInfo = new AliTRDclusterInfo;
 451       fClResolution->Add(clInfo);
 452       clInfo->SetCluster(c);
 453       clInfo->SetMC(pdg, label);
 454       clInfo->SetGlobalPosition(yt, zt, dydx, dzdx);
 455       clInfo->SetResolution(dy);
 456       clInfo->SetAnisochronity(d);
 457       clInfo->SetDriftLength(dx);
 458       clInfo->SetTilt(tilt);
 459
 460       // Fill Debug Tree
 461       if(fDebugLevel>=2){
 462         //clInfo->Print();
 463         (*fDebugStream) << "MCcluster"
 464           <<"clInfo.=" << clInfo
 465           << "\n";
 466       }
 467     }
 468   }
 469   return h;
 470 }
 471
 472
 473 //________________________________________________________
 474 Bool_t AliTRDtrackingResolution::GetRefFigure(Int_t ifig)
 475 {
 476   Float_t y[2] = {0., 0.};
 477   TBox *b = 0x0;
 478   TAxis *ax = 0x0;
 479   TGraphErrors *g = 0x0;
 480   switch(ifig){
 481   case kCluster:
 482     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 483     g->Draw("apl");
 484     ax = g->GetHistogram()->GetYaxis();
 485     y[0] = -0.5; y[1] = 2.5;
 486     ax->SetRangeUser(y[0], y[1]);
 487     ax->SetTitle("Cluster-Track Pools #sigma/#mu [mm]");
 488     ax = g->GetHistogram()->GetXaxis();
 489     ax->SetTitle("tg(#phi)");
 490     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 491     g->Draw("pl");
 492     b = new TBox(-.15, y[0], .15, y[1]);
 493     b->SetFillStyle(3002);b->SetFillColor(kGreen);
 494     b->SetLineColor(0); b->Draw();
 495     return kTRUE;
 496   case kTrackletY:
 497     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 498     g->Draw("apl");
 499     ax = g->GetHistogram()->GetYaxis();
 500     ax->SetRangeUser(-.5, 3.);
 501     ax->SetTitle("Tracklet-Track Y-Pools #sigma/#mu [mm]");
 502     ax = g->GetHistogram()->GetXaxis();
 503     ax->SetTitle("tg(#phi)");
 504     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 505     g->Draw("pl");
 506     b = new TBox(-.15, y[0], .15, y[1]);
 507     b->SetFillStyle(3002);b->SetFillColor(kGreen);
 508     b->SetLineColor(0); b->Draw();
 509     return kTRUE;
 510   case kTrackletPhi:
 511     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 512     g->Draw("apl");
 513     ax = g->GetHistogram()->GetYaxis();
 514     ax->SetRangeUser(-.5, 2.);
 515     ax->SetTitle("Tracklet-Track Phi-Residuals #sigma/#mu [rad]");
 516     ax = g->GetHistogram()->GetXaxis();
 517     ax->SetTitle("tg(#phi)");
 518     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 519     g->Draw("pl");
 520     b = new TBox(-.15, y[0], .15, y[1]);
 521     b->SetFillStyle(3002);b->SetFillColor(kGreen);
 522     b->SetLineColor(0); b->Draw();
 523     return kTRUE;
 524   case kMCcluster:
 525     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 526     ax = g->GetHistogram()->GetYaxis();
 527     y[0] = -.1; y[1] = 1.5;
 528     ax->SetRangeUser(y[0], y[1]);
 529     ax->SetTitle("Y_{cluster} #sigma/#mu [mm]");
 530     ax = g->GetHistogram()->GetXaxis();
 531     ax->SetTitle("tg(#phi)");
 532     g->Draw("apl");
 533     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 534     g->Draw("pl");
 535     b = new TBox(-.15, y[0], .15, y[1]);
 536     b->SetFillStyle(3002);b->SetFillColor(kBlue);
 537     b->SetLineColor(0); b->Draw();
 538     return kTRUE;
 539   case kMCtrackletY:
 540     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 541     ax = g->GetHistogram()->GetYaxis();
 542     y[0] = -.05; y[1] = 0.3;
 543     ax->SetRangeUser(y[0], y[1]);
 544     ax->SetTitle("Y_{tracklet} #sigma/#mu [mm]");
 545     ax = g->GetHistogram()->GetXaxis();
 546     ax->SetTitle("tg(#phi)");
 547     g->Draw("apl");
 548     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 549     g->Draw("pl");
 550     b = new TBox(-.15, y[0], .15, y[1]);
 551     b->SetFillStyle(3002);b->SetFillColor(kBlue);
 552     b->SetLineColor(0); b->Draw();
 553     return kTRUE;
 554   case kMCtrackletZ:
 555     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 556     ax = g->GetHistogram()->GetYaxis();
 557     ax->SetRangeUser(-.5, 1.);
 558     ax->SetTitle("Z_{tracklet} #sigma/#mu [mm]");
 559     ax = g->GetHistogram()->GetXaxis();
 560     ax->SetTitle("tg(#theta)");
 561     g->Draw("apl");
 562     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 563     g->Draw("pl");
 564     return kTRUE;
 565   case kMCtrackletPhi:
 566     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 567     ax = g->GetHistogram()->GetYaxis();
 568     y[0] = -.05; y[1] = .2;
 569     ax->SetRangeUser(y[0], y[1]);
 570     ax->SetTitle("#Phi_{tracklet} #sigma/#mu [deg]");
 571     ax = g->GetHistogram()->GetXaxis();
 572     ax->SetTitle("tg(#phi)");
 573     g->Draw("apl");
 574     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 575     g->Draw("pl");
 576     return kTRUE;
 577   case kMCtrackY:
 578     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 579     ax = g->GetHistogram()->GetYaxis();
 580     y[0] = -.05; y[1] = 0.2;
 581     ax->SetRangeUser(y[0], y[1]);
 582     ax->SetTitle("Y_{track} #sigma/#mu [mm]");
 583     ax = g->GetHistogram()->GetXaxis();
 584     ax->SetTitle("tg(#phi)");
 585     g->Draw("apl");
 586     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 587     g->Draw("pl");
 588     b = new TBox(-.15, y[0], .15, y[1]);
 589     b->SetFillStyle(3002);b->SetFillColor(kBlue);
 590     b->SetLineColor(0); b->Draw();
 591     return kTRUE;
 592   case kMCtrackZ:
 593     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 594     ax = g->GetHistogram()->GetYaxis();
 595     ax->SetRangeUser(-.5, 2.);
 596     ax->SetTitle("Z_{track} #sigma/#mu [mm]");
 597     ax = g->GetHistogram()->GetXaxis();
 598     ax->SetTitle("tg(#theta)");
 599     g->Draw("apl");
 600     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 601     g->Draw("pl");
 602     return kTRUE;
 603   case kMCtrackPt:
 604     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
 605     ax = g->GetHistogram()->GetYaxis();
 606     ax->SetRangeUser(-.5, 2.);
 607     ax->SetTitle("#epsilon_{P_{t}}^{track} / #mu [%]");
 608     ax = g->GetHistogram()->GetXaxis();
 609     ax->SetTitle("1/p_{t}");
 610     g->Draw("apl");
 611     if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
 612     g->Draw("pl");
 613     return kTRUE;
 614   }
 615   AliInfo(Form("Reference plot [%d] missing result", ifig));
 616   return kFALSE;
 617 }
 618
 619
 620 //________________________________________________________
 621 Bool_t AliTRDtrackingResolution::PostProcess()
 622 {
 623   //fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
 624   if (!fContainer) {
 625     Printf("ERROR: list not available");
 626     return kFALSE;
 627   }
 628   fNRefFigures = fContainer->GetEntriesFast();
 629   TGraphErrors *gm = 0x0, *gs = 0x0;
 630   if(!fGraphS){
 631     fGraphS = new TObjArray(fNRefFigures);
 632     fGraphS->SetOwner();
 633     for(Int_t ig=0; ig<fNRefFigures; ig++){
 634       gs = new TGraphErrors();
 635       gs->SetLineColor(kRed);
 636       gs->SetMarkerStyle(23);
 637       gs->SetMarkerColor(kRed);
 638       gs->SetNameTitle(Form("s_%d", ig), "");
 639       fGraphS->AddAt(gs, ig);
 640     }
 641   }
 642   if(!fGraphM){
 643     fGraphM = new TObjArray(fNRefFigures);
 644     fGraphM->SetOwner();
 645     for(Int_t ig=0; ig<fNRefFigures; ig++){
 646       gm = new TGraphErrors();
 647       gm->SetLineColor(kBlue);
 648       gm->SetMarkerStyle(7);
 649       gm->SetMarkerColor(kBlue);
 650       gm->SetNameTitle(Form("m_%d", ig), "");
 651       fGraphM->AddAt(gm, ig);
 652     }
 653   }
 654
 655   TH2I *h2 = 0x0;
 656   TH1D *h = 0x0;
 657
 658   // define models
 659   TF1 f("f1", "gaus", -.5, .5);
 660
 661   TF1 fb("fb", "[0]*exp(-0.5*((x-[1])/[2])**2)+[3]", -.5, .5);
 662
 663   TF1 fc("fc", "[0]*exp(-0.5*((x-[1])/[2])**2)+[3]*exp(-0.5*((x-[4])/[5])**2)", -.5, .5);
 664
 665   TCanvas *c = 0x0;
 666   if(IsVisual()) c = new TCanvas("c", Form("%s Visual", GetName()), 500, 500);
 667   char opt[5];
 668   sprintf(opt, "%c%c", IsVerbose() ? ' ' : 'Q', IsVisual() ? ' ': 'N');
 669
 670
 671   //PROCESS RESIDUAL DISTRIBUTIONS
 672
 673   // Clusters residuals
 674   h2 = (TH2I *)(fContainer->At(kCluster));
 675   gm = (TGraphErrors*)fGraphM->At(kCluster);
 676   gs = (TGraphErrors*)fGraphS->At(kCluster);
 677   for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
 678     Double_t phi = h2->GetXaxis()->GetBinCenter(ibin);
 679     h = h2->ProjectionY("py", ibin, ibin);
 680     if(h->GetEntries()<100) continue;
 681     AdjustF1(h, &f);
 682
 683     if(IsVisual()){c->cd(); c->SetLogy();}
 684     h->Fit(&f, opt, "", -0.5, 0.5);
 685     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 686
 687     Int_t ip = gm->GetN();
 688     gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
 689     gm->SetPointError(ip, 0., 10.*f.GetParError(1));
 690     gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
 691     gs->SetPointError(ip, 0., 10.*f.GetParError(2));
 692   }
 693
 694   // Tracklet y residuals
 695   h2 = (TH2I *)(fContainer->At(kTrackletY));
 696   gm = (TGraphErrors*)fGraphM->At(kTrackletY);
 697   gs = (TGraphErrors*)fGraphS->At(kTrackletY);
 698   for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
 699     Double_t phi = h2->GetXaxis()->GetBinCenter(ibin);
 700     h = h2->ProjectionY("py", ibin, ibin);
 701     if(h->GetEntries()<100) continue;
 702     AdjustF1(h, &f);
 703
 704     if(IsVisual()){c->cd(); c->SetLogy();}
 705     h->Fit(&f, opt, "", -0.5, 0.5);
 706     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 707
 708     Int_t ip = gm->GetN();
 709     gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
 710     gm->SetPointError(ip, 0., 10.*f.GetParError(1));
 711     gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
 712     gs->SetPointError(ip, 0., 10.*f.GetParError(2));
 713   }
 714
 715   // Tracklet phi residuals
 716   h2 = (TH2I *)(fContainer->At(kTrackletPhi));
 717   gm = (TGraphErrors*)fGraphM->At(kTrackletPhi);
 718   gs = (TGraphErrors*)fGraphS->At(kTrackletPhi);
 719   for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
 720     Double_t phi = h2->GetXaxis()->GetBinCenter(ibin);
 721     h = h2->ProjectionY("py", ibin, ibin);
 722     if(h->GetEntries()<100) continue;
 723     AdjustF1(h, &f);
 724
 725     if(IsVisual()){c->cd(); c->SetLogy();}
 726     h->Fit(&f, opt, "", -0.5, 0.5);
 727     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 728
 729     Int_t ip = gm->GetN();
 730     gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
 731     gm->SetPointError(ip, 0., 10.*f.GetParError(1));
 732     gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
 733     gs->SetPointError(ip, 0., 10.*f.GetParError(2));
 734   }
 735
 736   if(!HasMCdata()){
 737     if(c) delete c;
 738     return kTRUE;
 739   }
 740
 741   //PROCESS MC RESIDUAL DISTRIBUTIONS
 742
 743   // cluster y resolution
 744   h2 = (TH2I*)fContainer->At(kMCcluster);
 745   gm = (TGraphErrors*)fGraphM->At(kMCcluster);
 746   gs = (TGraphErrors*)fGraphS->At(kMCcluster);
 747   for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
 748     h = h2->ProjectionY("py", iphi, iphi);
 749     if(h->GetEntries()<100) continue;
 750     AdjustF1(h, &f);
 751
 752     if(IsVisual()){c->cd(); c->SetLogy();}
 753     h->Fit(&f, opt, "", -0.5, 0.5);
 754     if(IsVerbose()){
 755       printf("phi[%d] mean[%e] sigma[%e]\n\n", iphi, 10.*f.GetParameter(1), 10.*f.GetParameter(2));
 756     }
 757     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 758
 759     Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
 760     Int_t ip = gm->GetN();
 761     gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
 762     gm->SetPointError(ip, 0., 10.*f.GetParError(1));
 763     gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
 764     gs->SetPointError(ip, 0., 10.*f.GetParError(2));
 765   }
 766
 767   // tracklet y resolution
 768   h2 = (TH2I*)fContainer->At(kMCtrackletY);
 769   gm = (TGraphErrors*)fGraphM->At(kMCtrackletY);
 770   gs = (TGraphErrors*)fGraphS->At(kMCtrackletY);
 771   for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
 772     h = h2->ProjectionY("py", iphi, iphi);
 773     if(h->GetEntries()<100) continue;
 774     AdjustF1(h, &f);
 775
 776     if(IsVisual()){c->cd(); c->SetLogy();}
 777     h->Fit(&f, opt, "", -0.5, 0.5);
 778     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 779
 780     Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
 781     Int_t ip = gm->GetN();
 782     gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
 783     gm->SetPointError(ip, 0., 10.*f.GetParError(1));
 784     gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
 785     gs->SetPointError(ip, 0., 10.*f.GetParError(2));
 786   }
 787
 788   // tracklet z resolution
 789   h2 = (TH2I*)fContainer->At(kMCtrackletZ);
 790   gm = (TGraphErrors*)fGraphM->At(kMCtrackletZ);
 791   gs = (TGraphErrors*)fGraphS->At(kMCtrackletZ);
 792   for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
 793     h = h2->ProjectionY("py", iphi, iphi);
 794     if(h->GetEntries()<100) continue;
 795     AdjustF1(h, &fb);
 796
 797     if(IsVisual()){c->cd(); c->SetLogy();}
 798     h->Fit(&fb, opt, "", -0.5, 0.5);
 799     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 800
 801     Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
 802     Int_t ip = gm->GetN();
 803     gm->SetPoint(ip, phi, 10.*fb.GetParameter(1));
 804     gm->SetPointError(ip, 0., 10.*fb.GetParError(1));
 805     gs->SetPoint(ip, phi, 10.*fb.GetParameter(2));
 806     gs->SetPointError(ip, 0., 10.*fb.GetParError(2));
 807   }
 808
 809   //tracklet phi resolution
 810   h2 = (TH2I*)fContainer->At(kMCtrackletPhi);
 811   gm = (TGraphErrors*)fGraphM->At(kMCtrackletPhi);
 812   gs = (TGraphErrors*)fGraphS->At(kMCtrackletPhi);
 813   for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
 814     h = h2->ProjectionY("py", iphi, iphi);
 815     if(h->GetEntries()<100) continue;
 816
 817     if(IsVisual()){c->cd(); c->SetLogy();}
 818     h->Fit(&f, opt, "", -0.5, 0.5);
 819     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 820
 821     Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
 822     Int_t ip = gm->GetN();
 823     gm->SetPoint(ip, phi, f.GetParameter(1));
 824     gm->SetPointError(ip, 0., f.GetParError(1));
 825     gs->SetPoint(ip, phi, f.GetParameter(2));
 826     gs->SetPointError(ip, 0., f.GetParError(2));
 827   }
 828
 829   // track y resolution
 830   h2 = (TH2I*)fContainer->At(kMCtrackY);
 831   gm = (TGraphErrors*)fGraphM->At(kMCtrackY);
 832   gs = (TGraphErrors*)fGraphS->At(kMCtrackY);
 833   for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
 834     h = h2->ProjectionY("py", iphi, iphi);
 835     if(h->GetEntries()<100) continue;
 836     AdjustF1(h, &f);
 837
 838     if(IsVisual()){c->cd(); c->SetLogy();}
 839     h->Fit(&f, opt, "", -0.5, 0.5);
 840     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 841
 842     Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
 843     Int_t ip = gm->GetN();
 844     gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
 845     gm->SetPointError(ip, 0., 10.*f.GetParError(1));
 846     gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
 847     gs->SetPointError(ip, 0., 10.*f.GetParError(2));
 848   }
 849
 850   // track z resolution
 851   h2 = (TH2I*)fContainer->At(kMCtrackZ);
 852   gm = (TGraphErrors*)fGraphM->At(kMCtrackZ);
 853   gs = (TGraphErrors*)fGraphS->At(kMCtrackZ);
 854   for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
 855     h = h2->ProjectionY("pz", iphi, iphi);
 856     if(h->GetEntries()<70) continue;
 857     AdjustF1(h, &f);
 858
 859     if(IsVisual()){c->cd(); c->SetLogy();}
 860     h->Fit(&f, opt, "", -0.5, 0.5);
 861     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 862
 863     Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
 864     Int_t ip = gm->GetN();
 865     gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
 866     gm->SetPointError(ip, 0., 10.*f.GetParError(1));
 867     gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
 868     gs->SetPointError(ip, 0., 10.*f.GetParError(2));
 869   }
 870
 871   // track Pt resolution
 872   h2 = (TH2I*)fContainer->At(kMCtrackPt);
 873   TAxis *ax = h2->GetXaxis();
 874   gm = (TGraphErrors*)fGraphM->At(kMCtrackPt);
 875   gs = (TGraphErrors*)fGraphS->At(kMCtrackPt);
 876   TF1 fg("fg", "gaus", -1.5, 1.5);
 877   TF1 fl("fl", "landau", -4., 15.);
 878   TF1 fgl("fgl", "gaus(0)+landau(3)", -5., 20.);
 879   for(Int_t ip=1; ip<=ax->GetNbins(); ip++){
 880     h = h2->ProjectionY("ppt", ip, ip);
 881     if(h->GetEntries()<70) continue;
 882
 883     h->Fit(&fg, "QN", "", -1.5, 1.5);
 884     fgl.SetParameter(0, fg.GetParameter(0));
 885     fgl.SetParameter(1, fg.GetParameter(1));
 886     fgl.SetParameter(2, fg.GetParameter(2));
 887     h->Fit(&fl, "QN", "", -4., 15.);
 888     fgl.SetParameter(3, fl.GetParameter(0));
 889     fgl.SetParameter(4, fl.GetParameter(1));
 890     fgl.SetParameter(5, fl.GetParameter(2));
 891
 892     if(IsVisual()){c->cd(); c->SetLogy();}
 893     h->Fit(&fgl, opt, "", -5., 20.);
 894     if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
 895
 896     Float_t invpt = ax->GetBinCenter(ip);
 897     Int_t ip = gm->GetN();
 898     gm->SetPoint(ip, invpt, fgl.GetParameter(1));
 899     gm->SetPointError(ip, 0., fgl.GetParError(1));
 900     gs->SetPoint(ip, invpt, fgl.GetParameter(2)*invpt);
 901     gs->SetPointError(ip, 0., fgl.GetParError(2));
 902     // fgl.GetParameter(4) // Landau MPV
 903     // fgl.GetParameter(5) // Landau Sigma
 904   }
 905
 906
 907   if(c) delete c;
 908   return kTRUE;
 909 }
 910
 911
 912 //________________________________________________________
 913 void AliTRDtrackingResolution::Terminate(Option_t *)
 914 {
 915   if(fDebugStream){
 916     delete fDebugStream;
 917     fDebugStream = 0x0;
 918     fDebugLevel = 0;
 919   }
 920   if(HasPostProcess()) PostProcess();
 921 }
 922
 923 //________________________________________________________
 924 void AliTRDtrackingResolution::AdjustF1(TH1 *h, TF1 *f)
 925 {
 926 // Helper function to avoid duplication of code
 927 // Make first guesses on the fit parameters
 928
 929   // find the intial parameters of the fit !! (thanks George)
 930   Int_t nbinsy = Int_t(.5*h->GetNbinsX());
 931   Double_t sum = 0.;
 932   for(Int_t jbin=nbinsy-4; jbin<=nbinsy+4; jbin++) sum+=h->GetBinContent(jbin); sum/=9.;
 933   f->SetParLimits(0, 0., 3.*sum);
 934   f->SetParameter(0, .9*sum);
 935
 936   f->SetParLimits(1, -.2, .2);
 937   f->SetParameter(1, -0.1);
 938
 939   f->SetParLimits(2, 0., 4.e-1);
 940   f->SetParameter(2, 2.e-2);
 941   if(f->GetNpar() <= 4) return;
 942
 943   f->SetParLimits(3, 0., sum);
 944   f->SetParameter(3, .1*sum);
 945
 946   f->SetParLimits(4, -.3, .3);
 947   f->SetParameter(4, 0.);
 948
 949   f->SetParLimits(5, 0., 1.e2);
 950   f->SetParameter(5, 2.e-1);
 951 }
 952
 953 //________________________________________________________
 954 TObjArray* AliTRDtrackingResolution::Histos()
 955 {
 956   if(fContainer) return fContainer;
 957
 958   fContainer  = new TObjArray(10);
 959   //fContainer->SetOwner(kTRUE);
 960
 961   TH1 *h = 0x0;
 962   // cluster to tracklet residuals [2]
 963   if(!(h = (TH2I*)gROOT->FindObject("hCl"))){
 964     h = new TH2I("hCl", "Clusters-Track Residuals", 21, -.33, .33, 100, -.5, .5);
 965     h->GetXaxis()->SetTitle("tg(#phi)");
 966     h->GetYaxis()->SetTitle("#Delta y [cm]");
 967     h->GetZaxis()->SetTitle("entries");
 968   } else h->Reset();
 969   fContainer->AddAt(h, kCluster);
 970
 971   // tracklet to track residuals [2]
 972   if(!(h = (TH2I*)gROOT->FindObject("hTrkltY"))){
 973     h = new TH2I("hTrkltY", "Tracklets-Track Residuals (Y)", 21, -.33, .33, 100, -.5, .5);
 974     h->GetXaxis()->SetTitle("tg(#phi)");
 975     h->GetYaxis()->SetTitle("#Delta y [cm]");
 976     h->GetZaxis()->SetTitle("entries");
 977   } else h->Reset();
 978   fContainer->AddAt(h, kTrackletY);
 979
 980   // tracklet to track residuals angular [2]
 981   if(!(h = (TH2I*)gROOT->FindObject("hTrkltPhi"))){
 982     h = new TH2I("hTrkltPhi", "Tracklets-Track Residuals (#Phi)", 21, -.33, .33, 100, -.5, .5);
 983     h->GetXaxis()->SetTitle("tg(#phi)");
 984     h->GetYaxis()->SetTitle("#Delta phi [#circ]");
 985     h->GetZaxis()->SetTitle("entries");
 986   } else h->Reset();
 987   fContainer->AddAt(h, kTrackletPhi);
 988
 989
 990   // Resolution histos
 991   if(!HasMCdata()) return fContainer;
 992
 993   // cluster y resolution [0]
 994   if(!(h = (TH2I*)gROOT->FindObject("hMCcl"))){
 995     h = new TH2I("hMCcl", "Cluster Resolution", 31, -.48, .48, 100, -.5, .5);
 996     h->GetXaxis()->SetTitle("tg(#phi)");
 997     h->GetYaxis()->SetTitle("#Delta y [cm]");
 998     h->GetZaxis()->SetTitle("entries");
 999   } else h->Reset();
1000   fContainer->AddAt(h, kMCcluster);
1001
1002   // tracklet y resolution [0]
1003   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltY"))){
1004     h = new TH2I("hMCtrkltY", "Tracklet Resolution (Y)", 31, -.48, .48, 100, -.5, .5);
1005     h->GetXaxis()->SetTitle("tg(#phi)");
1006     h->GetYaxis()->SetTitle("#Delta y [cm]");
1007     h->GetZaxis()->SetTitle("entries");
1008   } else h->Reset();
1009   fContainer->AddAt(h, kMCtrackletY);
1010
1011   // tracklet y resolution [0]
1012   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltZ"))){
1013     h = new TH2I("hMCtrkltZ", "Tracklet Resolution (Z)", 31, -.48, .48, 100, -.5, .5);
1014     h->GetXaxis()->SetTitle("tg(#theta)");
1015     h->GetYaxis()->SetTitle("#Delta z [cm]");
1016     h->GetZaxis()->SetTitle("entries");
1017   } else h->Reset();
1018   fContainer->AddAt(h, kMCtrackletZ);
1019
1020   // tracklet angular resolution [1]
1021   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltPhi"))){
1022     h = new TH2I("hMCtrkltPhi", "Tracklet Resolution (#Phi)", 31, -.48, .48, 100, -10., 10.);
1023     h->GetXaxis()->SetTitle("tg(#phi)");
1024     h->GetYaxis()->SetTitle("#Delta #phi [deg]");
1025     h->GetZaxis()->SetTitle("entries");
1026   } else h->Reset();
1027   fContainer->AddAt(h, kMCtrackletPhi);
1028
1029   // Kalman track y resolution
1030   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkY"))){
1031     h = new TH2I("hMCtrkY", "Kalman Track Resolution (Y)", 31, -.48, .48, 100, -.5, .5);
1032     h->GetXaxis()->SetTitle("tg(#phi)");
1033     h->GetYaxis()->SetTitle("#Delta y [cm]");
1034     h->GetZaxis()->SetTitle("entries");
1035   } else h->Reset();
1036   fContainer->AddAt(h, kMCtrackY);
1037
1038   // Kalman track Z resolution
1039   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZ"))){
1040     h = new TH2I("hMCtrkZ", "Kalman Track Resolution (Z)", 20, -1., 1., 100, -1.5, 1.5);
1041     h->GetXaxis()->SetTitle("tg(#theta)");
1042     h->GetYaxis()->SetTitle("#Delta z [cm]");
1043     h->GetZaxis()->SetTitle("entries");
1044   } else h->Reset();
1045   fContainer->AddAt(h, kMCtrackZ);
1046
1047   // Kalman track Pt resolution
1048   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkPt"))){
1049     h = new TH2I("hMCtrkPt", "Kalman Track Resolution (Pt)", 100, 0., 2., 150, -5., 20.);
1050     h->GetXaxis()->SetTitle("1/p_{t}");
1051     h->GetYaxis()->SetTitle("#Delta p_{t} [GeV/c]");
1052     h->GetZaxis()->SetTitle("entries");
1053   } else h->Reset();
1054   fContainer->AddAt(h, kMCtrackPt);
1055
1056   return fContainer;
1057 }
1058
1059
1060 //________________________________________________________
1061 void AliTRDtrackingResolution::SetRecoParam(AliTRDrecoParam *r)
1062 {
1063
1064   fReconstructor->SetRecoParam(r);
1065 }