PWGGA/GammaConv/AliAnalysisTaskdPhi.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                                                                                        *
   4  * Authors: Svein Lindal                                                  *
   5  * Version 1.0                                                                    *
   6  *                                                                                                                                        *
   7  *                                                                                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is                  *
  14  * provided "as is" without express or implied warranty.                                  *
  15  **************************************************************************/
  16
  17 ////////////////////////////////////////////////
  18 //---------------------------------------------
  19 // Class doing conversion gamma dPhi correlations
  20 // Gamma Conversion analysis
  21 //---------------------------------------------
  22 ////////////////////////////////////////////////
  23
  24 #include "AliAnalysisTaskdPhi.h"
  25
  26 #include <TH2I.h>
  27 #include <TList.h>
  28 #include <TChain.h>
  29 #include <TFile.h>
  30 #include <AliAnalysisManager.h>
  31 #include <AliInputEventHandler.h>
  32 #include <AliESDInputHandler.h>
  33 #include <AliAODInputHandler.h>
  34 #include <TGeoGlobalMagField.h>
  35
  36 #include "AliConversionTrackCuts.h"
  37 #include "AliConversionCuts.h"
  38 #include "AliConversionMesonCuts.h"
  39 #include "AliAODConversionPhoton.h"
  40 #include "AliAODConversionMother.h"
  41 #include "TGrid.h"
  42 // #include "AliAnaConvCorrPhoton.h"
  43 // #include "AliAnaConvCorrPion.h"
  44 // #include "AliAnaConvIsolation.h"
  45 #include "AliV0ReaderV1.h"
  46 // Author Svein Lindal <slindal@fys.uio.no>
  47 using namespace std;
  48
  49 ClassImp(AliAnalysisTaskdPhi)
  50
  51
  52 //________________________________________________________________________
  53 AliAnalysisTaskdPhi::AliAnalysisTaskdPhi(const char *name) : AliAnalysisTaskSE(name),
  54   fHistograms(NULL),
  55   fCorrSparse(NULL),
  56   fTrigSparse(NULL),
  57   fTrackSparse(NULL),
  58   fMassSparse(NULL),
  59   fV0Reader(NULL),
  60   fSaveReaderHists(kFALSE),
  61   fV0Filter(NULL),
  62   fV0Filters(),
  63   fPhotonFilter(NULL),
  64   fMesonFilter(NULL),
  65   fMesonFilters(),
  66   fTrackFilter(NULL),
  67   fTrackFilters(),
  68   fGammas(),
  69   fTracks(),
  70   hMEvents(NULL),
  71   hTrackCent(NULL),
  72   hTrigPt(NULL),
  73   hTrackPt(NULL),
  74   hTrigPhi(NULL),
  75   fDeltaAODBranchName("AliAODGammaConversion_gamma"),
  76   fAxistPt(),
  77   fAxiscPt(),
  78   fAxisdEta(),
  79   fAxisTrigEta(),
  80   fAxisAssEta(),
  81   fAxisdPhi(),
  82   fAxisCent(),
  83   fAxisZ(),
  84   fAxisPiM(),
  85   fAxisTrackFilters(),
  86   fAxisV0Filters(),
  87   fAxisMesonFilters(),
  88   fkTrackAxis(kFALSE),
  89   fkV0Axis(kFALSE),
  90   fkPionAxis(kFALSE),
  91   fAxesList(),
  92   fTrigAxesList(),
  93   fTrackAxesList(),
  94   fMassAxesList(),
  95   fDoPhoton(kFALSE),
  96   fCorrectionMap(NULL)
  97  {
  98    //constructor
  99    SetUpBins();
 100
 101    DefineInput(0, TChain::Class());
 102    DefineOutput(1, TList::Class());
 103
 104    fGammas.SetOwner(kTRUE);
 105    fTracks.SetOwner(kTRUE);
 106
 107
 108
 109
 110  }
 111
 112
 113
 114  //________________________________________________________________________
 115  AliAnalysisTaskdPhi::~AliAnalysisTaskdPhi(){
 116    //destructor
 117
 118    if(fV0Filter)
 119          delete fV0Filter;
 120    fV0Filter = NULL;
 121
 122    if(fMesonFilter)
 123          delete fMesonFilter;
 124    fMesonFilter = NULL;
 125
 126    if(fPhotonFilter)
 127          delete fPhotonFilter;
 128    fPhotonFilter = NULL;
 129
 130    if(fHistograms)
 131          delete fHistograms;
 132    fHistograms = NULL;
 133
 134    if(fTrackFilter)
 135      delete fTrackFilter;
 136    fTrackFilter = NULL;
 137
 138    fGammas.Delete();
 139    fTracks.Delete();
 140
 141  }
 142
 143  ///________________________________________________________________________
 144 void AliAnalysisTaskdPhi::SetUpBins() {
 145
 146    fAxisTrigEta.SetNameTitle("tEta", "Eta");
 147    fAxisTrigEta.Set(320, -0.8, 0.8);
 148
 149    fAxisAssEta.SetNameTitle("aEta", "Eta");
 150    fAxisAssEta.Set(360, -0.9, 0.9);
 151
 152    fAxisdEta.SetNameTitle("dEta", "Eta");
 153    fAxisdEta.Set(34, -1.7, 1.7);
 154
 155    fAxisdPhi.SetNameTitle("dPhi", "delta Phi");
 156    fAxisdPhi.Set(32, -TMath::PiOver2(), 3*TMath::PiOver2());
 157
 158    Double_t tptbins[10] = {0.1, 3.0, 4.0, 5.0, 6.0, 8.0, 10.0, 15, 50, 100};                                    fAxistPt.SetNameTitle("tPt", "tPt");
 159    fAxistPt.Set(9, tptbins);
 160
 161    Double_t cptbins[13] = {0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 4.0, 5.0, 8.0, 10.0, 25, 50, 100};                     fAxiscPt.SetNameTitle("cPt", "cPt");
 162    fAxiscPt.Set(12, cptbins);
 163
 164    fAxisZ.SetNameTitle("vertexz", "Z");
 165    fAxisZ.Set(4, -10, 10);
 166
 167    Double_t centbins[5] = {0, 10, 30, 60, 100.1};
 168    fAxisCent.SetNameTitle("centrality", "Cent");
 169    fAxisCent.Set(4, centbins);
 170
 171    fAxisPiM.SetNameTitle("InvMassPi0", "Invariant mass");
 172    Double_t mbins[7] = {0.1, 0.11, 0.12, 0.15, 0.16, 0.18, 0.2};
 173    fAxisPiM.Set(6, mbins);
 174
 175 }
 176
 177
 178  ///________________________________________________________________________
 179 // void AliAnalysisTaskdPhi::SetUpCorrObjects() {
 180    // //Set up corr objects
 181    // AliDebug(AliLog::kDebug + 5, "Set Up corr objects");
 182
 183    // if(fDoPhoton) {
 184    //   fPhotonCorr = new AliAnaConvCorrPhoton("PhotonCorr","photon %s");
 185    //   SetUpCorrAxes(fPhotonCorr);
 186    //   fPhotonCorr->CreateHistograms();
 187    //   fHistograms->Add(fPhotonCorr->GetHistograms());
 188    // }
 189
 190    // fPionCorr = new AliAnaConvCorrPion("PionCorr", "pion");
 191    // SetUpCorrAxes(fPionCorr);
 192    // fPionCorr->GetAxisM().Set(fAxisPiM.GetNbins(), fAxisPiM.GetXbins()->GetArray());
 193    // fPionCorr->CreateHistograms();
 194    // fHistograms->Add(fPionCorr->GetHistograms());
 195    //}
 196
 197  ///________________________________________________________________________
 198 // void AliAnalysisTaskdPhi::SetUpCorrAxes(AliAnaConvCorrBase * corr) {
 199    ///Set up axes in corr object
 200    // corr->GetAxisCent().Set(fAxisCent.GetNbins(), fAxisCent.GetXbins()->GetArray());
 201    // const Double_t * zbins = fAxisZ.GetXbins()->GetArray();
 202    // if(zbins) {
 203    //   corr->GetAxisZ().Set(fAxisZ.GetNbins(), fAxisZ.GetXbins()->GetArray());
 204    // } else {
 205    //   corr->GetAxisZ().Set(fAxisZ.GetNbins(), fAxisZ.GetBinLowEdge(1), fAxisZ.GetBinUpEdge(fAxisZ.GetNbins()));
 206    // }
 207
 208    // corr->GetAxistPt().Set(fAxistPt.GetNbins(), fAxistPt.GetXbins()->GetArray());
 209    // corr->GetAxiscPt().Set(fAxiscPt.GetNbins(), fAxiscPt.GetXbins()->GetArray());
 210    //   corr->GetAxisdEta().Set(fAxisdEta.GetNbins(), fAxisdEta.GetBinLowEdge(1), fAxisdEta.GetBinUpEdge(fAxisdEta.GetNbins()));
 211    // corr->GetAxisTrigEta().Set(fAxisTrigEta.GetNbins(), fAxisTrigEta.GetBinLowEdge(1), fAxisTrigEta.GetBinUpEdge(fAxisTrigEta.GetNbins()));
 212    // corr->GetAxisAssEta().Set(fAxisAssEta.GetNbins(), fAxisAssEta.GetBinLowEdge(1), fAxisAssEta.GetBinUpEdge(fAxisAssEta.GetNbins()));
 213 // }
 214
 215
 216 //________________________________________________________________________
 217 void AliAnalysisTaskdPhi::UserCreateOutputObjects() {
 218   // Create histograms
 219   // TGrid::Connect("alien://",0,0,"t");
 220   // if(!gGrid) AliWarning("no GGrid");
 221   // TFile *tfile = TFile::Open("alien:///alice/cern.ch/user/s/slindal/trackMap.root", "READ");
 222   // if(tfile) {
 223   //   THnF * corrmap = dynamic_cast<THnF*>(tfile->Get("hTrackCorr"));
 224   //   if (corrmap) {
 225   //     fCorrectionMap = dynamic_cast<THnF*>(THn::CreateHn("corr", "corr", corrmap));
 226   //     for(Int_t i = 0; i < fCorrectionMap->GetNdimensions(); i++) {
 227   //    TAxis * axis = fCorrectionMap->GetAxis(i);
 228   //    axis->SetRange(1, axis->GetNbins());
 229   //     }
 230
 231   //     cout << "yessssssssssssssssssssssssssssssssssssssssssssssssss"<<endl;
 232   //   } else {
 233   //     cout << "xxxxxxxxxxxxxxxxx xxxxxxxxxxxxxxx xxxxxxxxxxxxxxxx"<<endl;
 234   //   }
 235   //   tfile->Close();
 236   // } else {
 237   //   cout << "no tfile shit shit shit "<<endl;
 238   //   AliFatal("file not ther!!!");
 239   // }
 240
 241
 242   fHistograms = new TList();
 243   fHistograms->SetName("dPhi_histograms");
 244   fHistograms->SetOwner(kTRUE);
 245
 246
 247   if(!fV0Reader){
 248     fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
 249   }
 250
 251   if(!fV0Reader){
 252     printf("Error: No V0 Reader");
 253   } // GetV0Reader
 254
 255
 256   if(fSaveReaderHists) {
 257     AliConversionCuts * v0cuts = fV0Reader->GetConversionCuts();
 258     if(v0cuts) {
 259       TList * histograms = v0cuts->GetCutHistograms();
 260       if(!histograms) {
 261         AliWarning("initializing v0 reader hists");
 262         v0cuts->InitCutHistograms("V0Reader", kTRUE);
 263       }
 264       histograms = v0cuts->GetCutHistograms();
 265       if(histograms) {
 266         fHistograms->Add(histograms);
 267       }
 268     }
 269   }
 270
 271   for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf ++){
 272     AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[0].At(igf));
 273     if(f) {
 274       TList * histograms = f->GetCutHistograms();
 275       if(histograms) fHistograms->Add(f->GetCutHistograms());
 276     }
 277   }
 278
 279   for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf ++){
 280     AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[1].At(igf));
 281     if(f) {
 282       TList * histograms = f->GetCutHistograms();
 283       if(histograms) fHistograms->Add(f->GetCutHistograms());
 284     }
 285   }
 286
 287   for(Int_t igf = 0; igf < fMesonFilters[0].GetEntriesFast(); igf ++){
 288     AliConversionMesonCuts * f = dynamic_cast<AliConversionMesonCuts*>(fMesonFilters[0].At(igf));
 289     if(f) {
 290       TList * histograms = f->GetCutHistograms();
 291       if(histograms) fHistograms->Add(f->GetCutHistograms());
 292     }
 293   }
 294
 295   for(Int_t igf = 0; igf < fMesonFilters[1].GetEntriesFast(); igf ++){
 296     AliConversionMesonCuts * f = dynamic_cast<AliConversionMesonCuts*>(fMesonFilters[1].At(igf));
 297     if(f) {
 298       TList * histograms = f->GetCutHistograms();
 299       if(histograms) fHistograms->Add(f->GetCutHistograms());
 300     }
 301   }
 302
 303   if(fV0Filter) {
 304     fV0Filter->InitCutHistograms("V0Filter", kFALSE);
 305     fHistograms->Add(fV0Filter->GetCutHistograms());
 306   }
 307   if(fMesonFilter) {
 308     fMesonFilter->InitCutHistograms("PionFilter", kFALSE);
 309     fHistograms->Add(fMesonFilter->GetCutHistograms());
 310   }
 311   if(fPhotonFilter) {
 312     fPhotonFilter->InitCutHistograms("PhotonFilter", kFALSE);
 313     fHistograms->Add(fPhotonFilter->GetCutHistograms());
 314   }
 315
 316
 317   AliConversionTrackCuts * tc = dynamic_cast<AliConversionTrackCuts*>(fTrackFilter);
 318   if(tc) {
 319     fHistograms->Add(tc->CreateHistograms());
 320   }
 321
 322
 323    // for(Int_t i = 0; i < fTrackFilters.GetEntriesFast(); i++) {
 324    //   AliConversionTrackCuts * tc = dynamic_cast<AliConversionTrackCuts*>(fTrackFilters.At(i));
 325    //   if(tc) fHistograms->Add(tc->CreateHistograms());
 326    // }
 327
 328   //SetUpCorrObjects();
 329
 330    ///Set up ME histograms
 331    TList * MEHistograms = new TList();
 332    MEHistograms->SetName("MEHistograms");
 333    MEHistograms->SetOwner(kTRUE);
 334    fHistograms->Add(MEHistograms);
 335
 336    hMEvents = new TH2I("hMEvents", "Nevents vs centrality vertexz",
 337                        fAxisZ.GetNbins(), fAxisZ.GetXbins()->GetArray(),
 338                        fAxisCent.GetNbins(), fAxisCent.GetXbins()->GetArray());
 339    MEHistograms->Add(hMEvents);
 340
 341    hTrackCent = new TH2I("hTrackCent", "N accepted tracks vs centrality",
 342                          fAxisCent.GetNbins() > 1 ? (int) (10*(fAxisCent.GetXmax() - fAxisCent.GetXmin()))  : 1,
 343                          fAxisCent.GetXmin(), fAxisCent.GetXmax(),
 344                          fAxisCent.GetNbins() > 1 ? 900 : 50,
 345                          0,
 346                          fAxisCent.GetNbins() > 1 ? 1800 : 50);
 347    MEHistograms->Add(hTrackCent);
 348
 349    hTrigPt = new TH3F("hTrigPt", "trigger pt", 100, 0., 10.,
 350                       10, 0., 50.,
 351                       5,  0.05, 0.2);
 352    MEHistograms->Add(hTrigPt);
 353    hTrackPt = new TH2F("hTrackPt", "track pt", 100, 0, 10, 10, 0, 50);//fAxisCent.GetNbins(), fAxisCent.GetXbins()->GetArray());
 354    MEHistograms->Add(hTrackPt);
 355    hTrigPhi = new TH1F("hTrigPhi", "trigger pt", 32, 0, 2*TMath::Pi());
 356    MEHistograms->Add(hTrigPhi);
 357
 358
 359
 360    Int_t ntrackfilters[2] = {fTrackFilters[0].GetEntriesFast(), fTrackFilters[1].GetEntriesFast()};
 361    fkTrackAxis = kTRUE;
 362    fAxisTrackFilters.SetNameTitle("trackCuts", "trackCuts");
 363    fAxisTrackFilters.Set(ntrackfilters[0] + ntrackfilters[1] + 1, -ntrackfilters[0] -0.5, ntrackfilters[1]  + 0.5);
 364
 365    Int_t nV0filters[2] = {fV0Filters[0].GetEntriesFast(), fV0Filters[1].GetEntriesFast()};
 366    fkV0Axis = kTRUE;
 367    fAxisV0Filters.SetNameTitle("V0Cuts", "V0Cuts");
 368    fAxisV0Filters.Set(nV0filters[0] + nV0filters[1] + 1, -nV0filters[0] -0.5, nV0filters[1]  + 0.5);
 369
 370    Int_t nmesonfilters[2] = {fMesonFilters[0].GetEntriesFast(), fMesonFilters[1].GetEntriesFast()};
 371    fkPionAxis = kTRUE;
 372    fAxisMesonFilters.SetNameTitle("mesonCuts", "mesonCuts");
 373    fAxisMesonFilters.Set(nmesonfilters[0] + nmesonfilters[1] + 1, -nmesonfilters[0] -0.5, nmesonfilters[1]  + 0.5);
 374
 375   fAxesList.AddAt(&fAxisdEta, 0);
 376   fAxesList.AddAt(&fAxisdPhi, 1);
 377   fAxesList.AddAt(&fAxistPt, 2);
 378   fAxesList.AddAt(&fAxiscPt, 3);
 379   fAxesList.AddAt(&fAxisCent, 4);
 380   fAxesList.AddAt(&fAxisZ, 5);
 381   fAxesList.AddAt(&fAxisPiM, 6);
 382   fAxesList.AddAt(&fAxisTrackFilters, 7);
 383   fAxesList.AddAt(&fAxisV0Filters, 8);
 384   fAxesList.AddAt(&fAxisMesonFilters, 9);
 385
 386   fTrackAxesList.AddAt(&fAxisAssEta, 0);
 387   fTrackAxesList.AddAt(&fAxistPt, 1);
 388   fTrackAxesList.AddAt(&fAxiscPt, 2);
 389   fTrackAxesList.AddAt(&fAxisCent, 3);
 390   fTrackAxesList.AddAt(&fAxisZ, 4);
 391   //fTrackAxesList.AddAt(&fAxisPiM, 5);
 392   fTrackAxesList.AddAt(&fAxisTrackFilters, 5);
 393   fTrackAxesList.AddAt(&fAxisV0Filters, 6);
 394   fTrackAxesList.AddAt(&fAxisMesonFilters, 7);
 395
 396   fTrigAxesList.AddAt(&fAxisTrigEta, 0);
 397   fTrigAxesList.AddAt(&fAxistPt, 1);
 398   fTrigAxesList.AddAt(&fAxisCent, 2);
 399   fTrigAxesList.AddAt(&fAxisZ, 3);
 400   fTrigAxesList.AddAt(&fAxisPiM, 4);
 401   fTrigAxesList.AddAt(&fAxisV0Filters, 5);
 402   fTrigAxesList.AddAt(&fAxisMesonFilters, 6);
 403
 404   TList masslist;
 405   TAxis massax;
 406   massax.SetNameTitle("mass", "mass");
 407   massax.Set(360, 0.04, 0.4); //hardcoded! change also in filling!
 408
 409   masslist.AddAt(&massax, 0);
 410   masslist.AddAt(&fAxistPt, 1);
 411   masslist.AddAt(&fAxisCent, 2);
 412   masslist.AddAt(&fAxisV0Filters, 3);
 413   masslist.AddAt(&fAxisMesonFilters, 4);
 414
 415   fCorrSparse = CreateSparse(TString("pionSparse"), TString("pionSparse"), &fAxesList);
 416   fTrackSparse = CreateSparse(TString("trackSparse"), TString("trackSparse"), &fTrackAxesList);
 417   fTrigSparse = CreateSparse(TString("trigSparse"), TString("trigSparse"), &fTrigAxesList);
 418   fMassSparse = CreateSparse("massSparse", "massSparse", &masslist);
 419
 420   fHistograms->Add(fCorrSparse);
 421   fHistograms->Add(fTrackSparse);
 422   fHistograms->Add(fTrigSparse);
 423   fHistograms->Add(fMassSparse);
 424
 425
 426   ///Add gamma and track containers:
 427   for(Int_t i = 0; i < fV0Filters[1].GetEntriesFast() + 1; i++) {
 428     fGammas.Add(new TObjArray());
 429   }
 430
 431   for(Int_t i = 0; i < fTrackFilters[1].GetEntriesFast() + 1; i++) {
 432     fTracks.Add(new TObjArray());
 433   }
 434
 435
 436
 437   PostData(1, fHistograms);
 438 }
 439
 440 ///________________________________________________________________________
 441 THnSparseF * AliAnalysisTaskdPhi::CreateSparse(TString nameString, TString titleString, TList * axesList) {
 442   //Create sparse
 443   const Int_t dim = axesList->GetSize();
 444
 445   TAxis * axes[dim];
 446   Int_t   bins[dim];
 447   Double_t min[dim];
 448   Double_t max[dim];
 449
 450   for(Int_t i = 0; i<dim; i++) {
 451     TAxis * axis = dynamic_cast<TAxis*>(axesList->At(i));
 452     if(axis) axes[i] = axis;
 453     else {
 454       cout << "AliAnalysisTaskdPhi::CreateSparse: Error error, all the axes are not present in axis list" << endl;
 455       return NULL;
 456     }
 457   }
 458
 459   for(Int_t i = 0; i<dim; i++) {
 460     //cout << axes[i]->GetTitle() << endl;
 461     bins[i] = axes[i]->GetNbins();
 462     min[i] = axes[i]->GetBinLowEdge(1);
 463     max[i] = axes[i]->GetBinUpEdge(axes[i]->GetNbins());
 464   }
 465
 466   THnSparseF * sparse = new THnSparseF(Form("%s", nameString.Data()),
 467                                        Form("%s", titleString.Data()),
 468                                        dim, bins, min, max);
 469
 470   for(Int_t i = 0; i<dim; i++) {
 471     sparse->GetAxis(i)->SetNameTitle(axes[i]->GetName(), axes[i]->GetTitle() );
 472     if(axes[i]->GetXbins()->GetSize() > 0) {
 473       sparse->SetBinEdges(i, axes[i]->GetXbins()->GetArray() );
 474     }
 475   }
 476   return sparse;
 477 }
 478
 479 //________________________________________________________________________
 480 void AliAnalysisTaskdPhi::UserExec(Option_t *) {
 481   ///User exec.
 482   ///This is a very ugly function, cut the complexity of the logic demands it.
 483
 484
 485   AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
 486   if(esdEvent) {
 487     if (!TGeoGlobalMagField::Instance()->GetField()) esdEvent->InitMagneticField();
 488   }
 489
 490
 491   //if(! fV0Filter->EventIsSelected(fInputEvent)) return;
 492   if(!fV0Reader){
 493     AliError("Error: No V0 Reader");
 494     return;
 495   } // GetV0Reader
 496
 497   if(!fV0Reader->IsEventSelected()) {
 498     return;
 499   }
 500
 501   AliDebug(5, "Processing event");
 502
 503
 504   for(Int_t i = 0; i < fGammas.GetEntriesFast(); i++) {
 505     static_cast<TObjArray*>(fGammas.At(i))->Clear();
 506   }
 507
 508   for(Int_t i = 0; i < fTracks.GetEntriesFast(); i++) {
 509     static_cast<TObjArray*>(fTracks.At(i))->Clear();
 510   }
 511
 512
 513
 514   AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 515   Bool_t isAOD=man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class();
 516
 517   AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 518   if (!inputHandler) {
 519     cout << "cout no input event handler"<<endl;
 520     return;
 521   }
 522
 523
 524
 525
 526   for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf++) {
 527     AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[0].At(igf));
 528     if ( f && !f->GetPIDResponse() ) {
 529       if ( inputHandler->GetPIDResponse() ){
 530         f->SetPIDResponse( inputHandler->GetPIDResponse() );
 531       } else {
 532         if (isAOD){
 533           if (!f->GetPIDResponse()){
 534             f->InitAODpidUtil(1);
 535           }
 536         }
 537       }
 538     } else {
 539       break;
 540     }
 541   }
 542
 543   for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf++) {
 544     AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[1].At(igf));
 545     if ( f && !f->GetPIDResponse() ) {
 546       if ( inputHandler->GetPIDResponse() ){
 547         f->SetPIDResponse( inputHandler->GetPIDResponse() );
 548       } else {
 549         if (isAOD){
 550           if (!f->GetPIDResponse()){
 551             f->InitAODpidUtil(1);
 552           }
 553         }
 554       }
 555     } else {
 556       break;
 557     }
 558   }
 559
 560
 561   if ( fV0Filter && !fV0Filter->GetPIDResponse() ) {
 562     if ( inputHandler->GetPIDResponse() ){
 563       fV0Filter->SetPIDResponse( inputHandler->GetPIDResponse() );
 564     } else {
 565
 566       //AOD case
 567       if (isAOD){
 568         if (!fV0Filter->GetPIDResponse()){
 569           fV0Filter->InitAODpidUtil(1);
 570         }
 571       }
 572     }
 573   }
 574
 575
 576   ///Initialize track cuts. Delete tracks that have been constrained to vertex (copies)
 577   AliConversionTrackCuts * tc = dynamic_cast<AliConversionTrackCuts*>(fTrackFilter);
 578   if(tc) {
 579     tc->SetEvent(fInputEvent);
 580     tc->DeleteTracks();
 581   }
 582
 583   for(Int_t i = 0; i < fTrackFilters[0].GetEntriesFast(); i++){
 584     AliConversionTrackCuts * tct = dynamic_cast<AliConversionTrackCuts*>(fTrackFilters[0].At(i));
 585     if(tct) {
 586       tct->SetEvent(fInputEvent);
 587       tct->DeleteTracks();
 588     }
 589   }
 590   for(Int_t i = 0; i < fTrackFilters[1].GetEntriesFast(); i++){
 591     AliConversionTrackCuts * tct = dynamic_cast<AliConversionTrackCuts*>(fTrackFilters[1].At(i));
 592     if(tct) {
 593       tct->SetEvent(fInputEvent);
 594       tct->DeleteTracks();
 595     }
 596   }
 597
 598
 599   Double_t centrality = 0.0;
 600   Double_t vertexz = fInputEvent->GetPrimaryVertex()->GetZ();
 601   if(isAOD) {
 602     AliAODHeader * header = static_cast<AliAODHeader*>(fInputEvent->GetHeader());
 603     centrality = header->GetCentrality();
 604   } else {
 605     centrality = static_cast<AliESDEvent*>(fInputEvent)->GetCentrality()->GetCentralityPercentile("V0M");
 606   }
 607
 608
 609   const Int_t centBin = GetBin(fAxisCent, centrality);
 610   const Int_t vertexBin = GetBin(fAxisZ, vertexz);
 611
 612
 613   if(DebugLevel () > 4) {
 614     cout << "centrality: " << centrality <<  " " << GetBin(fAxisCent, centrality) << endl;
 615     cout << "vertexz: " << vertexz <<  " " << GetBin(fAxisZ, vertexz) << endl;
 616   }
 617
 618
 619
 620   if(centBin < 0 || vertexBin < 0) {
 621     //  AliError("bin out of range");
 622     return;
 623   }
 624
 625
 626
 627   //TClonesArray * aodGammas = GetConversionGammas(isAOD);
 628   TClonesArray * aodGammas = fV0Reader->GetReconstructedGammas();
 629   if(!aodGammas) {
 630     AliError("no aod gammas found!");
 631     return;
 632   }
 633
 634   TObjArray * ggammas = static_cast<TObjArray*>(fGammas.At(0));
 635
 636   ///Fill arrays of accepted gammas
 637   if(DebugLevel() > 1) printf("Number of conversion gammas %d \n", aodGammas->GetEntriesFast());
 638   for(Int_t ig = 0; ig < aodGammas->GetEntriesFast(); ig++) {
 639     AliAODConversionPhoton * photon = dynamic_cast<AliAODConversionPhoton*>(aodGammas->At(ig));
 640     if(!photon) continue;
 641     if(!fV0Filter || fV0Filter->PhotonIsSelected(photon, fInputEvent)) {
 642       ggammas->Add(photon);
 643     } else {
 644       for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf++) {
 645         AliConversionCuts * gfilter = dynamic_cast<AliConversionCuts*>(fV0Filters[1].At(igf));
 646         if(gfilter && gfilter->PhotonIsSelected(photon, fInputEvent)) {
 647           static_cast<TObjArray*>(fGammas.At(igf+1))->Add(photon);
 648         }
 649       }
 650     }
 651   }
 652
 653   Bool_t lowgmap[fV0Filters[0].GetEntriesFast()][ggammas->GetEntriesFast()];
 654
 655   for(Int_t ig = 0; ig < ggammas->GetEntriesFast(); ig++ ) {
 656     AliAODConversionPhoton * gamma = static_cast<AliAODConversionPhoton*>(ggammas->At(ig));
 657     for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf++) {
 658
 659       AliConversionCuts * v0cuts = static_cast<AliConversionCuts*>(fV0Filters[0].At(igf));
 660       if(!(v0cuts->PhotonIsSelected(gamma, fInputEvent))) {
 661         lowgmap[igf][ig] = kTRUE;
 662       } else {
 663         lowgmap[igf][ig] = kFALSE;
 664       }
 665     }
 666   }
 667
 668
 669   if(DebugLevel() > 4) printf("Number of accepted gammas %d \n", ggammas->GetEntriesFast());
 670   hMEvents->Fill(vertexz, centrality);
 671
 672   ///create track array
 673   const Int_t ntrackfilters[2] = { fTrackFilters[0].GetEntriesFast(), fTrackFilters[1].GetEntriesFast()};
 674
 675   TObjArray * ttracks = static_cast<TObjArray*>(fTracks.At(0));
 676   const Double_t aetalim[2] = { fAxisAssEta.GetXmin(), fAxisAssEta.GetXmax()};
 677   const Double_t aptlim[2] = { fAxiscPt.GetXmin(), fAxiscPt.GetXmax()};
 678   for(Int_t iTrack = 0; iTrack < fInputEvent->GetNumberOfTracks(); iTrack++) {
 679     AliVTrack * track = static_cast<AliVTrack*>(fInputEvent->GetTrack(iTrack));
 680     if(track->Pt() < aptlim[0] || track->Pt() > aptlim[1]) continue;
 681     if(track->Eta() < aetalim[0] || track->Eta() > aetalim[1]) continue;
 682     if(fTrackFilter->IsSelected(track)) {
 683       hTrackPt->Fill(track->Pt(), centrality);
 684       ttracks->Add(track);
 685     } else {
 686       ///upside cuts
 687       for(Int_t itf = 1; itf < ntrackfilters[1] + 1; itf++) {
 688         AliAnalysisCuts * trackCuts = static_cast<AliAnalysisCuts*>(fTrackFilters[1].At(itf -1));
 689         if(trackCuts->IsSelected(track)) {
 690           static_cast<TObjArray*>(fTracks.At(itf))->Add(track);
 691         }
 692       }
 693     }
 694   }
 695
 696
 697   Bool_t lowtrackmap[ntrackfilters[0]][ttracks->GetEntriesFast()];
 698   ///Check lowside cuts
 699   for(Int_t iTrack = 0; iTrack < ttracks->GetEntriesFast(); iTrack++ ) {
 700     AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(iTrack));
 701     for(Int_t itf = 0; itf < ntrackfilters[0]; itf++) {
 702       AliAnalysisCuts * trackCuts = static_cast<AliAnalysisCuts*>(fTrackFilters[0].At(itf));
 703       if(!trackCuts->IsSelected(track)) {
 704         lowtrackmap[itf][iTrack] = kTRUE;
 705       } else {
 706         lowtrackmap[itf][iTrack] = kFALSE;
 707       }
 708     }
 709   }
 710
 711   hTrackCent->Fill(centrality, ttracks->GetEntriesFast());
 712
 713   const Double_t etalim[2] = { fAxisTrigEta.GetXmin(), fAxisTrigEta.GetXmax()};
 714   if(DebugLevel() > 4) printf("Number of accepted gammas, tracks %d  %d \n", ggammas->GetEntriesFast(), ttracks->GetEntriesFast());
 715
 716   //AliAnaConvCorrBase * gCorr = fPhotonCorr; //GetCorrObject(vertexBin, centBin, fPhotonCorr);
 717   //  AliAnaConvCorrPion * piCorr = fPionCorr; //static_cast<AliAnaConvCorrPion*>(GetCorrObject(vertexBin, centBin, fPionCorr));
 718   // if(!piCorr) {
 719   //   AliError("corr object missing");
 720   //   return;
 721   // }
 722
 723   TObjArray pions;
 724
 725   ///corr values
 726   Double_t dphivalues[fAxesList.GetSize()];
 727   dphivalues[4] = centrality;
 728   dphivalues[5] = vertexz;
 729
 730   ///Trigger me counters
 731   Double_t trigValues[7];
 732   trigValues[2] = centrality;
 733   trigValues[3] = vertexz;
 734
 735   ///Mass histogram
 736   Double_t massval[5];
 737   massval[2] = centrality;
 738
 739   ///Set up track me counters and initialize
 740   const Int_t nbins = fAxistPt.GetNbins();
 741   Bool_t tmap[fAxistPt.GetNbins()];
 742   Bool_t lv0tmap[fAxistPt.GetNbins()][fV0Filters[0].GetEntriesFast()];
 743   Bool_t uv0tmap[fAxistPt.GetNbins()][fV0Filters[1].GetEntriesFast()];
 744
 745   Bool_t lpitmap[fAxistPt.GetNbins()][fMesonFilters[0].GetEntriesFast()];
 746   Bool_t upitmap[fAxistPt.GetNbins()][fMesonFilters[1].GetEntriesFast()];
 747
 748   for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf++) {
 749     for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
 750       lv0tmap[ptbin][igf] = kFALSE;
 751     }
 752   }
 753
 754   for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf++) {
 755     for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
 756       uv0tmap[ptbin][igf] = kFALSE;
 757     }
 758   }
 759
 760   for(Int_t igf = 0; igf < fMesonFilters[0].GetEntriesFast(); igf++) {
 761     for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
 762       lpitmap[ptbin][igf] = kFALSE;
 763     }
 764   }
 765
 766   for(Int_t igf = 0; igf < fMesonFilters[1].GetEntriesFast(); igf++) {
 767     for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
 768       upitmap[ptbin][igf] = kFALSE;
 769     }
 770   }
 771
 772   for(Int_t ptbin = 0; ptbin < nbins; ptbin++){
 773     tmap[ptbin] = kFALSE;
 774   }
 775
 776   /////////////////////////////////////////
 777
 778
 779   for(Int_t igf1 = 0; igf1 < fV0Filters[1].GetEntriesFast() + 1; igf1++) {
 780     TObjArray * gamm1 = static_cast<TObjArray*>(fGammas.At(igf1));
 781     for(Int_t i1 = 0; i1 < gamm1->GetEntriesFast(); i1++) {
 782       AliAODConversionPhoton * ph1 = static_cast<AliAODConversionPhoton*>(gamm1->UncheckedAt(i1));
 783       Int_t tIDs[4] = {ph1->GetLabel(0), ph1->GetLabel(1), -1, -1};
 784
 785       ///Combine gamma into pions
 786       Int_t igmax = 0;
 787       for(Int_t igf2 = 0; igf2 <= igf1; igf2++) {
 788         TObjArray * gamm2 = NULL;
 789         if(igf2 == igf1) {
 790           gamm2 = gamm1;
 791           igmax = i1;
 792         } else {
 793           gamm2 = static_cast<TObjArray*>(fGammas.At(igf2));
 794           igmax = gamm2->GetEntriesFast();
 795         }
 796
 797         for(Int_t i2 = 0; i2 < igmax; i2++) {
 798           AliAODConversionPhoton * ph2 = static_cast<AliAODConversionPhoton*>(gamm2->UncheckedAt(i2));
 799
 800           if( ph2->GetTrackLabelPositive()==ph1->GetTrackLabelPositive()
 801               || ph2->GetTrackLabelNegative()==ph1->GetTrackLabelNegative()
 802               || ph2->GetTrackLabelNegative()==ph1->GetTrackLabelPositive()
 803               || ph2->GetTrackLabelPositive()==ph1->GetTrackLabelNegative()) {
 804             continue;
 805           }
 806
 807           AliAODConversionMother * pion = new AliAODConversionMother(ph1, ph2);
 808           if(pion->Eta() < etalim[0] || pion->Eta() > etalim[1]) continue;
 809           pion->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
 810
 811
 812           tIDs[2] = ph2->GetLabel(0);
 813           tIDs[3] = ph2->GetLabel(1);
 814
 815           massval[0] = pion->M();
 816           massval[1] = pion->Pt();
 817           massval[3] = igf1;
 818           massval[4] = 0;
 819
 820           dphivalues[2] = pion->Pt();
 821           dphivalues[6] = pion->M();
 822           dphivalues[8] = igf1;
 823           dphivalues[9] = 0;
 824
 825           trigValues[0] = pion->Eta();
 826           trigValues[1] = pion->Pt();
 827           trigValues[4] = pion->M();
 828           trigValues[5] = igf1;
 829           trigValues[6] = 0;
 830
 831
 832           ///Check that particle is in histo phase space
 833           if(pion->Pt() > fAxistPt.GetBinLowEdge(1) && pion->Pt() < fAxistPt.GetXmax() &&
 834              pion->Eta() > etalim[0] && pion->Eta() < etalim[1] &&
 835              pion->M() > 0.04 && pion->M() < 0.4
 836              ) {
 837
 838
 839             if(fMesonFilter->MesonIsSelected(pion, kTRUE)) {
 840
 841               fMassSparse->Fill(massval);
 842
 843               Bool_t lpimap[fMesonFilters[0].GetEntriesFast()];
 844               ///See if it passes lowside cuts
 845               if(igf1 == 0 && igf2 == 0) {
 846
 847                 ///Low side pion
 848                 massval[3] = 0;
 849                 for(Int_t ilpf = 0; ilpf < fMesonFilters[0].GetEntriesFast(); ilpf++) {
 850                   if(!(static_cast<AliConversionMesonCuts*>(fMesonFilters[0].At(ilpf))->MesonIsSelected(pion, kTRUE))) {
 851                     lpimap[ilpf] = kTRUE;
 852                     massval[4] = -(ilpf + 1);
 853                     fMassSparse->Fill(massval);
 854                   } else {
 855                     lpimap[ilpf] = kFALSE;
 856                   }
 857                   massval[4] = 0;
 858                 }
 859                 ///Lowside v0
 860                 for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
 861                   if(lowgmap[iglf][i1] || lowgmap[iglf][i2]){
 862                     massval[3] = -(iglf+1);
 863                     fMassSparse->Fill(massval);
 864                   }
 865                 }
 866                 massval[3] = 0;
 867               }  /// End lowside mass histo fillers
 868
 869
 870               ///Check that particle is in histo phase space
 871               if(pion->Pt() > fAxistPt.GetBinLowEdge(1) && pion->Pt() < fAxistPt.GetXmax() &&
 872                  pion->M() > fAxisPiM.GetXmin() && pion->M() < fAxisPiM.GetXmax() &&
 873                  pion->Eta() > etalim[0] && pion->Eta() < etalim[1]) {
 874
 875
 876                 const Int_t tbin = fAxistPt.FindFixBin(pion->Pt());
 877
 878                 ///Fill standard triggers including upside v0 filters
 879                 fTrigSparse->Fill(trigValues);
 880
 881
 882                 if(igf1 == 0 && igf2 == 0) {
 883                   // piCorr->FillTriggerCounters(pion);
 884                   // piCorr->CorrelateWithTracks(pion, &tracks[0], tIDs, centrality, vertexz);
 885
 886
 887                   ////Only mix events with pion in signal region
 888                   hTrigPt->Fill(pion->Pt(), centrality, pion->M());
 889                   if(pion->M() > 0.1 && pion->M() < 0.15) {
 890                     hTrigPhi->Fill(pion->Phi());
 891
 892                     ///Check trigger bin
 893                     if (tbin > 0 && tbin < (nbins + 1)) {
 894                       tmap[tbin-1] = kTRUE;
 895                     }
 896
 897                     ///Check if trigger also in low side (both gamma present in low side!)
 898                     for(Int_t ilgf = 0; ilgf < fV0Filters[0].GetEntriesFast(); ilgf++) {
 899                       if(!lowgmap[ilgf][i1] || !lowgmap[ilgf][i2]) {
 900                         lv0tmap[tbin-1][ilgf] = kTRUE;
 901                       }
 902                     }
 903
 904                     ///See if the lowside pion filter also passes this, if not
 905                     for(Int_t ilpf = 0; ilpf < fMesonFilters[0].GetEntriesFast(); ilpf++) {
 906                       if(!lpimap[ilpf]) {
 907                         lpitmap[tbin-1][ilpf] = kTRUE;
 908                       }
 909                     }
 910                   }
 911
 912                 } else {
 913                   if(pion->M() > 0.1 && pion->M() < 0.15) {
 914                     uv0tmap[tbin-1][igf1 - 1] = kTRUE;
 915                   }
 916                 }
 917
 918
 919                 ///Fill the triggers not selected in lowside filters only if passsing standard v0 filter
 920                 if(igf1 == 0 && igf2 == 0) {
 921
 922                   ///Lowside v0 filters
 923                   for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
 924                     if(lowgmap[iglf][i1] || lowgmap[iglf][i2]){
 925                       trigValues[5] = -(iglf+1);
 926                       fTrigSparse->Fill(trigValues);
 927                     }
 928                   }
 929
 930                   ////Low side pion filters
 931                   trigValues[5] = 0;
 932                   for(Int_t iplf = 0; iplf < fMesonFilters[0].GetEntriesFast(); iplf ++) {
 933                     if(lpimap[iplf]) {
 934                       trigValues[6] = -(iplf + 1);
 935                       fTrigSparse->Fill(trigValues);
 936                     }
 937                   }
 938                 } // ifg1 == 0
 939
 940                 trigValues[5] = igf1;
 941                 trigValues[6] = 0;
 942
 943                 ///////////////////////////////////////////////
 944                 /// Correlate with tracks
 945                 ///////////////////////////////////////////////
 946
 947                 Int_t ntf = 1;
 948                 if(igf1 == 0 && igf2 == 0) {
 949                   ntf = fTrackFilters[1].GetEntriesFast() + 1;
 950                 }
 951
 952                 for(Int_t itf = 0; itf < ntf; itf++) {
 953                   TObjArray * tracks = static_cast<TObjArray*>(fTracks.At(itf));
 954                   for(int ij = 0; ij < tracks->GetEntriesFast(); ij++) {
 955                     AliVTrack * track = static_cast<AliVTrack*>(tracks->At(ij));
 956                     Int_t tid = track->GetID();
 957
 958
 959                     if(tid == tIDs[0] || tid == tIDs[1] || tid == tIDs[2] || tid == tIDs[3])  {
 960                       continue;
 961                     }
 962
 963                     if(tid < 0) {
 964                       if(-tid-1 == tIDs[0]+1 || -tid-1 == tIDs[1]+1 || -tid-1 == tIDs[2]+1 || -tid-1 == tIDs[3]+1)  {
 965                         continue;
 966                       }
 967                     }
 968
 969
 970                     dphivalues[0] = pion->Eta() - track->Eta();
 971                     dphivalues[1] = GetDPhi(pion->Phi() - track->Phi());
 972                     dphivalues[3] = track->Pt();
 973                     dphivalues[7] = itf;
 974                     dphivalues[8] = igf1;
 975                     dphivalues[9] = 0;
 976                     fCorrSparse->Fill(dphivalues, GetTrackCorrection(vertexz, track));
 977
 978                     if(itf == 0 && igf1 == 0 && igf2 == 0) {
 979                       ///Fill the low side track filters
 980                       for(Int_t itlf = 0; itlf < fTrackFilters[0].GetEntriesFast(); itlf++) {
 981                         if(lowtrackmap[itlf][ij]){
 982                           dphivalues[7] = -(itlf+1);
 983                           fCorrSparse->Fill(dphivalues, GetTrackCorrection(vertexz, track));
 984                         }
 985                       }
 986                       ///Fill the low side v0 filters
 987                       dphivalues[7] = 0;
 988                       for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
 989                         if(lowgmap[iglf][i1] || lowgmap[iglf][i2]){
 990                           dphivalues[8] = -(iglf+1);
 991                           fCorrSparse->Fill(dphivalues, GetTrackCorrection(vertexz, track));
 992                         }
 993                       }
 994
 995                       ///Fill the low side pi filter
 996                       dphivalues[7] = 0;
 997                       dphivalues[8] = 0;
 998                       for(Int_t iplf = 0; iplf < fMesonFilters[0].GetEntriesFast(); iplf ++) {
 999                         if(lpimap[iplf]) {
1000                           dphivalues[9] = -(iplf + 1);
1001                           fCorrSparse->Fill(dphivalues, GetTrackCorrection(vertexz, track));
1002                         }
1003                       }
1004                     }  /// end non standard filters track corr
1005
1006                   } // end for tracks
1007                 } // end trackfilters loop
1008               } //end check pion in histogram range to prevent overflow
1009             } else {
1010               /////////////////////////////
1011               //// Not passing standard meson cuts, check upside filters
1012               ////////////////////////////
1013
1014               ///Only check the pions from standard v0 filter
1015               if(igf1 == 0 && igf2 == 0) {
1016                 for(Int_t ipuf = 0; ipuf < fMesonFilters[1].GetEntriesFast(); ipuf++) {
1017                   if(static_cast<AliConversionMesonCuts*>(fMesonFilters[1].At(ipuf))->MesonIsSelected(pion, kTRUE)) {
1018                     ///Fill invariant mass hist
1019                     massval[4] = (ipuf + 1);
1020                     fMassSparse->Fill(massval);
1021                     ///Check that particle is in histo phase space --- redundant!
1022
1023
1024
1025                     if(pion->Pt() > fAxistPt.GetBinLowEdge(1) && pion->Pt() < fAxistPt.GetXmax() &&
1026                        pion->M() > fAxisPiM.GetXmin() && pion->M() < fAxisPiM.GetXmax() &&
1027                        pion->Eta() > etalim[0] && pion->Eta() < etalim[1]) {
1028
1029
1030                       ////Only mix events with pion in signal region
1031                       if(pion->M() > 0.1 && pion->M() < 0.15) {
1032                         const Int_t tbin = fAxistPt.FindFixBin(pion->Pt());
1033                         upitmap[tbin-1][ipuf] = kTRUE;
1034                       }
1035
1036                       ///Fill trigger counters
1037                       trigValues[6] = (ipuf + 1);
1038                       fTrigSparse->Fill(trigValues);
1039
1040                       ///Correlate with standard tracks
1041                       for(int ij = 0; ij < ttracks->GetEntriesFast(); ij++) {
1042                         AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(ij));
1043                         Int_t tid = track->GetID();
1044
1045                         if(tid == tIDs[0] || tid == tIDs[1] || tid == tIDs[2] || tid == tIDs[3] ) {
1046                           continue;
1047                         }
1048
1049                         if(tid < 0) {
1050                           if(-tid-1 == tIDs[0]+1 || -tid-1 == tIDs[1]+1 || -tid-1 == tIDs[2]+1 || -tid-1 == tIDs[3]+1)  {
1051                             continue;
1052                           }
1053                         }
1054
1055                         dphivalues[0] = pion->Eta() - track->Eta();
1056                         dphivalues[1] = GetDPhi(pion->Phi() - track->Phi());
1057                         dphivalues[3] = track->Pt();
1058                         dphivalues[7] = 0; // track filter
1059                         dphivalues[8] = 0; // v0 filter
1060                         dphivalues[9] = ipuf + 1; // pion filter
1061                         fCorrSparse->Fill(dphivalues, GetTrackCorrection(vertexz, track));
1062                       } /// end track corr
1063                     }
1064                   } // MesonIsSelected
1065                 }
1066               }
1067             } /// end else ..  end upside meson filters
1068             /////////////////////////////////////////////
1069           } ///Etalim && pt limits
1070
1071
1072         } // i2 second gamma
1073       }
1074     } // i1 first gamma
1075   }
1076   //FillCounters(&pions, tracks, ntrackfilters, centrality, vertexz);
1077
1078   ///////Fill track counters after entire event has been passed through
1079   ////
1080
1081   Double_t trackValues[fTrackAxesList.GetSize()];
1082   trackValues[3] = centrality;
1083   trackValues[4] = vertexz;
1084   //trackValues[5] = particle->M(); remove !!!
1085
1086   for(Int_t tbin = 0; tbin < fAxistPt.GetNbins(); tbin++) {
1087     trackValues[1] = fAxistPt.GetBinCenter(tbin+1);
1088
1089     if(tmap[tbin]) {
1090
1091       for(Int_t itf = 0; itf < fTrackFilters[1].GetEntriesFast() + 1; itf++) {
1092         TObjArray * tracks = static_cast<TObjArray*>(fTracks.At(itf));
1093         for(Int_t iTrack = 0; iTrack < tracks->GetEntriesFast(); iTrack++) {
1094           AliVTrack * track = static_cast<AliVTrack*>(tracks->At(iTrack));
1095           trackValues[0] = track->Eta();
1096           trackValues[2] = track->Pt();
1097           trackValues[5] = itf;
1098           trackValues[6] = 0;  ///v0 filter
1099           trackValues[7] = 0; ////Pi filter
1100           fTrackSparse->Fill(trackValues, GetTrackCorrection(vertexz, track));
1101
1102           if(itf == 0) {
1103             for(Int_t itlf = 0; itlf < fTrackFilters[0].GetEntriesFast(); itlf++) {
1104               if(lowtrackmap[itlf][iTrack]) {
1105                 trackValues[5] = -(itlf + 1);
1106                 fTrackSparse->Fill(trackValues, GetTrackCorrection(vertexz, track) );
1107               }
1108             }
1109             trackValues[5] = 0;
1110
1111             ///Check lowside gamma
1112             for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
1113               if(!lv0tmap[tbin][iglf]) {
1114                 trackValues[6] = -(iglf + 1);
1115                 fTrackSparse->Fill(trackValues, GetTrackCorrection(vertexz, track));
1116               }
1117             }
1118             trackValues[6] = 0;
1119
1120             ////Check lowside pion
1121             for(Int_t iplf = 0; iplf < fMesonFilters[0].GetEntriesFast(); iplf++) {
1122               if(!lpitmap[tbin][iplf]) {
1123                 trackValues[7] = -(iplf + 1);
1124                 fTrackSparse->Fill(trackValues, GetTrackCorrection(vertexz, track));
1125               }
1126             }
1127
1128
1129
1130           } // itf == 00
1131         }
1132       }
1133     } else {
1134       ///If not in main, see if in upside filters
1135       ///Do upside v0 filters
1136       for(Int_t iguf = 0; iguf < fV0Filters[1].GetEntriesFast(); iguf++) {
1137         if (uv0tmap[tbin][iguf] ) {
1138
1139           //cout << "c vtx " <<  centrality << vertexz << endl;
1140
1141           for(Int_t iTrack = 0; iTrack < ttracks->GetEntriesFast(); iTrack++) {
1142             AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(iTrack));
1143             trackValues[0] = track->Eta();
1144             trackValues[2] = track->Pt();
1145             trackValues[5] = 0;
1146             trackValues[6] = iguf+1;  ///v0 filter
1147             trackValues[7] = 0; ////Pi filter
1148             fTrackSparse->Fill(trackValues, GetTrackCorrection(vertexz, track));
1149           }
1150         }
1151       }
1152       ///Do upside pi filter
1153       for(Int_t ipuf = 0; ipuf < fMesonFilters[1].GetEntriesFast(); ipuf++) {
1154         if (upitmap[tbin][ipuf] ) {
1155
1156           //cout << "c2 vtx2 " <<  centrality << vertexz << endl;
1157
1158
1159           for(Int_t iTrack = 0; iTrack < ttracks->GetEntriesFast(); iTrack++) {
1160             AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(iTrack));
1161             trackValues[0] = track->Eta();
1162             trackValues[2] = track->Pt();
1163             trackValues[5] = 0;
1164             trackValues[6] = 0;  ///v0 filter
1165             trackValues[7] = ipuf+1; ////Pi filter
1166             fTrackSparse->Fill(trackValues, GetTrackCorrection(vertexz, track));
1167           }
1168         }
1169       }
1170     }
1171   }
1172
1173
1174
1175   //////
1176   //////
1177
1178
1179
1180   //piCorr->FillCounters(&pions, tracks, centrality, vertexz);
1181
1182   PostData(1, fHistograms);
1183
1184 }
1185
1186 //_______________________________________________________________________________
1187 // void AliAnalysisTaskdPhi::FillCounters(TObjArray * particles, TObjArray tracks[], Int_t ntrackfilters, Float_t cent, Float_t vtxz) {
1188
1189 //   return;
1190
1191
1192 //  //Fill ME Counters
1193 //   const Int_t nbins = fAxistPt.GetNbins();
1194 //   Bool_t tmap[nbins];
1195 //   for(Int_t ptbin = 0; ptbin < nbins; ptbin++){
1196 //       tmap[ptbin] = kFALSE;
1197 //     }
1198 //   }
1199
1200
1201 //   Double_t trackValues[fTrackAxesList.GetSize()];
1202 //   trackValues[3] = cent;
1203 //   trackValues[4] = vtxz;
1204
1205 //   for(Int_t ip = 0; ip < particles->GetEntriesFast(); ip++){
1206 //     AliAODConversionParticle * particle = static_cast<AliAODConversionParticle*>(particles->At(ip));
1207
1208 //     Int_t tbin = fAxistPt.FindFixBin(particle->Pt());
1209 //     if (tbin > 0 && tbin < nbins + 1) {
1210 //       if(tmap[tbin - 1] == kTRUE) {
1211 //      continue;
1212 //       } else {
1213 //      tmap[tbin -1 ] = kTRUE;
1214
1215 //      trackValues[5] = particle->M();
1216 //      trackValues[1] = particle->Pt();
1217
1218 //      for(Int_t itf = 0; itf < ntrackfilters; itf++) {
1219 //        if(fkTrackAxis) trackValues[6] = itf;
1220 //        for(int ij = 0; ij < tracks[itf].GetEntriesFast(); ij++) {
1221 //          AliVTrack * track = static_cast<AliVTrack*>(tracks->UncheckedAt(ij));
1222 //          trackValues[0] = track->Eta();
1223 //          trackValues[2] = track->Pt();
1224 //          fTrackSparse->Fill(trackValues);
1225 //        }
1226 //      }
1227 //       }
1228 //     }
1229 //   }
1230 // }
1231
1232
1233
1234
1235 // //________________________________________________________________
1236 // void AliAnalysisTaskdPhi::CorrelateWithTracks(AliAODConversionParticle * particle, TObjArray tracks[], Int_t ntrackfilters, Bool_t ** lowtrackmap, Int_t nltf, Int_t const tIDs[4], Double_t dphivalues[]) {
1237 //   //Correlate particle with tracks
1238 //   ///Correlate with tracks
1239
1240 // }
1241
1242 //________________________________________________________________________
1243 void AliAnalysisTaskdPhi::Terminate(Option_t *) {
1244
1245   // Draw result to the screen
1246   // Called once at the end of the query
1247 }
1248
1249 //________________________________________________________________________
1250 TClonesArray * AliAnalysisTaskdPhi::GetConversionGammas(Bool_t isAOD) {
1251   if(isAOD) {
1252
1253     TClonesArray * gammas = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fDeltaAODBranchName.Data()));
1254     if(gammas) {
1255       return gammas;
1256     }
1257
1258     FindDeltaAODBranchName(fInputEvent);
1259     gammas = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fDeltaAODBranchName.Data()));
1260     return gammas;
1261
1262   } else {
1263     TClonesArray * gammas = dynamic_cast<TClonesArray*>(GetInputData(1));
1264     return gammas;
1265   }
1266
1267 }
1268
1269 //________________________________________________________________________
1270 void AliAnalysisTaskdPhi::FindDeltaAODBranchName(AliVEvent * event){
1271   ///Find aod branch
1272   TList *list=event->GetList();
1273   for(Int_t ii=0;ii<list->GetEntries();ii++){
1274     TString name((list->At(ii))->GetName());
1275     if(name.BeginsWith("GammaConv")&&name.EndsWith("gamma")){
1276       fDeltaAODBranchName=name;
1277       AliDebug(AliLog::kDebug + 5, Form("Set DeltaAOD BranchName to: %s",fDeltaAODBranchName.Data()));
1278       return;
1279     }
1280   }
1281 }
1282
1283
1284 //________________________________________________________________________
1285 Double_t AliAnalysisTaskdPhi::GetTrackCorrection(Double_t vtxz, AliVTrack * track) {
1286   ////Get track correction from map
1287   Int_t coord[4] = {-1, -1, -1, -1};
1288   if(fCorrectionMap) {
1289     Double_t values[4] = { vtxz, track->Pt(), track->Eta(), track->Phi() };
1290     Double_t correction = fCorrectionMap->GetBinContent(fCorrectionMap->GetBin(values, kFALSE), coord);
1291     if (fCorrectionMap->IsInRange(coord)) {
1292       return correction;
1293     }
1294   }
1295   return 1.0;
1296 }
1297
1298