PWG1/TPC/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,90,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,90,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->GetEM()==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->GetEM()==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->GetEM()==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 = (AliMCParticle*) 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->GetEM()==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   Double_t field[3]; track->GetBxByBz(field);
 541   Bool_t isOK = track->PropagateBxByBz(alpha,trRadius,field);
 542   if(!isOK) return;
 543
 544   Float_t mceta =  -TMath::Log(TMath::Tan(0.5 * ref0->Theta()));
 545   Float_t mcphi =  ref0->Phi();
 546   if(mcphi<0) mcphi += 2.*TMath::Pi();
 547   Float_t mcpt = ref0->Pt();
 548   Float_t mcsnp = TMath::Sin(TMath::ATan2(ref0->Py(),ref0->Px()));
 549   Float_t mctgl = TMath::Tan(TMath::ATan2(ref0->Pz(),ref0->Pt()));
 550
 551   if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
 552   if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
 553
 554   Float_t deltaPtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
 555   Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
 556
 557   // select primaries
 558   if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
 559   {
 560     if(mcpt == 0) return;
 561
 562     deltaYTPC= track->GetY(); // already rotated
 563     deltaZTPC = track->GetZ()-ref0->Z();
 564     deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(ref0->Pz(),ref0->Pt());
 565     deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(ref0->Py(),ref0->Px());
 566     //delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
 567     deltaPtTPC = (track->Pt()-mcpt) / mcpt;
 568
 569     pullYTPC= track->GetY() / TMath::Sqrt(track->GetSigmaY2());
 570     pullZTPC = (track->GetZ()-ref0->Z()) / TMath::Sqrt(track->GetSigmaZ2());
 571
 572     //Double_t sigma_lambda = 1./(1.+track->GetTgl()*track->GetTgl()) * TMath::Sqrt(track->GetSigmaTgl2());
 573     //Double_t sigma_phi = 1./TMath::Sqrt(1-track->GetSnp()*track->GetSnp()) * TMath::Sqrt(track->GetSigmaSnp2());
 574     pullPhiTPC = (track->GetSnp() - mcsnp) / TMath::Sqrt(track->GetSigmaSnp2());
 575     pullLambdaTPC = (track->GetTgl() - mctgl) / TMath::Sqrt(track->GetSigmaTgl2());
 576
 577     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
 578     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
 579     if (mcpt) pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
 580     else pull1PtTPC = 0.;
 581
 582     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,deltaPtTPC,ref0->Y(),ref0->Z(),mcphi,mceta,mcpt};
 583     fResolHisto->Fill(vResolHisto);
 584
 585     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,ref0->Y(),ref0->Z(),mcsnp,mctgl,1./mcpt};
 586     fPullHisto->Fill(vPullHisto);
 587   }
 588
 589   if(track) delete track;
 590 }
 591
 592 //_____________________________________________________________________________
 593 void AliPerformanceRes::ProcessOuterTPC(AliMCEvent *const mcEvent, AliESDtrack *const esdTrack, AliESDfriendTrack *const friendTrack)
 594 {
 595   //
 596   // Fill resolution comparison information (outer params at TPC reference point)
 597   //
 598   if(!friendTrack) return;
 599
 600   const AliExternalTrackParam * outerParam = friendTrack->GetTPCOut();
 601   if(!outerParam) return;
 602
 603   // create new AliExternalTrackParam
 604   AliExternalTrackParam *track = new AliExternalTrackParam(*outerParam);
 605   if(!track) return;
 606
 607   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 608   esdTrack->GetImpactParametersTPC(dca,cov);
 609
 610   //
 611   // Fill rec vs MC information
 612   //
 613   if(!mcEvent) return;
 614
 615   Int_t label = TMath::Abs(esdTrack->GetLabel());
 616   AliMCParticle *mcParticle = (AliMCParticle*) mcEvent->GetTrack(label);
 617   if(!mcParticle) return;
 618
 619   // get the last TPC track reference
 620   AliTrackReference *ref0 = GetLastTPCTrackRef(mcParticle);
 621   if(!ref0) return;
 622
 623   // Only 5 charged particle species (e,mu,pi,K,p)
 624   TParticle *particle = mcParticle->Particle();
 625   if(!particle) return;
 626
 627   if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) return;
 628
 629   // exclude electrons
 630   if (fCutsMC->GetEM()==TMath::Abs(particle->GetPdgCode())) return;
 631
 632   // calculate alpha angle
 633   Double_t xyz[3] = {ref0->X(),ref0->Y(),ref0->Z()};
 634   Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
 635   //if (alpha<0) alpha += TMath::TwoPi();
 636
 637   // rotate outer track to local coordinate system
 638   // and propagate to the radius of the last track reference of TPC
 639   Double_t trRadius = TMath::Sqrt(xyz[1] * xyz[1] + xyz[0] * xyz[0]);
 640   //Bool_t isOK = track->Propagate(alpha,trRadius,AliTracker::GetBz());
 641   Double_t field[3]; track->GetBxByBz(field);
 642   Bool_t isOK = track->PropagateBxByBz(alpha,trRadius,field);
 643   if(!isOK) return;
 644
 645   Float_t mceta =  -TMath::Log(TMath::Tan(0.5 * ref0->Theta()));
 646   Float_t mcphi =  ref0->Phi();
 647   if(mcphi<0) mcphi += 2.*TMath::Pi();
 648   Float_t mcpt = ref0->Pt();
 649   Float_t mcsnp = TMath::Sin(TMath::ATan2(ref0->Py(),ref0->Px()));
 650   Float_t mctgl = TMath::Tan(TMath::ATan2(ref0->Pz(),ref0->Pt()));
 651
 652   if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
 653   if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
 654
 655   Float_t deltaPtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
 656   Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
 657
 658   // select primaries
 659   if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
 660   {
 661     if(mcpt == 0) return;
 662
 663     deltaYTPC= track->GetY(); // already rotated
 664     deltaZTPC = track->GetZ()-ref0->Z();
 665     deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(ref0->Pz(),ref0->Pt());
 666     deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(ref0->Py(),ref0->Px());
 667     //delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
 668     deltaPtTPC = (track->Pt()-mcpt) / mcpt;
 669
 670     pullYTPC= track->GetY() / TMath::Sqrt(track->GetSigmaY2());
 671     pullZTPC = (track->GetZ()-ref0->Z()) / TMath::Sqrt(track->GetSigmaZ2());
 672
 673     //Double_t sigma_lambda = 1./(1.+track->GetTgl()*track->GetTgl()) * TMath::Sqrt(track->GetSigmaTgl2());
 674     //Double_t sigma_phi = 1./TMath::Sqrt(1-track->GetSnp()*track->GetSnp()) * TMath::Sqrt(track->GetSigmaSnp2());
 675     pullPhiTPC = (track->GetSnp() - mcsnp) / TMath::Sqrt(track->GetSigmaSnp2());
 676     pullLambdaTPC = (track->GetTgl() - mctgl) / TMath::Sqrt(track->GetSigmaTgl2());
 677
 678     //pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
 679     //pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
 680     if (mcpt) pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
 681     else pull1PtTPC = 0.;
 682
 683     Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,deltaPtTPC,ref0->Y(),ref0->Z(),mcphi,mceta,mcpt};
 684     fResolHisto->Fill(vResolHisto);
 685
 686     Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,ref0->Y(),ref0->Z(),mcsnp,mctgl,1./mcpt};
 687     fPullHisto->Fill(vPullHisto);
 688   }
 689
 690   if(track) delete track;
 691 }
 692
 693 //_____________________________________________________________________________
 694 AliTrackReference * AliPerformanceRes::GetFirstTPCTrackRef(AliMCParticle *mcParticle)
 695 {
 696   // return first TPC track reference
 697   if(!mcParticle) return 0;
 698
 699   // find first track reference
 700   // check direction to select proper reference point for looping tracks
 701   Int_t nTrackRef = mcParticle->GetNumberOfTrackReferences();
 702   AliTrackReference *ref = 0;
 703   AliTrackReference *refIn = 0;
 704   for (Int_t iref = 0; iref < nTrackRef; iref++) {
 705     ref = mcParticle->GetTrackReference(iref);
 706     if(ref && (ref->DetectorId()==AliTrackReference::kTPC))
 707     {
 708       Float_t dir = ref->X()*ref->Px()+ref->Y()*ref->Py();
 709       if(dir < 0.) break;
 710
 711       refIn = ref;
 712       break;
 713     }
 714   }
 715
 716 return refIn;
 717 }
 718
 719 //_____________________________________________________________________________
 720 AliTrackReference * AliPerformanceRes::GetLastTPCTrackRef(AliMCParticle *mcParticle)
 721 {
 722   // return last TPC track reference
 723   if(!mcParticle) return 0;
 724
 725   // find last track reference
 726   // check direction to select proper reference point for looping tracks
 727   Int_t nTrackRef = mcParticle->GetNumberOfTrackReferences();
 728   AliTrackReference *ref = 0;
 729   AliTrackReference *refOut = 0;
 730   for (Int_t iref = 0; iref < nTrackRef; iref++) {
 731     ref = mcParticle->GetTrackReference(iref);
 732     if(ref && (ref->DetectorId()==AliTrackReference::kTPC)) {
 733       Float_t dir=ref->X()*ref->Px()+ref->Y()*ref->Py();
 734       if(dir< 0.0) break;
 735       refOut = ref;
 736     }
 737   }
 738
 739 return refOut;
 740 }
 741
 742 //_____________________________________________________________________________
 743 void AliPerformanceRes::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const esdFriend, const Bool_t bUseMC, const Bool_t bUseESDfriend)
 744 {
 745   // Process comparison information
 746   //
 747   if(!esdEvent)
 748   {
 749     Error("Exec","esdEvent not available");
 750     return;
 751   }
 752   AliHeader* header = 0;
 753   AliGenEventHeader* genHeader = 0;
 754   AliStack* stack = 0;
 755   TArrayF vtxMC(3);
 756
 757   if(bUseMC)
 758   {
 759     if(!mcEvent) {
 760       Error("Exec","mcEvent not available");
 761       return;
 762     }
 763     // get MC event header
 764     header = mcEvent->Header();
 765     if (!header) {
 766       Error("Exec","Header not available");
 767       return;
 768     }
 769     // MC particle stack
 770     stack = mcEvent->Stack();
 771     if (!stack) {
 772       Error("Exec","Stack not available");
 773       return;
 774     }
 775     // get MC vertex
 776     genHeader = header->GenEventHeader();
 777     if (!genHeader) {
 778       Error("Exec","Could not retrieve genHeader from Header");
 779       return;
 780     }
 781     genHeader->PrimaryVertex(vtxMC);
 782   }
 783   else {
 784     Error("Exec","MC information required!");
 785     return;
 786   }
 787
 788   // use ESD friends
 789   if(bUseESDfriend) {
 790     if(!esdFriend) {
 791       Error("Exec","esdFriend not available");
 792       return;
 793     }
 794   }
 795
 796   //  Process events
 797   for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
 798   {
 799     AliESDtrack *track = esdEvent->GetTrack(iTrack);
 800     if(!track) continue;
 801
 802     AliESDfriendTrack *friendTrack=0;
 803     if(bUseESDfriend) {
 804       friendTrack=esdFriend->GetTrack(iTrack);
 805       if(!friendTrack) continue;
 806     }
 807
 808     if(GetAnalysisMode() == 0) ProcessTPC(stack,track);
 809     else if(GetAnalysisMode() == 1) ProcessTPCITS(stack,track);
 810     else if(GetAnalysisMode() == 2) ProcessConstrained(stack,track);
 811     else if(GetAnalysisMode() == 3) ProcessInnerTPC(mcEvent,track);
 812     else if(GetAnalysisMode() == 4) ProcessOuterTPC(mcEvent,track,friendTrack);
 813     else {
 814       printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
 815       return;
 816     }
 817   }
 818 }
 819
 820 //_____________________________________________________________________________
 821 TH1F* AliPerformanceRes::MakeResol(TH2F * his, Int_t integ, Bool_t type, Int_t cut){
 822   // Create resolution histograms
 823
 824   TH1F *hisr, *hism;
 825   if (!gPad) new TCanvas;
 826   hisr = AliTreeDraw::CreateResHistoII(his,&hism,integ,kTRUE,cut);
 827   if (type) return hism;
 828   else
 829     return hisr;
 830 }
 831
 832 //_____________________________________________________________________________
 833 void AliPerformanceRes::Analyse() {
 834   // Analyse comparison information and store output histograms
 835   // in the folder "folderRes"
 836   //
 837   TH1::AddDirectory(kFALSE);
 838   TH1F *h=0;
 839   TH2F *h2D=0;
 840   TObjArray *aFolderObj = new TObjArray;
 841
 842   // write results in the folder
 843   TCanvas * c = new TCanvas("Phi resol Tan","Phi resol Tan");
 844   c->cd();
 845
 846   char name[256];
 847   char title[256];
 848
 849   for(Int_t i=0; i<5; i++)
 850   {
 851     for(Int_t j=5; j<10; j++)
 852     {
 853       //if(j!=8) fResolHisto->GetAxis(8)->SetRangeUser(-0.9,0.89); // eta window
 854       if(j!=8) fResolHisto->GetAxis(8)->SetRangeUser(0.0,0.89); // eta window
 855       else fResolHisto->GetAxis(8)->SetRangeUser(-1.5,1.49);
 856       fResolHisto->GetAxis(9)->SetRangeUser(0.16,100.); // pt threshold
 857       if(GetAnalysisMode() == 3) fResolHisto->GetAxis(5)->SetRangeUser(-80.,80.); // y range
 858
 859       h2D = (TH2F*)fResolHisto->Projection(i,j);
 860
 861       h = AliPerformanceRes::MakeResol(h2D,1,0,100);
 862       sprintf(name,"h_res_%d_vs_%d",i,j);
 863       h->SetName(name);
 864
 865       h->GetXaxis()->SetTitle(fResolHisto->GetAxis(j)->GetTitle());
 866       sprintf(title,"%s %s",fResolHisto->GetAxis(i)->GetTitle(),"(resolution)");
 867       h->GetYaxis()->SetTitle(title);
 868       sprintf(title,"%s vs %s",title,fResolHisto->GetAxis(j)->GetTitle());
 869       h->SetTitle(title);
 870
 871       if(j==9) h->SetBit(TH1::kLogX);
 872       aFolderObj->Add(h);
 873
 874       h = AliPerformanceRes::MakeResol(h2D,1,1,100);
 875       //h = (TH1F*)arr->At(1);
 876       sprintf(name,"h_mean_res_%d_vs_%d",i,j);
 877       h->SetName(name);
 878
 879       h->GetXaxis()->SetTitle(fResolHisto->GetAxis(j)->GetTitle());
 880       sprintf(title,"%s %s",fResolHisto->GetAxis(i)->GetTitle(),"(mean)");
 881       h->GetYaxis()->SetTitle(title);
 882
 883       sprintf(title,"%s vs %s",title,fResolHisto->GetAxis(j)->GetTitle());
 884       h->SetTitle(title);
 885
 886       if(j==9) h->SetBit(TH1::kLogX);
 887       aFolderObj->Add(h);
 888
 889       fResolHisto->GetAxis(8)->SetRangeUser(-1.5,1.5);
 890       fResolHisto->GetAxis(9)->SetRangeUser(0.0,10.);
 891
 892       //
 893       //if(j!=8) fPullHisto->GetAxis(8)->SetRangeUser(-0.9,0.89); // eta window
 894       if(j!=8) fPullHisto->GetAxis(8)->SetRangeUser(0.0,0.89);    // eta window
 895       else  fPullHisto->GetAxis(8)->SetRangeUser(-1.5,1.49);      // eta window
 896       fPullHisto->GetAxis(9)->SetRangeUser(0.16,100.);            // pt threshold
 897       if(GetAnalysisMode() == 3) fPullHisto->GetAxis(5)->SetRangeUser(-80.,80.); // y range
 898
 899       h2D = (TH2F*)fPullHisto->Projection(i,j);
 900
 901       h = AliPerformanceRes::MakeResol(h2D,1,0,100);
 902       sprintf(name,"h_pull_%d_vs_%d",i,j);
 903       h->SetName(name);
 904
 905       h->GetXaxis()->SetTitle(fPullHisto->GetAxis(j)->GetTitle());
 906       sprintf(title,"%s %s",fPullHisto->GetAxis(i)->GetTitle(),"(resolution)");
 907       h->GetYaxis()->SetTitle(title);
 908       sprintf(title,"%s vs %s",title,fPullHisto->GetAxis(j)->GetTitle());
 909       h->SetTitle(title);
 910
 911       //if(j==9) h->SetBit(TH1::kLogX);
 912       aFolderObj->Add(h);
 913
 914       h = AliPerformanceRes::MakeResol(h2D,1,1,100);
 915       sprintf(name,"h_mean_pull_%d_vs_%d",i,j);
 916       h->SetName(name);
 917
 918       h->GetXaxis()->SetTitle(fPullHisto->GetAxis(j)->GetTitle());
 919       sprintf(title,"%s %s",fPullHisto->GetAxis(i)->GetTitle(),"(mean)");
 920       h->GetYaxis()->SetTitle(title);
 921       sprintf(title,"%s vs %s",title,fPullHisto->GetAxis(j)->GetTitle());
 922       h->SetTitle(title);
 923
 924       //if(j==9) h->SetBit(TH1::kLogX);
 925       aFolderObj->Add(h);
 926     }
 927   }
 928
 929   // export objects to analysis folder
 930   fAnalysisFolder = ExportToFolder(aFolderObj);
 931
 932   // delete only TObjArray
 933   if(aFolderObj) delete aFolderObj;
 934 }
 935
 936 //_____________________________________________________________________________
 937 TFolder* AliPerformanceRes::ExportToFolder(TObjArray * array)
 938 {
 939   // recreate folder avery time and export objects to new one
 940   //
 941   AliPerformanceRes * comp=this;
 942   TFolder *folder = comp->GetAnalysisFolder();
 943
 944   TString name, title;
 945   TFolder *newFolder = 0;
 946   Int_t i = 0;
 947   Int_t size = array->GetSize();
 948
 949   if(folder) {
 950      // get name and title from old folder
 951      name = folder->GetName();
 952      title = folder->GetTitle();
 953
 954          // delete old one
 955      delete folder;
 956
 957          // create new one
 958      newFolder = CreateFolder(name.Data(),title.Data());
 959      newFolder->SetOwner();
 960
 961          // add objects to folder
 962      while(i < size) {
 963            newFolder->Add(array->At(i));
 964            i++;
 965          }
 966   }
 967
 968 return newFolder;
 969 }
 970
 971 //_____________________________________________________________________________
 972 Long64_t AliPerformanceRes::Merge(TCollection* const list)
 973 {
 974   // Merge list of objects (needed by PROOF)
 975
 976   if (!list)
 977   return 0;
 978
 979   if (list->IsEmpty())
 980   return 1;
 981
 982   TIterator* iter = list->MakeIterator();
 983   TObject* obj = 0;
 984
 985   // collection of generated histograms
 986   Int_t count=0;
 987   while((obj = iter->Next()) != 0)
 988   {
 989   AliPerformanceRes* entry = dynamic_cast<AliPerformanceRes*>(obj);
 990   if (entry == 0) continue;
 991
 992   fResolHisto->Add(entry->fResolHisto);
 993   fPullHisto->Add(entry->fPullHisto);
 994
 995   count++;
 996   }
 997
 998 return count;
 999 }
1000
1001 //_____________________________________________________________________________
1002 TFolder* AliPerformanceRes::CreateFolder(TString name,TString title) {
1003 // create folder for analysed histograms
1004 //
1005 TFolder *folder = 0;
1006   folder = new TFolder(name.Data(),title.Data());
1007
1008   return folder;
1009 }