PWGLF/SPECTRA/PiKaPr/TestAOD/AliAnalysisTaskV2AllChAOD.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2009, 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-commercial 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 //-----------------------------------------------------------------
  17 //         AliAnalysisTaskV2AllChAOD class
  18 //-----------------------------------------------------------------
  19
  20 #include "TChain.h"
  21 #include "TTree.h"
  22 #include "TLegend.h"
  23 #include "TH1F.h"
  24 #include "TH1D.h"
  25 #include "TH2F.h"
  26 #include "THnSparse.h"
  27 #include "TProfile.h"
  28 #include "TCanvas.h"
  29 #include "AliAnalysisTask.h"
  30 #include "AliAODTrack.h"
  31 #include "AliAODMCParticle.h"
  32 #include "AliVParticle.h"
  33 #include "AliAODEvent.h"
  34 #include "AliAODInputHandler.h"
  35 #include "AliAnalysisTaskV2AllChAOD.h"
  36 #include "AliAnalysisTaskESDfilter.h"
  37 #include "AliAnalysisDataContainer.h"
  38 #include "AliSpectraAODTrackCuts.h"
  39 #include "AliSpectraAODEventCuts.h"
  40 #include "AliPIDCombined.h"
  41 #include "AliCentrality.h"
  42 #include "TProof.h"
  43 #include "AliVEvent.h"
  44 #include "AliStack.h"
  45 #include <TMCProcess.h>
  46 #include <TRandom.h>
  47
  48 #include <iostream>
  49
  50 using namespace std;
  51
  52 ClassImp(AliAnalysisTaskV2AllChAOD)
  53
  54 //________________________________________________________________________
  55 AliAnalysisTaskV2AllChAOD::AliAnalysisTaskV2AllChAOD(const char *name) : AliAnalysisTaskSE(name),
  56   fAOD(0x0),
  57   fTrackCuts(0x0),
  58   fEventCuts(0x0),
  59   fIsMC(0),
  60   fCharge(0),
  61   fVZEROside(0),
  62   fOutput(0x0),
  63   fOutput_lq(0x0),
  64   fOutput_sq(0x0),
  65   fnCentBins(20),
  66   fnQvecBins(100),
  67   fQvecUpperLim(100),
  68   fCutLargeQperc(90.),
  69   fCutSmallQperc(10.),
  70   fEtaGapMin(-0.5),
  71   fEtaGapMax(0.5),
  72   fTrkBit(128),
  73   fEtaCut(0.8),
  74   fMinPt(0),
  75   fMaxPt(20.0),
  76   fMinTPCNcls(70),
  77   fFillTHn(kFALSE),
  78   fCentrality(0),
  79   fResSP(0),
  80   fResSP_vs_Cent(0),
  81   f2partCumQA_vs_Cent(0),
  82   f2partCumQB_vs_Cent(0),
  83   fEta_vs_Phi_bef(0),
  84   fEta_vs_PhiA(0),
  85   fEta_vs_PhiB(0),
  86   fResSP_lq(0),
  87   fResSP_vs_Cent_lq(0),
  88   f2partCumQA_vs_Cent_lq(0),
  89   f2partCumQB_vs_Cent_lq(0),
  90   fResSP_sq(0),
  91   fResSP_vs_Cent_sq(0),
  92   f2partCumQA_vs_Cent_sq(0),
  93   f2partCumQB_vs_Cent_sq(0),
  94   fResSP_inclusive(0),
  95   fv2SPGap1A_inclusive_mb(0),
  96   fv2SPGap1B_inclusive_mb(0),
  97   fv2SPGap1A_inclusive_lq(0),
  98   fv2SPGap1B_inclusive_lq(0),
  99   fv2SPGap1A_inclusive_sq(0),
 100   fv2SPGap1B_inclusive_sq(0),
 101   fResSPmc_inclusive(0),
 102   fv2SPGap1Amc_inclusive_mb(0),
 103   fv2SPGap1Bmc_inclusive_mb(0),
 104   fv2SPGap1Amc_inclusive_lq(0),
 105   fv2SPGap1Bmc_inclusive_lq(0),
 106   fv2SPGap1Amc_inclusive_sq(0),
 107   fv2SPGap1Bmc_inclusive_sq(0),
 108   fIsRecoEff(0),
 109   fRecoEffList(0),
 110   fQvecGen(0),
 111   fQgenType(0),
 112   fnNchBins(400),
 113   fDoCentrSystCentrality(0)
 114 {
 115
 116   for (Int_t i = 0; i< 9; i++){
 117
 118     fResSP_vs_Qvec[i] = 0;
 119
 120     fv2SPGap1A[i] = 0;
 121     fv2SPGap1B[i] = 0;
 122
 123     fSinGap1Aq[i] = 0;
 124     fCosGap1Aq[i] = 0;
 125     fSinGap1Bq[i] = 0;
 126     fCosGap1Bq[i] = 0;
 127
 128     fSinGap1A[i] = 0;
 129     fCosGap1A[i] = 0;
 130     fSinGap1B[i] = 0;
 131     fCosGap1B[i] = 0;
 132
 133     //large q
 134     fv2SPGap1A_lq[i] = 0;
 135     fv2SPGap1B_lq[i] = 0;
 136
 137     fSinGap1Aq_lq[i] = 0;
 138     fCosGap1Aq_lq[i] = 0;
 139     fSinGap1Bq_lq[i] = 0;
 140     fCosGap1Bq_lq[i] = 0;
 141
 142     fSinGap1A_lq[i] = 0;
 143     fCosGap1A_lq[i] = 0;
 144     fSinGap1B_lq[i] = 0;
 145     fCosGap1B_lq[i] = 0;
 146
 147     //small q
 148     fv2SPGap1A_sq[i] = 0;
 149     fv2SPGap1B_sq[i] = 0;
 150
 151     fSinGap1Aq_sq[i] = 0;
 152     fCosGap1Aq_sq[i] = 0;
 153     fSinGap1Bq_sq[i] = 0;
 154     fCosGap1Bq_sq[i] = 0;
 155
 156     fSinGap1A_sq[i] = 0;
 157     fCosGap1A_sq[i] = 0;
 158     fSinGap1B_sq[i] = 0;
 159     fCosGap1B_sq[i] = 0;
 160
 161   }
 162
 163   fRecoEffList=new TList();
 164   fRecoEffList->SetOwner();
 165   fRecoEffList->SetName("fRecoEffList");
 166
 167   // Default constructor
 168   DefineInput(0, TChain::Class());
 169   DefineOutput(1, TList::Class());
 170   DefineOutput(2, AliSpectraAODEventCuts::Class());
 171   DefineOutput(3, AliSpectraAODTrackCuts::Class());
 172   DefineOutput(4, TList::Class());
 173   DefineOutput(5, TList::Class());
 174 }
 175
 176 //________________________________________________________________________
 177 void AliAnalysisTaskV2AllChAOD::UserCreateOutputObjects()
 178 {
 179   // create output objects
 180   fOutput=new TList();
 181   fOutput->SetOwner();
 182   fOutput->SetName("fOutput");
 183
 184   fOutput_lq=new TList();
 185   fOutput_lq->SetOwner();
 186   fOutput_lq->SetName("fOutput_lq");
 187
 188   fOutput_sq=new TList();
 189   fOutput_sq->SetOwner();
 190   fOutput_sq->SetName("fOutput_sq");
 191
 192   if (!fTrackCuts) AliFatal("Track Cuts should be set in the steering macro");
 193   if (!fEventCuts) AliFatal("Event Cuts should be set in the steering macro");
 194
 195   if( fFillTHn ){
 196     //dimensions of THnSparse for Q vector checks
 197     const Int_t nvarev=6;
 198     //                                             cent         q-rec_perc        qvec-rec      q-gen_tracks   qvec-gen_vzero          Nch
 199     Int_t    binsHistRealEv[nvarev] = {     fnCentBins,              100,        fnQvecBins,     fnQvecBins,       fnQvecBins,     fnNchBins};
 200     Double_t xminHistRealEv[nvarev] = {             0.,               0.,                0.,             0.,               0.,            0.};
 201     Double_t xmaxHistRealEv[nvarev] = {           100.,             100.,     fQvecUpperLim,  fQvecUpperLim,    fQvecUpperLim,         2000.};
 202
 203     THnSparseF* NSparseHistEv = new THnSparseF("NSparseHistEv","NSparseHistEv",nvarev,binsHistRealEv,xminHistRealEv,xmaxHistRealEv);
 204     NSparseHistEv->GetAxis(0)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
 205     NSparseHistEv->GetAxis(0)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
 206
 207     NSparseHistEv->GetAxis(1)->SetTitle("q-vec rec percentile");
 208     NSparseHistEv->GetAxis(1)->SetName("Qrec_perc");
 209
 210     NSparseHistEv->GetAxis(2)->SetTitle("q-vec rec");
 211     NSparseHistEv->GetAxis(2)->SetName("Qrec");
 212
 213     NSparseHistEv->GetAxis(3)->SetTitle("q-vec gen tracks");
 214     NSparseHistEv->GetAxis(3)->SetName("Qgen_tracks");
 215
 216     NSparseHistEv->GetAxis(4)->SetTitle("q-vec gen vzero");
 217     NSparseHistEv->GetAxis(4)->SetName("Qgen_vzero");
 218
 219     NSparseHistEv->GetAxis(5)->SetTitle("Ncharged");
 220     NSparseHistEv->GetAxis(5)->SetName("Nch");
 221     fOutput->Add(NSparseHistEv);
 222   }
 223
 224   fCentrality = new TH1D("fCentrality", "centrality distribution; centrality", 200, 0., 100);
 225   fOutput->Add(fCentrality);
 226
 227   // binning common to all the THn
 228   //change it according to your needs + move it to global variables -> setter/getter
 229 //   Double_t ptBins[] = {0., 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.4, 2.8, 3.2, 3.6, 4.0, 4.5, 5.0, 5.5, 6.0, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 20.0};
 230 //   const Int_t nptBins = 31;
 231   Double_t ptBins[] = {0., 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.2, 2.4, 2.6, 3.0, 3.4, 3.8, 4.2, 4.6, 5.0, 5.5, 6.0, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 20.0};
 232   const Int_t nptBins = 33;
 233
 234   fResSP = new TProfile("fResSP", "Resolution; centrality; Resolution", 9, -0.5, 8.5);
 235   fOutput->Add(fResSP);
 236
 237   fResSP_vs_Cent = new TProfile("fResSP_vs_Cent", "Resolution; centrality; Resolution", 20., 0., 100.);
 238   fOutput->Add(fResSP_vs_Cent);
 239
 240   f2partCumQA_vs_Cent = new TProfile("f2partCumQA_vs_Cent", "Resolution; centrality; Resolution", 100., 0., 100.);
 241   fOutput->Add(f2partCumQA_vs_Cent);
 242
 243   f2partCumQB_vs_Cent = new TProfile("f2partCumQB_vs_Cent", "Resolution; centrality; Resolution", 100., 0., 100.);
 244   fOutput->Add(f2partCumQB_vs_Cent);
 245
 246   fEta_vs_Phi_bef = new TH2D("fEta_vs_Phi_bef","eta vs phi distribution before eta gap;#eta;#phi",200.,-1.,1.,175.,0.,7.);
 247   fOutput->Add(fEta_vs_Phi_bef);
 248
 249   fEta_vs_PhiA = new TH2D("fEta_vs_PhiA","eta vs phi distribution;#eta;#phi",200.,-1.,1.,175.,0.,7.);
 250   fOutput->Add(fEta_vs_PhiA);
 251
 252   fEta_vs_PhiB = new TH2D("fEta_vs_PhiB","eta vs phi distribution;#eta;#phi",200.,-1.,1.,175.,0.,7.);
 253   fOutput->Add(fEta_vs_PhiB);
 254
 255     // MC closure test
 256     fResSP_inclusive = new TProfile("fResSP_inclusive", "Resolution; ese; Resolution", 3, 0., 3.);
 257     fOutput->Add(fResSP_inclusive);
 258
 259     fv2SPGap1A_inclusive_mb = new TProfile("fv2SPGap1A_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 260     fOutput->Add(fv2SPGap1A_inclusive_mb);
 261
 262     fv2SPGap1B_inclusive_mb = new TProfile("fv2SPGap1B_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 263     fOutput->Add(fv2SPGap1B_inclusive_mb);
 264
 265
 266   //large q
 267   fResSP_lq = new TProfile("fResSP_lq", "Resolution; centrality; Resolution", 9, -0.5, 8.5);
 268   fOutput_lq->Add(fResSP_lq);
 269
 270   fResSP_vs_Cent_lq = new TProfile("fResSP_vs_Cent_lq", "Resolution; centrality; Resolution", 20., 0., 100.);
 271   fOutput_lq->Add(fResSP_vs_Cent_lq);
 272
 273   f2partCumQA_vs_Cent_lq = new TProfile("f2partCumQA_vs_Cent_lq", "Resolution; centrality; Resolution", 100., 0., 100.);
 274   fOutput_lq->Add(f2partCumQA_vs_Cent_lq);
 275
 276   f2partCumQB_vs_Cent_lq = new TProfile("f2partCumQB_vs_Cent_lq", "Resolution; centrality; Resolution", 100., 0., 100.);
 277   fOutput_lq->Add(f2partCumQB_vs_Cent_lq);
 278
 279     // MC closure test
 280     fv2SPGap1A_inclusive_lq = new TProfile("fv2SPGap1A_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 281     fOutput_lq->Add(fv2SPGap1A_inclusive_lq);
 282
 283     fv2SPGap1B_inclusive_lq = new TProfile("fv2SPGap1B_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 284     fOutput_lq->Add(fv2SPGap1B_inclusive_lq);
 285
 286   //small q resolution
 287   fResSP_sq = new TProfile("fResSP_sq", "Resolution; centrality; Resolution", 9, -0.5, 8.5);
 288   fOutput_sq->Add(fResSP_sq);
 289
 290   fResSP_vs_Cent_sq = new TProfile("fResSP_vs_Cent_sq", "Resolution; centrality; Resolution", 20., 0., 100.);
 291   fOutput_sq->Add(fResSP_vs_Cent_sq);
 292
 293   f2partCumQA_vs_Cent_sq = new TProfile("f2partCumQA_vs_Cent_sq", "Resolution; centrality; Resolution", 100., 0., 100.);
 294   fOutput_sq->Add(f2partCumQA_vs_Cent_sq);
 295
 296   f2partCumQB_vs_Cent_sq = new TProfile("f2partCumQB_vs_Cent_sq", "Resolution; centrality; Resolution", 100., 0., 100.);
 297   fOutput_sq->Add(f2partCumQB_vs_Cent_sq);
 298
 299     // MC closure test
 300     fv2SPGap1A_inclusive_sq = new TProfile("fv2SPGap1A_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 301     fOutput_sq->Add(fv2SPGap1A_inclusive_sq);
 302
 303     fv2SPGap1B_inclusive_sq = new TProfile("fv2SPGap1B_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 304     fOutput_sq->Add(fv2SPGap1B_inclusive_sq);
 305
 306   for (Int_t iC = 0; iC < 9; iC++){
 307
 308     fResSP_vs_Qvec[iC] = new TProfile(Form("fResSP_vs_Qvec_%d", iC), "Resolution; Qvec (V0A); Resolution", 10., 0., 100.);
 309     fOutput->Add(fResSP_vs_Qvec[iC]);
 310
 311     fv2SPGap1A[iC] = new TProfile(Form("fv2SPGap1A_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 312     fOutput->Add(fv2SPGap1A[iC]);
 313
 314     fv2SPGap1B[iC] = new TProfile(Form("fv2SPGap1B_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 315     fOutput->Add(fv2SPGap1B[iC]);
 316
 317     fSinGap1Aq[iC] = new TProfile(Form("fSinGap1Aq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 318     fOutput->Add(fSinGap1Aq[iC]);
 319
 320     fCosGap1Aq[iC] = new TProfile(Form("fCosGap1Aq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 321     fOutput->Add(fCosGap1Aq[iC]);
 322
 323     fSinGap1Bq[iC] = new TProfile(Form("fSinGap1Bq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 324     fOutput->Add(fSinGap1Bq[iC]);
 325
 326     fCosGap1Bq[iC] = new TProfile(Form("fCosGap1Bq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 327     fOutput->Add(fCosGap1Bq[iC]);
 328
 329     fSinGap1A[iC] = new TProfile(Form("fSinGap1A_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 330     fOutput->Add(fSinGap1A[iC]);
 331
 332     fCosGap1A[iC] = new TProfile(Form("fCosGap1A_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 333     fOutput->Add(fCosGap1A[iC]);
 334
 335     fSinGap1B[iC] = new TProfile(Form("fSinGap1B_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 336     fOutput->Add(fSinGap1B[iC]);
 337
 338     fCosGap1B[iC] = new TProfile(Form("fCosGap1B_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 339     fOutput->Add(fCosGap1B[iC]);
 340
 341     //large q
 342     fv2SPGap1A_lq[iC] = new TProfile(Form("fv2SPGap1A_lq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 343     fOutput_lq->Add(fv2SPGap1A_lq[iC]);
 344
 345     fv2SPGap1B_lq[iC] = new TProfile(Form("fv2SPGap1B_lq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 346     fOutput_lq->Add(fv2SPGap1B_lq[iC]);
 347
 348     fSinGap1Aq_lq[iC] = new TProfile(Form("fSinGap1Aq_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 349     fOutput_lq->Add(fSinGap1Aq_lq[iC]);
 350
 351     fCosGap1Aq_lq[iC] = new TProfile(Form("fCosGap1Aq_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 352     fOutput_lq->Add(fCosGap1Aq_lq[iC]);
 353
 354     fSinGap1Bq_lq[iC] = new TProfile(Form("fSinGap1Bq_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 355     fOutput_lq->Add(fSinGap1Bq_lq[iC]);
 356
 357     fCosGap1Bq_lq[iC] = new TProfile(Form("fCosGap1Bq_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 358     fOutput_lq->Add(fCosGap1Bq_lq[iC]);
 359
 360     fSinGap1A_lq[iC] = new TProfile(Form("fSinGap1A_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 361     fOutput_lq->Add(fSinGap1A_lq[iC]);
 362
 363     fCosGap1A_lq[iC] = new TProfile(Form("fCosGap1A_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 364     fOutput_lq->Add(fCosGap1A_lq[iC]);
 365
 366     fSinGap1B_lq[iC] = new TProfile(Form("fSinGap1B_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 367     fOutput_lq->Add(fSinGap1B_lq[iC]);
 368
 369     fCosGap1B_lq[iC] = new TProfile(Form("fCosGap1B_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 370     fOutput_lq->Add(fCosGap1B_lq[iC]);
 371
 372     //small q
 373     fv2SPGap1A_sq[iC] = new TProfile(Form("fv2SPGap1A_sq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 374     fOutput_sq->Add(fv2SPGap1A_sq[iC]);
 375
 376     fv2SPGap1B_sq[iC] = new TProfile(Form("fv2SPGap1B_sq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 377     fOutput_sq->Add(fv2SPGap1B_sq[iC]);
 378
 379     fSinGap1Aq_sq[iC] = new TProfile(Form("fSinGap1Aq_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 380     fOutput_sq->Add(fSinGap1Aq_sq[iC]);
 381
 382     fCosGap1Aq_sq[iC] = new TProfile(Form("fCosGap1Aq_sq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 383     fOutput_sq->Add(fCosGap1Aq_sq[iC]);
 384
 385     fSinGap1Bq_sq[iC] = new TProfile(Form("fSinGap1Bq_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 386     fOutput_sq->Add(fSinGap1Bq_sq[iC]);
 387
 388     fCosGap1Bq_sq[iC] = new TProfile(Form("fCosGap1Bq_sq_%d", iC), "p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 389     fOutput_sq->Add(fCosGap1Bq_sq[iC]);
 390
 391     fSinGap1A_sq[iC] = new TProfile(Form("fSinGap1A_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 392     fOutput_sq->Add(fSinGap1A_sq[iC]);
 393
 394     fCosGap1A_sq[iC] = new TProfile(Form("fCosGap1A_sq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 395     fOutput_sq->Add(fCosGap1A_sq[iC]);
 396
 397     fSinGap1B_sq[iC] = new TProfile(Form("fSinGap1B_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
 398     fOutput_sq->Add(fSinGap1B_sq[iC]);
 399
 400     fCosGap1B_sq[iC] = new TProfile(Form("fCosGap1B_sq_%d", iC), "p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
 401     fOutput_sq->Add(fCosGap1B_sq[iC]);
 402   };
 403
 404   if(fIsMC){
 405     fResSPmc_inclusive = new TProfile("fResSPmc_inclusive", "Resolution; ese; Resolution", 3, 0., 3.);
 406     fOutput->Add(fResSPmc_inclusive);
 407
 408     fv2SPGap1Amc_inclusive_mb = new TProfile("fv2SPGap1Amc_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 409     fOutput->Add(fv2SPGap1Amc_inclusive_mb);
 410
 411     fv2SPGap1Bmc_inclusive_mb = new TProfile("fv2SPGap1Bmc_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 412     fOutput->Add(fv2SPGap1Bmc_inclusive_mb);
 413
 414     //large-q
 415     fv2SPGap1Amc_inclusive_lq = new TProfile("fv2SPGap1Amc_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 416     fOutput_lq->Add(fv2SPGap1Amc_inclusive_lq);
 417
 418     fv2SPGap1Bmc_inclusive_lq = new TProfile("fv2SPGap1Bmc_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 419     fOutput_lq->Add(fv2SPGap1Bmc_inclusive_lq);
 420
 421     //small-q
 422     fv2SPGap1Amc_inclusive_sq = new TProfile("fv2SPGap1Amc_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 423     fOutput_sq->Add(fv2SPGap1Amc_inclusive_sq);
 424
 425     fv2SPGap1Bmc_inclusive_sq = new TProfile("fv2SPGap1Bmc_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
 426     fOutput_sq->Add(fv2SPGap1Bmc_inclusive_sq);
 427   }
 428
 429
 430   PostData(1, fOutput  );
 431   PostData(2, fEventCuts);
 432   PostData(3, fTrackCuts);
 433   PostData(4, fOutput_lq  );
 434   PostData(5, fOutput_sq  );
 435 }
 436
 437 //________________________________________________________________________
 438
 439 void AliAnalysisTaskV2AllChAOD::UserExec(Option_t *)
 440 {
 441   //Printf("An event");
 442   // main event loop
 443   fAOD = dynamic_cast<AliAODEvent*>(fInputEvent);
 444   if (!fAOD) {
 445     AliWarning("ERROR: AliAODEvent not available \n");
 446     return;
 447   }
 448
 449   if (strcmp(fAOD->ClassName(), "AliAODEvent"))
 450     {
 451       AliFatal("Not processing AODs");
 452     }
 453
 454   if(!fEventCuts->IsSelected(fAOD,fTrackCuts))return;//event selection
 455
 456   //Get q-vector percentile.
 457   Double_t Qvec=0.;
 458   if(fIsMC && fQvecGen) Qvec = fEventCuts->GetQvecPercentileMC(fVZEROside, fQgenType);
 459   else Qvec = fEventCuts->GetQvecPercentile(fVZEROside);
 460
 461
 462   Double_t Cent=(fDoCentrSystCentrality)?1.01*fEventCuts->GetCent():fEventCuts->GetCent();
 463   fCentrality->Fill(Cent);
 464
 465   Int_t centV0 = -1;
 466   if ((Cent > 0) && (Cent <= 5.0))
 467       centV0 = 0;
 468     else if ((Cent > 5.0) && (Cent <= 10.0))
 469       centV0 = 1;
 470     else if ((Cent > 10.0) && (Cent <= 20.0))
 471       centV0 = 2;
 472     else if ((Cent > 20.0) && (Cent <= 30.0))
 473       centV0 = 3;
 474     else if ((Cent > 30.0) && (Cent <= 40.0))
 475       centV0 = 4;
 476     else if ((Cent > 40.0) && (Cent <= 50.0))
 477       centV0 = 5;
 478     else if ((Cent > 50.0) && (Cent <= 60.0))
 479       centV0 = 6;
 480     else if ((Cent > 60.0) && (Cent <= 70.0))
 481       centV0 = 7;
 482     else if ((Cent > 70.0) && (Cent <= 80.0))
 483       centV0 = 8;
 484
 485   if(fIsMC) MCclosure(Qvec); // fill mc histograms for montecarlo closure
 486
 487   Double_t QxGap1A = 0., QyGap1A = 0.;
 488   Double_t QxGap1B = 0., QyGap1B = 0.;
 489   Int_t multGap1A = 0, multGap1B = 0;
 490
 491   for (Int_t loop = 0; loop < 2; loop++){
 492
 493     //main loop on tracks
 494     for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
 495       AliAODTrack* track = fAOD->GetTrack(iTracks);
 496       if(fCharge != 0 && track->Charge() != fCharge) continue;//if fCharge != 0 only select fCharge
 497       if (!fTrackCuts->IsSelected(track,kTRUE)) continue; //track selection (rapidity selection NOT in the standard cuts)
 498
 499       fEta_vs_Phi_bef->Fill( track->Eta(), track->Phi() );
 500
 501       if (fIsRecoEff){
 502
 503         // 2) reject randomly tracks at high pT until the reconstruction efficiency becomes flat (add the following before the loop == 0 part): (mail by Alexandru)
 504
 505         Double_t recoEff = GetRecoEff(track->Pt(), centV0);
 506         if (recoEff < 0){
 507           cout<<"No reconstruction efficiency!"<<endl;
 508           continue;
 509         }
 510
 511         Double_t rndPt = gRandom->Rndm();
 512 //         cout<<"rndPt: "<<rndPt<<endl;
 513
 514         Double_t minRecPt = GetRecoEff(0.200001, centV0);
 515 //         cout<<"minRecPt: "<<minRecPt<<endl;
 516
 517         if (rndPt > minRecPt/recoEff){
 518 //        cout<<"Track rejected: "<<iTracks<<"  from "<<fAOD->GetNumberOfTracks()<<endl;
 519           continue;
 520         }
 521
 522       } // end fIsRecoEff
 523
 524       if (loop == 0) {
 525
 526         if (track->Eta() > fEtaGapMax){
 527           QxGap1A += TMath::Cos(2.*track->Phi());
 528           QyGap1A += TMath::Sin(2.*track->Phi());
 529           multGap1A++;
 530
 531           fSinGap1Aq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 532           fCosGap1Aq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 533
 534           fEta_vs_PhiA->Fill( track->Eta(), track->Phi() );
 535
 536           if (Qvec > fCutLargeQperc && Qvec < 100.){
 537             fSinGap1Aq_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 538             fCosGap1Aq_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 539           }
 540
 541           if (Qvec > 0. && Qvec < fCutSmallQperc){
 542             fSinGap1Aq_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 543             fCosGap1Aq_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 544           }
 545
 546         }
 547
 548       if (track->Eta() < fEtaGapMin){
 549         QxGap1B += TMath::Cos(2.*track->Phi());
 550         QyGap1B += TMath::Sin(2.*track->Phi());
 551         multGap1B++;
 552
 553         fCosGap1Bq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 554         fSinGap1Bq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 555
 556         fEta_vs_PhiB->Fill( track->Eta(), track->Phi() );
 557
 558         if (Qvec > fCutLargeQperc && Qvec < 100.){
 559           fSinGap1Bq_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 560           fCosGap1Bq_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 561         }
 562
 563         if (Qvec > 0. && Qvec < fCutSmallQperc){
 564           fSinGap1Bq_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 565           fCosGap1Bq_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 566         }
 567
 568       }
 569
 570     } else {
 571
 572         //eval v2 scalar product
 573         if (track->Eta() < fEtaGapMin && multGap1A > 0){
 574           Double_t v2SPGap1A = (TMath::Cos(2.*track->Phi())*QxGap1A + TMath::Sin(2.*track->Phi())*QyGap1A)/(Double_t)multGap1A;
 575           fv2SPGap1A[centV0]->Fill(track->Pt(), v2SPGap1A);
 576
 577           fv2SPGap1A_inclusive_mb->Fill(track->Pt(), v2SPGap1A); //mb v2 for mc closure
 578
 579           fSinGap1A[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 580           fCosGap1A[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 581
 582           if (Qvec > fCutLargeQperc && Qvec < 100.){
 583             fv2SPGap1A_lq[centV0]->Fill(track->Pt(), v2SPGap1A);
 584             fSinGap1A_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 585             fCosGap1A_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 586
 587             fv2SPGap1A_inclusive_lq->Fill(track->Pt(), v2SPGap1A); //lq v2 for mc closure
 588           }
 589
 590           if (Qvec > 0. && Qvec < fCutSmallQperc){
 591             fv2SPGap1A_sq[centV0]->Fill(track->Pt(), v2SPGap1A);
 592             fSinGap1A_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 593             fCosGap1A_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 594
 595             fv2SPGap1A_inclusive_sq->Fill(track->Pt(), v2SPGap1A); //sq v2 for mc closure
 596           }
 597
 598         }
 599
 600         if (track->Eta() > fEtaGapMax && multGap1B > 0){
 601           Double_t v2SPGap1B = (TMath::Cos(2.*track->Phi())*QxGap1B + TMath::Sin(2.*track->Phi())*QyGap1B)/(Double_t)multGap1B;
 602           fv2SPGap1B[centV0]->Fill(track->Pt(), v2SPGap1B);
 603
 604           fv2SPGap1B_inclusive_mb->Fill(track->Pt(), v2SPGap1B); //mb v2 for mc closure
 605
 606           fCosGap1B[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 607           fSinGap1B[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 608
 609           if (Qvec > fCutLargeQperc && Qvec < 100.){
 610             fv2SPGap1B_lq[centV0]->Fill(track->Pt(), v2SPGap1B);
 611             fSinGap1B_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 612             fCosGap1B_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 613
 614             fv2SPGap1B_inclusive_lq->Fill(track->Pt(), v2SPGap1B); //lq v2 for mc closure
 615           }
 616
 617           if (Qvec > 0. && Qvec < fCutSmallQperc){
 618             fv2SPGap1B_sq[centV0]->Fill(track->Pt(), v2SPGap1B);
 619             fSinGap1B_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
 620             fCosGap1B_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
 621
 622             fv2SPGap1B_inclusive_sq->Fill(track->Pt(), v2SPGap1B); //sq v2 for mc closure
 623           }
 624
 625         }
 626       }// end else
 627     } // end loop on tracks
 628   } // end loop
 629
 630
 631   if (multGap1A > 0 && multGap1B > 0){
 632     Double_t res = (QxGap1A*QxGap1B + QyGap1A*QyGap1B)/(Double_t)multGap1A/(Double_t)multGap1B;
 633     fResSP->Fill((Double_t)centV0, res);
 634
 635     fResSP_inclusive->Fill(0., res); //mb v2 for mc closure
 636     fResSP_vs_Cent->Fill(Cent, res);
 637     fResSP_vs_Qvec[centV0]->Fill(Qvec,res);
 638
 639     Double_t f2partCumQA = -999.;
 640     if(multGap1A>1)
 641       f2partCumQA = ( ( (QxGap1A*QxGap1A + QyGap1A*QyGap1A) - (Double_t)multGap1A ) / ((Double_t)multGap1A*((Double_t)multGap1A-1)) );
 642     if(f2partCumQA>0)f2partCumQA_vs_Cent->Fill((Double_t)Cent,f2partCumQA);
 643
 644     Double_t f2partCumQB = -999.;
 645     if(multGap1B>1)
 646       f2partCumQB = ( ( (QxGap1B*QxGap1B + QyGap1B*QyGap1B) - (Double_t)multGap1B ) / ((Double_t)multGap1B*((Double_t)multGap1B-1)) );
 647     if(f2partCumQB>0)f2partCumQB_vs_Cent->Fill((Double_t)Cent,f2partCumQB);
 648
 649     if (Qvec > fCutLargeQperc && Qvec < 100.){
 650       fResSP_lq->Fill((Double_t)centV0, res);
 651       fResSP_vs_Cent_lq->Fill(Cent, res);
 652       if(f2partCumQA>0)f2partCumQA_vs_Cent_lq->Fill((Double_t)Cent,f2partCumQA);
 653       if(f2partCumQB>0)f2partCumQB_vs_Cent_lq->Fill((Double_t)Cent,f2partCumQB);
 654
 655       fResSP_inclusive->Fill(1., res); //lq v2 for mc closure
 656     }
 657
 658     if (Qvec > 0. && Qvec < fCutSmallQperc){
 659       fResSP_sq->Fill((Double_t)centV0, res);
 660       fResSP_vs_Cent_sq->Fill(Cent, res);
 661       if(f2partCumQA>0)f2partCumQA_vs_Cent_sq->Fill((Double_t)Cent,f2partCumQA);
 662       if(f2partCumQB>0)f2partCumQB_vs_Cent_sq->Fill((Double_t)Cent,f2partCumQB);
 663
 664       fResSP_inclusive->Fill(2., res); //sq v2 for mc closure
 665     }
 666   }// end multiplicity if
 667
 668   if( fFillTHn ){
 669
 670
 671     Double_t varEv[6];
 672     varEv[0]=Cent;
 673     varEv[1]=(Double_t)Qvec; // qvec_rec_perc
 674
 675     Double_t qvzero = 0.;
 676     if(fVZEROside==0)qvzero=(Double_t)fEventCuts->GetqV0A();
 677     else if (fVZEROside==1)qvzero=(Double_t)fEventCuts->GetqV0C(); // qvec_rec
 678     varEv[2]=(Double_t)qvzero; // qvec from VZERO
 679
 680     Double_t qgen_tracks = (Double_t)fEventCuts->CalculateQVectorMC(fVZEROside, 0);
 681     varEv[3]= (Double_t)qgen_tracks;
 682
 683     Double_t qgen_vzero = (Double_t)fEventCuts->CalculateQVectorMC(fVZEROside, 1);
 684     varEv[4]= (Double_t)qgen_vzero;
 685
 686     varEv[5]=(Double_t)fEventCuts->GetNch(); // Nch
 687
 688     ((THnSparseF*)fOutput->FindObject("NSparseHistEv"))->Fill(varEv);//event loop
 689
 690   }
 691
 692
 693   PostData(1, fOutput  );
 694   PostData(2, fEventCuts);
 695   PostData(3, fTrackCuts);
 696   PostData(4, fOutput_lq  );
 697   PostData(5, fOutput_sq  );
 698 }
 699
 700 //_________________________________________________________________
 701 Bool_t  AliAnalysisTaskV2AllChAOD::GetDCA(const AliAODTrack* trk, Double_t * p){
 702
 703   //AliAODTrack::DCA(): for newest AOD fTrack->DCA() always gives -999. This should fix.
 704   //FIXME should update EventCuts?
 705   //FIXME add track->GetXYZ(p) method
 706
 707   double xyz[3],cov[3];
 708
 709   if (!trk->GetXYZ(xyz)) { // dca is not stored
 710     AliExternalTrackParam etp;
 711     etp.CopyFromVTrack(trk);
 712     AliVEvent* ev = (AliVEvent*)trk->GetEvent();
 713     if (!ev) {/*printf("Event is not connected to the track\n");*/ return kFALSE;}
 714     if (!etp.PropagateToDCA(ev->GetPrimaryVertex(), ev->GetMagneticField(),999,xyz,cov)) return kFALSE; // failed, track is too far from vertex
 715   }
 716   p[0] = xyz[0];
 717   p[1] = xyz[1];
 718   return kTRUE;
 719
 720 }
 721
 722 //_________________________________________________________________
 723 void  AliAnalysisTaskV2AllChAOD::MCclosure(Double_t qvec){
 724   // First do MC to fill up the MC particle array
 725
 726   TClonesArray *arrayMC = 0;
 727   if (fIsMC)
 728     {
 729       arrayMC = (TClonesArray*) fAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 730       if (!arrayMC) {
 731         AliFatal("Error: MC particles branch not found!\n");
 732       }
 733
 734       Double_t QxGap1Amc = 0., QyGap1Amc = 0.;
 735       Double_t QxGap1Bmc = 0., QyGap1Bmc = 0.;
 736       Int_t multGap1Amc = 0, multGap1Bmc = 0;
 737
 738       for (Int_t loop = 0; loop < 2; loop++){
 739
 740         Int_t nMC = arrayMC->GetEntries();
 741
 742         for (Int_t iMC = 0; iMC < nMC; iMC++)
 743           {
 744             AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(iMC);
 745             if(!partMC->Charge()) continue;//Skip neutrals
 746             if(fCharge != 0 && partMC->Charge()*fCharge < 0.) continue;//if fCharge != 0 only select fCharge
 747
 748             if (!(partMC->IsPhysicalPrimary()))
 749                 continue;
 750
 751             if(partMC->Eta()<fTrackCuts->GetEtaMin() || partMC->Eta()>fTrackCuts->GetEtaMax()) continue;
 752
 753             //Printf("a particle");
 754
 755
 756            if (loop == 0) {
 757
 758              if (partMC->Eta() > fEtaGapMax){
 759                QxGap1Amc += TMath::Cos(2.*partMC->Phi());
 760                QyGap1Amc += TMath::Sin(2.*partMC->Phi());
 761                multGap1Amc++;
 762             }
 763
 764             if (partMC->Eta() < fEtaGapMin){
 765               QxGap1Bmc += TMath::Cos(2.*partMC->Phi());
 766               QyGap1Bmc += TMath::Sin(2.*partMC->Phi());
 767               multGap1Bmc++;
 768             }
 769
 770           } else {
 771
 772             //eval v2 scalar product
 773             if (partMC->Eta() < fEtaGapMin && multGap1Amc > 0){
 774               Double_t v2SPGap1Amc = (TMath::Cos(2.*partMC->Phi())*QxGap1Amc + TMath::Sin(2.*partMC->Phi())*QyGap1Amc)/(Double_t)multGap1Amc;
 775               fv2SPGap1Amc_inclusive_mb->Fill(partMC->Pt(), v2SPGap1Amc);
 776
 777               if (qvec > fCutLargeQperc && qvec < 100.){
 778                 fv2SPGap1Amc_inclusive_lq->Fill(partMC->Pt(), v2SPGap1Amc);
 779               }
 780
 781               if (qvec > 0. && qvec < fCutSmallQperc){
 782                 fv2SPGap1Amc_inclusive_sq->Fill(partMC->Pt(), v2SPGap1Amc);
 783               }
 784
 785             }
 786
 787             if (partMC->Eta() > fEtaGapMax && multGap1Bmc > 0){
 788               Double_t v2SPGap1Bmc = (TMath::Cos(2.*partMC->Phi())*QxGap1Bmc + TMath::Sin(2.*partMC->Phi())*QyGap1Bmc)/(Double_t)multGap1Bmc;
 789               fv2SPGap1Bmc_inclusive_mb->Fill(partMC->Pt(), v2SPGap1Bmc);
 790
 791               if (qvec > fCutLargeQperc && qvec < 100.){
 792                 fv2SPGap1Bmc_inclusive_lq->Fill(partMC->Pt(), v2SPGap1Bmc);
 793               }
 794
 795               if (qvec > 0. && qvec < fCutSmallQperc){
 796                 fv2SPGap1Bmc_inclusive_sq->Fill(partMC->Pt(), v2SPGap1Bmc);
 797               }
 798
 799             }
 800
 801
 802           }// end else
 803         } // end loop on partMCs
 804       } // end loop
 805
 806       if (multGap1Amc > 0 && multGap1Bmc > 0){
 807         Double_t resmc = (QxGap1Amc*QxGap1Bmc + QyGap1Amc*QyGap1Bmc)/(Double_t)multGap1Amc/(Double_t)multGap1Bmc;
 808         fResSPmc_inclusive->Fill(0.,resmc);
 809
 810         if (qvec > fCutLargeQperc && qvec < 100.){
 811           fResSPmc_inclusive->Fill(1.,resmc);
 812         }
 813
 814         if (qvec > 0. && qvec < fCutSmallQperc){
 815           fResSPmc_inclusive->Fill(2.,resmc);
 816         }
 817
 818       }
 819
 820       }// end if MC
 821 }
 822
 823 //_________________________________________________________________
 824 Double_t AliAnalysisTaskV2AllChAOD::GetRecoEff(Double_t pt, Int_t iC){
 825
 826   if(iC>8) return 1.;
 827
 828   if(pt<0.2 || pt>100.) return 1.;
 829
 830 // // //   //spectra ese binning
 831 // // //   //const Double_t ptBins[] = {0.20,0.30,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.4,2.8,3.2,3.6,4.0,5.0,6.0,7.0,8.0,9.0,10.,12.,15.,20.,25.,30.,35.,40.,50.,75.,100.};
 832 // // //   //const Int_t nptBins=34;
 833
 834   const Double_t fEpsilon=0.000001;
 835
 836   TH1F *h = (TH1F*)fRecoEffList->At(iC);
 837
 838   Int_t bin = h->FindBin(pt);
 839
 840   Double_t lowlim = h->GetBinLowEdge(bin);
 841   Double_t uplim = h->GetBinLowEdge(bin) + h->GetBinWidth(bin);
 842
 843   if( pt>lowlim && pt<uplim ) return h->GetBinContent(bin);
 844   if( pt == lowlim ) return h->GetBinContent( h->FindBin(pt+fEpsilon) );
 845   if( pt == uplim ) return h->GetBinContent( h->FindBin(pt-fEpsilon) );
 846
 847   else return 1.;
 848
 849 }
 850 //_________________________________________________________________
 851 void   AliAnalysisTaskV2AllChAOD::Terminate(Option_t *)
 852 {
 853   // Terminate
 854 }