]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/Correlations/DPhi/AliAnalysisTaskPhiCorrelations.cxx
update from Markus
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / AliAnalysisTaskPhiCorrelations.cxx
CommitLineData
e0331fd9 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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/* $Id:$ */
17
18#include <TROOT.h>
1fb21a38 19#include <TInterpreter.h>
e0331fd9 20#include <TChain.h>
21#include <TFile.h>
22#include <TList.h>
23#include <TMath.h>
24#include <TTree.h>
25#include <TH2F.h>
ebda4319 26#include <TH3F.h>
e0331fd9 27#include <TRandom.h>
dc44bc1f 28#include <TParameter.h>
e0331fd9 29
30#include "AliAnalysisTaskPhiCorrelations.h"
31#include "AliAnalyseLeadingTrackUE.h"
32#include "AliUEHistograms.h"
33#include "AliUEHist.h"
34
e0331fd9 35#include "AliAnalysisManager.h"
36#include "AliAODHandler.h"
37#include "AliAODInputHandler.h"
38#include "AliAODMCParticle.h"
e0331fd9 39#include "AliInputEventHandler.h"
40#include "AliLog.h"
41#include "AliMCEventHandler.h"
42#include "AliVParticle.h"
2a910c25 43#include "AliCFContainer.h"
e0331fd9 44
45#include "AliESDEvent.h"
46#include "AliESDInputHandler.h"
47#include "AliMultiplicity.h"
2a910c25 48#include "AliCentrality.h"
49#include "AliStack.h"
7dd4dec4 50#include "AliAODMCHeader.h"
1ccd8a0a 51#include "AliGenCocktailEventHeader.h"
52#include "AliGenEventHeader.h"
3e09dd3e 53#include "AliCollisionGeometry.h"
e0331fd9 54
2a910c25 55#include "AliEventPoolManager.h"
e0331fd9 56
d728f490 57#include "AliESDZDC.h"
3bbad7c1 58#include "AliESDtrackCuts.h"
d728f490 59
32e49607 60#include "AliHelperPID.h"
6e84bd0a 61#include "AliAnalysisUtils.h"
dc44bc1f 62#include "TMap.h"
e0331fd9 63
64////////////////////////////////////////////////////////////////////////
65//
66// Analysis class for azimuthal correlation studies
67// Based on the UE task from Sara Vallero and Jan Fiete
68//
69// This class needs input AODs.
70// The output is a list of analysis-specific containers.
71//
72// The AOD can be either connected to the InputEventHandler
73// for a chain of AOD files
74// or
75// to the OutputEventHandler
76// for a chain of ESD files,
77// in this case the class should be in the train after the jet-finder
78//
79// Authors:
80// Jan Fiete Grosse-Oetringhaus
81//
82////////////////////////////////////////////////////////////////////////
83
84
85ClassImp( AliAnalysisTaskPhiCorrelations )
a1c31636 86ClassImp( AliDPhiBasicParticle )
e0331fd9 87
88//____________________________________________________________________
89AliAnalysisTaskPhiCorrelations:: AliAnalysisTaskPhiCorrelations(const char* name):
90AliAnalysisTask(name,""),
91// general configuration
92fDebug(0),
93fMode(0),
94fReduceMemoryFootprint(kFALSE),
eed401dc 95fFillMixed(kTRUE),
ac647b0f 96fMixingTracks(50000),
1bba939a 97fTwoTrackEfficiencyStudy(kFALSE),
d4b3dbfc 98fTwoTrackEfficiencyCut(0),
7e9608f2 99fTwoTrackCutMinRadius(0.8),
44af28f9 100fUseVtxAxis(kFALSE),
9894bedd 101fCourseCentralityBinning(kFALSE),
04af8d15 102fSkipTrigger(kFALSE),
1ccd8a0a 103fInjectedSignals(kFALSE),
f89c31a7 104fRandomizeReactionPlane(kFALSE),
32e49607 105fHelperPID(0x0),
6e84bd0a 106fAnalysisUtils(0x0),
dc44bc1f 107fMap(0x0),
108// pointers to UE classes
e0331fd9 109fAnalyseUE(0x0),
110fHistos(0x0),
111fHistosMixed(0),
418b56c5 112fEfficiencyCorrectionTriggers(0),
113fEfficiencyCorrectionAssociated(0),
24d8278b 114fCentralityWeights(0),
e0331fd9 115// handlers and events
85bfac17 116fAOD(0x0),
117fESD(0x0),
e0331fd9 118fArrayMC(0x0),
119fInputHandler(0x0),
120fMcEvent(0x0),
121fMcHandler(0x0),
2a910c25 122fPoolMgr(0x0),
e0331fd9 123// histogram settings
124fListOfHistos(0x0),
125// event QA
126fnTracksVertex(1), // QA tracks pointing to principal vertex (= 3 default)
85bfac17 127fZVertex(7.),
576bfde8 128fAcceptOnlyMuEvents(kFALSE),
2a910c25 129fCentralityMethod("V0M"),
e0331fd9 130// track cuts
131fTrackEtaCut(0.8),
97dc2b4e 132fTrackEtaCutMin(-1.),
9da2f080 133fOnlyOneEtaSide(0),
e0331fd9 134fPtMin(0.5),
97270ee2 135fDCAXYCut(0),
7d756130 136fSharedClusterCut(-1),
132b6963 137fCrossedRowsCut(-1),
138fFoundFractionCut(-1),
e0331fd9 139fFilterBit(0xFF),
418b56c5 140fTrackStatus(0),
5cabd2cd 141fSelectBit(AliVEvent::kMB|AliVEvent::kUserDefined),
2a910c25 142fUseChargeHadrons(kFALSE),
5c9b9fa6 143fParticleSpeciesTrigger(-1),
144fParticleSpeciesAssociated(-1),
9e35c487 145fCheckMotherPDG(kTRUE),
46ec3872 146fTrackletDphiCut(9999999.),
c05ff6be 147fSelectCharge(0),
7a77d480 148fTriggerSelectCharge(0),
15b0fdd0 149fAssociatedSelectCharge(0),
d38fa455 150fTriggerRestrictEta(-1),
00b6f3c6 151fEtaOrdering(kFALSE),
b0d56b29 152fCutConversions(kFALSE),
153fCutResonances(kFALSE),
97270ee2 154fRejectResonanceDaughters(-1),
a26093ba 155fFillOnlyStep0(kFALSE),
156fSkipStep6(kFALSE),
5e053cad 157fRejectCentralityOutliers(kFALSE),
ff9f18e4 158fRejectZeroTrackEvents(kFALSE),
f613255f 159fRemoveWeakDecays(kFALSE),
d6a8903f 160fRemoveDuplicates(kFALSE),
51d0a028 161fSkipFastCluster(kFALSE),
8a368fc2 162fWeightPerEvent(kFALSE),
3f3f12d9 163fCustomBinning(),
640b9425 164fPtOrder(kTRUE),
2b27d6f4 165fTriggersFromDetector(0),
677dc831 166fAssociatedFromDetector(0),
2313a5d0 167fMCUseUncheckedCentrality(kFALSE),
c05ff6be 168fFillpT(kFALSE)
e0331fd9 169{
170 // Default constructor
171 // Define input and output slots here
172 // Input slot #0 works with a TChain
173 DefineInput(0, TChain::Class());
174 // Output slot #0 writes into a TList container
175 DefineOutput(0, TList::Class());
e0331fd9 176}
177
178AliAnalysisTaskPhiCorrelations::~AliAnalysisTaskPhiCorrelations()
179{
180 // destructor
181
182 if (fListOfHistos && !AliAnalysisManager::GetAnalysisManager()->IsProofMode())
183 delete fListOfHistos;
184}
185
186//____________________________________________________________________
187void AliAnalysisTaskPhiCorrelations::ConnectInputData(Option_t* /*option*/)
188{
189
190 // Connect the input data
191 if (fDebug > 1) AliInfo("ConnectInputData() ");
192
193 // Since AODs can either be connected to the InputEventHandler
194 // or to the OutputEventHandler ( the AOD is created by a previus task in the train )
195 // we need to get the pointer to the AODEvent correctly.
196
197 // Delta AODs are also accepted.
198
199 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
200
85bfac17 201 if( handler && handler->InheritsFrom("AliAODInputHandler") )
202 { // input AOD
203 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
204 if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODInputHandler");
205 }
206 else
207 { //output AOD
208 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
209 if (handler && handler->InheritsFrom("AliAODHandler") )
210 {
211 fAOD = ((AliAODHandler*)handler)->GetAOD();
212 if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODHandler");
213 }
214 else
215 { // no AOD
216 AliWarning("I can't get any AOD Event Handler");
217 }
218 }
219
220 if (handler && handler->InheritsFrom("AliESDInputHandler") )
221 { // input ESD
222 // pointer received per event in ::Exec
223 if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliESDInputHandler");
224 }
e0331fd9 225
226 // Initialize common pointers
227 Initialize();
e0331fd9 228}
229
230//____________________________________________________________________
231void AliAnalysisTaskPhiCorrelations::CreateOutputObjects()
232{
233 // Create the output container
234
235 if (fDebug > 1) AliInfo("CreateOutputObjects()");
236
237 // Initialize class with main algorithms, event and track selection.
238 fAnalyseUE = new AliAnalyseLeadingTrackUE();
97dc2b4e 239 fAnalyseUE->SetParticleSelectionCriteria(fFilterBit, fUseChargeHadrons, fTrackEtaCut, fTrackEtaCutMin, fPtMin);
97270ee2 240 fAnalyseUE->SetDCAXYCut(fDCAXYCut);
7d756130 241 fAnalyseUE->SetSharedClusterCut(fSharedClusterCut);
132b6963 242 fAnalyseUE->SetCrossedRowsCut(fCrossedRowsCut);
243 fAnalyseUE->SetFoundFractionCut(fFoundFractionCut);
418b56c5 244 fAnalyseUE->SetTrackStatus(fTrackStatus);
9e35c487 245 fAnalyseUE->SetCheckMotherPDG(fCheckMotherPDG);
e0331fd9 246 fAnalyseUE->SetDebug(fDebug);
85bfac17 247 fAnalyseUE->DefineESDCuts(fFilterBit);
5cabd2cd 248 fAnalyseUE->SetEventSelection(fSelectBit);
32e49607 249 fAnalyseUE->SetHelperPID(fHelperPID);
250 if ((fParticleSpeciesTrigger != -1 || fParticleSpeciesAssociated != -1) && !fHelperPID)
251 AliFatal("HelperPID object should be set in the steering macro");
e0331fd9 252
253 // Initialize output list of containers
254 if (fListOfHistos != NULL){
255 delete fListOfHistos;
256 fListOfHistos = NULL;
257 }
258 if (!fListOfHistos){
259 fListOfHistos = new TList();
260 fListOfHistos->SetOwner(kTRUE);
261 }
262
263 // Initialize class to handle histograms
9894bedd 264 TString histType = "4R";
3bbad7c1 265 if (fUseVtxAxis == 1)
b752706a 266 histType = "5R";
3bbad7c1 267 else if (fUseVtxAxis == 2)
268 histType = "6R";
9894bedd 269 if (fCourseCentralityBinning)
270 histType += "C";
3f3f12d9 271 fHistos = new AliUEHistograms("AliUEHistogramsSame", histType, fCustomBinning);
272 fHistosMixed = new AliUEHistograms("AliUEHistogramsMixed", histType, fCustomBinning);
e0331fd9 273
2a910c25 274 fHistos->SetSelectCharge(fSelectCharge);
275 fHistosMixed->SetSelectCharge(fSelectCharge);
276
afe13b94 277 fHistos->SetSelectTriggerCharge(fTriggerSelectCharge);
278 fHistosMixed->SetSelectTriggerCharge(fTriggerSelectCharge);
15b0fdd0 279
280 fHistos->SetSelectAssociatedCharge(fAssociatedSelectCharge);
281 fHistosMixed->SetSelectAssociatedCharge(fAssociatedSelectCharge);
7a77d480 282
d38fa455 283 fHistos->SetTriggerRestrictEta(fTriggerRestrictEta);
284 fHistosMixed->SetTriggerRestrictEta(fTriggerRestrictEta);
285
9da2f080 286 fHistos->SetOnlyOneEtaSide(fOnlyOneEtaSide);
287 fHistosMixed->SetOnlyOneEtaSide(fOnlyOneEtaSide);
288
00b6f3c6 289 fHistos->SetEtaOrdering(fEtaOrdering);
290 fHistosMixed->SetEtaOrdering(fEtaOrdering);
291
b0d56b29 292 fHistos->SetPairCuts(fCutConversions, fCutResonances);
293 fHistosMixed->SetPairCuts(fCutConversions, fCutResonances);
294
97270ee2 295 fHistos->SetRejectResonanceDaughters(fRejectResonanceDaughters);
296 fHistosMixed->SetRejectResonanceDaughters(fRejectResonanceDaughters);
297
defad170 298 fHistos->SetTrackEtaCut(fTrackEtaCut);
299 fHistosMixed->SetTrackEtaCut(fTrackEtaCut);
300
8a368fc2 301 fHistos->SetWeightPerEvent(fWeightPerEvent);
302 fHistosMixed->SetWeightPerEvent(fWeightPerEvent);
418b56c5 303
640b9425 304 fHistos->SetPtOrder(fPtOrder);
305 fHistosMixed->SetPtOrder(fPtOrder);
306
7e9608f2 307 fHistos->SetTwoTrackCutMinRadius(fTwoTrackCutMinRadius);
308 fHistosMixed->SetTwoTrackCutMinRadius(fTwoTrackCutMinRadius);
309
418b56c5 310 if (fEfficiencyCorrectionTriggers)
311 {
312 fHistos->SetEfficiencyCorrectionTriggers(fEfficiencyCorrectionTriggers);
313 fHistosMixed->SetEfficiencyCorrectionTriggers((THnF*) fEfficiencyCorrectionTriggers->Clone());
314 }
315 if (fEfficiencyCorrectionAssociated)
408d1ac9 316 {
418b56c5 317 fHistos->SetEfficiencyCorrectionAssociated(fEfficiencyCorrectionAssociated);
318 fHistosMixed->SetEfficiencyCorrectionAssociated((THnF*) fEfficiencyCorrectionAssociated->Clone());
408d1ac9 319 }
13a404ed 320
e0331fd9 321 // add histograms to list
322 fListOfHistos->Add(fHistos);
323 fListOfHistos->Add(fHistosMixed);
32e49607 324 // add HelperPID to list
325 if (fHelperPID)
326 fListOfHistos->Add(fHelperPID);
dc44bc1f 327 // add TMap to list
328 if (fMap)
329 fListOfHistos->Add(fMap);
e0331fd9 330
2a910c25 331 fListOfHistos->Add(new TH2F("processIDs", ";#Delta#phi;process id", 100, -0.5 * TMath::Pi(), 1.5 * TMath::Pi(), kPNoProcess + 1, -0.5, kPNoProcess + 0.5));
c066889a 332 fListOfHistos->Add(new TH1F("eventStat", ";;events", 4, -0.5, 3.5));
5614e02f 333 fListOfHistos->Add(new TH2F("mixedDist", ";centrality;tracks;events", 101, 0, 101, 200, 0, fMixingTracks * 1.5));
f0a25b1d 334 fListOfHistos->Add(new TH2F("referenceMultiplicity", ";centrality;tracks;events", 101, 0, 101, 200, 0, 200));
2b27d6f4 335 if (fCentralityMethod == "V0A_MANUAL")
336 {
337 fListOfHistos->Add(new TH2F("V0AMult", "V0A multiplicity;V0A multiplicity;V0A multiplicity (scaled)", 1000, -.5, 999.5, 1000, -.5, 999.5));
338 fListOfHistos->Add(new TH2F("V0AMultCorrelation", "V0A multiplicity;V0A multiplicity;SPD tracklets", 1000, -.5, 999.5, 1000, -.5, 999.5));
339 }
8cf60d03 340 if (fTriggersFromDetector == 1 || fTriggersFromDetector == 2 || fAssociatedFromDetector == 1 || fAssociatedFromDetector == 2)
2b27d6f4 341 fListOfHistos->Add(new TH1F("V0SingleCells", "V0 single cell multiplicity;multiplicity;events", 100, -0.5, 99.5));
8cf60d03 342 if (fTriggersFromDetector == 3 || fAssociatedFromDetector == 3)
46ec3872 343 fListOfHistos->Add(new TH1F("DphiTrklets", "tracklets Dphi;#Delta#phi,trklets (mrad);entries", 100, -100, 100));
e0331fd9 344
1a32ecc0 345 PostData(0,fListOfHistos);
346
e0331fd9 347 // Add task configuration to output list
348 AddSettingsTree();
349
2a910c25 350 // event mixing
3f3f12d9 351 Int_t poolsize = 1000; // Maximum number of events, ignored in the present implemention of AliEventPoolManager
c066889a 352
2a910c25 353 Int_t nCentralityBins = fHistos->GetUEHist(2)->GetEventHist()->GetNBins(1);
354 Double_t* centralityBins = (Double_t*) fHistos->GetUEHist(2)->GetEventHist()->GetAxis(1, 0)->GetXbins()->GetArray();
355
5f54dac4 356 const Int_t kNZvtxBins = 10+(1+10)*4;
408d1ac9 357 // bins for further buffers are shifted by 100 cm
5f54dac4 358 Double_t vertexBins[kNZvtxBins+1] = { -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10,
359 90, 92, 94, 96, 98, 100, 102, 104, 106, 108, 110,
360 190, 192, 194, 196, 198, 200, 202, 204, 206, 208, 210,
361 290, 292, 294, 296, 298, 300, 302, 304, 306, 308, 310,
362 390, 392, 394, 396, 398, 400, 402, 404, 406, 408, 410 };
363
364 Int_t nZvtxBins = kNZvtxBins;
eb88bdfe 365 Double_t* zvtxbin = vertexBins;
c066889a 366
408d1ac9 367 if (fMode == 0 && fHistos->GetUEHist(2)->GetEventHist()->GetNVar() > 2)
eb88bdfe 368 {
369 nZvtxBins = fHistos->GetUEHist(2)->GetEventHist()->GetNBins(2);
370 zvtxbin = (Double_t*) fHistos->GetUEHist(2)->GetEventHist()->GetAxis(2, 0)->GetXbins()->GetArray();
371 }
85bfac17 372
3f3f12d9 373 fPoolMgr = new AliEventPoolManager(poolsize, fMixingTracks, nCentralityBins, centralityBins, nZvtxBins, zvtxbin);
f0df9076 374 fPoolMgr->SetTargetValues(fMixingTracks, 0.1, 5);
e0331fd9 375}
376
377//____________________________________________________________________
378void AliAnalysisTaskPhiCorrelations::Exec(Option_t */*option*/)
379{
416af868 380 // exec (per event)
381 fAnalyseUE->NextEvent();
382
85bfac17 383 // receive ESD pointer if we are not running AOD analysis
384 if (!fAOD)
385 {
386 AliVEventHandler* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
387 if (handler && handler->InheritsFrom("AliESDInputHandler"))
388 fESD = ((AliESDInputHandler*)handler)->GetEvent();
389 }
390
7dd4dec4 391 if (fMode)
392 {
393 // correction mode
394
395 if (fMcHandler)
396 fMcEvent = fMcHandler->MCEvent();
397
398 if (fAOD)
399 {
400 // array of MC particles
401 fArrayMC = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
402 if (!fArrayMC)
403 AliFatal("No array of MC particles found !!!");
404 }
405
406 AnalyseCorrectionMode();
407 }
e0331fd9 408 else AnalyseDataMode();
e0331fd9 409}
410
411/******************** ANALYSIS METHODS *****************************/
412
413//____________________________________________________________________
414void AliAnalysisTaskPhiCorrelations::AddSettingsTree()
415{
416 //Write settings to output list
417 TTree *settingsTree = new TTree("UEAnalysisSettings","Analysis Settings in UE estimation");
c05ff6be 418 settingsTree->Branch("fnTracksVertex", &fnTracksVertex,"nTracksVertex/I");
419 settingsTree->Branch("fZVertex", &fZVertex,"ZVertex/D");
576bfde8 420 settingsTree->Branch("fAcceptOnlyMuEvents", &fAcceptOnlyMuEvents,"AcceptOnlyMuEvents/O");
c05ff6be 421 //settingsTree->Branch("fCentralityMethod", fCentralityMethod.Data(),"CentralityMethod/C");
e0331fd9 422 settingsTree->Branch("fTrackEtaCut", &fTrackEtaCut, "TrackEtaCut/D");
97dc2b4e 423 settingsTree->Branch("fTrackEtaCutMin", &fTrackEtaCutMin, "TrackEtaCutMin/D");
9da2f080 424 settingsTree->Branch("fOnlyOneEtaSide", &fOnlyOneEtaSide,"OnlyOneEtaSide/I");
e0331fd9 425 settingsTree->Branch("fPtMin", &fPtMin, "PtMin/D");
c05ff6be 426 settingsTree->Branch("fFilterBit", &fFilterBit,"FilterBit/I");
7d756130 427 settingsTree->Branch("fSharedClusterCut", &fSharedClusterCut,"SharedClusterCut/D");
132b6963 428 settingsTree->Branch("fCrossedRowsCut", &fCrossedRowsCut,"CrossedRowsCut/I");
429 settingsTree->Branch("fFoundFractionCut", &fFoundFractionCut,"FoundFractionCut/D");
418b56c5 430 settingsTree->Branch("fTrackStatus", &fTrackStatus,"TrackStatus/I");
c05ff6be 431 settingsTree->Branch("fSelectBit", &fSelectBit,"EventSelectionBit/I");
e0331fd9 432 settingsTree->Branch("fUseChargeHadrons", &fUseChargeHadrons,"UseChHadrons/O");
5c9b9fa6 433 settingsTree->Branch("fParticleSpeciesTrigger", &fParticleSpeciesTrigger,"ParticleSpeciesTrigger/I");
434 settingsTree->Branch("fParticleSpeciesAssociated", &fParticleSpeciesAssociated,"ParticleSpeciesAssociated/I");
9e35c487 435 settingsTree->Branch("fCheckMotherPDG", &fCheckMotherPDG,"CheckMotherPDG/I");
c05ff6be 436 settingsTree->Branch("fSelectCharge", &fSelectCharge,"SelectCharge/I");
7a77d480 437 settingsTree->Branch("fTriggerSelectCharge", &fTriggerSelectCharge,"TriggerSelectCharge/I");
15b0fdd0 438 settingsTree->Branch("fAssociatedSelectCharge", &fAssociatedSelectCharge,"fAssociatedSelectCharge/I");
b0d56b29 439 settingsTree->Branch("fTriggerRestrictEta", &fTriggerRestrictEta,"TriggerRestrictEta/D");
00b6f3c6 440 settingsTree->Branch("fEtaOrdering", &fEtaOrdering,"EtaOrdering/O");
b0d56b29 441 settingsTree->Branch("fCutConversions", &fCutConversions,"CutConversions/O");
442 settingsTree->Branch("fCutResonances", &fCutResonances,"CutResonances/O");
97270ee2 443 settingsTree->Branch("fRejectResonanceDaughters", &fRejectResonanceDaughters,"RejectResonanceDaughters/I");
c05ff6be 444 settingsTree->Branch("fFillpT", &fFillpT,"FillpT/O");
ac647b0f 445 settingsTree->Branch("fMixingTracks", &fMixingTracks,"MixingTracks/I");
04af8d15 446 settingsTree->Branch("fSkipTrigger", &fSkipTrigger,"SkipTrigger/O");
f89c31a7
JFGO
447 settingsTree->Branch("fInjectedSignals", &fInjectedSignals,"InjectedSignals/O");
448 settingsTree->Branch("fRandomizeReactionPlane", &fRandomizeReactionPlane,"RandomizeReactionPlane/O");
5e053cad 449 settingsTree->Branch("fRejectCentralityOutliers", &fRejectCentralityOutliers,"RejectCentralityOutliers/O");
ff9f18e4 450 settingsTree->Branch("fRejectZeroTrackEvents", &fRejectZeroTrackEvents,"RejectZeroTrackEvents/O");
f613255f 451 settingsTree->Branch("fRemoveWeakDecays", &fRemoveWeakDecays,"RemoveWeakDecays/O");
d6a8903f 452 settingsTree->Branch("fRemoveDuplicates", &fRemoveDuplicates,"RemoveDuplicates/O");
51d0a028 453 settingsTree->Branch("fSkipFastCluster", &fSkipFastCluster,"SkipFastCluster/O");
8a368fc2 454 settingsTree->Branch("fWeightPerEvent", &fWeightPerEvent,"WeightPerEvent/O");
640b9425 455 settingsTree->Branch("fPtOrder", &fPtOrder,"PtOrder/O");
2b27d6f4 456 settingsTree->Branch("fTriggersFromDetector", &fTriggersFromDetector,"TriggersFromDetector/I");
677dc831 457 settingsTree->Branch("fAssociatedFromDetector", &fAssociatedFromDetector,"AssociatedFromDetector/I");
2313a5d0 458 settingsTree->Branch("fMCUseUncheckedCentrality", &fMCUseUncheckedCentrality,"MCUseUncheckedCentrality/O");
7e9608f2 459 settingsTree->Branch("fTwoTrackEfficiencyCut", &fTwoTrackEfficiencyCut,"TwoTrackEfficiencyCut/D");
460 settingsTree->Branch("fTwoTrackCutMinRadius", &fTwoTrackCutMinRadius,"TwoTrackCutMinRadius/D");
2b27d6f4 461
3f3f12d9 462 //fCustomBinning
04af8d15 463
e0331fd9 464 settingsTree->Fill();
465 fListOfHistos->Add(settingsTree);
466}
467
468//____________________________________________________________________
469void AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
470{
471 // Run the analysis on MC to get the correction maps
472 //
473
2a910c25 474 if ( fDebug > 3 ) AliInfo( " Processing event in Corrections mode ..." );
475
476 Double_t centrality = 0;
477
478 if (fCentralityMethod.Length() > 0)
479 {
3e09dd3e 480 if (fCentralityMethod == "MC_b")
2a910c25 481 {
3e09dd3e 482 AliGenEventHeader* eventHeader = GetFirstHeader();
483 if (!eventHeader)
484 {
485 // We avoid AliFatal here, because the AOD productions sometimes have events where the MC header is missing
486 // (due to unreadable Kinematics) and we don't want to loose the whole job because of a few events
487 AliError("Event header not found. Skipping this event.");
488 fHistos->FillEvent(0, AliUEHist::kCFStepAnaTopology);
489 return;
490 }
491
492 AliCollisionGeometry* collGeometry = dynamic_cast<AliCollisionGeometry*> (eventHeader);
493 if (!collGeometry)
494 {
495 eventHeader->Dump();
496 AliFatal("Asking for MC_b centrality, but event header has no collision geometry information");
497 }
2313a5d0 498
3e09dd3e 499 centrality = collGeometry->ImpactParameter();
2a910c25 500 }
018f2edc 501 else if (fCentralityMethod == "nano")
502 {
503 centrality = (Float_t) gROOT->ProcessLine(Form("100.0 + 100.0 * ((AliNanoAODHeader*) %p)->GetCentrality(\"%s\")", fAOD->GetHeader(), fCentralityMethod.Data())) / 100 - 1.0;
504 }
2a910c25 505 else
506 {
3e09dd3e 507 AliCentrality *centralityObj = 0;
508 if (fAOD)
509 centralityObj = fAOD->GetHeader()->GetCentralityP();
510 else if (fESD)
511 centralityObj = fESD->GetCentrality();
512
513 if (centralityObj)
514 {
515 if (fMCUseUncheckedCentrality)
516 centrality = centralityObj->GetCentralityPercentileUnchecked(fCentralityMethod);
517 else
518 centrality = centralityObj->GetCentralityPercentile(fCentralityMethod);
519 }
520 else
521 {
522 Printf("WARNING: Centrality object is 0");
523 centrality = -1;
524 }
525 }
526
527 AliInfo(Form("Centrality is %f", centrality));
2a910c25 528 }
529
85bfac17 530 // Support for ESD and AOD based analysis
531 AliVEvent* inputEvent = fAOD;
532 if (!inputEvent)
533 inputEvent = fESD;
d4b3dbfc 534
932155c7 535 Float_t bSign = 0;
536
d4b3dbfc 537 if (inputEvent)
932155c7 538 {
d4b3dbfc 539 fHistos->SetRunNumber(inputEvent->GetRunNumber());
932155c7 540 bSign = (inputEvent->GetMagneticField() > 0) ? 1 : -1;
541 }
542
85bfac17 543 TObject* mc = fArrayMC;
544 if (!mc)
545 mc = fMcEvent;
546
2a910c25 547 // count all events
548 fHistos->FillEvent(centrality, -1);
549
72b62cb8 550 if (centrality < 0)
551 return;
552
2a910c25 553 // Only consider MC events within the vtx-z region used also as cut on the reconstructed vertex
7dd4dec4 554 TObject* vertexSupplier = fMcEvent;
555 if (fAOD) // AOD
556 vertexSupplier = fAOD->FindListObject(AliAODMCHeader::StdBranchName());
557
558 if (!fAnalyseUE->VertexSelection(vertexSupplier, 0, fZVertex))
2a910c25 559 return;
576bfde8 560
561 Float_t zVtx = 0;
7dd4dec4 562 if (fAOD)
563 zVtx = ((AliAODMCHeader*) vertexSupplier)->GetVtxZ();
564 else
565 zVtx = fMcEvent->GetPrimaryVertex()->GetZ();
c05ff6be 566 Float_t weight = 1;
567 if (fFillpT)
568 weight = -1;
1ccd8a0a 569
570 // For productions with injected signals, figure out above which label to skip particles/tracks
571 Int_t skipParticlesAbove = 0;
572 if (fInjectedSignals)
573 {
574 AliGenEventHeader* eventHeader = 0;
575 Int_t headers = 0;
576
577 if (fMcEvent)
578 {
579 // ESD
580 AliHeader* header = (AliHeader*) fMcEvent->Header();
581 if (!header)
582 AliFatal("fInjectedSignals set but no MC header found");
583
584 AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
585 if (!cocktailHeader)
586 {
587 header->Dump();
588 AliFatal("fInjectedSignals set but no MC cocktail header found");
589 }
590
591 headers = cocktailHeader->GetHeaders()->GetEntries();
592 eventHeader = dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->First());
888a53cf 593
594 if (fDebug > 4)
595 {
596 for (Int_t i=0; i<cocktailHeader->GetHeaders()->GetEntries(); i++)
597 {
24d8278b 598 AliGenEventHeader* headerTmp = dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->At(i));
599 if (headerTmp)
600 Printf("%d particles in header:", headerTmp->NProduced());
888a53cf 601 cocktailHeader->GetHeaders()->At(i)->Dump();
602 }
603 }
1ccd8a0a 604 }
605 else
606 {
607 // AOD
608 AliAODMCHeader* header = (AliAODMCHeader*) fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName());
609 if (!header)
610 AliFatal("fInjectedSignals set but no MC header found");
611
612 headers = header->GetNCocktailHeaders();
613 eventHeader = header->GetCocktailHeader(0);
614 }
2a910c25 615
1ccd8a0a 616 if (!eventHeader)
46848f0b 617 {
618 // We avoid AliFatal here, because the AOD productions sometimes have events where the MC header is missing
619 // (due to unreadable Kinematics) and we don't want to loose the whole job because of a few events
620 AliError("First event header not found. Skipping this event.");
621 fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
622 return;
623 }
1ccd8a0a 624
625 skipParticlesAbove = eventHeader->NProduced();
888a53cf 626 AliInfo(Form("Injected signals in this event (%d headers). Keeping particles/tracks of %s. Will skip particles/tracks above %d.", headers, eventHeader->ClassName(), skipParticlesAbove));
1ccd8a0a 627 }
628
f89c31a7
JFGO
629 if (fCentralityWeights && !AcceptEventCentralityWeight(centrality))
630 {
631 AliInfo(Form("Rejecting event because of centrality weighting: %f", centrality));
632 fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
633 return;
634 }
635
2a910c25 636 // Get MC primaries
3f3f12d9 637 // triggers
5c9b9fa6 638 TObjArray* tmpList = fAnalyseUE->GetAcceptedParticles(mc, 0, kTRUE, fParticleSpeciesTrigger, kTRUE);
3f3f12d9 639 CleanUp(tmpList, mc, skipParticlesAbove);
762206c7 640 TObjArray* tracksMC = CloneAndReduceTrackList(tmpList);
641 delete tmpList;
2a910c25 642
3f3f12d9 643 // associated
644 TObjArray* tracksCorrelateMC = tracksMC;
5c9b9fa6 645 if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
3f3f12d9 646 {
b54ab0b1 647 tmpList = fAnalyseUE->GetAcceptedParticles(mc, 0, kTRUE, fParticleSpeciesAssociated, kTRUE);
648 CleanUp(tmpList, mc, skipParticlesAbove);
649 tracksCorrelateMC = CloneAndReduceTrackList(tmpList);
650 delete tmpList;
3f3f12d9 651 }
652
f89c31a7
JFGO
653 if (fRandomizeReactionPlane)
654 {
655 Double_t centralityDigits = centrality*1000. - (Int_t)(centrality*1000.);
656 Double_t angle = TMath::TwoPi() * centralityDigits;
657 AliInfo(Form("Shifting phi of all tracks by %f (digits %f)", angle, centralityDigits));
658 ShiftTracks(tracksMC, angle);
659 if (tracksCorrelateMC != tracksMC)
660 ShiftTracks(tracksCorrelateMC, angle);
661 }
662
9375b5e2 663 if (fFillOnlyStep0)
664 zVtx = 0;
665
0cff229b 666 // Event selection based on number of number of MC particles
667 if (fRejectZeroTrackEvents && tracksMC->GetEntriesFast() == 0)
668 {
669 AliInfo(Form("Rejecting event due to kinematic selection: %f %d", centrality, tracksMC->GetEntriesFast()));
670 fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
671 if (tracksMC != tracksCorrelateMC)
672 delete tracksCorrelateMC;
673 delete tracksMC;
674 return;
675 }
676
2a910c25 677 // (MC-true all particles)
678 // STEP 0
3f3f12d9 679 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepAll, tracksMC, tracksCorrelateMC, weight);
2a910c25 680
7a028750 681 // mixed event
8781a433 682 if (fFillMixed)
683 {
684 AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx);
888a53cf 685 if (fFillOnlyStep0)
686 ((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());
f0df9076 687 if (pool->IsReady())
8781a433 688 for (Int_t jMix=0; jMix<pool->GetCurrentNEvents(); jMix++)
689 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepAll, tracksMC, pool->GetEvent(jMix), 1.0 / pool->GetCurrentNEvents(), (jMix == 0));
3f3f12d9 690 pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelateMC));
8781a433 691 }
7a028750 692
04af8d15 693// Printf("trigger: %d", ((AliInputEventHandler*)fInputHandler)->IsEventSelected());
694
2a910c25 695 // Trigger selection ************************************************
a26093ba 696 if (!fFillOnlyStep0 && (fSkipTrigger || fAnalyseUE->TriggerSelection(fInputHandler)))
2a910c25 697 {
c32a0ca9 698 // (MC-true all particles)
699 // STEP 1
700 if (!fReduceMemoryFootprint)
3f3f12d9 701 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTriggered, tracksMC, tracksCorrelateMC, weight);
c32a0ca9 702 else
703 fHistos->FillEvent(centrality, AliUEHist::kCFStepTriggered);
704
872dc651 705 if (!inputEvent) {
c32a0ca9 706 AliFatal("UNEXPECTED: inputEvent is 0. Trigger selection should have failed");
872dc651 707 return;
708 }
2a910c25 709
710 // Vertex selection *************************************************
85bfac17 711 if (fAnalyseUE->VertexSelection(inputEvent, fnTracksVertex, fZVertex))
2a910c25 712 {
713 // fill here for tracking efficiency
714 // loop over particle species
715
716 for (Int_t particleSpecies = 0; particleSpecies < 4; particleSpecies++)
717 {
85bfac17 718 TObjArray* primMCParticles = fAnalyseUE->GetAcceptedParticles(mc, 0x0, kTRUE, particleSpecies, kTRUE);
59375a30 719 TObjArray* primRecoTracksMatched = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, particleSpecies, kTRUE, kFALSE);
720 TObjArray* allRecoTracksMatched = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, particleSpecies, kTRUE, kFALSE);
721 TObjArray* primRecoTracksMatchedPID = 0;
722 TObjArray* allRecoTracksMatchedPID = 0;
723
724 if (fHelperPID)
725 {
726 primRecoTracksMatchedPID = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, particleSpecies, kTRUE, kTRUE);
727 allRecoTracksMatchedPID = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, particleSpecies, kTRUE, kTRUE);
728 }
1ccd8a0a 729
3f3f12d9 730 CleanUp(primMCParticles, mc, skipParticlesAbove);
731 CleanUp(primRecoTracksMatched, mc, skipParticlesAbove);
732 CleanUp(allRecoTracksMatched, mc, skipParticlesAbove);
59375a30 733 CleanUp(primRecoTracksMatchedPID, mc, skipParticlesAbove);
734 CleanUp(allRecoTracksMatchedPID, mc, skipParticlesAbove);
29cfdb5f 735
736 // select charges
737 if (fTriggerSelectCharge != 0)
738 {
739 SelectCharge(primMCParticles);
740 SelectCharge(primRecoTracksMatched);
741 SelectCharge(allRecoTracksMatched);
59375a30 742 SelectCharge(primRecoTracksMatchedPID);
743 SelectCharge(allRecoTracksMatchedPID);
29cfdb5f 744 }
2a910c25 745
59375a30 746 fHistos->FillTrackingEfficiency(primMCParticles, primRecoTracksMatched, allRecoTracksMatched, primRecoTracksMatchedPID, allRecoTracksMatchedPID, 0, particleSpecies, centrality, zVtx);
2a910c25 747
eed401dc 748// Printf("%d --> %d %d %d", particleSpecies, primMCParticles->GetEntries(), primRecoTracksMatched->GetEntries(), allRecoTracksMatched->GetEntries());
749
750 delete primMCParticles;
2a910c25 751 delete primRecoTracksMatched;
752 delete allRecoTracksMatched;
753 }
3f3f12d9 754
b591fb9c 755 TObjArray* fakeParticles = fAnalyseUE->GetFakeParticles(inputEvent, mc, kFALSE, -1, kTRUE);
3f3f12d9 756 CleanUp((TObjArray*) fakeParticles->At(0), mc, skipParticlesAbove);
757 CleanUp((TObjArray*) fakeParticles->At(1), mc, skipParticlesAbove);
758
59375a30 759 fHistos->FillTrackingEfficiency(0, 0, 0, 0, 0, (TObjArray*) fakeParticles->At(2), 0, centrality, zVtx);
b591fb9c 760 fHistos->FillFakePt(fakeParticles, centrality);
d728f490 761// Printf(">>>>> %d %d %d fakes", ((TObjArray*) fakeParticles->At(0))->GetEntriesFast(), ((TObjArray*) fakeParticles->At(1))->GetEntriesFast(), ((TObjArray*) fakeParticles->At(2))->GetEntriesFast());
b591fb9c 762 delete fakeParticles;
2a910c25 763
764 // (MC-true all particles)
765 // STEP 2
85bfac17 766 if (!fReduceMemoryFootprint)
3f3f12d9 767 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepVertex, tracksMC, tracksCorrelateMC, weight);
c32a0ca9 768 else
769 fHistos->FillEvent(centrality, AliUEHist::kCFStepVertex);
2a910c25 770
771 // Get MC primaries that match reconstructed track
3f3f12d9 772 // triggers
b54ab0b1 773 tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, fParticleSpeciesTrigger, kTRUE);
774 CleanUp(tmpList, mc, skipParticlesAbove);
775 TObjArray* tracksRecoMatchedPrim = CloneAndReduceTrackList(tmpList);
776 delete tmpList;
777
3f3f12d9 778 // associated
779 TObjArray* tracksCorrelateRecoMatchedPrim = tracksRecoMatchedPrim;
5c9b9fa6 780 if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
3f3f12d9 781 {
b54ab0b1 782 tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, fParticleSpeciesAssociated, kTRUE);
783 CleanUp(tmpList, mc, skipParticlesAbove);
784 tracksCorrelateRecoMatchedPrim = CloneAndReduceTrackList(tmpList);
785 delete tmpList;
3f3f12d9 786 }
787
2a910c25 788 // (RECO-matched (quantities from MC particle) primary particles)
789 // STEP 4
3f3f12d9 790 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTrackedOnlyPrim, tracksRecoMatchedPrim, tracksCorrelateRecoMatchedPrim, weight);
408d1ac9 791
792 // mixed event
793 if (fFillMixed)
794 {
795 AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx + 200);
f0df9076 796 if (pool->IsReady())
408d1ac9 797 for (Int_t jMix=0; jMix<pool->GetCurrentNEvents(); jMix++)
798 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTrackedOnlyPrim, tracksRecoMatchedPrim, pool->GetEvent(jMix), 1.0 / pool->GetCurrentNEvents(), (jMix == 0));
3f3f12d9 799 pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelateRecoMatchedPrim));
408d1ac9 800 }
2a910c25 801
802 // Get MC primaries + secondaries that match reconstructed track
3f3f12d9 803 // triggers
b54ab0b1 804 tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, fParticleSpeciesTrigger, kTRUE);
805 CleanUp(tmpList, mc, skipParticlesAbove);
806 TObjArray* tracksRecoMatchedAll = CloneAndReduceTrackList(tmpList);
807 delete tmpList;
808
3f3f12d9 809 // associated
810 TObjArray* tracksCorrelateRecoMatchedAll = tracksRecoMatchedAll;
5c9b9fa6 811 if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
3f3f12d9 812 {
b54ab0b1 813 tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, fParticleSpeciesAssociated, kTRUE);
814 CleanUp(tmpList, mc, skipParticlesAbove);
815 tracksCorrelateRecoMatchedAll = CloneAndReduceTrackList(tmpList);
816 delete tmpList;
3f3f12d9 817 }
818
2a910c25 819 // (RECO-matched (quantities from MC particle) all particles)
820 // STEP 5
3f3f12d9 821 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTracked, tracksRecoMatchedAll, tracksCorrelateRecoMatchedAll, weight);
2a910c25 822
408d1ac9 823 // mixed event
824 if (fFillMixed)
825 {
826 AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx + 300);
f0df9076 827 if (pool->IsReady())
408d1ac9 828 for (Int_t jMix=0; jMix<pool->GetCurrentNEvents(); jMix++)
829 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTracked, tracksRecoMatchedAll, pool->GetEvent(jMix), 1.0 / pool->GetCurrentNEvents(), (jMix == 0));
3f3f12d9 830 pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelateRecoMatchedAll));
408d1ac9 831 }
832
2a910c25 833 // Get RECO tracks
3f3f12d9 834 // triggers
b54ab0b1 835 tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesTrigger, kTRUE);
836 CleanUp(tmpList, mc, skipParticlesAbove);
837 TObjArray* tracks = CloneAndReduceTrackList(tmpList);
838 delete tmpList;
839
3f3f12d9 840 // associated
841 TObjArray* tracksCorrelate = tracks;
5c9b9fa6 842 if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
3f3f12d9 843 {
b54ab0b1 844 tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesAssociated, kTRUE);
845 CleanUp(tmpList, mc, skipParticlesAbove);
846 tracksCorrelate = CloneAndReduceTrackList(tmpList);
847 delete tmpList;
3f3f12d9 848 }
1ccd8a0a 849
2a910c25 850 // (RECO all tracks)
851 // STEP 6
a26093ba 852 if (!fSkipStep6)
3f3f12d9 853 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracks, tracksCorrelate, weight);
2a910c25 854
932155c7 855 // two track cut, STEP 8
856 if (fTwoTrackEfficiencyCut > 0)
3f3f12d9 857 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracks, tracksCorrelate, weight, kTRUE, kTRUE, bSign, fTwoTrackEfficiencyCut);
932155c7 858
408d1ac9 859 // apply correction efficiency, STEP 10
418b56c5 860 if (fEfficiencyCorrectionTriggers || fEfficiencyCorrectionAssociated)
408d1ac9 861 {
418b56c5 862 // with or without two track efficiency depending on if fTwoTrackEfficiencyCut is set
408d1ac9 863 Bool_t twoTrackCut = (fTwoTrackEfficiencyCut > 0);
864
3f3f12d9 865 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepCorrected, tracks, tracksCorrelate, weight, kTRUE, twoTrackCut, bSign, fTwoTrackEfficiencyCut, kTRUE);
408d1ac9 866 }
418b56c5 867
7a028750 868 // mixed event
408d1ac9 869 if (fFillMixed)
8781a433 870 {
871 AliEventPool* pool2 = fPoolMgr->GetEventPool(centrality, zVtx + 100);
888a53cf 872 ((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool2->NTracksInPool());
f0df9076 873 if (pool2->IsReady())
408d1ac9 874 {
932155c7 875 for (Int_t jMix=0; jMix<pool2->GetCurrentNEvents(); jMix++)
876 {
877 // STEP 6
408d1ac9 878 if (!fSkipStep6)
879 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracks, pool2->GetEvent(jMix), 1.0 / pool2->GetCurrentNEvents(), (jMix == 0));
932155c7 880
881 // two track cut, STEP 8
882 if (fTwoTrackEfficiencyCut > 0)
883 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracks, pool2->GetEvent(jMix), 1.0 / pool2->GetCurrentNEvents(), (jMix == 0), kTRUE, bSign, fTwoTrackEfficiencyCut);
408d1ac9 884
885 // apply correction efficiency, STEP 10
418b56c5 886 if (fEfficiencyCorrectionTriggers || fEfficiencyCorrectionAssociated)
408d1ac9 887 {
888 // with or without two track efficiency depending on if fTwoTrackEfficiencyCut is set
889 Bool_t twoTrackCut = (fTwoTrackEfficiencyCut > 0);
890
891 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepCorrected, tracks, pool2->GetEvent(jMix), 1.0 / pool2->GetCurrentNEvents(), (jMix == 0), twoTrackCut, bSign, fTwoTrackEfficiencyCut, kTRUE);
892 }
932155c7 893 }
408d1ac9 894 }
3f3f12d9 895 pool2->UpdatePool(CloneAndReduceTrackList(tracksCorrelate));
8781a433 896 }
408d1ac9 897
85bfac17 898 if (0 && !fReduceMemoryFootprint)
2a910c25 899 {
900 // make list of secondaries (matched with MC)
901 TObjArray* tracksRecoMatchedSecondaries = new TObjArray;
b0d56b29 902 for (Int_t i=0; i<tracksRecoMatchedAll->GetEntriesFast(); i++)
2a910c25 903 if (((AliAODMCParticle*)tracksRecoMatchedAll->At(i))->IsPhysicalPrimary() == kFALSE)
904 tracksRecoMatchedSecondaries->Add(tracksRecoMatchedAll->At(i));
905
906 // Study: Use only secondaries as trigger particles and plot the correlation vs. all particles; store in step 9
c05ff6be 907 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy2, tracksRecoMatchedSecondaries, tracksRecoMatchedAll, weight);
2a910c25 908
909 // Study: Use only primaries as trigger particles and plot the correlation vs. secondaries; store in step 8
c05ff6be 910 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracksRecoMatchedPrim, tracksRecoMatchedSecondaries, weight);
2a910c25 911
912 // plot delta phi vs process id of secondaries
913 // trigger particles: primaries in 4 < pT < 10
914 // associated particles: secondaries in 1 < pT < 10
915
b0d56b29 916 for (Int_t i=0; i<tracksRecoMatchedPrim->GetEntriesFast(); i++)
2a910c25 917 {
918 AliVParticle* triggerParticle = (AliVParticle*) tracksRecoMatchedPrim->At(i);
919
920 if (triggerParticle->Pt() < 4 || triggerParticle->Pt() > 10)
921 continue;
922
b0d56b29 923 for (Int_t j=0; j<tracksRecoMatchedSecondaries->GetEntriesFast(); j++)
2a910c25 924 {
925 AliAODMCParticle* particle = (AliAODMCParticle*) tracksRecoMatchedSecondaries->At(j);
926
927 if (particle->Pt() < 1 || particle->Pt() > 10)
928 continue;
929
930 if (particle->Pt() > triggerParticle->Pt())
931 continue;
932
933 Double_t deltaPhi = triggerParticle->Phi() - particle->Phi();
934 if (deltaPhi > 1.5 * TMath::Pi())
935 deltaPhi -= TMath::TwoPi();
936 if (deltaPhi < -0.5 * TMath::Pi())
937 deltaPhi += TMath::TwoPi();
938
939 Int_t processID = fMcEvent->Stack()->Particle(particle->GetLabel())->GetUniqueID();
940
941 ((TH2F*) fListOfHistos->FindObject("processIDs"))->Fill(deltaPhi, processID);
942 }
943 }
944
945 delete tracksRecoMatchedSecondaries;
946 }
947
3f3f12d9 948 if (tracksCorrelateRecoMatchedPrim != tracksRecoMatchedPrim)
949 delete tracksCorrelateRecoMatchedPrim;
2a910c25 950 delete tracksRecoMatchedPrim;
3f3f12d9 951
952 if (tracksCorrelateRecoMatchedAll != tracksRecoMatchedAll)
953 delete tracksCorrelateRecoMatchedAll;
2a910c25 954 delete tracksRecoMatchedAll;
3f3f12d9 955
956 if (tracksCorrelate != tracks)
957 delete tracksCorrelate;
2a910c25 958 delete tracks;
959 }
960 }
961
3f3f12d9 962 if (tracksMC != tracksCorrelateMC)
963 delete tracksCorrelateMC;
2a910c25 964 delete tracksMC;
e0331fd9 965}
966
3e09dd3e 967//____________________________________________________________________
968AliGenEventHeader* AliAnalysisTaskPhiCorrelations::GetFirstHeader()
969{
970 // get first MC header from either ESD/AOD (including cocktail header if available)
971
972 if (fMcEvent)
973 {
974 // ESD
975 AliHeader* header = (AliHeader*) fMcEvent->Header();
976 if (!header)
977 return 0;
978
979 AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
980 if (cocktailHeader)
981 return dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->First());
982
16e7d490 983 return dynamic_cast<AliGenEventHeader*> (header->GenEventHeader());
3e09dd3e 984 }
985 else
986 {
987 // AOD
988 AliAODMCHeader* header = (AliAODMCHeader*) fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName());
989 if (!header)
990 return 0;
991
992 return header->GetCocktailHeader(0);
993 }
994}
995
e0331fd9 996//____________________________________________________________________
997void AliAnalysisTaskPhiCorrelations::AnalyseDataMode()
998{
999
1000 // Run the analysis on DATA or MC to get raw distributions
1001
e0331fd9 1002 if ( fDebug > 3 ) AliInfo( " Processing event in Data mode ..." );
e0331fd9 1003
c066889a 1004 ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(0);
1005
2a910c25 1006 if (!fInputHandler)
1007 return;
1008
5e053cad 1009 // skip not selected events here (the AOD is not updated for those)
d728f490 1010 if (!fSkipTrigger && !(fInputHandler->IsEventSelected() & fSelectBit))
2a910c25 1011 return;
51d0a028 1012
1013 // skip fast cluster events here if requested
8a368fc2 1014 if (fSkipFastCluster && (fInputHandler->IsEventSelected() & AliVEvent::kFastOnly))
51d0a028 1015 return;
6e84bd0a 1016
72b62cb8 1017 // Support for ESD and AOD based analysis
1018 AliVEvent* inputEvent = fAOD;
1019 if (!inputEvent)
1020 inputEvent = fESD;
1021
2a910c25 1022 Double_t centrality = 0;
1023
c3294f09 1024 AliCentrality *centralityObj = 0;
2a910c25 1025 if (fCentralityMethod.Length() > 0)
1026 {
d728f490 1027 if (fCentralityMethod == "ZNA_MANUAL")
1028 {
9894bedd 1029 Bool_t zna = kFALSE;
1030 for(Int_t j = 0; j < 4; ++j) {
1031 if (fESD->GetZDCData()->GetZDCTDCData(12,j) != 0) {
1032 zna = kTRUE;
1033 }
1034 }
1035
1036// Printf("%d %f", zna, fZNAtower[0]);
1037 if (zna)
1038 {
1039 // code from Chiara O (23.10.12)
1040 const Double_t *fZNAtower = fESD->GetZDCData()->GetZN2TowerEnergy();
72b62cb8 1041 Float_t znacut[4] = {681., 563., 413., 191.};
1042
9894bedd 1043 if(fZNAtower[0]>znacut[0]) centrality = 1;
72b62cb8 1044 else if(fZNAtower[0]>znacut[1]) centrality = 21;
1045 else if(fZNAtower[0]>znacut[2]) centrality = 41;
1046 else if(fZNAtower[0]>znacut[3]) centrality = 61;
1047 else centrality = 81;
9894bedd 1048 }
1049 else
1050 centrality = -1;
d728f490 1051 }
72b62cb8 1052 else if (fCentralityMethod == "TRACKS_MANUAL")
1053 {
1054 // for pp
1055 TObjArray* tracks = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, -1, kTRUE);
1056 centrality = tracks->GetEntriesFast();
1057 if (centrality > 40)
1058 centrality = 41;
1059// Printf("%d %f", tracks->GetEntriesFast(), centrality);
1060 delete tracks;
1061 }
17f8b959 1062 else if (fCentralityMethod == "V0A_MANUAL")
1063 {
ba874719 1064 // for pp
1065
dc44bc1f 1066 //Total multiplicity in the VZERO A detector
1067 Float_t MV0A=inputEvent->GetVZEROData()->GetMTotV0A();
1068 Float_t MV0AScaled=0.;
1069 if (fMap){
1070 TParameter<float>* sf=(TParameter<float>*)fMap->GetValue(Form("%d",inputEvent->GetRunNumber()));
1071 if(sf)MV0AScaled=MV0A*sf->GetVal();
1072 }
dc44bc1f 1073
ba874719 1074 if (MV0AScaled > 0)
1075 centrality = MV0AScaled;
1076 else
1077 centrality = -1;
17f8b959 1078 }
018f2edc 1079 else if (fCentralityMethod == "nano")
1080 {
1081// fAOD->GetHeader()->Dump();
1082// Printf("%p %p %d", dynamic_cast<AliNanoAODHeader*> (fAOD->GetHeader()), dynamic_cast<AliNanoAODHeader*> ((TObject*) (fAOD->GetHeader())), fAOD->GetHeader()->InheritsFrom("AliNanoAODHeader"));
1083
1fb21a38 1084 Int_t error = 0;
1085 centrality = (Float_t) gROOT->ProcessLine(Form("100.0 + 100.0 * ((AliNanoAODHeader*) %p)->GetCentrality(\"%s\")", fAOD->GetHeader(), fCentralityMethod.Data()), &error) / 100 - 1.0;
1086 if (error != TInterpreter::kNoError)
1087 centrality = -1;
018f2edc 1088 }
2a910c25 1089 else
a6d82f4e 1090 {
d728f490 1091 if (fAOD)
1092 centralityObj = fAOD->GetHeader()->GetCentralityP();
1093 else if (fESD)
1094 centralityObj = fESD->GetCentrality();
1095
1096 if (centralityObj)
1097 centrality = centralityObj->GetCentralityPercentile(fCentralityMethod);
1098 //centrality = centralityObj->GetCentralityPercentileUnchecked(fCentralityMethod);
1099 else
1100 centrality = -1;
0df1f11b 1101
d728f490 1102 if (fAOD)
a6d82f4e 1103 {
d728f490 1104 // remove outliers
1105 if (centrality == 0)
a6d82f4e 1106 {
d728f490 1107 if (fAOD->GetVZEROData())
a6d82f4e 1108 {
d728f490 1109 Float_t multV0 = 0;
1110 for (Int_t i=0; i<64; i++)
1111 multV0 += fAOD->GetVZEROData()->GetMultiplicity(i);
1112 if (multV0 < 19500)
1113 {
1114 centrality = -1;
1115 AliInfo("Rejecting event due to too small V0 multiplicity");
1116 }
a6d82f4e 1117 }
1118 }
1119 }
1120 }
1121
a8583211 1122 AliInfo(Form("Centrality is %f", centrality));
2a910c25 1123 }
1124
04af8d15 1125 Float_t bSign = (inputEvent->GetMagneticField() > 0) ? 1 : -1;
1126
85bfac17 1127 fHistos->SetRunNumber(inputEvent->GetRunNumber());
1128
2a910c25 1129 // Fill the "event-counting-container", it is needed to get the number of events remaining after each event-selection cut
1130 fHistos->FillEvent(centrality, AliUEHist::kCFStepAll);
1131
e0331fd9 1132 // Trigger selection ************************************************
d728f490 1133 if (!fSkipTrigger && !fAnalyseUE->TriggerSelection(fInputHandler)) return;
6e84bd0a 1134
e0331fd9 1135 // Fill the "event-counting-container", it is needed to get the number of events remaining after each event-selection cut
2a910c25 1136 fHistos->FillEvent(centrality, AliUEHist::kCFStepTriggered);
e0331fd9 1137
6e84bd0a 1138 // Pileup selection ************************************************
1139 if (fAnalysisUtils && fAnalysisUtils->IsPileUpEvent(inputEvent))
1140 {
1141 // count the removed events
1142 fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
1143
1144 return;
1145 }
1146
576bfde8 1147 // Reject events without a muon in the muon arm ************************************************
1148 if(fAcceptOnlyMuEvents && !IsMuEvent())return;
1149
e0331fd9 1150 // Vertex selection *************************************************
85bfac17 1151 if(!fAnalyseUE->VertexSelection(inputEvent, fnTracksVertex, fZVertex)) return;
e0331fd9 1152
1153 // Fill the "event-counting-container", it is needed to get the number of events remaining after each event-selection cut
2a910c25 1154 fHistos->FillEvent(centrality, AliUEHist::kCFStepVertex);
e0331fd9 1155
b101208f 1156 // fill V0 control histograms
1157 if (fCentralityMethod == "V0A_MANUAL")
1158 {
1159 ((TH2F*) fListOfHistos->FindObject("V0AMult"))->Fill(inputEvent->GetVZEROData()->GetMTotV0A(), centrality);
1160 if (fAOD)
1161 ((TH2F*) fListOfHistos->FindObject("V0AMultCorrelation"))->Fill(inputEvent->GetVZEROData()->GetMTotV0A(), fAOD->GetTracklets()->GetNumberOfTracklets());
1162 }
1163
c3294f09 1164 // optimization
1d53fb34 1165 if (centrality < 0)
2a910c25 1166 return;
e0331fd9 1167
2b27d6f4 1168 TObjArray* tracks = 0;
1169
1170 if (fTriggersFromDetector == 0)
1171 tracks = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesTrigger, kTRUE);
1db82805 1172 else if (fTriggersFromDetector <= 4)
677dc831 1173 tracks=GetParticlesFromDetector(inputEvent,fTriggersFromDetector);
2b27d6f4 1174 else
1175 AliFatal(Form("Invalid setting for fTriggersFromDetector: %d", fTriggersFromDetector));
677dc831 1176
636a6c1d 1177 //Printf("Accepted %d tracks", tracks->GetEntries());
5e053cad 1178
1179 // check for outlier in centrality vs number of tracks (rough constants extracted from correlation histgram)
1180 Bool_t reject = kFALSE;
1181 if (fRejectCentralityOutliers)
1182 {
1183 if (centrality > 40 && centrality <= 50 && tracks->GetEntriesFast() > 1160)
1184 reject = kTRUE;
1185 if (centrality > 50 && centrality <= 60 && tracks->GetEntriesFast() > 650)
1186 reject = kTRUE;
1187 if (centrality > 60 && centrality <= 70 && tracks->GetEntriesFast() > 370)
1188 reject = kTRUE;
1189 if (centrality > 70 && centrality <= 80 && tracks->GetEntriesFast() > 220)
1190 reject = kTRUE;
1191 if (centrality > 80 && centrality <= 90 && tracks->GetEntriesFast() > 130)
1192 reject = kTRUE;
1193 if (centrality > 90 && tracks->GetEntriesFast() > 75)
1194 reject = kTRUE;
1195 }
1196
1197 if (reject)
1198 {
1199 AliInfo(Form("Rejecting event due to centrality vs tracks correlation: %f %d", centrality, tracks->GetEntriesFast()));
1200 fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
ff9f18e4 1201 delete tracks;
1202 return;
1203 }
1204
1205 if (fRejectZeroTrackEvents && tracks->GetEntriesFast() == 0)
1206 {
1207 AliInfo(Form("Rejecting event because it has no tracks: %f %d", centrality, tracks->GetEntriesFast()));
1208 fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
5e053cad 1209 delete tracks;
1210 return;
1211 }
24d8278b 1212
1213 if (fCentralityWeights && !AcceptEventCentralityWeight(centrality))
1214 {
1215 AliInfo(Form("Rejecting event because of centrality weighting: %f", centrality));
1216 fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
1217 delete tracks;
1218 return;
1219 }
1220
636a6c1d 1221 // correlate particles with...
1222 TObjArray* tracksCorrelate = 0;
677dc831 1223 if(fAssociatedFromDetector==0){
1224 if (fParticleSpeciesAssociated != fParticleSpeciesTrigger || fTriggersFromDetector > 0 )
1225 tracksCorrelate = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesAssociated, kTRUE);
1226 }
1db82805 1227 else if (fAssociatedFromDetector <= 4){
677dc831 1228 if(fAssociatedFromDetector != fTriggersFromDetector)
1229 tracksCorrelate = GetParticlesFromDetector(inputEvent,fAssociatedFromDetector);
1230 }
1231 else
1232 AliFatal(Form("Invalid setting for fAssociatedFromDetector: %d", fAssociatedFromDetector));
636a6c1d 1233
3bbad7c1 1234 // reference multiplicity
1235 Int_t referenceMultiplicity = -1;
1236 if (fESD)
1237 referenceMultiplicity = AliESDtrackCuts::GetReferenceMultiplicity(fESD);
3f3f12d9 1238 else if (fAOD)
1239 referenceMultiplicity = tracks->GetEntriesFast(); // TODO to be replaced by the estimator once available in the AOD
418b56c5 1240// referenceMultiplicity = fAOD->GetHeader()->GetRefMultiplicityComb05();
5e053cad 1241
3bbad7c1 1242 ((TH2F*) fListOfHistos->FindObject("referenceMultiplicity"))->Fill(centrality, referenceMultiplicity);
5e053cad 1243
85bfac17 1244 const AliVVertex* vertex = inputEvent->GetPrimaryVertex();
2a910c25 1245 Double_t zVtx = vertex->GetZ();
1246
c05ff6be 1247 Float_t weight = 1;
1248 if (fFillpT)
1249 weight = -1;
1250
c3294f09 1251 // Fill containers at STEP 6 (reconstructed)
1252 if (centrality >= 0)
c066889a 1253 {
a26093ba 1254 if (!fSkipStep6)
636a6c1d 1255 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracks, tracksCorrelate, weight, kTRUE, kFALSE, 0, 0.02, kTRUE);
a26093ba 1256
c066889a 1257 ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(1);
04af8d15 1258
d4b3dbfc 1259 if (fTwoTrackEfficiencyCut > 0)
636a6c1d 1260 fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracks, tracksCorrelate, weight, kTRUE, kTRUE, bSign, fTwoTrackEfficiencyCut, kTRUE);
c066889a 1261 }
c3294f09 1262
636a6c1d 1263 // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool
1264 TObjArray* tracksClone = CloneAndReduceTrackList(tracks);
1265 delete tracks;
1266
eed401dc 1267 if (fFillMixed)
e0331fd9 1268 {
eed401dc 1269 // event mixing
2a910c25 1270
eed401dc 1271 // 1. First get an event pool corresponding in mult (cent) and
1272 // zvertex to the current event. Once initialized, the pool
1273 // should contain nMix (reduced) events. This routine does not
1274 // pre-scan the chain. The first several events of every chain
1275 // will be skipped until the needed pools are filled to the
1276 // specified depth. If the pool categories are not too rare, this
1277 // should not be a problem. If they are rare, you could lose
1278 // statistics.
1279
1280 // 2. Collect the whole pool's content of tracks into one TObjArray
1281 // (bgTracks), which is effectively a single background super-event.
1282
1283 // 3. The reduced and bgTracks arrays must both be passed into
1284 // FillCorrelations(). Also nMix should be passed in, so a weight
1285 // of 1./nMix can be applied.
1286
1287 AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx);
1288
1289 if (!pool)
1290 AliFatal(Form("No pool found for centrality = %f, zVtx = %f", centrality, zVtx));
1291
f0df9076 1292// pool->SetDebug(1);
2d7827bf 1293
f0df9076 1294 if (pool->IsReady())
eed401dc 1295 {
1296
1297 Int_t nMix = pool->GetCurrentNEvents();
a1c31636 1298// cout << "nMix = " << nMix << " tracks in pool = " << pool->NTracksInPool() << endl;
c066889a 1299
1300 ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(2);
5614e02f 1301 ((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());
c066889a 1302 if (pool->IsReady())
1303 ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(3);
eed401dc 1304
1305 // Fill mixed-event histos here
c05ff6be 1306 for (Int_t jMix=0; jMix<nMix; jMix++)
1307 {
eed401dc 1308 TObjArray* bgTracks = pool->GetEvent(jMix);
c5c840c5 1309
a26093ba 1310 if (!fSkipStep6)
408d1ac9 1311 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracksClone, bgTracks, 1.0 / nMix, (jMix == 0), kFALSE, 0, 0.02, kTRUE);
c5c840c5 1312
d4b3dbfc 1313 if (fTwoTrackEfficiencyCut > 0)
408d1ac9 1314 fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracksClone, bgTracks, 1.0 / nMix, (jMix == 0), kTRUE, bSign, fTwoTrackEfficiencyCut, kTRUE);
eed401dc 1315 }
e0331fd9 1316 }
eed401dc 1317
04af8d15 1318 // ownership is with the pool now
636a6c1d 1319 if (tracksCorrelate)
1320 {
1321 pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelate));
1322 delete tracksClone;
1323 }
1324 else
1325 pool->UpdatePool(tracksClone);
eed401dc 1326 //pool->PrintInfo();
e0331fd9 1327 }
04af8d15 1328 else
636a6c1d 1329 {
04af8d15 1330 delete tracksClone;
636a6c1d 1331 if (tracksCorrelate)
1332 delete tracksCorrelate;
1333 }
e0331fd9 1334}
1335
7a028750 1336TObjArray* AliAnalysisTaskPhiCorrelations::CloneAndReduceTrackList(TObjArray* tracks)
1337{
1338 // clones a track list by using AliDPhiBasicParticle which uses much less memory (used for event mixing)
1339
1340 TObjArray* tracksClone = new TObjArray;
1341 tracksClone->SetOwner(kTRUE);
1342
b0d56b29 1343 for (Int_t i=0; i<tracks->GetEntriesFast(); i++)
7a028750 1344 {
1d53fb34 1345 AliVParticle* particle = (AliVParticle*) tracks->UncheckedAt(i);
1346 AliDPhiBasicParticle* copy = new AliDPhiBasicParticle(particle->Eta(), particle->Phi(), particle->Pt(), particle->Charge());
1347 copy->SetUniqueID(particle->GetUniqueID());
1348 tracksClone->Add(copy);
7a028750 1349 }
1350
1351 return tracksClone;
1352}
1353
e0331fd9 1354//____________________________________________________________________
1355void AliAnalysisTaskPhiCorrelations::Initialize()
1356{
1357 // input handler
1358 fInputHandler = (AliInputEventHandler*)
1359 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
1360 // MC handler
d8ee8cd0 1361 fMcHandler = dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
e0331fd9 1362}
d6a8903f 1363
1364//____________________________________________________________________
1365void AliAnalysisTaskPhiCorrelations::RemoveDuplicates(TObjArray* tracks)
1366{
1367 // remove particles with the same label
1368
1369 Int_t before = tracks->GetEntriesFast();
1370
1371 for (Int_t i=0; i<before; ++i)
1372 {
1373 AliVParticle* part = (AliVParticle*) tracks->At(i);
1374
1375 for (Int_t j=i+1; j<before; ++j)
1376 {
1377 AliVParticle* part2 = (AliVParticle*) tracks->At(j);
1378
1379 if (part->GetLabel() == part2->GetLabel())
1380 {
1381 Printf("Removing %d with label %d (duplicated in %d)", i, part->GetLabel(), j); part->Dump(); part2->Dump();
1382 TObject* object = tracks->RemoveAt(i);
1383 if (tracks->IsOwner())
1384 delete object;
1385 break;
1386 }
1387 }
1388 }
1389
1390 tracks->Compress();
1391
1392 if (before > tracks->GetEntriesFast())
1393 AliInfo(Form("Reduced from %d to %d", before, tracks->GetEntriesFast()));
1394}
3f3f12d9 1395
1396void AliAnalysisTaskPhiCorrelations::CleanUp(TObjArray* tracks, TObject* mcObj, Int_t maxLabel)
1397{
1398 // calls RemoveInjectedSignals, RemoveWeakDecays and RemoveDuplicates
1399
59375a30 1400 if (!tracks)
1401 return;
1402
3f3f12d9 1403 if (fInjectedSignals)
1404 fAnalyseUE->RemoveInjectedSignals(tracks, mcObj, maxLabel);
1405 if (fRemoveWeakDecays)
1406 fAnalyseUE->RemoveWeakDecays(tracks, mcObj);
1407 if (fRemoveDuplicates)
1408 RemoveDuplicates(tracks);
1409}
29cfdb5f 1410
1411//____________________________________________________________________
1412void AliAnalysisTaskPhiCorrelations::SelectCharge(TObjArray* tracks)
1413{
1414 // remove particles with charge not selected (depending on fTriggerSelectCharge)
1415
59375a30 1416 if (!tracks)
1417 return;
1418
29cfdb5f 1419 Int_t before = tracks->GetEntriesFast();
1420
1421 for (Int_t i=0; i<before; ++i)
1422 {
1423 AliVParticle* part = (AliVParticle*) tracks->At(i);
1424
1425 if (part->Charge() * fTriggerSelectCharge < -1)
1426 {
1427// Printf("Removing %d with charge %d", i, part->Charge());
1428 TObject* object = tracks->RemoveAt(i);
1429 if (tracks->IsOwner())
1430 delete object;
1431 }
1432 }
1433
1434 tracks->Compress();
1435
1436 if (before > tracks->GetEntriesFast())
1437 AliInfo(Form("Reduced from %d to %d", before, tracks->GetEntriesFast()));
1438}
24d8278b 1439
1440//____________________________________________________________________
1441Bool_t AliAnalysisTaskPhiCorrelations::AcceptEventCentralityWeight(Double_t centrality)
1442{
1443 // rejects "randomly" events such that the centrality gets flat
1444 // uses fCentralityWeights histogram
1445
1446 // TODO code taken and adapted from AliRDHFCuts; waiting for general class AliCentralityFlattening
1447
1448 Double_t weight = fCentralityWeights->GetBinContent(fCentralityWeights->FindBin(centrality));
e245fe0a 1449 Double_t centralityDigits = centrality*100. - (Int_t)(centrality*100.);
24d8278b 1450
1451 Bool_t result = kFALSE;
1452 if (centralityDigits < weight)
1453 result = kTRUE;
1454
1455 AliInfo(Form("Centrality: %f; Digits: %f; Weight: %f; Result: %d", centrality, centralityDigits, weight, result));
1456
1457 return result;
1458}
f89c31a7
JFGO
1459
1460//____________________________________________________________________
1461void AliAnalysisTaskPhiCorrelations::ShiftTracks(TObjArray* tracks, Double_t angle)
1462{
1463 // shifts the phi angle of all tracks by angle
1464 // 0 <= angle <= 2pi
1465
1466 for (Int_t i=0; i<tracks->GetEntriesFast(); ++i)
1467 {
1468 AliDPhiBasicParticle* part = (AliDPhiBasicParticle*) tracks->At(i);
1469 Double_t newAngle = part->Phi() + angle;
1470 if (newAngle >= TMath::TwoPi())
1471 newAngle -= TMath::TwoPi();
1472
1473 part->SetPhi(newAngle);
1474 }
1475}
677dc831 1476
1477//____________________________________________________________________
1478TObjArray* AliAnalysisTaskPhiCorrelations::GetParticlesFromDetector(AliVEvent* inputEvent, Int_t idet)
1479{
1480 //1 = VZERO_A; 2 = VZERO_C; 3 = SPD tracklets
1481 TObjArray* obj = new TObjArray;
1482 obj->SetOwner(kTRUE);
1483
1484 if (idet == 1 || idet == 2)
1485 {
1486 AliVVZERO* vZero = inputEvent->GetVZEROData();
1487
1488 const Int_t vZeroStart = (idet == 1) ? 32 : 0;
1489
1490 TH1F* singleCells = (TH1F*) fListOfHistos->FindObject("V0SingleCells");
1491 for (Int_t i=vZeroStart; i<vZeroStart+32; i++)
1492 {
1493 Float_t weight = vZero->GetMultiplicity(i);
1494 singleCells->Fill(weight);
1495
1496 // rough estimate of multiplicity
1497 for (Int_t j=0; j<TMath::Nint(weight); j++)
1498 {
1499 AliDPhiBasicParticle* particle = new AliDPhiBasicParticle((AliVVZERO::GetVZEROEtaMax(i) + AliVVZERO::GetVZEROEtaMin(i)) / 2, AliVVZERO::GetVZEROAvgPhi(i), 1.1, 0); // fit pT = 1.1 and charge = 0
1500 particle->SetUniqueID(fAnalyseUE->GetEventCounter()* 100000 + j + i * 1000);
1501
1502 obj->Add(particle);
1503 }
1504 }
1505 }
1506 else if (idet == 3)
1507 {
1508 if(!fAOD)
1509 AliFatal("Tracklets only available on AOD");
1510 AliAODTracklets* trklets=(AliAODTracklets*)fAOD->GetTracklets();
1511 if(!trklets)
1512 AliFatal("AliAODTracklets not found");
1513 for(Int_t itrklets=0;itrklets<trklets->GetNumberOfTracklets();itrklets++)
1514 {
b2cf76f9 1515 Float_t eta=-TMath::Log(TMath::Tan(trklets->GetTheta(itrklets)/2));
677dc831 1516 if(TMath::Abs(eta)>fTrackEtaCut)continue;
b2cf76f9 1517 Float_t pT=1000*TMath::Abs(trklets->GetDeltaPhi(itrklets));//in mrad
46ec3872 1518 if(pT>fTrackletDphiCut)continue;
f334f670 1519 TH1F* DphiTrklets = (TH1F*)fListOfHistos->FindObject("DphiTrklets");
46ec3872 1520 DphiTrklets->Fill(1000*trklets->GetDeltaPhi(itrklets)); //in mrad
b2cf76f9 1521 Float_t phi=trklets->GetPhi(itrklets);
1522 phi+=trklets->GetDeltaPhi(itrklets)*39./34.; //correction dphi*39./34. (Dphi in rad)
8bdd8d01 1523 if (phi<0) phi+=TMath::TwoPi();
1524 if (phi>TMath::TwoPi()) phi-=TMath::TwoPi();
b2cf76f9 1525
1526 AliDPhiBasicParticle* particle = new AliDPhiBasicParticle(eta,phi, pT, 0); // pT = TMath::Abs(trklets->GetDeltaPhi(itrklets)) in mrad and charge = 0
677dc831 1527 particle->SetUniqueID(fAnalyseUE->GetEventCounter()* 100000 + itrklets);
1528
1529 obj->Add(particle);
1530 }
1531 }
1db82805 1532 else if (idet == 4)
1533 {
1534 if(!fAOD)
1535 AliFatal("Muon selection only implemented on AOD");//FIXME to be implemented also for ESDs as in AliAnalyseLeadingTrackUE::GetAcceptedPArticles
1536 for (Int_t iTrack = 0; iTrack < fAOD->GetNTracks(); iTrack++) {
1537 AliAODTrack* track = fAOD->GetTrack(iTrack);
1538 if (!track->IsMuonTrack()) continue;
1539 //Float_t dca = track->DCA();
1540 //Float_t chi2 = track->Chi2perNDF();
1db82805 1541 Float_t rabs = track->GetRAtAbsorberEnd();
46ec3872 1542 Float_t eta = track->Eta();
1543 Int_t matching = track->GetMatchTrigger();
1db82805 1544 if (rabs < 17.6 || rabs > 89.5) continue;
1545 if (eta < -4 || eta > -2.5) continue;
46ec3872 1546 if (matching < 2) continue;
1db82805 1547
1548 AliDPhiBasicParticle* particle = new AliDPhiBasicParticle(eta,track->Phi(), track->Pt(), track->Charge());
1549 particle->SetUniqueID(fAnalyseUE->GetEventCounter() * 100000 + iTrack);
1550
1551 obj->Add(particle);
1552 }
1553 }
677dc831 1554 else
1555 AliFatal(Form("GetParticlesFromDetector: Invalid idet value: %d", idet));
1556
1557 return obj;
1558}
1559
576bfde8 1560//____________________________________________________________________
1561Bool_t AliAnalysisTaskPhiCorrelations::IsMuEvent(){
1562
1563 if(!fAOD)
1564 AliFatal("Muon selection only implemented on AOD");//FIXME to be implemented also for ESDs as in AliAnalyseLeadingTrackUE::GetAcceptedPArticles
1565 for (Int_t iTrack = 0; iTrack < fAOD->GetNTracks(); iTrack++) {
1566 AliAODTrack* track = fAOD->GetTrack(iTrack);
1567 if (!track->IsMuonTrack()) continue;
1568 //Float_t dca = track->DCA();
1569 //Float_t chi2 = track->Chi2perNDF();
1570 Float_t rabs = track->GetRAtAbsorberEnd();
1571 Float_t eta = track->Eta();
1572 Int_t matching = track->GetMatchTrigger();
1573 if (rabs < 17.6 || rabs > 89.5) continue;
1574 if (eta < -4 || eta > -2.5) continue;
1575 if (matching < 2) continue;
1576 return kTRUE;
1577 }
1578 return kFALSE;
1579
1580}