1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Authors: Svein Lindal *
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 **************************************************************************/
17 ////////////////////////////////////////////////
18 //---------------------------------------------
19 // Class doing conversion gamma dPhi correlations
20 // Gamma Conversion analysis
21 //---------------------------------------------
22 ////////////////////////////////////////////////
24 #include "AliAnalysisTaskdPhi.h"
30 #include <AliAnalysisManager.h>
31 #include <AliInputEventHandler.h>
32 #include <AliESDInputHandler.h>
33 #include <AliAODInputHandler.h>
34 #include <TGeoGlobalMagField.h>
36 #include "AliConversionTrackCuts.h"
37 #include "AliConversionCuts.h"
38 #include "AliConversionMesonCuts.h"
39 #include "AliAODConversionPhoton.h"
40 #include "AliAODConversionMother.h"
41 // #include "AliAnaConvCorrPhoton.h"
42 // #include "AliAnaConvCorrPion.h"
43 // #include "AliAnaConvIsolation.h"
44 #include "AliV0ReaderV1.h"
45 // Author Svein Lindal <slindal@fys.uio.no>
48 ClassImp(AliAnalysisTaskdPhi)
51 //________________________________________________________________________
52 AliAnalysisTaskdPhi::AliAnalysisTaskdPhi(const char *name) : AliAnalysisTaskSE(name),
59 fSaveReaderHists(kFALSE),
73 fDeltaAODBranchName("AliAODGammaConversion_gamma"),
98 DefineInput(0, TChain::Class());
99 DefineOutput(1, TList::Class());
101 fGammas.SetOwner(kTRUE);
102 fTracks.SetOwner(kTRUE);
108 //________________________________________________________________________
109 AliAnalysisTaskdPhi::~AliAnalysisTaskdPhi(){
121 delete fPhotonFilter;
122 fPhotonFilter = NULL;
137 ///________________________________________________________________________
138 void AliAnalysisTaskdPhi::SetUpBins() {
140 fAxisTrigEta.SetNameTitle("tEta", "Eta");
141 fAxisTrigEta.Set(320, -0.8, 0.8);
143 fAxisAssEta.SetNameTitle("aEta", "Eta");
144 fAxisAssEta.Set(360, -0.9, 0.9);
146 fAxisdEta.SetNameTitle("dEta", "Eta");
147 fAxisdEta.Set(34, -1.7, 1.7);
149 fAxisdPhi.SetNameTitle("dPhi", "delta Phi");
150 fAxisdPhi.Set(32, -TMath::PiOver2(), 3*TMath::PiOver2());
152 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");
153 fAxistPt.Set(9, tptbins);
155 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");
156 fAxiscPt.Set(12, cptbins);
158 fAxisZ.SetNameTitle("vertexz", "Z");
159 fAxisZ.Set(4, -10, 10);
161 Double_t centbins[5] = {0, 10, 30, 60, 100.1};
162 fAxisCent.SetNameTitle("centrality", "Cent");
163 fAxisCent.Set(4, centbins);
165 fAxisPiM.SetNameTitle("InvMassPi0", "Invariant mass");
166 Double_t mbins[7] = {0.1, 0.11, 0.12, 0.15, 0.16, 0.18, 0.2};
167 fAxisPiM.Set(6, mbins);
172 ///________________________________________________________________________
173 // void AliAnalysisTaskdPhi::SetUpCorrObjects() {
174 // //Set up corr objects
175 // AliDebug(AliLog::kDebug + 5, "Set Up corr objects");
178 // fPhotonCorr = new AliAnaConvCorrPhoton("PhotonCorr","photon %s");
179 // SetUpCorrAxes(fPhotonCorr);
180 // fPhotonCorr->CreateHistograms();
181 // fHistograms->Add(fPhotonCorr->GetHistograms());
184 // fPionCorr = new AliAnaConvCorrPion("PionCorr", "pion");
185 // SetUpCorrAxes(fPionCorr);
186 // fPionCorr->GetAxisM().Set(fAxisPiM.GetNbins(), fAxisPiM.GetXbins()->GetArray());
187 // fPionCorr->CreateHistograms();
188 // fHistograms->Add(fPionCorr->GetHistograms());
191 ///________________________________________________________________________
192 // void AliAnalysisTaskdPhi::SetUpCorrAxes(AliAnaConvCorrBase * corr) {
193 ///Set up axes in corr object
194 // corr->GetAxisCent().Set(fAxisCent.GetNbins(), fAxisCent.GetXbins()->GetArray());
195 // const Double_t * zbins = fAxisZ.GetXbins()->GetArray();
197 // corr->GetAxisZ().Set(fAxisZ.GetNbins(), fAxisZ.GetXbins()->GetArray());
199 // corr->GetAxisZ().Set(fAxisZ.GetNbins(), fAxisZ.GetBinLowEdge(1), fAxisZ.GetBinUpEdge(fAxisZ.GetNbins()));
202 // corr->GetAxistPt().Set(fAxistPt.GetNbins(), fAxistPt.GetXbins()->GetArray());
203 // corr->GetAxiscPt().Set(fAxiscPt.GetNbins(), fAxiscPt.GetXbins()->GetArray());
204 // corr->GetAxisdEta().Set(fAxisdEta.GetNbins(), fAxisdEta.GetBinLowEdge(1), fAxisdEta.GetBinUpEdge(fAxisdEta.GetNbins()));
205 // corr->GetAxisTrigEta().Set(fAxisTrigEta.GetNbins(), fAxisTrigEta.GetBinLowEdge(1), fAxisTrigEta.GetBinUpEdge(fAxisTrigEta.GetNbins()));
206 // corr->GetAxisAssEta().Set(fAxisAssEta.GetNbins(), fAxisAssEta.GetBinLowEdge(1), fAxisAssEta.GetBinUpEdge(fAxisAssEta.GetNbins()));
210 //________________________________________________________________________
211 void AliAnalysisTaskdPhi::UserCreateOutputObjects() {
214 fHistograms = new TList();
215 fHistograms->SetName("dPhi_histograms");
216 fHistograms->SetOwner(kTRUE);
220 fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
224 printf("Error: No V0 Reader");
228 if(fSaveReaderHists) {
229 AliConversionCuts * v0cuts = fV0Reader->GetConversionCuts();
231 TList * histograms = v0cuts->GetCutHistograms();
233 AliWarning("initializing v0 reader hists");
234 v0cuts->InitCutHistograms("V0Reader", kTRUE);
236 histograms = v0cuts->GetCutHistograms();
238 fHistograms->Add(histograms);
243 for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf ++){
244 AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[0].At(igf));
246 f->InitCutHistograms(Form("V0Filter_%d", -(igf+1)), kFALSE);
247 fHistograms->Add(f->GetCutHistograms());
251 for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf ++){
252 AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[1].At(igf));
254 f->InitCutHistograms(Form("V0Filter_%d", igf+1), kFALSE);
255 fHistograms->Add(f->GetCutHistograms());
259 for(Int_t igf = 0; igf < fMesonFilters[0].GetEntriesFast(); igf ++){
260 AliConversionMesonCuts * f = dynamic_cast<AliConversionMesonCuts*>(fMesonFilters[0].At(igf));
262 f->InitCutHistograms(Form("PionFilter_%d", -(igf+1)), kFALSE);
263 fHistograms->Add(f->GetCutHistograms());
267 for(Int_t igf = 0; igf < fMesonFilters[1].GetEntriesFast(); igf ++){
268 AliConversionMesonCuts * f = dynamic_cast<AliConversionMesonCuts*>(fMesonFilters[1].At(igf));
270 f->InitCutHistograms(Form("PionFilter_%d", igf+1), kFALSE);
271 fHistograms->Add(f->GetCutHistograms());
276 fV0Filter->InitCutHistograms("V0Filter", kFALSE);
277 fHistograms->Add(fV0Filter->GetCutHistograms());
280 fMesonFilter->InitCutHistograms("PionFilter", kFALSE);
281 fHistograms->Add(fMesonFilter->GetCutHistograms());
284 fPhotonFilter->InitCutHistograms("PhotonFilter", kFALSE);
285 fHistograms->Add(fPhotonFilter->GetCutHistograms());
289 AliConversionTrackCuts * tc = dynamic_cast<AliConversionTrackCuts*>(fTrackFilter);
291 fHistograms->Add(tc->CreateHistograms());
295 // for(Int_t i = 0; i < fTrackFilters.GetEntriesFast(); i++) {
296 // AliConversionTrackCuts * tc = dynamic_cast<AliConversionTrackCuts*>(fTrackFilters.At(i));
297 // if(tc) fHistograms->Add(tc->CreateHistograms());
300 //SetUpCorrObjects();
302 ///Set up ME histograms
303 TList * MEHistograms = new TList();
304 MEHistograms->SetName("MEHistograms");
305 MEHistograms->SetOwner(kTRUE);
306 fHistograms->Add(MEHistograms);
308 hMEvents = new TH2I("hMEvents", "Nevents vs centrality vertexz",
309 fAxisZ.GetNbins(), fAxisZ.GetXbins()->GetArray(),
310 fAxisCent.GetNbins(), fAxisCent.GetXbins()->GetArray());
311 MEHistograms->Add(hMEvents);
313 hTrackCent = new TH2I("hTrackCent", "N accepted tracks vs centrality",
314 fAxisCent.GetNbins() > 2 ? 100 : 1, fAxisCent.GetBinLowEdge(1), fAxisCent.GetBinUpEdge(fAxisCent.GetNbins()),
316 MEHistograms->Add(hTrackCent);
318 Int_t ntrackfilters[2] = {fTrackFilters[0].GetEntriesFast(), fTrackFilters[1].GetEntriesFast()};
320 fAxisTrackFilters.SetNameTitle("trackCuts", "trackCuts");
321 fAxisTrackFilters.Set(ntrackfilters[0] + ntrackfilters[1] + 1, -ntrackfilters[0] -0.5, ntrackfilters[1] + 0.5);
323 Int_t nV0filters[2] = {fV0Filters[0].GetEntriesFast(), fV0Filters[1].GetEntriesFast()};
325 fAxisV0Filters.SetNameTitle("V0Cuts", "V0Cuts");
326 fAxisV0Filters.Set(nV0filters[0] + nV0filters[1] + 1, -nV0filters[0] -0.5, nV0filters[1] + 0.5);
328 Int_t nmesonfilters[2] = {fMesonFilters[0].GetEntriesFast(), fMesonFilters[1].GetEntriesFast()};
330 fAxisMesonFilters.SetNameTitle("mesonCuts", "mesonCuts");
331 fAxisMesonFilters.Set(nmesonfilters[0] + nmesonfilters[1] + 1, -nmesonfilters[0] -0.5, nmesonfilters[1] + 0.5);
333 fAxesList.AddAt(&fAxisdEta, 0);
334 fAxesList.AddAt(&fAxisdPhi, 1);
335 fAxesList.AddAt(&fAxistPt, 2);
336 fAxesList.AddAt(&fAxiscPt, 3);
337 fAxesList.AddAt(&fAxisCent, 4);
338 fAxesList.AddAt(&fAxisZ, 5);
339 fAxesList.AddAt(&fAxisPiM, 6);
340 fAxesList.AddAt(&fAxisTrackFilters, 7);
341 fAxesList.AddAt(&fAxisV0Filters, 8);
342 fAxesList.AddAt(&fAxisMesonFilters, 9);
344 fTrackAxesList.AddAt(&fAxisAssEta, 0);
345 fTrackAxesList.AddAt(&fAxistPt, 1);
346 fTrackAxesList.AddAt(&fAxiscPt, 2);
347 fTrackAxesList.AddAt(&fAxisCent, 3);
348 fTrackAxesList.AddAt(&fAxisZ, 4);
349 //fTrackAxesList.AddAt(&fAxisPiM, 5);
350 fTrackAxesList.AddAt(&fAxisTrackFilters, 5);
351 fTrackAxesList.AddAt(&fAxisV0Filters, 6);
352 fTrackAxesList.AddAt(&fAxisMesonFilters, 7);
354 fTrigAxesList.AddAt(&fAxisTrigEta, 0);
355 fTrigAxesList.AddAt(&fAxistPt, 1);
356 fTrigAxesList.AddAt(&fAxisCent, 2);
357 fTrigAxesList.AddAt(&fAxisZ, 3);
358 fTrigAxesList.AddAt(&fAxisPiM, 4);
359 fTrigAxesList.AddAt(&fAxisV0Filters, 5);
360 fTrigAxesList.AddAt(&fAxisMesonFilters, 6);
364 massax.SetNameTitle("mass", "mass");
365 massax.Set(360, 0.04, 0.4); //hardcoded! change also in filling!
367 masslist.AddAt(&massax, 0);
368 masslist.AddAt(&fAxistPt, 1);
369 masslist.AddAt(&fAxisCent, 2);
370 masslist.AddAt(&fAxisV0Filters, 3);
371 masslist.AddAt(&fAxisMesonFilters, 4);
373 fCorrSparse = CreateSparse(TString("pionSparse"), TString("pionSparse"), &fAxesList);
374 fTrackSparse = CreateSparse(TString("trackSparse"), TString("trackSparse"), &fTrackAxesList);
375 fTrigSparse = CreateSparse(TString("trigSparse"), TString("trigSparse"), &fTrigAxesList);
376 fMassSparse = CreateSparse("massSparse", "massSparse", &masslist);
378 fHistograms->Add(fCorrSparse);
379 fHistograms->Add(fTrackSparse);
380 fHistograms->Add(fTrigSparse);
381 fHistograms->Add(fMassSparse);
384 ///Add gamma and track containers:
385 for(Int_t i = 0; i < fV0Filters[1].GetEntriesFast() + 1; i++) {
386 fGammas.Add(new TObjArray());
389 for(Int_t i = 0; i < fTrackFilters[1].GetEntriesFast() + 1; i++) {
390 fTracks.Add(new TObjArray());
395 PostData(1, fHistograms);
398 ///________________________________________________________________________
399 THnSparseF * AliAnalysisTaskdPhi::CreateSparse(TString nameString, TString titleString, TList * axesList) {
401 const Int_t dim = axesList->GetSize();
408 for(Int_t i = 0; i<dim; i++) {
409 TAxis * axis = dynamic_cast<TAxis*>(axesList->At(i));
410 if(axis) axes[i] = axis;
412 cout << "AliAnalysisTaskdPhi::CreateSparse: Error error, all the axes are not present in axis list" << endl;
417 for(Int_t i = 0; i<dim; i++) {
418 //cout << axes[i]->GetTitle() << endl;
419 bins[i] = axes[i]->GetNbins();
420 min[i] = axes[i]->GetBinLowEdge(1);
421 max[i] = axes[i]->GetBinUpEdge(axes[i]->GetNbins());
424 THnSparseF * sparse = new THnSparseF(Form("%s", nameString.Data()),
425 Form("%s", titleString.Data()),
426 dim, bins, min, max);
428 for(Int_t i = 0; i<dim; i++) {
429 sparse->GetAxis(i)->SetNameTitle(axes[i]->GetName(), axes[i]->GetTitle() );
430 if(axes[i]->GetXbins()->GetSize() > 0) {
431 sparse->SetBinEdges(i, axes[i]->GetXbins()->GetArray() );
437 //________________________________________________________________________
438 void AliAnalysisTaskdPhi::UserExec(Option_t *) {
440 ///This is a very ugly function, cut the complexity of the logic demands it.
443 AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
445 if (!TGeoGlobalMagField::Instance()->GetField()) esdEvent->InitMagneticField();
449 //if(! fV0Filter->EventIsSelected(fInputEvent)) return;
451 AliError("Error: No V0 Reader");
455 if(!fV0Reader->IsEventSelected()) {
458 AliDebug(5, "Processing event");
461 for(Int_t i = 0; i < fGammas.GetEntriesFast(); i++) {
462 static_cast<TObjArray*>(fGammas.At(i))->Clear();
465 for(Int_t i = 0; i < fTracks.GetEntriesFast(); i++) {
466 static_cast<TObjArray*>(fTracks.At(i))->Clear();
471 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
472 Bool_t isAOD=man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class();
474 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
476 cout << "cout no input event handler"<<endl;
483 for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf++) {
484 AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[0].At(igf));
485 if ( f && !f->GetPIDResponse() ) {
486 if ( inputHandler->GetPIDResponse() ){
487 f->SetPIDResponse( inputHandler->GetPIDResponse() );
490 if (!f->GetPIDResponse()){
491 f->InitAODpidUtil(1);
500 for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf++) {
501 AliConversionCuts * f = dynamic_cast<AliConversionCuts*>(fV0Filters[1].At(igf));
502 if ( f && !f->GetPIDResponse() ) {
503 if ( inputHandler->GetPIDResponse() ){
504 f->SetPIDResponse( inputHandler->GetPIDResponse() );
507 if (!f->GetPIDResponse()){
508 f->InitAODpidUtil(1);
518 if ( fV0Filter && !fV0Filter->GetPIDResponse() ) {
519 if ( inputHandler->GetPIDResponse() ){
520 fV0Filter->SetPIDResponse( inputHandler->GetPIDResponse() );
525 if (!fV0Filter->GetPIDResponse()){
526 fV0Filter->InitAODpidUtil(1);
532 AliConversionTrackCuts * tc = dynamic_cast<AliConversionTrackCuts*>(fTrackFilter);
534 tc->SetEvent(fInputEvent);
538 Double_t centrality = 0.0;
539 Double_t eventPlane = 0.0;
540 Double_t vertexz = fInputEvent->GetPrimaryVertex()->GetZ();
542 AliAODHeader * header = static_cast<AliAODHeader*>(fInputEvent->GetHeader());
543 centrality = header->GetCentrality();
544 eventPlane = header->GetEventplane();
546 centrality = static_cast<AliESDEvent*>(fInputEvent)->GetCentrality()->GetCentralityPercentile("V0M");
547 eventPlane = fInputEvent->GetEventplane()->GetEventplane("Q");
551 const Int_t centBin = GetBin(fAxisCent, centrality);
552 const Int_t vertexBin = GetBin(fAxisZ, vertexz);
555 if(DebugLevel () > 4) {
556 cout << "centrality: " << centrality << " " << GetBin(fAxisCent, centrality) << endl;
557 cout << "vertexz: " << vertexz << " " << GetBin(fAxisZ, vertexz) << endl;
558 cout << "eventPlane: " << eventPlane << " " << endl;
563 if(centBin < 0 || vertexBin < 0) {
564 // AliError("bin out of range");
570 //TClonesArray * aodGammas = GetConversionGammas(isAOD);
571 TClonesArray * aodGammas = fV0Reader->GetReconstructedGammas();
573 AliError("no aod gammas found!");
577 TObjArray * ggammas = static_cast<TObjArray*>(fGammas.At(0));
579 ///Fill arrays of accepted gammas
580 if(DebugLevel() > 1) printf("Number of conversion gammas %d \n", aodGammas->GetEntriesFast());
581 for(Int_t ig = 0; ig < aodGammas->GetEntriesFast(); ig++) {
582 AliAODConversionPhoton * photon = dynamic_cast<AliAODConversionPhoton*>(aodGammas->At(ig));
583 if(!photon) continue;
584 if(!fV0Filter || fV0Filter->PhotonIsSelected(photon, fInputEvent)) {
585 ggammas->Add(photon);
587 for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf++) {
588 AliConversionCuts * gfilter = dynamic_cast<AliConversionCuts*>(fV0Filters[1].At(igf));
589 if(gfilter && gfilter->PhotonIsSelected(photon, fInputEvent)) {
590 static_cast<TObjArray*>(fGammas.At(igf+1))->Add(photon);
596 Bool_t lowgmap[fV0Filters[0].GetEntriesFast()][ggammas->GetEntriesFast()];
598 for(Int_t ig = 0; ig < ggammas->GetEntriesFast(); ig++ ) {
599 AliAODConversionPhoton * gamma = static_cast<AliAODConversionPhoton*>(ggammas->At(ig));
600 for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf++) {
602 AliConversionCuts * v0cuts = static_cast<AliConversionCuts*>(fV0Filters[0].At(igf));
603 if(!(v0cuts->PhotonIsSelected(gamma, fInputEvent))) {
604 lowgmap[igf][ig] = kTRUE;
606 lowgmap[igf][ig] = kFALSE;
612 if(DebugLevel() > 4) printf("Number of accepted gammas %d \n", ggammas->GetEntriesFast());
613 hMEvents->Fill(vertexz, centrality);
615 ///create track array
616 const Int_t ntrackfilters[2] = { fTrackFilters[0].GetEntriesFast(), fTrackFilters[1].GetEntriesFast()};
619 TObjArray * ttracks = static_cast<TObjArray*>(fTracks.At(0));
620 const Double_t aetalim[2] = { fAxisAssEta.GetXmin(), fAxisAssEta.GetXmax()};
621 const Double_t aptlim[2] = { fAxiscPt.GetXmin(), fAxiscPt.GetXmax()};
622 for(Int_t iTrack = 0; iTrack < fInputEvent->GetNumberOfTracks(); iTrack++) {
623 AliVTrack * track = static_cast<AliVTrack*>(fInputEvent->GetTrack(iTrack));
624 if(track->Pt() < aptlim[0] || track->Pt() > aptlim[1]) continue;
625 if(track->Eta() < aetalim[0] || track->Eta() > aetalim[1]) continue;
626 if(fTrackFilter->IsSelected(track)) {
630 for(Int_t itf = 1; itf < ntrackfilters[1] + 1; itf++) {
631 AliAnalysisCuts * trackCuts = static_cast<AliAnalysisCuts*>(fTrackFilters[1].At(itf -1));
632 if(trackCuts->IsSelected(track)) {
633 static_cast<TObjArray*>(fTracks.At(itf))->Add(track);
640 Bool_t lowtrackmap[ntrackfilters[0]][ttracks->GetEntriesFast()];
641 ///Check lowside cuts
642 for(Int_t iTrack = 0; iTrack < ttracks->GetEntriesFast(); iTrack++ ) {
643 AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(iTrack));
644 for(Int_t itf = 0; itf < ntrackfilters[0]; itf++) {
645 AliAnalysisCuts * trackCuts = static_cast<AliAnalysisCuts*>(fTrackFilters[0].At(itf));
646 if(!trackCuts->IsSelected(track)) {
647 lowtrackmap[itf][iTrack] = kTRUE;
649 lowtrackmap[itf][iTrack] = kFALSE;
654 hTrackCent->Fill(centrality, ttracks->GetEntriesFast());
656 const Double_t etalim[2] = { fAxisTrigEta.GetXmin(), fAxisTrigEta.GetXmax()};
657 if(DebugLevel() > 4) printf("Number of accepted gammas, tracks %d %d \n", ggammas->GetEntriesFast(), ttracks->GetEntriesFast());
659 //AliAnaConvCorrBase * gCorr = fPhotonCorr; //GetCorrObject(vertexBin, centBin, fPhotonCorr);
660 // AliAnaConvCorrPion * piCorr = fPionCorr; //static_cast<AliAnaConvCorrPion*>(GetCorrObject(vertexBin, centBin, fPionCorr));
663 // AliError("corr object missing");
670 Double_t dphivalues[fAxesList.GetSize()];
671 dphivalues[4] = centrality;
672 dphivalues[5] = vertexz;
674 ///Trigger me counters
675 Double_t trigValues[7];
676 trigValues[2] = centrality;
677 trigValues[3] = vertexz;
681 massval[2] = centrality;
683 ///Set up track me counters and initialize
684 const Int_t nbins = fAxistPt.GetNbins();
685 Bool_t tmap[fAxistPt.GetNbins()];
686 Bool_t lv0tmap[fAxistPt.GetNbins()][fV0Filters[0].GetEntriesFast()];
687 Bool_t uv0tmap[fAxistPt.GetNbins()][fV0Filters[1].GetEntriesFast()];
689 Bool_t lpitmap[fAxistPt.GetNbins()][fMesonFilters[0].GetEntriesFast()];
690 Bool_t upitmap[fAxistPt.GetNbins()][fMesonFilters[1].GetEntriesFast()];
694 for(Int_t igf = 0; igf < fV0Filters[0].GetEntriesFast(); igf++) {
695 for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
696 lv0tmap[ptbin][igf] = kFALSE;
700 for(Int_t igf = 0; igf < fV0Filters[1].GetEntriesFast(); igf++) {
701 for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
702 uv0tmap[ptbin][igf] = kFALSE;
706 for(Int_t igf = 0; igf < fMesonFilters[0].GetEntriesFast(); igf++) {
707 for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
708 lpitmap[ptbin][igf] = kFALSE;
712 for(Int_t igf = 0; igf < fMesonFilters[1].GetEntriesFast(); igf++) {
713 for(Int_t ptbin = 0; ptbin < fAxistPt.GetNbins(); ptbin++) {
714 upitmap[ptbin][igf] = kFALSE;
718 for(Int_t ptbin = 0; ptbin < nbins; ptbin++){
719 tmap[ptbin] = kFALSE;
722 /////////////////////////////////////////
725 for(Int_t igf1 = 0; igf1 < fV0Filters[1].GetEntriesFast() + 1; igf1++) {
726 TObjArray * gamm1 = static_cast<TObjArray*>(fGammas.At(igf1));
727 for(Int_t i1 = 0; i1 < gamm1->GetEntriesFast(); i1++) {
728 AliAODConversionPhoton * ph1 = static_cast<AliAODConversionPhoton*>(gamm1->UncheckedAt(i1));
729 Int_t tIDs[4] = {ph1->GetLabel(0), ph1->GetLabel(1), -1, -1};
731 ///Combine gamma into pions
733 for(Int_t igf2 = 0; igf2 <= igf1; igf2++) {
734 TObjArray * gamm2 = NULL;
739 gamm2 = static_cast<TObjArray*>(fGammas.At(igf2));
740 igmax = gamm2->GetEntriesFast();
743 for(Int_t i2 = 0; i2 < igmax; i2++) {
744 AliAODConversionPhoton * ph2 = static_cast<AliAODConversionPhoton*>(gamm2->UncheckedAt(i2));
746 if( ph2->GetTrackLabelPositive()==ph1->GetTrackLabelPositive()
747 || ph2->GetTrackLabelNegative()==ph1->GetTrackLabelNegative()
748 || ph2->GetTrackLabelNegative()==ph1->GetTrackLabelPositive()
749 || ph2->GetTrackLabelPositive()==ph1->GetTrackLabelNegative()) {
753 AliAODConversionMother * pion = new AliAODConversionMother(ph1, ph2);
754 if(pion->Eta() < etalim[0] || pion->Eta() > etalim[1]) continue;
755 pion->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
758 tIDs[2] = ph2->GetLabel(0);
759 tIDs[3] = ph2->GetLabel(1);
761 massval[0] = pion->M();
762 massval[1] = pion->Pt();
766 dphivalues[2] = pion->Pt();
767 dphivalues[6] = pion->M();
768 dphivalues[8] = igf1;
771 trigValues[0] = pion->Eta();
772 trigValues[1] = pion->Pt();
773 trigValues[4] = pion->M();
774 trigValues[5] = igf1;
778 ///Check that particle is in histo phase space
779 if(pion->Pt() > fAxistPt.GetBinLowEdge(1) && pion->Pt() < fAxistPt.GetXmax() &&
780 pion->Eta() > etalim[0] && pion->Eta() < etalim[1] &&
781 pion->M() > 0.04 && pion->M() < 0.4
785 if(fMesonFilter->MesonIsSelected(pion, kTRUE)) {
787 fMassSparse->Fill(massval);
789 Bool_t lpimap[fMesonFilters[0].GetEntriesFast()];
790 ///See if it passes lowside cuts
791 if(igf1 == 0 && igf2 == 0) {
795 for(Int_t ilpf = 0; ilpf < fMesonFilters[0].GetEntriesFast(); ilpf++) {
796 if(!(static_cast<AliConversionMesonCuts*>(fMesonFilters[0].At(ilpf))->MesonIsSelected(pion, kTRUE))) {
797 lpimap[ilpf] = kTRUE;
798 massval[4] = -(ilpf + 1);
799 fMassSparse->Fill(massval);
801 lpimap[ilpf] = kFALSE;
806 for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
807 if(lowgmap[iglf][i1] || lowgmap[iglf][i2]){
808 massval[3] = -(iglf+1);
809 fMassSparse->Fill(massval);
813 } /// End lowside mass histo fillers
816 ///Check that particle is in histo phase space
817 if(pion->Pt() > fAxistPt.GetBinLowEdge(1) && pion->Pt() < fAxistPt.GetXmax() &&
818 pion->M() > fAxisPiM.GetXmin() && pion->M() < fAxisPiM.GetXmax() &&
819 pion->Eta() > etalim[0] && pion->Eta() < etalim[1]) {
822 const Int_t tbin = fAxistPt.FindFixBin(pion->Pt());
824 ///Fill standard triggers including upside v0 filters
825 fTrigSparse->Fill(trigValues);
828 if(igf1 == 0 && igf2 == 0) {
829 // piCorr->FillTriggerCounters(pion);
830 // piCorr->CorrelateWithTracks(pion, &tracks[0], tIDs, centrality, vertexz);
833 ////Only mix events with pion in signal region
834 if(pion->M() > 0.1 && pion->M() < 0.15) {
838 if (tbin > 0 && tbin < (nbins + 1)) {
839 tmap[tbin-1] = kTRUE;
842 ///Check if trigger also in low side (both gamma present in low side!)
843 for(Int_t ilgf = 0; ilgf < fV0Filters[0].GetEntriesFast(); ilgf++) {
844 if(!lowgmap[ilgf][i1] && !lowgmap[ilgf][i2]) {
845 lv0tmap[tbin-1][ilgf] = kTRUE;
849 ///See if the lowside pion filter also passes this, if not
850 for(Int_t ilpf = 0; ilpf < fMesonFilters[0].GetEntriesFast(); ilpf++) {
852 lpitmap[tbin-1][ilpf] = kTRUE;
858 if(pion->M() > 0.1 && pion->M() < 0.15) {
859 uv0tmap[tbin-1][igf1 - 1] = kTRUE;
864 ///Fill the triggers not selected in lowside filters only if passsing standard v0 filter
865 if(igf1 == 0 && igf2 == 0) {
867 ///Lowside v0 filters
868 for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
869 if(lowgmap[iglf][i1] || lowgmap[iglf][i2]){
870 trigValues[5] = -(iglf+1);
871 fTrigSparse->Fill(trigValues);
875 ////Low side pion filters
877 for(Int_t iplf = 0; iplf < fMesonFilters[0].GetEntriesFast(); iplf ++) {
879 trigValues[6] = -(iplf + 1);
880 fTrigSparse->Fill(trigValues);
885 trigValues[5] = igf1;
888 ///////////////////////////////////////////////
889 /// Correlate with tracks
890 ///////////////////////////////////////////////
893 if(igf1 == 0 && igf2 == 0) {
894 ntf = fTrackFilters[1].GetEntriesFast() + 1;
897 for(Int_t itf = 0; itf < ntf; itf++) {
898 TObjArray * tracks = static_cast<TObjArray*>(fTracks.At(itf));
899 for(int ij = 0; ij < tracks->GetEntriesFast(); ij++) {
900 AliVTrack * track = static_cast<AliVTrack*>(tracks->At(ij));
901 Int_t tid = track->GetID();
903 if((tid > 0) && (tid == tIDs[0] || tid == tIDs[1] || tid == tIDs[2] || tid == tIDs[3]) ) {
907 dphivalues[0] = pion->Eta() - track->Eta();
908 dphivalues[1] = GetDPhi(pion->Phi() - track->Phi());
909 dphivalues[3] = track->Pt();
911 dphivalues[8] = igf1;
913 fCorrSparse->Fill(dphivalues);
916 if(itf == 0 && igf1 == 0 && igf2 == 0) {
917 ///Fill the low side track filters
918 for(Int_t itlf = 0; itlf < fTrackFilters[0].GetEntriesFast(); itlf++) {
919 if(lowtrackmap[itlf][ij]){
920 dphivalues[7] = -(itlf+1);
921 fCorrSparse->Fill(dphivalues);
924 ///Fill the low side v0 filters
926 for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
927 if(lowgmap[iglf][i1] || lowgmap[iglf][i2]){
928 dphivalues[8] = -(iglf+1);
929 fCorrSparse->Fill(dphivalues);
933 ///Fill the low side pi filter
936 for(Int_t iplf = 0; iplf < fMesonFilters[0].GetEntriesFast(); iplf ++) {
938 dphivalues[9] = -(iplf + 1);
939 fCorrSparse->Fill(dphivalues);
942 } /// end non standard filters track corr
945 } // end trackfilters loop
946 } //end check pion in histogram range to prevent overflow
948 /////////////////////////////
949 //// Not passing standard meson cuts, check upside filters
950 ////////////////////////////
952 ///Only check the pions from standard v0 filter
953 if(igf1 == 0 && igf2 == 0) {
954 for(Int_t ipuf = 0; ipuf < fMesonFilters[1].GetEntriesFast(); ipuf++) {
955 if(static_cast<AliConversionMesonCuts*>(fMesonFilters[1].At(ipuf))->MesonIsSelected(pion, kTRUE)) {
956 ///Fill invariant mass hist
957 massval[4] = (ipuf + 1);
958 fMassSparse->Fill(massval);
959 ///Check that particle is in histo phase space --- redundant!
963 if(pion->Pt() > fAxistPt.GetBinLowEdge(1) && pion->Pt() < fAxistPt.GetXmax() &&
964 pion->M() > fAxisPiM.GetXmin() && pion->M() < fAxisPiM.GetXmax() &&
965 pion->Eta() > etalim[0] && pion->Eta() < etalim[1]) {
968 ////Only mix events with pion in signal region
969 if(pion->M() > 0.1 && pion->M() < 0.15) {
970 const Int_t tbin = fAxistPt.FindFixBin(pion->Pt());
971 upitmap[tbin-1][ipuf] = kTRUE;
974 ///Fill trigger counters
975 trigValues[6] = (ipuf + 1);
976 fTrigSparse->Fill(trigValues);
978 ///Correlate with standard tracks
979 for(int ij = 0; ij < ttracks->GetEntriesFast(); ij++) {
980 AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(ij));
981 Int_t tid = track->GetID();
983 if((tid > 0) && (tid == tIDs[0] || tid == tIDs[1] || tid == tIDs[2] || tid == tIDs[3]) ) {
987 dphivalues[0] = pion->Eta() - track->Eta();
988 dphivalues[1] = GetDPhi(pion->Phi() - track->Phi());
989 dphivalues[3] = track->Pt();
990 dphivalues[7] = 0; // track filter
991 dphivalues[8] = 0; // v0 filter
992 dphivalues[9] = ipuf + 1; // pion filter
993 fCorrSparse->Fill(dphivalues);
999 } /// end else .. end upside meson filters
1000 /////////////////////////////////////////////
1001 } ///Etalim && pt limits
1004 } // i2 second gamma
1008 //FillCounters(&pions, tracks, ntrackfilters, centrality, vertexz);
1010 ///////Fill track counters after entire event has been passed through
1013 Double_t trackValues[fTrackAxesList.GetSize()];
1014 trackValues[3] = centrality;
1015 trackValues[4] = vertexz;
1016 //trackValues[5] = particle->M(); remove !!!
1018 for(Int_t tbin = 0; tbin < fAxistPt.GetNbins(); tbin++) {
1019 trackValues[1] = fAxistPt.GetBinCenter(tbin+1);
1023 for(Int_t itf = 0; itf < fTrackFilters[1].GetEntriesFast() + 1; itf++) {
1024 TObjArray * tracks = static_cast<TObjArray*>(fTracks.At(itf));
1025 for(Int_t iTrack = 0; iTrack < tracks->GetEntriesFast(); iTrack++) {
1026 AliVTrack * track = static_cast<AliVTrack*>(tracks->At(iTrack));
1027 trackValues[0] = track->Eta();
1028 trackValues[2] = track->Pt();
1029 trackValues[5] = itf;
1030 trackValues[6] = 0; ///v0 filter
1031 trackValues[7] = 0; ////Pi filter
1032 fTrackSparse->Fill(trackValues);
1035 for(Int_t itlf = 0; itlf < fTrackFilters[0].GetEntriesFast(); itlf++) {
1036 if(lowtrackmap[itlf][iTrack]) {
1037 trackValues[5] = -(itlf + 1);
1038 fTrackSparse->Fill(trackValues);
1043 ///Check lowside gamma
1044 for(Int_t iglf = 0; iglf < fV0Filters[0].GetEntriesFast(); iglf++) {
1045 if(!lv0tmap[tbin][iglf]) {
1046 trackValues[6] = -(iglf + 1);
1047 fTrackSparse->Fill(trackValues);
1052 ////Check lowside pion
1053 for(Int_t iplf = 0; iplf < fMesonFilters[0].GetEntriesFast(); iplf++) {
1054 if(!lpitmap[tbin][iplf]) {
1055 trackValues[7] = -(iplf + 1);
1056 fTrackSparse->Fill(trackValues);
1066 ///If not in main, see if in upside filters
1067 ///Do upside v0 filters
1068 for(Int_t iguf = 0; iguf < fV0Filters[1].GetEntriesFast(); iguf++) {
1069 if (uv0tmap[tbin][iguf] ) {
1071 cout << "c vtx " << centrality << vertexz << endl;
1073 for(Int_t iTrack = 0; iTrack < ttracks->GetEntriesFast(); iTrack++) {
1074 AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(iTrack));
1075 trackValues[0] = track->Eta();
1076 trackValues[2] = track->Pt();
1078 trackValues[6] = iguf+1; ///v0 filter
1079 trackValues[7] = 0; ////Pi filter
1080 fTrackSparse->Fill(trackValues);
1084 ///Do upside pi filter
1085 for(Int_t ipuf = 0; ipuf < fMesonFilters[1].GetEntriesFast(); ipuf++) {
1086 if (upitmap[tbin][ipuf] ) {
1088 cout << "c2 vtx2 " << centrality << vertexz << endl;
1091 for(Int_t iTrack = 0; iTrack < ttracks->GetEntriesFast(); iTrack++) {
1092 AliVTrack * track = static_cast<AliVTrack*>(ttracks->At(iTrack));
1093 trackValues[0] = track->Eta();
1094 trackValues[2] = track->Pt();
1096 trackValues[6] = 0; ///v0 filter
1097 trackValues[7] = ipuf+1; ////Pi filter
1098 fTrackSparse->Fill(trackValues);
1113 //piCorr->FillCounters(&pions, tracks, centrality, vertexz);
1115 PostData(1, fHistograms);
1119 //_______________________________________________________________________________
1120 // void AliAnalysisTaskdPhi::FillCounters(TObjArray * particles, TObjArray tracks[], Int_t ntrackfilters, Float_t cent, Float_t vtxz) {
1125 // //Fill ME Counters
1126 // const Int_t nbins = fAxistPt.GetNbins();
1127 // Bool_t tmap[nbins];
1128 // for(Int_t ptbin = 0; ptbin < nbins; ptbin++){
1129 // tmap[ptbin] = kFALSE;
1134 // Double_t trackValues[fTrackAxesList.GetSize()];
1135 // trackValues[3] = cent;
1136 // trackValues[4] = vtxz;
1138 // for(Int_t ip = 0; ip < particles->GetEntriesFast(); ip++){
1139 // AliAODConversionParticle * particle = static_cast<AliAODConversionParticle*>(particles->At(ip));
1141 // Int_t tbin = fAxistPt.FindFixBin(particle->Pt());
1142 // if (tbin > 0 && tbin < nbins + 1) {
1143 // if(tmap[tbin - 1] == kTRUE) {
1146 // tmap[tbin -1 ] = kTRUE;
1148 // trackValues[5] = particle->M();
1149 // trackValues[1] = particle->Pt();
1151 // for(Int_t itf = 0; itf < ntrackfilters; itf++) {
1152 // if(fkTrackAxis) trackValues[6] = itf;
1153 // for(int ij = 0; ij < tracks[itf].GetEntriesFast(); ij++) {
1154 // AliVTrack * track = static_cast<AliVTrack*>(tracks->UncheckedAt(ij));
1155 // trackValues[0] = track->Eta();
1156 // trackValues[2] = track->Pt();
1157 // fTrackSparse->Fill(trackValues);
1168 // //________________________________________________________________
1169 // void AliAnalysisTaskdPhi::CorrelateWithTracks(AliAODConversionParticle * particle, TObjArray tracks[], Int_t ntrackfilters, Bool_t ** lowtrackmap, Int_t nltf, Int_t const tIDs[4], Double_t dphivalues[]) {
1170 // //Correlate particle with tracks
1171 // ///Correlate with tracks
1175 //________________________________________________________________________
1176 void AliAnalysisTaskdPhi::Terminate(Option_t *) {
1178 // Draw result to the screen
1179 // Called once at the end of the query
1182 //________________________________________________________________________
1183 TClonesArray * AliAnalysisTaskdPhi::GetConversionGammas(Bool_t isAOD) {
1186 TClonesArray * gammas = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fDeltaAODBranchName.Data()));
1191 FindDeltaAODBranchName(fInputEvent);
1192 gammas = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fDeltaAODBranchName.Data()));
1196 TClonesArray * gammas = dynamic_cast<TClonesArray*>(GetInputData(1));
1202 //________________________________________________________________________
1203 void AliAnalysisTaskdPhi::FindDeltaAODBranchName(AliVEvent * event){
1205 TList *list=event->GetList();
1206 for(Int_t ii=0;ii<list->GetEntries();ii++){
1207 TString name((list->At(ii))->GetName());
1208 if(name.BeginsWith("GammaConv")&&name.EndsWith("gamma")){
1209 fDeltaAODBranchName=name;
1210 AliDebug(AliLog::kDebug + 5, Form("Set DeltaAOD BranchName to: %s",fDeltaAODBranchName.Data()));