PWG1/AliPerformanceRes.cxx

   1 //------------------------------------------------------------------------------
   2 // Implementation of AliPerformanceRes class. It keeps information from
   3 // comparison of reconstructed and MC particle tracks. In addtion,
   4 // it keeps selection cuts used during comparison. The comparison
   5 // information is stored in the ROOT histograms. Analysis of these
   6 // histograms can be done by using Analyse() class function. The result of
   7 // the analysis (histograms/graphs) are stored in the folder which is
   8 // a data member of AliPerformanceRes.
   9 //
  10 // Author: J.Otwinowski 04/02/2008
  11 //------------------------------------------------------------------------------
  12
  13 /*
  14
  15   // after running comparison task, read the file, and get component
  16   gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
  17   LoadMyLibs();
  18
  19   TFile f("Output.root");
  20   AliPerformanceRes * compObj = (AliPerformanceRes*)coutput->FindObject("AliPerformanceRes");
  21
  22   // analyse comparison data
  23   compObj->Analyse();
  24
  25   // the output histograms/graphs will be stored in the folder "folderRes"
  26   compObj->GetAnalysisFolder()->ls("*");
  27
  28   // user can save whole comparison object (or only folder with anlysed histograms)
  29   // in the seperate output file (e.g.)
  30   TFile fout("Analysed_Res.root","recreate");
  31   compObj->Write(); // compObj->GetAnalysisFolder()->Write();
  32   fout.Close();
  33
  34 */
  35
  36 #include "TCanvas.h"
  37 #include "TH1.h"
  38 #include "TH2.h"
  39 #include "TAxis.h"
  40 #include "TF1.h"
  41
  42 #include "AliPerformanceRes.h"
  43 #include "AliESDEvent.h"
  44 #include "AliESDVertex.h"
  45 #include "AliESDtrack.h"
  46 #include "AliESDfriendTrack.h"
  47 #include "AliESDfriend.h"
  48 #include "AliLog.h"
  49 #include "AliMCEvent.h"
  50 #include "AliMCParticle.h"
  51 #include "AliHeader.h"
  52 #include "AliGenEventHeader.h"
  53 #include "AliStack.h"
  54 #include "AliMCInfoCuts.h"
  55 #include "AliRecInfoCuts.h"
  56 #include "AliTracker.h"
  57 #include "AliTreeDraw.h"
  58
  59 using namespace std;
  60
  61 ClassImp(AliPerformanceRes)
  62
  63 //_____________________________________________________________________________
  64 AliPerformanceRes::AliPerformanceRes():
  65   AliPerformanceObject("AliPerformanceRes"),
  66   fResolHisto(0),
  67   fPullHisto(0),
  68
  69   // Cuts
  70   fCutsRC(0),
  71   fCutsMC(0),
  72
  73   // histogram folder
  74   fAnalysisFolder(0)
  75 {
  76   Init();
  77 }
  78
  79 //_____________________________________________________________________________
  80 AliPerformanceRes::AliPerformanceRes(Char_t* name="AliPerformanceRes", Char_t* title="AliPerformanceRes",Int_t analysisMode=0,Bool_t hptGenerator=kFALSE):
  81   AliPerformanceObject(name,title),
  82   fResolHisto(0),
  83   fPullHisto(0),
  84
  85   // Cuts
  86   fCutsRC(0),
  87   fCutsMC(0),
  88
  89   // histogram folder
  90   fAnalysisFolder(0)
  91 {
  92   // named constructor
  93   //
  94   SetAnalysisMode(analysisMode);
  95   SetHptGenerator(hptGenerator);
  96
  97   Init();
  98 }
  99
 100 //_____________________________________________________________________________
 101 AliPerformanceRes::~AliPerformanceRes()
 102 {
 103   // destructor
 104
 105   if(fResolHisto) delete fResolHisto; fResolHisto=0;
 106   if(fPullHisto)  delete fPullHisto;  fPullHisto=0;
 107
 108   if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
 109 }
 110
 111 //_____________________________________________________________________________
 112 void AliPerformanceRes::Init(){
 113
 114   //
 115   // histogram bining
 116   //
 117
 118   // set pt bins
 119   Int_t nPtBins = 50;
 120   Double_t ptMin = 1.e-2, ptMax = 10.;
 121
 122   Double_t *binsPt = 0;
 123   if (IsHptGenerator())  {
 124     nPtBins = 100; ptMax = 100.;
 125     binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
 126   } else {
 127     binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
 128   }
 129
 130   Double_t yMin = -0.02, yMax = 0.02;
 131   Double_t zMin = -12.0, zMax = 12.0;
 132   if(GetAnalysisMode() == 3) { // TrackRef coordinate system
 133     yMin = -100.; yMax = 100.;
 134     zMin = -100.; zMax = 100.;
 135   }
 136
 137   // res_y:res_z:res_phi,res_lambda:res_pt:y:z:eta:phi:pt
 138   Int_t binsResolHisto[10]={100,100,100,100,100,25,50,144,30,nPtBins};
 139   Double_t minResolHisto[10]={-1.,-1.,-0.03,-0.03,-0.2, yMin, zMin, 0., -1.5, ptMin};
 140   Double_t maxResolHisto[10]={ 1., 1., 0.03, 0.03, 0.2, yMax, zMax, 2.*TMath::Pi(), 1.5, ptMax};
 141
 142   fResolHisto = new THnSparseF("fResolHisto","res_y:res_z:res_phi:res_lambda:res_pt:y:z:phi:eta:pt",10,binsResolHisto,minResolHisto,maxResolHisto);
 143   fResolHisto->SetBinEdges(9,binsPt);
 144
 145   fResolHisto->GetAxis(0)->SetTitle("y-y_{mc} (cm)");
 146   fResolHisto->GetAxis(1)->SetTitle("z-z_{mc} (cm)");
 147   fResolHisto->GetAxis(2)->SetTitle("#phi-#phi_{mc} (rad)");
 148   fResolHisto->GetAxis(3)->SetTitle("#lambda-#lambda_{mc} (rad)");
 149   fResolHisto->GetAxis(4)->SetTitle("(p_{T}/p_{Tmc}-1)");
 150   fResolHisto->GetAxis(5)->SetTitle("y_{mc} (cm)");
 151   fResolHisto->GetAxis(6)->SetTitle("z_{mc} (cm)");
 152   fResolHisto->GetAxis(7)->SetTitle("#phi_{mc} (rad)");
 153   fResolHisto->GetAxis(8)->SetTitle("#eta_{mc}");
 154   fResolHisto->GetAxis(9)->SetTitle("p_{Tmc} (GeV/c)");
 155   fResolHisto->Sumw2();
 156
 157   ////pull_y:pull_z:pull_phi:pull_lambda:pull_1pt:y:z:eta:phi:pt
 158   //Int_t binsPullHisto[10]={100,100,100,100,100,50,50,30,144,nPtBins};
 159   //Double_t minPullHisto[10]={-5.,-5.,-5.,-5.,-5.,yMin, zMin,-1.5, 0., ptMin};
 160   //Double_t maxPullHisto[10]={ 5., 5., 5., 5., 5., yMax, zMax, 1.5, 2.*TMath::Pi(),ptMax};
 161   //fPullHisto = new THnSparseF("fPullHisto","pull_y:pull_z:pull_phi:pull_lambda:pull_1pt:y:z:eta:phi:pt",10,binsPullHisto,minPullHisto,maxPullHisto);
 162
 163   //pull_y:pull_z:pull_snp:pull_tgl:pull_1pt:y:z:snp:tgl:1pt
 164   Int_t binsPullHisto[10]={100,100,100,100,100,50,50,50,50,nPtBins};
 165   Double_t minPullHisto[10]={-5.,-5.,-5.,-5.,-5.,yMin, zMin,-1., -2.0, ptMin};
 166   Double_t maxPullHisto[10]={ 5., 5., 5., 5., 5., yMax, zMax, 1., 2.0, ptMax};
 167   fPullHisto = new THnSparseF("fPullHisto","pull_y:pull_z:pull_y:pull_z:pull_snp:pull_tgl:pull_1pt:y:z:snp:tgl:1pt",10,binsPullHisto,minPullHisto,maxPullHisto);
 168
 169   /*
 170   if(!IsHptGenerator()) fPullHisto->SetBinEdges(9,bins1Pt);
 171   fPullHisto->GetAxis(0)->SetTitle("(y-y_{mc})/#sigma");
 172   fPullHisto->GetAxis(1)->SetTitle("(z-z_{mc})/#sigma");
 173   fPullHisto->GetAxis(2)->SetTitle("(#phi-#phi_{mc})/#sigma");
 174   fPullHisto->GetAxis(3)->SetTitle("(#lambda-#lambda_{mc})/#sigma");
 175   fPullHisto->GetAxis(4)->SetTitle("(p_{Tmc}/p_{T}-1)/#sigma");
 176   fPullHisto->GetAxis(5)->SetTitle("y_{mc} (cm)");
 177   fPullHisto->GetAxis(6)->SetTitle("z_{mc} (cm)");
 178   fPullHisto->GetAxis(7)->SetTitle("#eta_{mc}");
 179   fPullHisto->GetAxis(8)->SetTitle("#phi_{mc} (rad)");
 180   fPullHisto->GetAxis(9)->SetTitle("p_{Tmc} (GeV/c)");
 181   fPullHisto->Sumw2();
 182   */
 183
 184   fPullHisto->GetAxis(0)->SetTitle("(y-y_{mc})/#sigma");
 185   fPullHisto->GetAxis(1)->SetTitle("(z-z_{mc})/#sigma");
 186   fPullHisto->GetAxis(2)->SetTitle("(sin#phi-sin#phi_{mc})/#sigma");
 187   fPullHisto->GetAxis(3)->SetTitle("(tan#lambda-tan#lambda_{mc})/#sigma");
 188   fPullHisto->GetAxis(4)->SetTitle("(p_{Tmc}/p_{T}-1)/#sigma");
 189   fPullHisto->GetAxis(5)->SetTitle("y_{mc} (cm)");
 190   fPullHisto->GetAxis(6)->SetTitle("z_{mc} (cm)");
 191   fPullHisto->GetAxis(7)->SetTitle("sin#phi_{mc}");
 192   fPullHisto->GetAxis(8)->SetTitle("tan#lambda_{mc}");
 193   fPullHisto->GetAxis(9)->SetTitle("1/p_{Tmc} (GeV/c)^{-1}");
 194   fPullHisto->Sumw2();
 195
 196   // Init cuts
 197   if(!fCutsMC)
 198     AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
 199   if(!fCutsRC)
 200     AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
 201
 202   // init folder
 203   fAnalysisFolder = CreateFolder("folderRes","Analysis Resolution Folder");
 204 }
 205
 206 //_____________________________________________________________________________
 207 void AliPerformanceRes::ProcessTPC(AliStack* const stack, AliESDtrack *const esdTrack)
 208 {
 209   if(!esdTrack) return;
 210
 211   // Fill TPC only resolution comparison information
 212   const AliExternalTrackParam *track = esdTrack->GetTPCInnerParam();
 213   if(!track) return;
 214
 215   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 216   esdTrack->GetImpactParametersTPC(dca,cov);
 217
 218   //
 219   // Fill rec vs MC information
 220   //
 221   if(!stack) return;
 222   Int_t label = TMath::Abs(esdTrack->GetLabel());
 223   TParticle* particle = stack->Particle(label);
 224   if(!particle) return;
 225   if(!particle->GetPDG()) return;
 226   if(particle->GetPDG()->Charge()==0) return;
 227   //printf("charge %d \n",particle->GetPDG()->Charge());
 228
 229   // Only 5 charged particle species (e,mu,pi,K,p)
 230   if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) return;
 231
 232   // exclude electrons
 233   if (fCutsMC->GetEP()==TMath::Abs(particle->GetPdgCode())) return;
 234
 235   Float_t deltaPtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
 236   Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
 237
 238   Float_t mceta =  particle->Eta();
 239   Float_t mcphi =  particle->Phi();
 240   if(mcphi<0) mcphi += 2.*TMath::Pi();
 241   Float_t mcpt = particle->Pt();
 242   Float_t mcsnp = TMath::Sin(TMath::ATan2(particle->Py(),particle->Px()));
 243   Float_t mctgl = TMath::Tan(TMath::ATan2(particle->Pz(),particle->Pt()));
 244
 245   // nb. TPC clusters cut
 246   if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
 247
 248   // select primaries
 249   if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
 250   {
 251     if(mcpt == 0) return;
 252
 253     deltaYTPC= track->GetY()-particle->Vy();
 254     deltaZTPC = track->GetZ()-particle->Vz();
 255     deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(particle->Pz(),particle->Pt());
 256     deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(particle->Py(),particle->Px());
 257     //delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
 258     deltaPtTPC = (track->Pt()-mcpt) / mcpt;
 259
 260     pullYTPC= (track->GetY()-particle->Vy()) / TMath::Sqrt(track->GetSigmaY2());
 261     pullZTPC = (track->GetZ()-particle->Vz()) / TMath::Sqrt(track->GetSigmaZ2());
 262
 263     //Double_t sigma_lambda = 1./(1.+track->GetTgl()*track->GetTgl()) * TMath::Sqrt(track->GetSigmaTgl2());
 264     //Double_t sigma_phi = 1./TMath::Sqrt(1-track->GetSnp()*track->GetSnp()) * TMath::Sqrt(track->GetSigmaSnp2());
 265     pullPhiTPC = (track->GetSnp() - mcsnp) / TMath::Sqrt(track->GetSigmaSnp2());
 266     pullLambdaTPC = (track->GetTgl() - mctgl) / TMath::Sqrt(track->GetSigmaTgl2());
 267
 268     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
 269     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
 270     if (mcpt) pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
 271     else pull1PtTPC = 0.;
 272
 273     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,deltaPtTPC,particle->Vy(),particle->Vz(),mcphi,mceta,mcpt};
 274     fResolHisto->Fill(vResolHisto);
 275
 276     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,particle->Vy(),particle->Vz(),mcsnp,mctgl,1./mcpt};
 277     fPullHisto->Fill(vPullHisto);
 278   }
 279 }
 280
 281 //_____________________________________________________________________________
 282 void AliPerformanceRes::ProcessTPCITS(AliStack* const stack, AliESDtrack *const esdTrack)
 283 {
 284   // Fill resolution comparison information (TPC+ITS)
 285   if(!esdTrack) return;
 286
 287   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 288   esdTrack->GetImpactParameters(dca,cov);
 289
 290   //
 291   // Fill rec vs MC information
 292   //
 293   if(!stack) return;
 294
 295   Int_t label = TMath::Abs(esdTrack->GetLabel());
 296   TParticle* particle = stack->Particle(label);
 297   if(!particle) return;
 298   if(particle->GetPDG()->Charge()==0) return;
 299   //printf("charge %d \n",particle->GetPDG()->Charge());
 300
 301
 302   // Only 5 charged particle species (e,mu,pi,K,p)
 303   if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) return;
 304
 305   // exclude electrons
 306   if (fCutsMC->GetEP()==TMath::Abs(particle->GetPdgCode())) return;
 307
 308   Float_t mceta =  particle->Eta();
 309   Float_t mcphi =  particle->Phi();
 310   if(mcphi<0) mcphi += 2.*TMath::Pi();
 311   Float_t mcpt = particle->Pt();
 312   Float_t mcsnp = TMath::Sin(TMath::ATan2(particle->Py(),particle->Px()));
 313   Float_t mctgl = TMath::Tan(TMath::ATan2(particle->Pz(),particle->Pt()));
 314
 315   if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
 316   if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
 317   Int_t clusterITS[200];
 318   if(esdTrack->GetITSclusters(clusterITS)<fCutsRC->GetMinNClustersITS()) return;  // min. nb. ITS clusters
 319
 320   Float_t deltaPtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
 321   Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
 322
 323   // select primaries
 324   if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
 325   {
 326     if(mcpt == 0) return;
 327
 328     deltaYTPC= esdTrack->GetY()-particle->Vy();
 329     deltaZTPC = esdTrack->GetZ()-particle->Vz();
 330     deltaLambdaTPC = TMath::ATan2(esdTrack->Pz(),esdTrack->Pt())-TMath::ATan2(particle->Pz(),particle->Pt());
 331     deltaPhiTPC = TMath::ATan2(esdTrack->Py(),esdTrack->Px())-TMath::ATan2(particle->Py(),particle->Px());
 332     //delta1PtTPC = (esdTrack->OneOverPt()-1./mcpt)*mcpt;
 333     deltaPtTPC = (esdTrack->Pt()-mcpt) / mcpt;
 334
 335     pullYTPC= (esdTrack->GetY()-particle->Vy()) / TMath::Sqrt(esdTrack->GetSigmaY2());
 336     pullZTPC = (esdTrack->GetZ()-particle->Vz()) / TMath::Sqrt(esdTrack->GetSigmaZ2());
 337
 338     //Double_t sigma_lambda = 1./(1.+esdTrack->GetTgl()*esdTrack->GetTgl()) * TMath::Sqrt(esdTrack->GetSigmaTgl2());
 339     //Double_t sigma_phi = 1./TMath::Sqrt(1-esdTrack->GetSnp()*esdTrack->GetSnp()) * TMath::Sqrt(esdTrack->GetSigmaSnp2());
 340     pullPhiTPC = (esdTrack->GetSnp() - mcsnp) / TMath::Sqrt(esdTrack->GetSigmaSnp2());
 341     pullLambdaTPC = (esdTrack->GetTgl() - mctgl) / TMath::Sqrt(esdTrack->GetSigmaTgl2());
 342
 343     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(esdTrack->GetSigmaTgl2());
 344     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(esdTrack->GetSigmaSnp2());
 345     if (mcpt) pull1PtTPC = (esdTrack->OneOverPt()-1./mcpt) / TMath::Sqrt(esdTrack->GetSigma1Pt2());
 346     else pull1PtTPC = 0.;
 347
 348     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,deltaPtTPC,particle->Vy(),particle->Vz(),mcphi,mceta,mcpt};
 349     fResolHisto->Fill(vResolHisto);
 350
 351     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,particle->Vy(),particle->Vz(),mcsnp,mctgl,1./mcpt};
 352     fPullHisto->Fill(vPullHisto);
 353
 354
 355     /*
 356     deltaYTPC= esdTrack->GetY()-particle->Vy();
 357     deltaZTPC = esdTrack->GetZ()-particle->Vz();
 358     deltaLambdaTPC = TMath::ATan2(esdTrack->Pz(),esdTrack->Pt())-TMath::ATan2(particle->Pz(),particle->Pt());
 359     deltaPhiTPC = TMath::ATan2(esdTrack->Py(),esdTrack->Px())-TMath::ATan2(particle->Py(),particle->Px());
 360     delta1PtTPC = (esdTrack->OneOverPt()-1./mcpt)*mcpt;
 361
 362     pullYTPC= (esdTrack->GetY()-particle->Vy()) / TMath::Sqrt(esdTrack->GetSigmaY2());
 363     pullZTPC = (esdTrack->GetZ()-particle->Vz()) / TMath::Sqrt(esdTrack->GetSigmaZ2());
 364
 365     Double_t sigma_lambda = 1./(1.+esdTrack->GetTgl()*esdTrack->GetTgl()) * TMath::Sqrt(esdTrack->GetSigmaTgl2());
 366     Double_t sigma_phi = 1./TMath::Sqrt(1-esdTrack->GetSnp()*esdTrack->GetSnp()) * TMath::Sqrt(esdTrack->GetSigmaSnp2());
 367     pullLambdaTPC = deltaLambdaTPC / sigma_lambda;
 368     pullPhiTPC = deltaPhiTPC / sigma_phi;
 369
 370     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(esdTrack->GetSigmaTgl2());
 371     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(esdTrack->GetSigmaSnp2());
 372     if (mcpt) pull1PtTPC = (esdTrack->OneOverPt()-1./mcpt) / TMath::Sqrt(esdTrack->GetSigma1Pt2());
 373     else pull1PtTPC = 0.;
 374
 375     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,delta1PtTPC,particle->Vy(),particle->Vz(),mceta,mcphi,mcpt};
 376     fResolHisto->Fill(vResolHisto);
 377
 378     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,particle->Vy(),particle->Vz(),mceta,mcphi,mcpt};
 379     fPullHisto->Fill(vPullHisto);
 380     */
 381   }
 382 }
 383
 384 //_____________________________________________________________________________
 385 void AliPerformanceRes::ProcessConstrained(AliStack* const stack, AliESDtrack *const esdTrack)
 386 {
 387   // Fill resolution comparison information (constarained parameters)
 388   //
 389   if(!esdTrack) return;
 390
 391   const AliExternalTrackParam * track = esdTrack->GetConstrainedParam();
 392   if(!track) return;
 393
 394   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 395   esdTrack->GetImpactParameters(dca,cov);
 396
 397   //
 398   // Fill rec vs MC information
 399   //
 400   if(!stack) return;
 401
 402   Int_t label = TMath::Abs(esdTrack->GetLabel());
 403   TParticle* particle = stack->Particle(label);
 404   if(!particle) return;
 405   if(particle->GetPDG()->Charge()==0) return;
 406   //printf("charge %d \n",particle->GetPDG()->Charge());
 407
 408   // Only 5 charged particle species (e,mu,pi,K,p)
 409   if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) return;
 410
 411   // exclude electrons
 412   if (fCutsMC->GetEP()==TMath::Abs(particle->GetPdgCode())) return;
 413
 414   Float_t mceta =  particle->Eta();
 415   Float_t mcphi =  particle->Phi();
 416   if(mcphi<0) mcphi += 2.*TMath::Pi();
 417   Float_t mcpt = particle->Pt();
 418   Float_t mcsnp = TMath::Sin(TMath::ATan2(particle->Py(),particle->Px()));
 419   Float_t mctgl = TMath::Tan(TMath::ATan2(particle->Pz(),particle->Pt()));
 420
 421   if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
 422   if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
 423
 424   Float_t deltaPtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
 425   Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
 426
 427   // select primaries
 428   if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
 429   {
 430
 431     if(mcpt == 0) return;
 432
 433     deltaYTPC= track->GetY()-particle->Vy();
 434     deltaZTPC = track->GetZ()-particle->Vz();
 435     deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(particle->Pz(),particle->Pt());
 436     deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(particle->Py(),particle->Px());
 437     //delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
 438     deltaPtTPC = (track->Pt()-mcpt) / mcpt;
 439
 440     pullYTPC= (track->GetY()-particle->Vy()) / TMath::Sqrt(track->GetSigmaY2());
 441     pullZTPC = (track->GetZ()-particle->Vz()) / TMath::Sqrt(track->GetSigmaZ2());
 442
 443     //Double_t sigma_lambda = 1./(1.+track->GetTgl()*track->GetTgl()) * TMath::Sqrt(track->GetSigmaTgl2());
 444     //Double_t sigma_phi = 1./TMath::Sqrt(1-track->GetSnp()*track->GetSnp()) * TMath::Sqrt(track->GetSigmaSnp2());
 445     pullPhiTPC = (track->GetSnp() - mcsnp) / TMath::Sqrt(track->GetSigmaSnp2());
 446     pullLambdaTPC = (track->GetTgl() - mctgl) / TMath::Sqrt(track->GetSigmaTgl2());
 447
 448     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
 449     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
 450     if (mcpt) pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
 451     else pull1PtTPC = 0.;
 452
 453     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,deltaPtTPC,particle->Vy(),particle->Vz(),mcphi,mceta,mcpt};
 454     fResolHisto->Fill(vResolHisto);
 455
 456     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,particle->Vy(),particle->Vz(),mcsnp,mctgl,1./mcpt};
 457     fPullHisto->Fill(vPullHisto);
 458
 459     /*
 460
 461     deltaYTPC= track->GetY()-particle->Vy();
 462     deltaZTPC = track->GetZ()-particle->Vz();
 463     deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(particle->Pz(),particle->Pt());
 464     deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(particle->Py(),particle->Px());
 465     delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
 466
 467     pullYTPC= (track->GetY()-particle->Vy()) / TMath::Sqrt(track->GetSigmaY2());
 468     pullZTPC = (track->GetZ()-particle->Vz()) / TMath::Sqrt(track->GetSigmaZ2());
 469
 470     Double_t sigma_lambda = 1./(1.+track->GetTgl()*track->GetTgl()) * TMath::Sqrt(track->GetSigmaTgl2());
 471     Double_t sigma_phi = 1./TMath::Sqrt(1-track->GetSnp()*track->GetSnp()) * TMath::Sqrt(track->GetSigmaSnp2());
 472     pullLambdaTPC = deltaLambdaTPC / sigma_lambda;
 473     pullPhiTPC = deltaPhiTPC / sigma_phi;
 474
 475     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
 476     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
 477
 478     if (mcpt) pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
 479     else pull1PtTPC = 0.;
 480
 481     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,delta1PtTPC,particle->Vy(),particle->Vz(),mceta,mcphi,mcpt};
 482     fResolHisto->Fill(vResolHisto);
 483
 484     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,particle->Vy(),particle->Vz(),mceta,mcphi,mcpt};
 485     fPullHisto->Fill(vPullHisto);
 486
 487     */
 488   }
 489 }
 490
 491 //_____________________________________________________________________________
 492 void AliPerformanceRes::ProcessInnerTPC(AliMCEvent *const mcEvent, AliESDtrack *const esdTrack)
 493 {
 494   //
 495   // Fill resolution comparison information (inner params at TPC reference point)
 496   //
 497   if(!esdTrack) return;
 498
 499   const AliExternalTrackParam * innerParam = esdTrack->GetInnerParam();
 500   if(!innerParam) return;
 501
 502   // create new AliExternalTrackParam
 503   AliExternalTrackParam *track = new AliExternalTrackParam(*innerParam);
 504   if(!track) return;
 505
 506   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 507   esdTrack->GetImpactParametersTPC(dca,cov);
 508
 509   //
 510   // Fill rec vs MC information
 511   //
 512   if(!mcEvent) return;
 513
 514   Int_t label = TMath::Abs(esdTrack->GetLabel());
 515   AliMCParticle *mcParticle = mcEvent->GetTrack(label);
 516   if(!mcParticle) return;
 517
 518   // get the first TPC track reference
 519   AliTrackReference *ref0 = GetFirstTPCTrackRef(mcParticle);
 520   if(!ref0) return;
 521
 522   // Only 5 charged particle species (e,mu,pi,K,p)
 523   TParticle *particle = mcParticle->Particle();
 524   if(!particle) return;
 525
 526   if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) return;
 527
 528   // exclude electrons
 529   if (fCutsMC->GetEP()==TMath::Abs(particle->GetPdgCode())) return;
 530
 531   // calculate alpha angle
 532   Double_t xyz[3] = {ref0->X(),ref0->Y(),ref0->Z()};
 533   Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
 534   //if (alpha<0) alpha += TMath::TwoPi();
 535
 536   // rotate inner track to local coordinate system
 537   // and propagate to the radius of the first track referenco of TPC
 538   Double_t trRadius = TMath::Sqrt(xyz[1] * xyz[1] + xyz[0] * xyz[0]);
 539   Bool_t isOK = track->Propagate(alpha,trRadius,AliTracker::GetBz());
 540   if(!isOK) return;
 541
 542   Float_t mceta =  -TMath::Log(TMath::Tan(0.5 * ref0->Theta()));
 543   Float_t mcphi =  ref0->Phi();
 544   if(mcphi<0) mcphi += 2.*TMath::Pi();
 545   Float_t mcpt = ref0->Pt();
 546   Float_t mcsnp = TMath::Sin(TMath::ATan2(ref0->Py(),ref0->Px()));
 547   Float_t mctgl = TMath::Tan(TMath::ATan2(ref0->Pz(),ref0->Pt()));
 548
 549   if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
 550   if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
 551
 552   Float_t deltaPtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
 553   Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
 554
 555   // select primaries
 556   if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
 557   {
 558     if(mcpt == 0) return;
 559
 560     deltaYTPC= track->GetY(); // already rotated
 561     deltaZTPC = track->GetZ()-ref0->Z();
 562     deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(ref0->Pz(),ref0->Pt());
 563     deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(ref0->Py(),ref0->Px());
 564     //delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
 565     deltaPtTPC = (track->Pt()-mcpt) / mcpt;
 566
 567     pullYTPC= track->GetY() / TMath::Sqrt(track->GetSigmaY2());
 568     pullZTPC = (track->GetZ()-ref0->Z()) / TMath::Sqrt(track->GetSigmaZ2());
 569
 570     //Double_t sigma_lambda = 1./(1.+track->GetTgl()*track->GetTgl()) * TMath::Sqrt(track->GetSigmaTgl2());
 571     //Double_t sigma_phi = 1./TMath::Sqrt(1-track->GetSnp()*track->GetSnp()) * TMath::Sqrt(track->GetSigmaSnp2());
 572     pullPhiTPC = (track->GetSnp() - mcsnp) / TMath::Sqrt(track->GetSigmaSnp2());
 573     pullLambdaTPC = (track->GetTgl() - mctgl) / TMath::Sqrt(track->GetSigmaTgl2());
 574
 575     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
 576     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
 577     if (mcpt) pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
 578     else pull1PtTPC = 0.;
 579
 580     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,deltaPtTPC,ref0->Y(),ref0->Z(),mcphi,mceta,mcpt};
 581     fResolHisto->Fill(vResolHisto);
 582
 583     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,ref0->Y(),ref0->Z(),mcsnp,mctgl,1./mcpt};
 584     fPullHisto->Fill(vPullHisto);
 585   }
 586
 587   if(track) delete track;
 588 }
 589
 590 //_____________________________________________________________________________
 591 void AliPerformanceRes::ProcessOuterTPC(AliMCEvent *const mcEvent, AliESDtrack *const esdTrack, AliESDfriendTrack *const friendTrack)
 592 {
 593   //
 594   // Fill resolution comparison information (outer params at TPC reference point)
 595   //
 596   if(!friendTrack) return;
 597
 598   const AliExternalTrackParam * outerParam = friendTrack->GetTPCOut();
 599   if(!outerParam) return;
 600
 601   // create new AliExternalTrackParam
 602   AliExternalTrackParam *track = new AliExternalTrackParam(*outerParam);
 603   if(!track) return;
 604
 605   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 606   esdTrack->GetImpactParametersTPC(dca,cov);
 607
 608   //
 609   // Fill rec vs MC information
 610   //
 611   if(!mcEvent) return;
 612
 613   Int_t label = TMath::Abs(esdTrack->GetLabel());
 614   AliMCParticle *mcParticle = mcEvent->GetTrack(label);
 615   if(!mcParticle) return;
 616
 617   // get the last TPC track reference
 618   AliTrackReference *ref0 = GetLastTPCTrackRef(mcParticle);
 619   if(!ref0) return;
 620
 621   // Only 5 charged particle species (e,mu,pi,K,p)
 622   TParticle *particle = mcParticle->Particle();
 623   if(!particle) return;
 624
 625   if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) return;
 626
 627   // exclude electrons
 628   if (fCutsMC->GetEP()==TMath::Abs(particle->GetPdgCode())) return;
 629
 630   // calculate alpha angle
 631   Double_t xyz[3] = {ref0->X(),ref0->Y(),ref0->Z()};
 632   Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
 633   //if (alpha<0) alpha += TMath::TwoPi();
 634
 635   // rotate outer track to local coordinate system
 636   // and propagate to the radius of the last track reference of TPC
 637   Double_t trRadius = TMath::Sqrt(xyz[1] * xyz[1] + xyz[0] * xyz[0]);
 638   Bool_t isOK = track->Propagate(alpha,trRadius,AliTracker::GetBz());
 639   if(!isOK) return;
 640
 641   Float_t mceta =  -TMath::Log(TMath::Tan(0.5 * ref0->Theta()));
 642   Float_t mcphi =  ref0->Phi();
 643   if(mcphi<0) mcphi += 2.*TMath::Pi();
 644   Float_t mcpt = ref0->Pt();
 645   Float_t mcsnp = TMath::Sin(TMath::ATan2(ref0->Py(),ref0->Px()));
 646   Float_t mctgl = TMath::Tan(TMath::ATan2(ref0->Pz(),ref0->Pt()));
 647
 648   if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
 649   if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
 650
 651   Float_t deltaPtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
 652   Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
 653
 654   // select primaries
 655   if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
 656   {
 657     if(mcpt == 0) return;
 658
 659     deltaYTPC= track->GetY(); // already rotated
 660     deltaZTPC = track->GetZ()-ref0->Z();
 661     deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(ref0->Pz(),ref0->Pt());
 662     deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(ref0->Py(),ref0->Px());
 663     //delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
 664     deltaPtTPC = (track->Pt()-mcpt) / mcpt;
 665
 666     pullYTPC= track->GetY() / TMath::Sqrt(track->GetSigmaY2());
 667     pullZTPC = (track->GetZ()-ref0->Z()) / TMath::Sqrt(track->GetSigmaZ2());
 668
 669     //Double_t sigma_lambda = 1./(1.+track->GetTgl()*track->GetTgl()) * TMath::Sqrt(track->GetSigmaTgl2());
 670     //Double_t sigma_phi = 1./TMath::Sqrt(1-track->GetSnp()*track->GetSnp()) * TMath::Sqrt(track->GetSigmaSnp2());
 671     pullPhiTPC = (track->GetSnp() - mcsnp) / TMath::Sqrt(track->GetSigmaSnp2());
 672     pullLambdaTPC = (track->GetTgl() - mctgl) / TMath::Sqrt(track->GetSigmaTgl2());
 673
 674     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
 675     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
 676     if (mcpt) pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
 677     else pull1PtTPC = 0.;
 678
 679     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,deltaPtTPC,ref0->Y(),ref0->Z(),mcphi,mceta,mcpt};
 680     fResolHisto->Fill(vResolHisto);
 681
 682     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,ref0->Y(),ref0->Z(),mcsnp,mctgl,1./mcpt};
 683     fPullHisto->Fill(vPullHisto);
 684   }
 685
 686   if(track) delete track;
 687 }
 688
 689 //_____________________________________________________________________________
 690 AliTrackReference * AliPerformanceRes::GetFirstTPCTrackRef(AliMCParticle *mcParticle)
 691 {
 692   // return first TPC track reference
 693   if(!mcParticle) return 0;
 694
 695   // find first track reference
 696   // check direction to select proper reference point for looping tracks
 697   Int_t nTrackRef = mcParticle->GetNumberOfTrackReferences();
 698   AliTrackReference *ref = 0;
 699   AliTrackReference *refIn = 0;
 700   for (Int_t iref = 0; iref < nTrackRef; iref++) {
 701     ref = mcParticle->GetTrackReference(iref);
 702     if(ref && (ref->DetectorId()==AliTrackReference::kTPC))
 703     {
 704       Float_t dir = ref->X()*ref->Px()+ref->Y()*ref->Py();
 705       if(dir < 0.) break;
 706
 707       refIn = ref;
 708       break;
 709     }
 710   }
 711
 712 return refIn;
 713 }
 714
 715 //_____________________________________________________________________________
 716 AliTrackReference * AliPerformanceRes::GetLastTPCTrackRef(AliMCParticle *mcParticle)
 717 {
 718   // return last TPC track reference
 719   if(!mcParticle) return 0;
 720
 721   // find last track reference
 722   // check direction to select proper reference point for looping tracks
 723   Int_t nTrackRef = mcParticle->GetNumberOfTrackReferences();
 724   AliTrackReference *ref = 0;
 725   AliTrackReference *refOut = 0;
 726   for (Int_t iref = 0; iref < nTrackRef; iref++) {
 727     ref = mcParticle->GetTrackReference(iref);
 728     if(ref && (ref->DetectorId()==AliTrackReference::kTPC)) {
 729       Float_t dir=ref->X()*ref->Px()+ref->Y()*ref->Py();
 730       if(dir< 0.0) break;
 731       refOut = ref;
 732     }
 733   }
 734
 735 return refOut;
 736 }
 737
 738 //_____________________________________________________________________________
 739 void AliPerformanceRes::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const esdFriend, const Bool_t bUseMC, const Bool_t bUseESDfriend)
 740 {
 741   // Process comparison information
 742   //
 743   if(!esdEvent)
 744   {
 745       AliDebug(AliLog::kError, "esdEvent not available");
 746       return;
 747   }
 748   AliHeader* header = 0;
 749   AliGenEventHeader* genHeader = 0;
 750   AliStack* stack = 0;
 751   TArrayF vtxMC(3);
 752
 753   if(bUseMC)
 754   {
 755     if(!mcEvent) {
 756       AliDebug(AliLog::kError, "mcEvent not available");
 757       return;
 758     }
 759
 760     // get MC event header
 761     header = mcEvent->Header();
 762     if (!header) {
 763       AliDebug(AliLog::kError, "Header not available");
 764       return;
 765     }
 766     // MC particle stack
 767     stack = mcEvent->Stack();
 768     if (!stack) {
 769       AliDebug(AliLog::kError, "Stack not available");
 770       return;
 771     }
 772
 773     // get MC vertex
 774     genHeader = header->GenEventHeader();
 775     if (!genHeader) {
 776       AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");
 777       return;
 778     }
 779     genHeader->PrimaryVertex(vtxMC);
 780
 781   } // end bUseMC
 782
 783   // use ESD friends
 784   if(bUseESDfriend) {
 785     if(!esdFriend) {
 786       AliDebug(AliLog::kError, "esdFriend not available");
 787       return;
 788     }
 789   }
 790
 791   //  Process events
 792   for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
 793   {
 794     AliESDtrack *track = esdEvent->GetTrack(iTrack);
 795     if(!track) continue;
 796
 797     AliESDfriendTrack *friendTrack=0;
 798     if(bUseESDfriend) {
 799       friendTrack=esdFriend->GetTrack(iTrack);
 800       if(!friendTrack) continue;
 801     }
 802
 803     if(GetAnalysisMode() == 0) ProcessTPC(stack,track);
 804     else if(GetAnalysisMode() == 1) ProcessTPCITS(stack,track);
 805     else if(GetAnalysisMode() == 2) ProcessConstrained(stack,track);
 806     else if(GetAnalysisMode() == 3) ProcessInnerTPC(mcEvent,track);
 807     else if(GetAnalysisMode() == 4) ProcessOuterTPC(mcEvent,track,friendTrack);
 808     else {
 809       printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
 810       return;
 811     }
 812   }
 813 }
 814
 815 //_____________________________________________________________________________
 816 TH1F* AliPerformanceRes::MakeResol(TH2F * his, Int_t integ, Bool_t type, Int_t cut){
 817   // Create resolution histograms
 818
 819   TH1F *hisr, *hism;
 820   if (!gPad) new TCanvas;
 821   hisr = AliTreeDraw::CreateResHistoII(his,&hism,integ,kTRUE,cut);
 822   if (type) return hism;
 823   else
 824     return hisr;
 825 }
 826
 827 //_____________________________________________________________________________
 828 void AliPerformanceRes::Analyse() {
 829   // Analyse comparison information and store output histograms
 830   // in the folder "folderRes"
 831   //
 832   TH1::AddDirectory(kFALSE);
 833   TH1F *h=0;
 834   TH2F *h2D=0;
 835   TObjArray *aFolderObj = new TObjArray;
 836
 837   // write results in the folder
 838   TCanvas * c = new TCanvas("Phi resol Tan","Phi resol Tan");
 839   c->cd();
 840
 841   char name[256];
 842   char title[256];
 843
 844   for(Int_t i=0; i<5; i++)
 845   {
 846     for(Int_t j=5; j<10; j++)
 847     {
 848       if(j!=8) fResolHisto->GetAxis(8)->SetRangeUser(-0.9,0.89); // eta window
 849       fResolHisto->GetAxis(9)->SetRangeUser(0.16,10.); // pt threshold
 850       //if(j!=8) fResolHisto->GetAxis(8)->SetRangeUser(0.,0.9); // eta window
 851       //fResolHisto->GetAxis(9)->SetRangeUser(0.16,3.); // pt window
 852       if(GetAnalysisMode() == 3) fResolHisto->GetAxis(5)->SetRangeUser(-80.,80.); // y range
 853
 854       h2D = (TH2F*)fResolHisto->Projection(i,j);
 855
 856       h = AliPerformanceRes::MakeResol(h2D,1,0,100);
 857       sprintf(name,"h_res_%d_vs_%d",i,j);
 858       h->SetName(name);
 859
 860       h->GetXaxis()->SetTitle(fResolHisto->GetAxis(j)->GetTitle());
 861       sprintf(title,"%s %s",fResolHisto->GetAxis(i)->GetTitle(),"(resolution)");
 862       h->GetYaxis()->SetTitle(title);
 863       sprintf(title,"%s vs %s",title,fResolHisto->GetAxis(j)->GetTitle());
 864       h->SetTitle(title);
 865
 866       if(j==9) h->SetBit(TH1::kLogX);
 867       aFolderObj->Add(h);
 868
 869       h = AliPerformanceRes::MakeResol(h2D,1,1,100);
 870       //h = (TH1F*)arr->At(1);
 871       sprintf(name,"h_mean_res_%d_vs_%d",i,j);
 872       h->SetName(name);
 873
 874       h->GetXaxis()->SetTitle(fResolHisto->GetAxis(j)->GetTitle());
 875       sprintf(title,"%s %s",fResolHisto->GetAxis(i)->GetTitle(),"(mean)");
 876       h->GetYaxis()->SetTitle(title);
 877
 878       sprintf(title,"%s vs %s",title,fResolHisto->GetAxis(j)->GetTitle());
 879       h->SetTitle(title);
 880
 881       if(j==9) h->SetBit(TH1::kLogX);
 882       aFolderObj->Add(h);
 883
 884       fResolHisto->GetAxis(8)->SetRangeUser(-1.5,1.5);
 885       fResolHisto->GetAxis(9)->SetRangeUser(0.0,10.);
 886
 887       //
 888       if(j!=8) fPullHisto->GetAxis(8)->SetRangeUser(-1.0,0.99); // theta window
 889       else  fPullHisto->GetAxis(8)->SetRangeUser(-1.5,1.5); // theta window
 890       fPullHisto->GetAxis(9)->SetRangeUser(0.0,10.);  // pt threshold
 891       if(GetAnalysisMode() == 3) fPullHisto->GetAxis(5)->SetRangeUser(-80.,80.); // y range
 892
 893       h2D = (TH2F*)fPullHisto->Projection(i,j);
 894
 895       h = AliPerformanceRes::MakeResol(h2D,1,0,100);
 896       sprintf(name,"h_pull_%d_vs_%d",i,j);
 897       h->SetName(name);
 898
 899       h->GetXaxis()->SetTitle(fPullHisto->GetAxis(j)->GetTitle());
 900       sprintf(title,"%s %s",fPullHisto->GetAxis(i)->GetTitle(),"(resolution)");
 901       h->GetYaxis()->SetTitle(title);
 902       sprintf(title,"%s vs %s",title,fPullHisto->GetAxis(j)->GetTitle());
 903       h->SetTitle(title);
 904
 905       //if(j==9) h->SetBit(TH1::kLogX);
 906       aFolderObj->Add(h);
 907
 908       h = AliPerformanceRes::MakeResol(h2D,1,1,100);
 909       sprintf(name,"h_mean_pull_%d_vs_%d",i,j);
 910       h->SetName(name);
 911
 912       h->GetXaxis()->SetTitle(fPullHisto->GetAxis(j)->GetTitle());
 913       sprintf(title,"%s %s",fPullHisto->GetAxis(i)->GetTitle(),"(mean)");
 914       h->GetYaxis()->SetTitle(title);
 915       sprintf(title,"%s vs %s",title,fPullHisto->GetAxis(j)->GetTitle());
 916       h->SetTitle(title);
 917
 918       //if(j==9) h->SetBit(TH1::kLogX);
 919       aFolderObj->Add(h);
 920     }
 921   }
 922
 923   // export objects to analysis folder
 924   fAnalysisFolder = ExportToFolder(aFolderObj);
 925
 926   // delete only TObjArray
 927   if(aFolderObj) delete aFolderObj;
 928 }
 929
 930 //_____________________________________________________________________________
 931 TFolder* AliPerformanceRes::ExportToFolder(TObjArray * array)
 932 {
 933   // recreate folder avery time and export objects to new one
 934   //
 935   AliPerformanceRes * comp=this;
 936   TFolder *folder = comp->GetAnalysisFolder();
 937
 938   TString name, title;
 939   TFolder *newFolder = 0;
 940   Int_t i = 0;
 941   Int_t size = array->GetSize();
 942
 943   if(folder) {
 944      // get name and title from old folder
 945      name = folder->GetName();
 946      title = folder->GetTitle();
 947
 948          // delete old one
 949      delete folder;
 950
 951          // create new one
 952      newFolder = CreateFolder(name.Data(),title.Data());
 953      newFolder->SetOwner();
 954
 955          // add objects to folder
 956      while(i < size) {
 957            newFolder->Add(array->At(i));
 958            i++;
 959          }
 960   }
 961
 962 return newFolder;
 963 }
 964
 965 //_____________________________________________________________________________
 966 Long64_t AliPerformanceRes::Merge(TCollection* const list)
 967 {
 968   // Merge list of objects (needed by PROOF)
 969
 970   if (!list)
 971   return 0;
 972
 973   if (list->IsEmpty())
 974   return 1;
 975
 976   TIterator* iter = list->MakeIterator();
 977   TObject* obj = 0;
 978
 979   // collection of generated histograms
 980   Int_t count=0;
 981   while((obj = iter->Next()) != 0)
 982   {
 983   AliPerformanceRes* entry = dynamic_cast<AliPerformanceRes*>(obj);
 984   if (entry == 0) continue;
 985
 986   fResolHisto->Add(entry->fResolHisto);
 987   fPullHisto->Add(entry->fPullHisto);
 988
 989   count++;
 990   }
 991
 992 return count;
 993 }
 994
 995 //_____________________________________________________________________________
 996 TFolder* AliPerformanceRes::CreateFolder(TString name,TString title) {
 997 // create folder for analysed histograms
 998 //
 999 TFolder *folder = 0;
1000   folder = new TFolder(name.Data(),title.Data());
1001
1002   return folder;
1003 }