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