]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/hfe/AliAnalysisTaskHFE.cxx
change in handling of V0 tagged tracks
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliAnalysisTaskHFE.cxx
CommitLineData
01745870 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// The analysis task:
17// Filling an AliCFContainer with the quantities pt, eta and phi
18// for tracks which survivied the particle cuts (MC resp. ESD tracks)
19// Track selection is done using the AliHFE package
20//
21// Author:
22// Raphaelle Bailhache <R.Bailhache@gsi.de>
23// Markus Fasel <M.Fasel@gsi.de>
24// Matus Kalisky <matus.kalisky@cern.ch>
25// MinJung Kweon <minjung@physi.uni-heidelberg.de>
26//
27#include <TAxis.h>
28#include <TBits.h>
29#include <TCanvas.h>
30#include <TChain.h>
31#include <TDirectory.h>
32#include <TFile.h>
33#include <TH3D.h>
34#include <TIterator.h>
35#include <TList.h>
36#include <TLegend.h>
37#include <TMath.h>
38#include <TObjArray.h>
39#include <TObjString.h>
40#include <TParticle.h>
41#include <TProfile.h>
42#include <TString.h>
43#include <TF1.h>
44#include <TTree.h>
45
46#include "AliESDtrackCuts.h"
47#include "AliAnalysisManager.h"
48#include "AliAnalysisUtils.h"
49#include "AliAODInputHandler.h"
50#include "AliAODMCParticle.h"
51#include "AliAODTrack.h"
52#include "AliAODVertex.h"
53#include "AliCentrality.h"
54#include "AliCFContainer.h"
55#include "AliCFManager.h"
56#include "AliESDEvent.h"
57#include "AliESDInputHandler.h"
58#include "AliESDtrack.h"
59#include "AliLog.h"
60#include "AliMCEvent.h"
61#include "AliMCEventHandler.h"
62#include "AliMCParticle.h"
63#include "AliMultiplicity.h"
64#include "AliPID.h"
65#include "AliPIDResponse.h"
66#include "AliOADBContainer.h"
67#include "AliStack.h"
68#include "AliTriggerAnalysis.h"
4437a0d2 69#include "AliTRDTriggerAnalysis.h"
01745870 70#include "AliVVertex.h"
71
72#include "AliHFEcollection.h"
73#include "AliHFEcontainer.h"
74#include "AliHFEcuts.h"
75#include "AliHFEelecbackground.h"
76#include "AliHFENonPhotonicElectron.h"
77#include "AliHFEmcQA.h"
78#include "AliHFEpairs.h"
79#include "AliHFEpid.h"
80#include "AliHFEpidQAmanager.h"
81#include "AliHFEsecVtxs.h"
82#include "AliHFEsecVtx.h"
83#include "AliHFEsignalCuts.h"
84#include "AliHFEtaggedTrackAnalysis.h"
85#include "AliHFEtools.h"
7986e54e 86#include "AliHFEV0taginfo.h"
01745870 87#include "AliHFEvarManager.h"
88#include "AliAnalysisTaskHFE.h"
89#include "AliAODMCHeader.h"
90#include "TClonesArray.h"
91
92ClassImp(AliAnalysisTaskHFE)
93
94//____________________________________________________________
95AliAnalysisTaskHFE::AliAnalysisTaskHFE():
96AliAnalysisTaskSE("PID efficiency Analysis")
97 , fAODMCHeader(NULL)
98 , fAODArrayMCInfo(NULL)
99 , fQAlevel(0)
100 , fPlugins(0)
101 , fCollisionSystem(3)
102 , fFillSignalOnly(kTRUE)
103 , fFillNoCuts(kFALSE)
4437a0d2 104 , fApplyCutAOD(kTRUE)
01745870 105 , fBackGroundFactorApply(kFALSE)
106 , fRemovePileUp(kFALSE)
107 , fIdentifiedAsPileUp(kFALSE)
108 , fIdentifiedAsOutInz(kFALSE)
109 , fPassTheEventCut(kFALSE)
110 , fRejectKinkMother(kTRUE)
111 , fisppMultiBin(kFALSE)
112 , fPbPbUserCentralityBinning(kFALSE)
113 , fRemoveFirstEvent(kFALSE)
114 , fisNonHFEsystematics(kFALSE)
115 , fSpecialTrigger(NULL)
116 , fCentralityF(-1)
117 , fCentralityPercent(-1)
118 , fCentralityEstimator("V0M")
119 , fContributors(0.5)
120 , fWeightBackGround(0.)
121 , fVz(0.0)
122 , fContainer(NULL)
123 , fVarManager(NULL)
124 , fSignalCuts(NULL)
125 , fCFM(NULL)
126 , fTriggerAnalysis(NULL)
127 , fPID(NULL)
128 , fPIDqa(NULL)
618038b6 129 , fTRDTriggerAnalysis(NULL)
01745870 130 , fPIDpreselect(NULL)
131 , fCuts(NULL)
132 , fTaggedTrackCuts(NULL)
133 , fCleanTaggedTrack(kFALSE)
134 , fVariablesTRDTaggedTrack(kFALSE)
135 , fAnalysisUtils(NULL)
136 , fCutspreselect(NULL)
137 , fSecVtx(NULL)
138 , fElecBackGround(NULL)
139 , fMCQA(NULL)
140 , fTaggedTrackAnalysis(NULL)
141 , fExtraCuts(NULL)
142 , fBackgroundSubtraction(NULL)
143 , fTRDTrigger(kFALSE)
144 , fWhichTRDTrigger(0)
7986e54e 145 , fV0Tagger(NULL)
01745870 146 , fQA(NULL)
147 , fOutput(NULL)
148 , fHistMCQA(NULL)
149 , fHistSECVTX(NULL)
150 , fHistELECBACKGROUND(NULL)
151 , fQACollection(NULL)
152{
153 //
154 // Dummy constructor
155 //
156 memset(fElecBackgroundFactor, 0, sizeof(Double_t) * kElecBgSpecies * kBgPtBins * kCentBins * kBgLevels);
157 memset(fkBackGroundFactorArray, 0, sizeof(TF1 *) * 12);
158 memset(fBinLimit, 0, sizeof(Double_t) * (kBgPtBins+1));
159 memset(&fisppMultiBin, kFALSE, sizeof(fisppMultiBin));
160 memset(fCentralityLimits, 0, sizeof(Float_t) * 12);
161
162 SetppAnalysis();
163}
164
165//____________________________________________________________
166AliAnalysisTaskHFE::AliAnalysisTaskHFE(const char * name):
167 AliAnalysisTaskSE(name)
168 , fAODMCHeader(NULL)
169 , fAODArrayMCInfo(NULL)
170 , fQAlevel(0)
171 , fPlugins(0)
172 , fCollisionSystem(3)
173 , fFillSignalOnly(kTRUE)
174 , fFillNoCuts(kFALSE)
4437a0d2 175 , fApplyCutAOD(kTRUE)
01745870 176 , fBackGroundFactorApply(kFALSE)
177 , fRemovePileUp(kFALSE)
178 , fIdentifiedAsPileUp(kFALSE)
179 , fIdentifiedAsOutInz(kFALSE)
180 , fPassTheEventCut(kFALSE)
181 , fRejectKinkMother(kTRUE)
182 , fisppMultiBin(kFALSE)
183 , fPbPbUserCentralityBinning(kFALSE)
184 , fRemoveFirstEvent(kFALSE)
185 , fisNonHFEsystematics(kFALSE)
186 , fSpecialTrigger(NULL)
187 , fCentralityF(-1)
188 , fCentralityPercent(-1)
189 , fCentralityEstimator("V0M")
190 , fContributors(0.5)
191 , fWeightBackGround(0.)
192 , fVz(0.0)
193 , fContainer(NULL)
194 , fVarManager(NULL)
195 , fSignalCuts(NULL)
196 , fCFM(NULL)
197 , fTriggerAnalysis(NULL)
198 , fPID(NULL)
199 , fPIDqa(NULL)
618038b6 200 , fTRDTriggerAnalysis(NULL)
01745870 201 , fPIDpreselect(NULL)
202 , fCuts(NULL)
203 , fTaggedTrackCuts(NULL)
204 , fCleanTaggedTrack(kFALSE)
205 , fVariablesTRDTaggedTrack(kFALSE)
206 , fAnalysisUtils(NULL)
207 , fCutspreselect(NULL)
208 , fSecVtx(NULL)
209 , fElecBackGround(NULL)
210 , fMCQA(NULL)
211 , fTaggedTrackAnalysis(NULL)
212 , fExtraCuts(NULL)
213 , fBackgroundSubtraction(NULL)
214 , fTRDTrigger(kFALSE)
215 , fWhichTRDTrigger(0)
7986e54e 216 , fV0Tagger(NULL)
01745870 217 , fQA(NULL)
218 , fOutput(NULL)
219 , fHistMCQA(NULL)
220 , fHistSECVTX(NULL)
221 , fHistELECBACKGROUND(NULL)
222 , fQACollection(0x0)
223{
224 //
225 // Default constructor
226 //
227 DefineOutput(1, TList::Class());
228 DefineOutput(2, TList::Class());
229
7986e54e 230 fV0Tagger = new AliHFEV0taginfo("Tagger");
01745870 231 fPID = new AliHFEpid("hfePid");
232 fPIDqa = new AliHFEpidQAmanager;
233 fVarManager = new AliHFEvarManager("hfeVarManager");
234 fAnalysisUtils = new AliAnalysisUtils;
618038b6 235 fTRDTriggerAnalysis = new AliTRDTriggerAnalysis();
01745870 236
237 memset(fElecBackgroundFactor, 0, sizeof(Double_t) * kElecBgSpecies * kBgPtBins * kCentBins * kBgLevels);
238 memset(fkBackGroundFactorArray, 0, sizeof(TF1 *) * 12);
239 memset(fBinLimit, 0, sizeof(Double_t) * (kBgPtBins+1));
240 memset(&fisppMultiBin, kFALSE, sizeof(fisppMultiBin));
241 memset(fCentralityLimits, 0, sizeof(Float_t) * 12);
242
243 SetppAnalysis();
244}
245
246//____________________________________________________________
247AliAnalysisTaskHFE::AliAnalysisTaskHFE(const AliAnalysisTaskHFE &ref):
248 AliAnalysisTaskSE(ref)
249 , fAODMCHeader(NULL)
250 , fAODArrayMCInfo(NULL)
251 , fQAlevel(0)
252 , fPlugins(0)
253 , fCollisionSystem(ref.fCollisionSystem)
254 , fFillSignalOnly(ref.fFillSignalOnly)
255 , fFillNoCuts(ref.fFillNoCuts)
256 , fApplyCutAOD(ref.fApplyCutAOD)
257 , fBackGroundFactorApply(ref.fBackGroundFactorApply)
258 , fRemovePileUp(ref.fRemovePileUp)
259 , fIdentifiedAsPileUp(ref.fIdentifiedAsPileUp)
260 , fIdentifiedAsOutInz(ref.fIdentifiedAsOutInz)
261 , fPassTheEventCut(ref.fPassTheEventCut)
262 , fRejectKinkMother(ref.fRejectKinkMother)
263 , fisppMultiBin(ref.fisppMultiBin)
264 , fPbPbUserCentralityBinning(ref.fPbPbUserCentralityBinning)
265 , fRemoveFirstEvent(ref.fRemoveFirstEvent)
266 , fisNonHFEsystematics(ref.fisNonHFEsystematics)
267 , fSpecialTrigger(ref.fSpecialTrigger)
268 , fCentralityF(ref.fCentralityF)
269 , fCentralityPercent(ref.fCentralityPercent)
270 , fCentralityEstimator(ref.fCentralityEstimator)
271 , fContributors(ref.fContributors)
272 , fWeightBackGround(ref.fWeightBackGround)
273 , fVz(ref.fVz)
274 , fContainer(NULL)
275 , fVarManager(NULL)
276 , fSignalCuts(NULL)
277 , fCFM(NULL)
278 , fTriggerAnalysis(NULL)
279 , fPID(NULL)
280 , fPIDqa(NULL)
618038b6 281 , fTRDTriggerAnalysis(NULL)
01745870 282 , fPIDpreselect(NULL)
283 , fCuts(NULL)
284 , fTaggedTrackCuts(NULL)
285 , fCleanTaggedTrack(ref.fCleanTaggedTrack)
286 , fVariablesTRDTaggedTrack(ref.fVariablesTRDTaggedTrack)
287 , fAnalysisUtils(NULL)
288 , fCutspreselect(NULL)
289 , fSecVtx(NULL)
290 , fElecBackGround(NULL)
291 , fMCQA(NULL)
292 , fTaggedTrackAnalysis(NULL)
293 , fExtraCuts(NULL)
294 , fBackgroundSubtraction(NULL)
295 , fTRDTrigger(ref.fTRDTrigger)
296 , fWhichTRDTrigger(ref.fWhichTRDTrigger)
7986e54e 297 , fV0Tagger(NULL)
01745870 298 , fQA(NULL)
299 , fOutput(NULL)
300 , fHistMCQA(NULL)
301 , fHistSECVTX(NULL)
302 , fHistELECBACKGROUND(NULL)
303 , fQACollection(NULL)
304{
305 //
306 // Copy Constructor
307 //
308 ref.Copy(*this);
309}
310
311//____________________________________________________________
312AliAnalysisTaskHFE &AliAnalysisTaskHFE::operator=(const AliAnalysisTaskHFE &ref){
313 //
314 // Assignment operator
315 //
316 if(this == &ref)
317 ref.Copy(*this);
318 return *this;
319}
320
321//____________________________________________________________
322void AliAnalysisTaskHFE::Copy(TObject &o) const {
323 //
324 // Copy into object o
325 //
326 AliAnalysisTaskHFE &target = dynamic_cast<AliAnalysisTaskHFE &>(o);
327 target.fAODMCHeader = fAODMCHeader;
328 target.fAODArrayMCInfo = fAODArrayMCInfo;
329 target.fQAlevel = fQAlevel;
330 target.fPlugins = fPlugins;
331 target.fCollisionSystem = fCollisionSystem;
332 target.fFillSignalOnly = fFillSignalOnly;
333 target.fFillNoCuts = fFillNoCuts;
334 target.fApplyCutAOD = fApplyCutAOD;
335 target.fBackGroundFactorApply = fBackGroundFactorApply;
336 target.fRemovePileUp = fRemovePileUp;
337 target.fIdentifiedAsPileUp = fIdentifiedAsPileUp;
338 target.fIdentifiedAsOutInz = fIdentifiedAsOutInz;
339 target.fPassTheEventCut = fPassTheEventCut;
340 target.fRejectKinkMother = fRejectKinkMother;
341 target.fisppMultiBin = fisppMultiBin;
342 target.fPbPbUserCentralityBinning = fPbPbUserCentralityBinning;
343 target.fRemoveFirstEvent = fRemoveFirstEvent;
344 target.fisNonHFEsystematics = fisNonHFEsystematics;
345 target.fSpecialTrigger = fSpecialTrigger;
346 target.fCentralityF = fCentralityF;
347 target.fCentralityPercent = fCentralityPercent;
348 target.fCentralityEstimator = fCentralityEstimator;
349 target.fContributors = fContributors;
350 target.fWeightBackGround = fWeightBackGround;
351 target.fVz = fVz;
352 target.fContainer = fContainer;
353 target.fVarManager = fVarManager;
354 target.fSignalCuts = fSignalCuts;
355 target.fCFM = fCFM;
356 target.fTriggerAnalysis = fTriggerAnalysis;
357 target.fPID = fPID;
358 target.fPIDqa = fPIDqa;
618038b6 359 target.fTRDTriggerAnalysis = fTRDTriggerAnalysis;
01745870 360 target.fPIDpreselect = fPIDpreselect;
361 target.fCuts = fCuts;
362 target.fTaggedTrackCuts = fTaggedTrackCuts;
363 target.fCleanTaggedTrack = fCleanTaggedTrack;
364 target.fVariablesTRDTaggedTrack = fVariablesTRDTaggedTrack;
365 target.fAnalysisUtils = fAnalysisUtils;
366 target.fCutspreselect = fCutspreselect;
367 target.fSecVtx = fSecVtx;
368 target.fElecBackGround = fElecBackGround;
369 target.fMCQA = fMCQA;
370 target.fTaggedTrackAnalysis = fTaggedTrackAnalysis;
371 target.fExtraCuts = fExtraCuts;
372 target.fBackgroundSubtraction = fBackgroundSubtraction;
373 target.fTRDTrigger = fTRDTrigger;
374 target.fWhichTRDTrigger = fWhichTRDTrigger;
7986e54e 375 target.fV0Tagger = fV0Tagger;
01745870 376 target.fQA = fQA;
377 target.fOutput = fOutput;
378 target.fHistMCQA = fHistMCQA;
379 target.fHistSECVTX = fHistSECVTX;
380 target.fHistELECBACKGROUND = fHistELECBACKGROUND;
381 target.fQACollection = fQACollection;
382}
383
384//____________________________________________________________
385AliAnalysisTaskHFE::~AliAnalysisTaskHFE(){
386 //
387 // Destructor
388 //
389 if(fPID) delete fPID;
390 if(fPIDpreselect) delete fPIDpreselect;
391 if(fVarManager) delete fVarManager;
618038b6 392 if(fTRDTriggerAnalysis) delete fTRDTriggerAnalysis;
01745870 393 if(fCFM) delete fCFM;
394 if(fTriggerAnalysis) delete fTriggerAnalysis;
395 if(fSignalCuts) delete fSignalCuts;
396 if(fSecVtx) delete fSecVtx;
397 if(fMCQA) delete fMCQA;
398 if(fElecBackGround) delete fElecBackGround;
399 if(fBackgroundSubtraction) delete fBackgroundSubtraction;
400 if(fSpecialTrigger) delete fSpecialTrigger;
401 if(fAnalysisUtils) delete fAnalysisUtils;
7986e54e 402 if(fV0Tagger) delete fV0Tagger;
01745870 403 // Delete output objects only if we are not running in PROOF mode because otherwise this produces a crash during merging
404 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
405 if(mgr && mgr->GetAnalysisType() != AliAnalysisManager::kProofAnalysis){
406 if(fPIDqa) delete fPIDqa;
407 if(fOutput) delete fOutput;
408 if(fQA) delete fQA;
409 }
410}
411
412//____________________________________________________________
413void AliAnalysisTaskHFE::UserCreateOutputObjects(){
414 //
415 // Creating output container and output objects
416 // Here we also Initialize the correction framework container and
417 // the objects for
418 // - PID
419 // - MC QA
420 // - SecVtx
421 // QA histograms are created if requested
422 // Called once per worker
423 //
424 AliDebug(3, "Creating Output Objects");
425
426 // Make lists for Output
427 if(!fQA) fQA = new TList;
428 fQA->SetOwner();
429 if(!fOutput) fOutput = new TList;
430 fOutput->SetOwner();
431
432 // Automatic determination of the analysis mode
433 AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
434 if(!TString(inputHandler->IsA()->GetName()).CompareTo("AliAODInputHandler")){
435 SetAODAnalysis();
436 } else {
437 SetESDAnalysis();
438 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())
439 SetHasMCData();
440 }
441 printf("Analysis Mode: %s Analysis\n", IsAODanalysis() ? "AOD" : "ESD");
442 printf("MC Data available %s\n", HasMCData() ? "Yes" : "No");
443
444 // Enable Trigger Analysis
445 fTriggerAnalysis = new AliTriggerAnalysis;
446 fTriggerAnalysis->EnableHistograms();
447 fTriggerAnalysis->SetAnalyzeMC(HasMCData());
448
449 // First Part: Make QA histograms
450 fQACollection = new AliHFEcollection("TaskQA", "QA histos from the Electron Task");
451 fQACollection->CreateTH1F("nElectronTracksEvent", "Number of Electron Candidates", 100, 0, 100);
452 fQACollection->CreateTH1F("nElectron", "Number of electrons", 100, 0, 100);
453 fQACollection->CreateTH2F("radius", "Production Vertex", 100, 0.0, 5.0, 100, 0.0, 5.0);
454 fQACollection->CreateTH1F("nTriggerBit", "Histo Trigger Bit", 22, 0, 22);
47faab26 455 fQACollection->CreateTH1F("Filterend", "AOD filter of tracks after all cuts", 21, -1, 20);
01745870 456
457 InitHistoITScluster();
458 InitContaminationQA();
459 fQA->Add(fQACollection);
460
461 // Initialize PID
462 fPID->SetHasMCData(HasMCData());
463 if(!fPID->GetNumberOfPIDdetectors()) fPID->AddDetector("TPC", 0);
464 if(IsQAOn(kPIDqa)){
465 AliInfo("PID QA switched on");
466 fPIDqa->Initialize(fPID);
467 fQA->Add(fPIDqa->MakeList("HFEpidQA"));
468 }
469 fPID->SortDetectors();
470
471 // Background subtraction-------------------------------------------------------------------
472 if (GetPlugin(kNonPhotonicElectron)) {
473 if(!fBackgroundSubtraction) fBackgroundSubtraction = new AliHFENonPhotonicElectron();
474 if(IsAODanalysis()) fBackgroundSubtraction->SetAOD(kTRUE);
475 fBackgroundSubtraction->Init();
476 fOutput->Add(fBackgroundSubtraction->GetListOutput());
477 }
478 //------------------------------------------------------------------------------------------
479
480
481 // Initialize correction Framework and Cuts
482 const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack + AliHFEcuts::kNcutStepsSecvtxTrack;
483 fCFM = new AliCFManager;
484 fCFM->SetNStepParticle(kNcutSteps);
485 MakeParticleContainer();
486 MakeEventContainer();
487 // Temporary fix: Initialize particle cuts with NULL
488 for(Int_t istep = 0; istep < kNcutSteps; istep++)
489 fCFM->SetParticleCutsList(istep, NULL);
490 if(!fCuts){
491 AliWarning("Cuts not available. Default cuts will be used");
492 fCuts = new AliHFEcuts;
493 fCuts->CreateStandardCuts();
494 }
495 if(IsAODanalysis()) fCuts->SetAOD();
496 // Make clone for V0 tagging step
497 fCuts->Initialize(fCFM);
498 if(fCuts->IsQAOn()) fQA->Add(fCuts->GetQAhistograms());
499 fSignalCuts = new AliHFEsignalCuts("HFEsignalCuts", "HFE MC Signal definition");
500 fVarManager->SetSignalCuts(fSignalCuts);
501
502 // add output objects to the List
503 fOutput->AddAt(fContainer, 0);
504 fOutput->AddAt(fCFM->GetEventContainer(), 1);
505
506 // mcQA----------------------------------
507 if (HasMCData() && IsQAOn(kMCqa)) {
508 AliInfo("MC QA on");
509 if(!fMCQA) fMCQA = new AliHFEmcQA;
510 if(!fHistMCQA) fHistMCQA = new TList();
511 fHistMCQA->SetOwner();
512 if(IsPbPb()) fMCQA->SetPbPb();
513 if(fisppMultiBin) fMCQA->SetPPMultiBin();
514 if(TestBit(kTreeStream)){
515 fMCQA->EnableDebugStreamer();
516 }
517 fMCQA->CreatDefaultHistograms(fHistMCQA);
518 fMCQA->SetBackgroundWeightFactor(fElecBackgroundFactor[0][0][0],fBinLimit);
519 fQA->Add(fHistMCQA);
520 }
521
522 // secvtx----------------------------------
523 if (GetPlugin(kSecVtx)) {
524 AliInfo("Secondary Vertex Analysis on");
525 if(!fSecVtx) fSecVtx = new AliHFEsecVtx;
526 fSecVtx->SetHasMCData(HasMCData());
527
528 if(!fHistSECVTX) fHistSECVTX = new TList();
529 fHistSECVTX->SetOwner();
530 fSecVtx->CreateHistograms(fHistSECVTX);
531 fOutput->Add(fHistSECVTX);
532 }
533
534 // background----------------------------------
535 if (GetPlugin(kIsElecBackGround)) {
536 AliInfo("Electron BackGround Analysis on");
537 if(!fElecBackGround){
538 AliWarning("ElecBackGround not available. Default elecbackground will be used");
539 fElecBackGround = new AliHFEelecbackground;
540 }
541 fElecBackGround->SetHasMCData(HasMCData());
542
543 if(!fHistELECBACKGROUND) fHistELECBACKGROUND = new TList();
544 fHistELECBACKGROUND->SetOwner();
545 fElecBackGround->CreateHistograms(fHistELECBACKGROUND);
546 fOutput->Add(fHistELECBACKGROUND);
547 }
548
549 // tagged tracks
550 if(GetPlugin(kTaggedTrackAnalysis)){
551 AliInfo("Analysis on V0-tagged tracks enabled");
552 fTaggedTrackAnalysis = new AliHFEtaggedTrackAnalysis(Form("taggedTrackAnalysis%s", GetName()));
553 fTaggedTrackAnalysis->SetCuts(fTaggedTrackCuts);
554 fTaggedTrackAnalysis->SetClean(fCleanTaggedTrack);
555 AliHFEvarManager *varManager = fTaggedTrackAnalysis->GetVarManager();
556 TObjArray *array = fVarManager->GetVariables();
557 Int_t nvars = array->GetEntriesFast();
558 TString namee;
559 for(Int_t v = 0; v < nvars; v++) {
560 AliHFEvarManager::AliHFEvariable *variable = (AliHFEvarManager::AliHFEvariable *) array->At(v);
561 if(!variable) continue;
562 TString name(((AliHFEvarManager::AliHFEvariable *)variable)->GetName());
563 if(!name.CompareTo("source")) namee = TString("species");
564 else namee = TString(name);
565 Int_t nbins = variable->GetNumberOfBins();
566 if(variable->HasUserDefinedBinning()){
567 varManager->AddVariable(namee, nbins, variable->GetBinning());
568 } else {
569 varManager->AddVariable(namee, nbins, variable->GetMinimum(), variable->GetMaximum());
570 }
571 //printf("For AliTaggedTrackAnalysis, had the variable %s and the one used %s\n",(const char*)variable->GetName(),(const char*) namee);
572 }
573 if(fPIDqa->HasHighResolutionHistos())
574 fTaggedTrackAnalysis->GetPIDqa()->SetHighResolutionHistos();
575 fTaggedTrackAnalysis->SetPID(fPID);
576 fTaggedTrackAnalysis->SetVariablesTRD(fVariablesTRDTaggedTrack);
577 fTaggedTrackAnalysis->InitContainer();
578 fOutput->Add(fTaggedTrackAnalysis->GetContainer());
579 fQA->Add(fTaggedTrackAnalysis->GetPIDQA());
580 fQA->Add(fTaggedTrackAnalysis->GetCutQA());
581 fQA->Add(fTaggedTrackAnalysis->GetQAcollection());
582 }
583
584 //fQA->Print();
585
586 PrintStatus();
587 // Done!!!
588 PostData(1, fOutput);
589 PostData(2, fQA);
590}
591
592//____________________________________________________________
593void AliAnalysisTaskHFE::UserExec(Option_t *){
594 //
595 // Run the analysis
596 //
597
598 //printf("test00\n");
599
600 AliDebug(3, "Starting Single Event Analysis");
601 if(!fInputEvent){
602 AliError("Reconstructed Event not available");
603 //printf("Reconstructed Event not available");
604 return;
605 }
606 if(HasMCData() && IsESDanalysis()){
607 AliDebug(4, Form("MC Event: %p", fMCEvent));
608 if(!fMCEvent){
609 AliError("No MC Event, but MC Data required");
610 //printf("No MC Event, but MC Data required");
611 return;
612 }
613 }
614 if(!fCuts){
615 AliError("HFE cuts not available");
616 //printf("HFE cuts not available");
617 return;
618 }
619 if(!fPID->IsInitialized()){
620 // Initialize PID with the given run number
621 fPID->InitializePID(fInputEvent->GetRunNumber());
622 }
623
624 if(fRemoveFirstEvent){
625 if(fAnalysisUtils->IsFirstEventInChunk(fInputEvent)) return;
626 }
627
628 AliESDEvent *ev = dynamic_cast<AliESDEvent *>(fInputEvent);
4437a0d2 629 if(ev && fTRDTrigger && (fWhichTRDTrigger<6))
01745870 630 {
4437a0d2 631 if(!CheckTRDTriggerESD(ev)) return;
632 }
633 if(fInputEvent && fTRDTrigger && (fWhichTRDTrigger>5))
634 {
635 if(!CheckTRDTrigger(fInputEvent)) return;
01745870 636 }
01745870 637
638 if(IsESDanalysis() && HasMCData()){
639 // Protect against missing MC trees
640 AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
641 if(!mcH){
642 AliError("No MC Event Handler available");
643 return;
644 }
645 if(!mcH->InitOk()) return;
646 if(!mcH->TreeK()) return;
647 if(!mcH->TreeTR()) return;
648
649 // Background subtraction-------------------------------------------------------------------
650 if(GetPlugin(kNonPhotonicElectron)) fBackgroundSubtraction->SetMCEvent(fMCEvent);
651 //------------------------------------------------------------------------------------------
652 }
653
654 if(IsAODanalysis() && HasMCData()){
655 // take MC info
656 AliAODEvent *aodE = dynamic_cast<AliAODEvent *>(fInputEvent);
657 if(!aodE){
658 AliError("No AOD Event");
659 return;
660 }
661 fAODMCHeader = dynamic_cast<AliAODMCHeader *>(fInputEvent->FindListObject(AliAODMCHeader::StdBranchName()));
662 if(!fAODMCHeader){
663 AliError("No AliAODMCHeader");
664 //printf("No AliAODMCHeader");
665 return;
666 }
667 fAODArrayMCInfo = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
668 if(!fAODArrayMCInfo){
669 AliError("No AOD MC particles");
670 //printf("No AOD MC particles");
671 return;
672 }
673 fSignalCuts->SetMCAODInfo(fAODArrayMCInfo);
674 // Background subtraction-------------------------------------------------------------------
675 if (GetPlugin(kNonPhotonicElectron)) fBackgroundSubtraction->SetAODArrayMCInfo(fAODArrayMCInfo);
676 //------------------------------------------------------------------------------------------
677 }
678
679 //printf("test2\n");
680
681 // need the centrality for everything (MC also)
682 fCentralityF = -1;
683 if(!ReadCentrality()) fCentralityF = -1;
684 //printf("pass centrality\n");
685 //printf("Reading fCentralityF %d\n",fCentralityF);
686
687 // See if pile up and z in the range
688 RejectionPileUpVertexRangeEventCut();
689
690 //printf("test3\n");
691
692 // Protect agains missing
693 if(HasMCData()){
694 //printf("Has MC data\n");
695 fSignalCuts->SetMCEvent(fMCEvent);
696 ProcessMC(); // Run the MC loop + MC QA in case MC Data are available
697 }
698
699 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
700 if(!pidResponse){
701 AliDebug(1, "Using default PID Response");
702 pidResponse = AliHFEtools::GetDefaultPID(HasMCData(), fInputEvent->IsA() == AliESDEvent::Class());
703 }
704 fPID->SetPIDResponse(pidResponse);
705 if(fPIDpreselect) fPIDpreselect->SetPIDResponse(pidResponse);
706
707 // Background subtraction-------------------------------------------------------------------
708 if(GetPlugin(kNonPhotonicElectron)) fBackgroundSubtraction->InitRun(fInputEvent,pidResponse);
709 //------------------------------------------------------------------------------------------
710
711 // Event loop
712 if(IsAODanalysis()){
713 //printf("test4\n");
714 ProcessAOD();
715 } else {
716 const char *specialTrigger = GetSpecialTrigger(fInputEvent->GetRunNumber());
717 // Check Trigger selection
718 if(specialTrigger){
719 AliDebug(2, Form("Special Trigger requested: %s", specialTrigger));
720 if(!(ev && ev->IsTriggerClassFired(specialTrigger))){
721 AliDebug(2, "Event not selected");
722 return;
723 } else AliDebug(2, "Event Selected");
724 } else AliDebug(2, "No Special Trigger requested");
725
726 ProcessESD();
727 }
728 // Done!!!
729 PostData(1, fOutput);
730 PostData(2, fQA);
731}
732
733//____________________________________________________________
734void AliAnalysisTaskHFE::Terminate(Option_t *){
735 //
736 // Terminate not implemented at the moment
737 //
738}
739
740//_______________________________________________________________
741Bool_t AliAnalysisTaskHFE::IsEventInBinZero() {
742 //
743 //
744 //
745
746 //printf("test in IsEventInBinZero\n");
747 if(!fInputEvent){
748 AliError("Reconstructed Event not available");
749 return kFALSE;
750 }
751
752 // check vertex
753 const AliVVertex *vertex = fInputEvent->GetPrimaryVertex();
754 if(!vertex) return kTRUE;
755 //if(vertex) return kTRUE;
756
757 // check tracks
758 if(fInputEvent->GetNumberOfTracks()<=0) return kTRUE;
759 //if(fInputEvent->GetNumberOfTracks()>0) return kTRUE;
760
761
762 return kFALSE;
763
764}
765//____________________________________________________________
766void AliAnalysisTaskHFE::ProcessMC(){
767 //
768 // Runs the MC Loop (filling the container for the MC Cut Steps with the observables pt, eta and phi)
769 // In case MC QA is on also MC QA loop is done
770 //
771 AliDebug(3, "Processing MC Information");
772 Double_t eventContainer [4] = {0., 0., 0., 0.};
773 if(IsESDanalysis()) eventContainer[0] = fMCEvent->GetPrimaryVertex()->GetZ();
774 else eventContainer[0] = fAODMCHeader->GetVtxZ();
775 eventContainer[2] = fCentralityF;
776 eventContainer[3] = fContributors;
777 fVz = eventContainer[0];
778 //printf("z position is %f\n",eventContainer[0]);
779 //if(fCFM->CheckEventCuts(AliHFEcuts::kEventStepGenerated, fMCEvent))
780 fCFM->GetEventContainer()->Fill(eventContainer,AliHFEcuts::kEventStepGenerated);
781 Int_t nElectrons = 0;
782 if(IsESDanalysis()){
783 if(!((fIdentifiedAsPileUp) || (TMath::Abs(fVz) > fCuts->GetVertexRange()) || (fCentralityF < 0))){ //kStepMCGeneratedZOutNoPileUpCentralityFine
784 if (HasMCData() && IsQAOn(kMCqa)) {
785 AliDebug(2, "Running MC QA");
786
787 if(fMCEvent->Stack()){
788 fMCQA->SetMCEvent(fMCEvent);
789 fMCQA->SetGenEventHeader(fMCEvent->GenEventHeader());
790 fMCQA->SetCentrality(fCentralityF);
791 fMCQA->SetPercentrality(static_cast<Int_t>(fCentralityPercent));
792 if(IsPbPb()) { fMCQA->SetPbPb();}
793 else
794 {
795 if(fisppMultiBin) fMCQA->SetPPMultiBin();
796 else fMCQA->SetPP();
797 }
798 fMCQA->Init();
799
800 fMCQA->GetMesonKine();
801
802 // loop over all tracks for decayed electrons
803 for (Int_t igen = 0; igen < fMCEvent->GetNumberOfTracks(); igen++){
804 TParticle* mcpart = fMCEvent->Stack()->Particle(igen);
805 if(!mcpart) continue;
806 fMCQA->GetQuarkKine(mcpart, igen, AliHFEmcQA::kCharm);
807 fMCQA->GetQuarkKine(mcpart, igen, AliHFEmcQA::kBeauty);
808 fMCQA->GetHadronKine(mcpart, AliHFEmcQA::kCharm);
809 fMCQA->GetHadronKine(mcpart, AliHFEmcQA::kBeauty);
810 fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kCharm, AliHFEmcQA::kElectronPDG); // no accept cut
811 fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kBeauty, AliHFEmcQA::kElectronPDG); // no accept cut
812 fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kOthers, AliHFEmcQA::kElectronPDG); // no accept cut
813 }
814 //fMCQA->EndOfEventAna(AliHFEmcQA::kCharm);
815 //fMCQA->EndOfEventAna(AliHFEmcQA::kBeauty);
816 }
817
818 } // end of MC QA loop
819 }
820 // -----------------------------------------------------------------
821 fCFM->SetMCEventInfo(fMCEvent);
822 // fCFM->CheckEventCuts(AliCFManager::kEvtRecCuts, fESD);
823 } else {
824 fCFM->SetMCEventInfo(fInputEvent);
825 }
826 // Run MC loop
827 AliVParticle *mctrack = NULL;
828 Int_t numberofmctracks = 0;
829 if(IsESDanalysis()){
830 numberofmctracks = fMCEvent->GetNumberOfTracks();
831 }
832 else {
833 numberofmctracks = fAODArrayMCInfo->GetEntriesFast();
834 }
835 AliDebug(3, Form("Number of Tracks: %d",numberofmctracks));
836 //printf("Number of MC track %d\n",numberofmctracks);
837 for(Int_t imc = 0; imc <numberofmctracks; imc++){
838 if(IsESDanalysis()) {
839 if(!(mctrack = fMCEvent->GetTrack(imc))) continue;
840 }
841 else {
842 if(!(mctrack = (AliVParticle *) fAODArrayMCInfo->At(imc))) continue;
843 }
844 //printf("Test in ProcessMC\n");
845 AliDebug(4, "Next MC Track");
846 if(ProcessMCtrack(mctrack)) nElectrons++;
847 }
848
849 // fCFM->CheckEventCuts(AliCFManager::kEvtRecCuts, fESD);
850 fQACollection->Fill("nElectron", nElectrons);
851}
852
853//____________________________________________________________
854void AliAnalysisTaskHFE::ProcessESD(){
855 //
856 // Run Analysis of reconstructed event in ESD Mode
857 // Loop over Tracks, filter according cut steps defined in AliHFEcuts
858 //
859 AliDebug(1, Form("Task %s", GetName()));
860 AliDebug(3, "Processing ESD Event");
861 AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
862 if(!fESD){
863 AliError("ESD Event required for ESD Analysis");
864 return;
865 }
866
7986e54e 867 // Tag all v0s in current event
868 if(fV0Tagger){
869 fV0Tagger->Reset();
870 fV0Tagger->TagV0Tracks(fESD);
871 }
01745870 872 // Set magnetic field if V0 task on
873 if(fTaggedTrackAnalysis) {
874 fTaggedTrackAnalysis->SetMagneticField(fESD->GetMagneticField());
875 fTaggedTrackAnalysis->SetCentrality(fCentralityF);
4437a0d2 876 if(IsPbPb()) fTaggedTrackAnalysis->SetPbPb();
877 else {
878 if(IspPb()) fTaggedTrackAnalysis->SetpPb();
879 else fTaggedTrackAnalysis->SetPP();
880 }
881
01745870 882 }
883
884 // Do event Normalization
885 Double_t eventContainer[4];
886 eventContainer[0] = 0.;
887 if(HasMCData()) eventContainer[0] = fVz;
888 else {
889 const AliESDVertex* vtxESD = fESD->GetPrimaryVertex();
890 if(vtxESD) eventContainer[0] = vtxESD->GetZ();
891 }
892 eventContainer[1] = 0.;
893 eventContainer[2] = fCentralityF;
894 eventContainer[3] = fContributors;
895 if(fTriggerAnalysis->IsOfflineTriggerFired(fESD, AliTriggerAnalysis::kV0AND))
896 eventContainer[1] = 1.;
897
898 //
899 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoCut);
900
901 //
902 if(fIdentifiedAsPileUp) return;
903 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoPileUp);
904
905 //
906 if(TMath::Abs(fCentralityF) < 0) return;
907 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecCentralityOk);
908 //printf("In ProcessESD %f\n",fCentralityF);
909
910 //
911 if(fIdentifiedAsOutInz) return;
912 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepZRange);
913
914 //
915 if(!fPassTheEventCut) return;
916 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepReconstructed);
917
918
919 fContainer->NewEvent();
920
921 if (GetPlugin(kIsElecBackGround)) {
922 fElecBackGround->SetEvent(fESD);
923 }
924 if (GetPlugin(kSecVtx)) {
925 fSecVtx->SetEvent(fESD);
926 fSecVtx->GetPrimaryCondition();
927 }
928
929 if(HasMCData()){
930 if (GetPlugin(kSecVtx)) {
931 fSecVtx->SetMCEvent(fMCEvent);
932 fSecVtx->SetMCQA(fMCQA);
933 }
934 if (GetPlugin(kIsElecBackGround)) {
935 fElecBackGround->SetMCEvent(fMCEvent);
936 }
937 }
938
939 Double_t container[10];
940 memset(container, 0, sizeof(Double_t) * 10);
941 // container for the output THnSparse
942 Double_t dataDca[6]; // [source, pT, dca, centrality]
943 Int_t nElectronCandidates = 0;
944 AliESDtrack *track = NULL, *htrack = NULL;
945 AliMCParticle *mctrack = NULL;
946 AliMCParticle *mctrackmother = NULL;
947
948 Bool_t signal = kTRUE;
949
950 fCFM->SetRecEventInfo(fESD);
951
952 // Get Number of contributors to the primary vertex for multiplicity-dependent correction
953 Int_t ncontribVtx = 0;
954 const AliESDVertex *priVtx = fESD->GetPrimaryVertexTracks();
955 if(priVtx){
956 ncontribVtx = priVtx->GetNContributors();
957 }
958
959 // minjung for IP QA(temporary ~ 2weeks)
960 if(!fExtraCuts){
961 fExtraCuts = new AliHFEextraCuts("hfetmpCuts","HFE tmp Cuts");
962 }
963 fExtraCuts->SetRecEventInfo(fESD);
964
965 // Electron background analysis
966 if (GetPlugin(kIsElecBackGround)) {
967
968 AliDebug(2, "Running BackGround Analysis");
969
970 fElecBackGround->Reset();
971
972 } // end of electron background analysis
973
974
975 // Background subtraction-------------------------------------------------------------------
976 if (GetPlugin(kNonPhotonicElectron)) fBackgroundSubtraction->FillPoolAssociatedTracks(fInputEvent, fCentralityF);
977 //------------------------------------------------------------------------------------------
978
979 //
980 // Loop ESD
981 //
982 AliDebug(3, Form("Number of Tracks: %d", fESD->GetNumberOfTracks()));
983 for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
984 AliDebug(4, "New ESD track");
985 track = fESD->GetTrack(itrack);
986 track->SetESDEvent(fESD);
4437a0d2 987
01745870 988 // fill counts of v0-identified particles
7986e54e 989 AliPID::EParticleType v0pid = fV0Tagger->GetV0Info(track->GetID());
01745870 990 // here the tagged track analysis will run
7986e54e 991 if(fTaggedTrackAnalysis && v0pid != AliPID::kUnknown){
01745870 992 AliDebug(1, Form("Track identified as %s", AliPID::ParticleName(v0pid)));
993 fTaggedTrackAnalysis->ProcessTrack(track, v0pid);
994 AliDebug(1, "V0 PID done");
995 }
996
997
998 //Fill non-HFE source containers at reconstructed events cut step
999 AliDebug(3, Form("Doing track %d, %p", itrack, track));
1000
1001
1002 //////////////////////////////////////
1003 // preselect
1004 //////////////////////////////////////
1005 if(fPIDpreselect || fCutspreselect) {
1006 if(!PreSelectTrack(track)) continue;
1007 }
1008
1009 signal = kTRUE;
1010
1011 // Fill step without any cut
1012
1013 if(HasMCData()){
1014 // Check if it is electrons near the vertex
1015 if(!(mctrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(track->GetLabel()))))) continue;
1016
1017 if(fFillSignalOnly && !fCFM->CheckParticleCuts(AliHFEcuts::kStepMCGenerated, mctrack)) signal = kFALSE;
1018 else AliDebug(3, "Signal Electron");
1019
1020 // Fill K pt for Ke3 contributions
1021 if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode())==321)) fQACollection->Fill("Kptspectra",mctrack->Pt());
1022 else if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode())==130)) fQACollection->Fill("K0Lptspectra",mctrack->Pt());
1023 }
1024 // Cache new Track information inside the var manager
1025 fVarManager->NewTrack(track, mctrack, fCentralityF, -1, signal);
1026
1027 if(fFillNoCuts) {
1028 if(signal || !fFillSignalOnly){
1029 fVarManager->FillContainer(fContainer, "recTrackContReco", AliHFEcuts::kStepRecNoCut, kFALSE);
1030 fVarManager->FillContainer(fContainer, "recTrackContMC", AliHFEcuts::kStepRecNoCut, kTRUE);
1031 }
1032 }
1033
1034 // RecKine: ITSTPC cuts
1035 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
1036
1037
1038 // RecPrim
1039 if(fRejectKinkMother) {
1040 if(track->GetKinkIndex(0) != 0) continue; } // Quick and dirty fix to reject both kink mothers and daughters
1041 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
1042
1043 // HFEcuts: ITS layers cuts
1044 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
1045
1046 // HFE cuts: TOF PID and mismatch flag
1047 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTOF, track)) continue;
1048
1049 // HFE cuts: TPC PID cleanup
1050 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
1051
1052 // HFEcuts: Nb of tracklets TRD0
1053 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTRD, track)) continue;
1054
1055 // Fill correlation maps before PID
1056 if(signal && fContainer->GetCorrelationMatrix("correlationstepbeforePID")) {
1057 //printf("Fill correlation maps before PID\n");
1058 fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepbeforePID"));
1059 }
1060
1061 if(HasMCData()){
1062 FillProductionVertex(track);
1063
1064 if(fMCQA && signal){
1065 fMCQA->SetCentrality(fCentralityF);
1066 if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) == 11)){
1067 Double_t weightElecBgV0[kBgLevels] = {0.,0.,0.};
1068 Double_t hfeimpactRtmp=0., hfeimpactnsigmaRtmp=0.;
1069 fExtraCuts->GetHFEImpactParameters(track, hfeimpactRtmp, hfeimpactnsigmaRtmp);
1070 UChar_t itsPixel = track->GetITSClusterMap();
1071 Double_t ilyrhit=0, ilyrstat=0;
1072 for(Int_t ilyr=0; ilyr<6; ilyr++){
1073 if(TESTBIT(itsPixel, ilyr)) ilyrhit += TMath::Power(2,ilyr);
1074 if(fExtraCuts->CheckITSstatus(fExtraCuts->GetITSstatus(track,ilyr))) ilyrstat += TMath::Power(2,ilyr);
1075 }
1076 fMCQA->SetITSInfo(ilyrhit,ilyrstat);
1077 fMCQA->SetHFEImpactParameters(hfeimpactRtmp, hfeimpactnsigmaRtmp);
1078 fMCQA->SetTrkKine(track->Pt(),track->Eta(), track->Phi());
1079 fMCQA->SetContainerStep(3);
1080 for(Int_t iLevel = 0; iLevel < kBgLevels; iLevel++){
1081 weightElecBgV0[iLevel] = fMCQA->GetWeightFactor(mctrack, iLevel); // positive:conversion e, negative: nonHFE
1082 if(!fisNonHFEsystematics || IsPbPb())break;
1083 }
1084
1085 if(fisNonHFEsystematics){
1086 //Fill additional containers for electron source distinction
1087 Int_t elecSource = 0;
1088 elecSource = fMCQA->GetElecSource(mctrack->Particle());
1089 const Char_t *sourceName[kElecBgSpecies]={"Pion","Eta","Omega","Phi","EtaPrime","Rho"};
1090 const Char_t *levelName[kBgLevels]={"Best","Lower","Upper"};
1091 Int_t iName = 0;
1092 for(Int_t iSource = AliHFEmcQA::kPi0; iSource <= AliHFEmcQA::kGammaRho0; iSource++){
1093 if((iSource == AliHFEmcQA::kElse)||(iSource == AliHFEmcQA::kMisID)) continue;
1094 if(elecSource == iSource){
1095 for(Int_t iLevel = 0; iLevel < kBgLevels; iLevel++){
1096 if(weightElecBgV0[iLevel]>0){
1097 fVarManager->FillContainer(fContainer, Form("conversionElecs%s%s",sourceName[iName], levelName[iLevel]), 3, kFALSE, weightElecBgV0[iLevel]);
1098 }
1099 else if(weightElecBgV0[iLevel]<0){
1100 fVarManager->FillContainer(fContainer, Form("mesonElecs%s%s",sourceName[iName], levelName[iLevel]), 3, kFALSE, -1*weightElecBgV0[iLevel]);
1101 }
1102 if(IsPbPb())break;
1103 }
1104 break;
1105 }
1106 iName++;
1107 if(iName == kElecBgSpecies)iName = 0;
1108 }
1109 }
1110 //else{
1111 if(weightElecBgV0[0]>0) {
1112 fVarManager->FillContainer(fContainer, "conversionElecs", 3, kFALSE, weightElecBgV0[0]);
1113 fVarManager->FillContainer(fContainer, "conversionElecs", 4, kTRUE, weightElecBgV0[0]);
1114 }
1115 else if(weightElecBgV0[0]<0) {
1116 fVarManager->FillContainer(fContainer, "mesonElecs", 3, kFALSE, -1*weightElecBgV0[0]);
1117 fVarManager->FillContainer(fContainer, "mesonElecs", 4, kTRUE, -1*weightElecBgV0[0]);
1118 }
1119 //}
1120 }
1121 }
1122
1123 Double_t hfeimpactR4all=0., hfeimpactnsigmaR4all=0.;
1124 Int_t sourceDca =-1;
1125 if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) == 211)){
1126 if(track->Pt()>4.){
1127 fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4all, hfeimpactnsigmaR4all);
1128 dataDca[0]=0; //pion
1129 dataDca[1]=track->Pt();
1130 dataDca[2]=hfeimpactR4all;
1131 dataDca[3]=fCentralityF;
1132 dataDca[4] = v0pid;
1133 dataDca[5] = double(track->Charge());
1134 fQACollection->Fill("Dca", dataDca);
1135 }
1136 }
1137 else if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) == 11)){ // to increas statistics for Martin
1138 if(signal){
1139 fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4all, hfeimpactnsigmaR4all);
1140 if(fSignalCuts->IsCharmElectron(track)){
1141 sourceDca=1;
1142 }
1143 else if(fSignalCuts->IsBeautyElectron(track)){
1144 sourceDca=2;
1145 }
1146 else if(fSignalCuts->IsGammaElectron(track)){
1147 sourceDca=3;
1148 }
1149 else if(fSignalCuts->IsNonHFElectron(track)){
1150 sourceDca=4;
1151 }
1152 else if(fSignalCuts->IsJpsiElectron(track)){
1153 sourceDca=5;
1154 }
1155 else {
1156 sourceDca=6;
1157 }
1158 dataDca[0]=sourceDca;
1159 dataDca[1]=track->Pt();
1160 dataDca[2]=hfeimpactR4all;
1161 dataDca[3]=fCentralityF;
1162 dataDca[4] = v0pid;
1163 dataDca[5] = double(track->Charge());
1164 if(signal) fQACollection->Fill("Dca", dataDca);
1165 }
1166 }
1167 }
1168
1169 AliHFEpidObject hfetrack;
1170 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1171 hfetrack.SetRecTrack(track);
1172 if(HasMCData()) hfetrack.SetMCTrack(mctrack);
1173 hfetrack.SetCentrality(fCentralityF);
1174 hfetrack.SetMulitplicity(ncontribVtx);
4437a0d2 1175 if(IsPbPb()) hfetrack.SetPbPb();
1176 else {
1177 if(IspPb()) hfetrack.SetpPb();
1178 else hfetrack.SetPP();
1179 }
01745870 1180 fPID->SetVarManager(fVarManager);
1181 if(!fPID->IsSelected(&hfetrack, fContainer, "recTrackCont", fPIDqa)) continue;
1182 nElectronCandidates++;
1183
1184 // Background subtraction------------------------------------------------------------------------------------------
1185 if (GetPlugin(kNonPhotonicElectron)) {
1186 Int_t indexmother = -1;
1187 Int_t mcsource = -1;
1188 if(HasMCData()){
1189 mcsource = fBackgroundSubtraction->FindMother(mctrack->GetLabel(),indexmother);
1190 }
1191 fBackgroundSubtraction->LookAtNonHFE(itrack, track, fInputEvent, 1, fCentralityF, -1, mcsource, indexmother);
1192 }
1193 //-----------------------------------------------------------------------------------------------------------------
1194
1195 // Temporary histogram for chi2/ITS cluster
1196 if(IsPbPb()) {
1197 TBits shared = track->GetTPCSharedMap();
1198 Int_t sharebit=0;
1199 if(shared.CountBits() >= 2) sharebit=1;
1200
1201 Double_t itschi2percluster = 0.0;
1202 Double_t itsnbcls = static_cast<Double_t>(track->GetNcls(0));
1203 if(itsnbcls > 0) itschi2percluster = track->GetITSchi2()/itsnbcls;
1204
1205 Double_t itsChi2[7] = {track->Pt(),track->Eta(), track->Phi(),
1206 static_cast<Double_t>(fCentralityF),static_cast<Double_t>(track->GetTPCsignalN()), static_cast<Double_t>(sharebit),itschi2percluster};
1207 fQACollection->Fill("fChi2perITScluster", itsChi2);
1208 }
1209 else{
1210
1211 Double_t itschi2percluster = 0.0;
1212 Double_t itsnbcls = static_cast<Double_t>(track->GetNcls(0));
1213 if(itsnbcls > 0) itschi2percluster = track->GetITSchi2()/itsnbcls;
1214
1215 Double_t itsChi2[3] = {track->Pt(), static_cast<Double_t>(fCentralityF), itschi2percluster};
1216 fQACollection->Fill("fChi2perITScluster", itsChi2);
1217 }
1218
1219 // Fill Histogram for Hadronic Background
1220 if(HasMCData()){
1221 if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11))
1222 fVarManager->FillContainer(fContainer, "hadronicBackground", UInt_t(0), kFALSE);
1223 else if(mctrack){
1224 // Fill Ke3 contributions
1225 Int_t glabel=TMath::Abs(mctrack->GetMother());
1226 if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
1227 if(TMath::Abs(mctrackmother->Particle()->GetPdgCode())==321)
1228 fQACollection->Fill("Ke3Kecorr",mctrack->Pt(),mctrackmother->Pt());
1229 else if(TMath::Abs(mctrackmother->Particle()->GetPdgCode())==130)
1230 fQACollection->Fill("Ke3K0Lecorr",mctrack->Pt(),mctrackmother->Pt());
1231 }
1232 }
1233 }
1234
1235 // Fill Containers
1236 if(signal) {
1237 // Apply weight for background contamination
1238 if(fBackGroundFactorApply) {
4437a0d2 1239 if(IsPbPb() && fCentralityF >= 0) fWeightBackGround = fkBackGroundFactorArray[fCentralityF >= 0 ? fCentralityF : 0]->Eval(TMath::Abs(track->P()));
1240 else fWeightBackGround = fkBackGroundFactorArray[0]->Eval(TMath::Abs(track->P())); // for pp and pPb
01745870 1241
1242 if(fWeightBackGround < 0.0) fWeightBackGround = 0.0;
1243 else if(fWeightBackGround > 1.0) fWeightBackGround = 1.0;
1244 // weightBackGround as special weight
1245 fVarManager->FillContainer(fContainer, "hadronicBackground", 1, kFALSE, fWeightBackGround);
1246 }
1247 fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepafterPID"));
1248 }
1249
1250 Bool_t bTagged=kFALSE;
1251 if(GetPlugin(kSecVtx)) {
1252 AliDebug(2, "Running Secondary Vertex Analysis");
1253 if(fSecVtx->Process(track) && signal) {
1254 fVarManager->FillContainer(fContainer, "recTrackContSecvtxReco", AliHFEcuts::kStepHFEcutsSecvtx, kFALSE);
1255 fVarManager->FillContainer(fContainer, "recTrackContSecvtxMC", AliHFEcuts::kStepHFEcutsSecvtx, kTRUE);
1256 bTagged=kTRUE;
1257 }
1258 }
1259
1260 // Electron background analysis
1261 if (GetPlugin(kIsElecBackGround)) {
1262
1263 AliDebug(2, "Running BackGround Analysis");
1264
1265 for(Int_t jtrack = 0; jtrack < fESD->GetNumberOfTracks(); jtrack++){
1266 htrack = fESD->GetTrack(jtrack);
1267 if ( itrack == jtrack ) continue;
1268 fElecBackGround->PairAnalysis(track, htrack);
1269 }
1270 } // end of electron background analysis
1271
1272 if (GetPlugin(kDEstep)) {
1273 Double_t weightElecBgV0[kBgLevels] = {0.,0.,0.,};
1274 Int_t elecSource = 0;
1275 Double_t hfeimpactR=0., hfeimpactnsigmaR=0.;
1276 fExtraCuts->GetHFEImpactParameters(track, hfeimpactR, hfeimpactnsigmaR);
1277 if(HasMCData())
1278 {
1279 if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
1280 fQACollection->Fill("hadronsBeforeIPcut",track->Pt());
1281 }
1282 if(fMCQA && signal) {
1283
1284 fMCQA->SetContainerStep(0);
1285 for(Int_t iLevel = 0; iLevel < kBgLevels; iLevel++){
1286 weightElecBgV0[iLevel] = fMCQA->GetWeightFactor(mctrack, iLevel); // positive:conversion e, negative: nonHFE
1287 if(!fisNonHFEsystematics || IsPbPb())break;
1288 }
1289
1290 if(fisNonHFEsystematics){
1291 //Fill additional containers for electron source distinction
1292 elecSource = fMCQA->GetElecSource(mctrack->Particle());
1293 const Char_t *sourceName[kElecBgSpecies]={"Pion","Eta","Omega","Phi","EtaPrime","Rho"};
1294 const Char_t *levelName[kBgLevels]={"Best","Lower","Upper"};
1295 Int_t iName = 0;
1296 for(Int_t iSource = AliHFEmcQA::kPi0; iSource <= AliHFEmcQA::kGammaRho0; iSource++){
1297 if((iSource == AliHFEmcQA::kElse)||(iSource == AliHFEmcQA::kMisID)) continue;
1298 if(elecSource == iSource){
1299 for(Int_t iLevel = 0; iLevel < kBgLevels; iLevel++){
1300 if(weightElecBgV0[iLevel]>0) fVarManager->FillContainer(fContainer, Form("conversionElecs%s%s",sourceName[iName], levelName[iLevel]), 0, kFALSE, weightElecBgV0[iLevel]);
1301 else if(weightElecBgV0[iLevel]<0) fVarManager->FillContainer(fContainer, Form("mesonElecs%s%s",sourceName[iName], levelName[iLevel]), 0, kFALSE, -1*weightElecBgV0[iLevel]);
1302 if(IsPbPb())break;
1303 }
1304 break;
1305 }
1306 iName++;
1307 if(iName == kElecBgSpecies)iName = 0;
1308 }
1309 }
1310 //else{
1311 if(weightElecBgV0[0]>0) {
1312 fVarManager->FillContainer(fContainer, "conversionElecs", 0, kFALSE, weightElecBgV0[0]);
1313 fVarManager->FillContainer(fContainer, "conversionElecs", 5, kTRUE, weightElecBgV0[0]);
1314 }
1315 else if(weightElecBgV0[0]<0) {
1316 fVarManager->FillContainer(fContainer, "mesonElecs", 0, kFALSE, -1*weightElecBgV0[0]);
1317 fVarManager->FillContainer(fContainer, "mesonElecs", 5, kTRUE, -1*weightElecBgV0[0]);
1318 }
1319 //}
1320 if(bTagged){ // bg estimation for the secondary vertex tagged signals
1321 if(weightElecBgV0[0]>0) fVarManager->FillContainer(fContainer, "conversionElecs", 2, kFALSE, weightElecBgV0[0]);
1322 else if(weightElecBgV0[0]<0) fVarManager->FillContainer(fContainer, "mesonElecs", 2, kFALSE, -1*weightElecBgV0[0]);
1323 }
1324 }
1325 } // end of MC
1326
1327 dataDca[0]=-1; //for data, don't know the origin
1328 dataDca[1]=track->Pt();
1329 dataDca[2]=hfeimpactR;
1330 dataDca[3]=fCentralityF;
1331 dataDca[4] = v0pid;
1332 dataDca[5] = double(track->Charge());
1333 if (!HasMCData()) fQACollection->Fill("Dca", dataDca);
1334
1335 // Fill Containers for impact parameter analysis
1336 if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsDca + AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack,track)) continue;
1337 if(signal) {
1338 // Apply weight for background contamination after ip cut
1339 if(fBackGroundFactorApply) {
1340 fWeightBackGround = fkBackGroundFactorArray[0]->Eval(TMath::Abs(track->P())); // pp case
1341 if(fWeightBackGround < 0.0) fWeightBackGround = 0.0;
1342 else if(fWeightBackGround > 1.0) fWeightBackGround = 1.0;
1343 // weightBackGround as special weight
1344 fVarManager->FillContainer(fContainer, "hadronicBackground", 2, kFALSE, fWeightBackGround);
1345 }
1346 }
1347
1348 if(HasMCData()){
1349 if(fMCQA && signal) {
1350 fMCQA->SetContainerStep(1);
1351 for(Int_t iLevel = 0; iLevel < kBgLevels; iLevel++){
1352 weightElecBgV0[iLevel] = fMCQA->GetWeightFactor(mctrack, iLevel); // positive:conversion e, negative: nonHFE
1353 if(!fisNonHFEsystematics || IsPbPb())break;
1354 }
1355 if(fisNonHFEsystematics){
1356 //Fill additional containers for electron source distinction
1357 elecSource = fMCQA->GetElecSource(mctrack->Particle());
1358 const Char_t *sourceName[kElecBgSpecies]={"Pion","Eta","Omega","Phi","EtaPrime","Rho"};
1359 const Char_t *levelName[kBgLevels]={"Best","Lower","Upper"};
1360 Int_t iName = 0;
1361 for(Int_t iSource = AliHFEmcQA::kPi0; iSource <= AliHFEmcQA::kGammaRho0; iSource++){
1362 if((iSource == AliHFEmcQA::kElse)||(iSource == AliHFEmcQA::kMisID)) continue;
1363 if(elecSource == iSource){
1364 for(Int_t iLevel = 0; iLevel < kBgLevels; iLevel++){
1365 if(weightElecBgV0[iLevel]>0) fVarManager->FillContainer(fContainer, Form("conversionElecs%s%s",sourceName[iName], levelName[iLevel]), 1, kFALSE, weightElecBgV0[iLevel]);
1366 else if(weightElecBgV0[iLevel]<0) fVarManager->FillContainer(fContainer, Form("mesonElecs%s%s",sourceName[iName], levelName[iLevel]), 1, kFALSE, -1*weightElecBgV0[iLevel]);
1367 if(IsPbPb())break;
1368 }
1369 break;
1370 }
1371 iName++;
1372 if(iName == kElecBgSpecies)iName = 0;
1373 }
1374 }
1375 // else{
1376 if(weightElecBgV0[0]>0) {
1377 fVarManager->FillContainer(fContainer, "conversionElecs", 1, kFALSE, weightElecBgV0[0]);
1378 fVarManager->FillContainer(fContainer, "conversionElecs", 6, kTRUE, weightElecBgV0[0]);
1379 }
1380 else if(weightElecBgV0[0]<0) {
1381 fVarManager->FillContainer(fContainer, "mesonElecs", 1, kFALSE, -1*weightElecBgV0[0]);
1382 fVarManager->FillContainer(fContainer, "mesonElecs", 6, kTRUE, -1*weightElecBgV0[0]);
1383 }
1384 //}
1385 }
1386 }
1387 if(signal) {
1388 fVarManager->FillContainer(fContainer, "recTrackContDEReco", AliHFEcuts::kStepHFEcutsDca, kFALSE);
1389 fVarManager->FillContainer(fContainer, "recTrackContDEMC", AliHFEcuts::kStepHFEcutsDca, kTRUE);
1390 fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepafterDE"));
1391 }
1392 if(HasMCData()){
1393 if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
1394 fQACollection->Fill("hadronsAfterIPcut",track->Pt());
1395 }
1396 }
1397 }
1398
1399 }
1400
1401 // Background subtraction-------------------------------------------------------------------
1402 if (GetPlugin(kNonPhotonicElectron)) fBackgroundSubtraction->CountPoolAssociated(fInputEvent, fCentralityF);
1403 //------------------------------------------------------------------------------------------
1404
1405 fQACollection->Fill("nElectronTracksEvent", nElectronCandidates);
1406}
1407
1408//____________________________________________________________
1409void AliAnalysisTaskHFE::ProcessAOD(){
1410 //
1411 // Run Analysis in AOD Mode
1412 // Function is still in development
1413 //
1414 //printf("Process AOD\n");
1415 AliDebug(3, "Processing AOD Event");
1416 Double_t eventContainer[4];
1417 eventContainer[0] = 0.0;
1418 if(HasMCData()) eventContainer[0] = fVz;
1419 else {
1420 if(fInputEvent->GetPrimaryVertex()) eventContainer[0] = fInputEvent->GetPrimaryVertex()->GetZ();
1421 }
1422 eventContainer[1] = 1.; // No Information available in AOD analysis, assume all events have V0AND
1423 eventContainer[2] = fCentralityF;
1424 eventContainer[3] = fContributors;
1425
1426 //printf("value event container %f, %f, %f, %f\n",eventContainer[0],eventContainer[1],eventContainer[2],eventContainer[3]);
1427
1428 AliAODEvent *fAOD = dynamic_cast<AliAODEvent *>(fInputEvent);
1429 if(!fAOD){
1430 AliError("AOD Event required for AOD Analysis");
1431 return;
1432 }
1433
1434 //printf("Will fill\n");
1435 //
1436 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoCut);
1437 //printf("Fill\n");
1438 //
1439 if(fIdentifiedAsPileUp) return;
1440 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoPileUp);
1441
1442 //
1443 if(fIdentifiedAsOutInz) return;
1444 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepZRange);
1445
1446 //
1447 if(!fPassTheEventCut) return;
1448 fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepReconstructed);
1449 //printf("pass\n");
1450
1451 fContainer->NewEvent();
1452
1453 fCFM->SetRecEventInfo(fAOD);
1454
1455 if(!fExtraCuts){
1456 fExtraCuts = new AliHFEextraCuts("hfeExtraCuts","HFE Extra Cuts");
1457 }
1458 fExtraCuts->SetRecEventInfo(fAOD);
1459
1460 // Get Number of contributors to the primary vertex for multiplicity-dependent correction
1461 Int_t ncontribVtx = 0;
1462 AliAODVertex *priVtx = fAOD->GetPrimaryVertex();
1463 if(priVtx){
1464 ncontribVtx = priVtx->GetNContributors();
1465 }
1466
1467 // Look for kink mother
1468 Int_t numberofvertices = fAOD->GetNumberOfVertices();
1469 Double_t listofmotherkink[numberofvertices];
1470 Int_t numberofmotherkink = 0;
1471 for(Int_t ivertex=0; ivertex < numberofvertices; ivertex++) {
1472 AliAODVertex *aodvertex = fAOD->GetVertex(ivertex);
1473 if(!aodvertex) continue;
1474 if(aodvertex->GetType()==AliAODVertex::kKink) {
1475 AliAODTrack *mother = (AliAODTrack *) aodvertex->GetParent();
1476 if(!mother) continue;
1477 Int_t idmother = mother->GetID();
1478 listofmotherkink[numberofmotherkink] = idmother;
1479 //printf("ID %d\n",idmother);
1480 numberofmotherkink++;
1481 }
1482 }
1483 //printf("Number of kink mother in the events %d\n",numberofmotherkink);
1484
1485 // Background subtraction-------------------------------------------------------------------
1486 if (GetPlugin(kNonPhotonicElectron)) fBackgroundSubtraction->FillPoolAssociatedTracks(fInputEvent, fCentralityF);
1487 //------------------------------------------------------------------------------------------
1488
1489 // Loop over tracks
1490 AliAODTrack *track = NULL;
1491 AliAODMCParticle *mctrack = NULL;
1492 Double_t dataDca[6]; // [source, pT, dca, centrality]
1493 Int_t nElectronCandidates = 0;
1494 Bool_t signal;
1495
1496 //printf("Number of track %d\n",(Int_t) fAOD->GetNumberOfTracks());
1497 for(Int_t itrack = 0; itrack < fAOD->GetNumberOfTracks(); itrack++){
1498 track = fAOD->GetTrack(itrack); mctrack = NULL;
1499 if(!track) continue;
1500
1501 signal = kTRUE;
1502 if(HasMCData()){
1503
1504 Int_t label = TMath::Abs(track->GetLabel());
1505 if(label && label < fAODArrayMCInfo->GetEntriesFast())
1506 mctrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(label));
1507 if(fFillSignalOnly && !fCFM->CheckParticleCuts(AliHFEcuts::kStepMCGenerated, mctrack)) signal = kFALSE;
1508 }
1509
1510 fVarManager->NewTrack(track, mctrack, fCentralityF, -1, signal);
1511
1512 if(fFillNoCuts) {
1513 if(signal || !fFillSignalOnly){
1514 fVarManager->FillContainer(fContainer, "recTrackContReco", AliHFEcuts::kStepRecNoCut, kFALSE);
1515 fVarManager->FillContainer(fContainer, "recTrackContMC", AliHFEcuts::kStepRecNoCut, kTRUE);
1516 }
1517 }
1518
1519 if(fApplyCutAOD) {
1520 //printf("Apply cuts\n");
1521 // RecKine: ITSTPC cuts
1522 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
1523
1524 // Reject kink mother
1525 if(fRejectKinkMother) {
4437a0d2 1526 Bool_t kinkmotherpass = kTRUE;
1527 for(Int_t kinkmother = 0; kinkmother < numberofmotherkink; kinkmother++) {
1528 if(track->GetID() == listofmotherkink[kinkmother]) {
1529 kinkmotherpass = kFALSE;
1530 continue;
1531 }
1532 }
1533 if(!kinkmotherpass) continue;
01745870 1534 }
1535
1536 // RecPrim
1537 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
1538
1539 // HFEcuts: ITS layers cuts
1540 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
1541
1542 // HFE cuts: TOF PID and mismatch flag
1543 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTOF, track)) continue;
1544
1545 // HFE cuts: TPC PID cleanup
1546 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
1547
1548 // HFEcuts: Nb of tracklets TRD0
1549 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTRD, track)) continue;
1550 }
1551
1552 // Fill correlation maps before PID
1553 if(signal && fContainer->GetCorrelationMatrix("correlationstepbeforePID")) {
1554 //printf("Fill correlation maps before PID\n");
1555 fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepbeforePID"));
1556 }
1557
1558 if(HasMCData()){
1559 Double_t hfeimpactR4all=0., hfeimpactnsigmaR4all=0.;
1560 Int_t sourceDca =-1;
1561 if(mctrack && (TMath::Abs(mctrack->GetPdgCode()) == 211)){
1562 if(track->Pt()>4.){
1563 fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4all, hfeimpactnsigmaR4all);
1564 dataDca[0]=0; //pion
1565 dataDca[1]=track->Pt();
1566 dataDca[2]=hfeimpactR4all;
1567 dataDca[3]=fCentralityF;
1568 dataDca[4] = -1; // not store V0 for the moment
1569 dataDca[5] = double(track->Charge());
1570 fQACollection->Fill("Dca", dataDca);
1571 }
1572 }
1573 else if(mctrack && (TMath::Abs(mctrack->GetPdgCode()) == 11)){ // to increas statistics for Martin
1574 if(signal){
1575 fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4all, hfeimpactnsigmaR4all);
1576 if(fSignalCuts->IsCharmElectron(track)){
1577 sourceDca=1;
1578 }
1579 else if(fSignalCuts->IsBeautyElectron(track)){
1580 sourceDca=2;
1581 }
1582 else if(fSignalCuts->IsGammaElectron(track)){
1583 sourceDca=3;
1584 }
1585 else if(fSignalCuts->IsNonHFElectron(track)){
1586 sourceDca=4;
1587 }
1588 else if(fSignalCuts->IsJpsiElectron(track)){
1589 sourceDca=5;
1590 }
1591 else {
1592 sourceDca=6;
1593 }
1594 dataDca[0]=sourceDca;
1595 dataDca[1]=track->Pt();
1596 dataDca[2]=hfeimpactR4all;
1597 dataDca[3]=fCentralityF;
1598 dataDca[4] = -1; // not store V0 for the moment
1599 dataDca[5] = double(track->Charge());
1600 if(signal) fQACollection->Fill("Dca", dataDca);
1601 }
1602 }
1603 }
1604
1605 //printf("Will process to PID\n");
1606
1607 // track accepted, do PID
1608 AliHFEpidObject hfetrack;
1609 hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);
1610 hfetrack.SetRecTrack(track);
1611 if(HasMCData()) hfetrack.SetMCTrack(mctrack);
1612 hfetrack.SetCentrality(fCentralityF);
1613 hfetrack.SetMulitplicity(ncontribVtx); // for correction
4437a0d2 1614 if(IsPbPb()) hfetrack.SetPbPb();
1615 else{
1616 if(IspPb()) hfetrack.SetpPb();
1617 else hfetrack.SetPP();
1618 }
01745870 1619 fPID->SetVarManager(fVarManager);
1620 if(!fPID->IsSelected(&hfetrack, fContainer, "recTrackCont", fPIDqa)) continue;
1621 // we will do PID here as soon as possible
1622
1623 // Background subtraction----------------------------------------------------------------------------------------------
1624 if (GetPlugin(kNonPhotonicElectron)) {
1625 Int_t indexmother = -1;
1626 Int_t mcsource = -1;
39427568 1627 if(HasMCData() && mctrack) mcsource = fBackgroundSubtraction->FindMother(TMath::Abs(track->GetLabel()),indexmother);
01745870 1628 fBackgroundSubtraction->LookAtNonHFE(itrack, track, fInputEvent, 1, fCentralityF, -1,mcsource, indexmother);
1629 }
1630 //---------------------------------------------------------------------------------------------------------------------
1631
1632 // end AOD QA
1633 fQACollection->Fill("Filterend", -1);
1634 for(Int_t k=0; k<20; k++) {
1635 Int_t u = 1<<k;
1636 if((track->TestFilterBit(u))) {
1637 fQACollection->Fill("Filterend", k);
1638 }
1639 }
1640
1641 // Apply weight for background contamination
1642 //Double_t weightBackGround = 1.0;
1643 if(signal) {
1644 // Apply weight for background contamination
1645 if(fBackGroundFactorApply) {
4437a0d2 1646 if(IsPbPb() && fCentralityF >= 0) fWeightBackGround = fkBackGroundFactorArray[fCentralityF >= 0 ? fCentralityF : 0]->Eval(TMath::Abs(track->P()));
01745870 1647 else fWeightBackGround = fkBackGroundFactorArray[0]->Eval(TMath::Abs(track->P())); // pp case
1648
1649 if(fWeightBackGround < 0.0) fWeightBackGround = 0.0;
1650 else if(fWeightBackGround > 1.0) fWeightBackGround = 1.0;
1651 // weightBackGround as special weight
1652 fVarManager->FillContainer(fContainer, "hadronicBackground", 1, kFALSE, fWeightBackGround);
1653 }
1654 fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepafterPID"));
1655 }
1656
1657 nElectronCandidates++;
1658
1659 if (GetPlugin(kDEstep)) {
1660 if (!HasMCData()){
1661 Double_t hfeimpactR=0., hfeimpactnsigmaR=0.;
1662 fExtraCuts->GetHFEImpactParameters(track, hfeimpactR, hfeimpactnsigmaR);
1663 dataDca[0]=-1; //for data, don't know the origin
1664 dataDca[1]=track->Pt();
1665 dataDca[2]=hfeimpactR;
1666 dataDca[3]=fCentralityF;
1667 dataDca[4] = -1; // not store V0 for the moment
1668 dataDca[5] = double(track->Charge());
1669 fQACollection->Fill("Dca", dataDca);
1670 }
1671
1672 // Fill Containers for impact parameter analysis
1673 if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsDca + AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack,track)) continue;
1674 if(signal) {
1675 // Apply weight for background contamination after ip cut
1676 if(fBackGroundFactorApply) {
1677 fWeightBackGround = fkBackGroundFactorArray[0]->Eval(TMath::Abs(track->P())); // pp case
1678 if(fWeightBackGround < 0.0) fWeightBackGround = 0.0;
1679 else if(fWeightBackGround > 1.0) fWeightBackGround = 1.0;
1680 // weightBackGround as special weight
1681 fVarManager->FillContainer(fContainer, "hadronicBackground", 2, kFALSE, fWeightBackGround);
1682 }
1683 }
1684 }
1685 }
1686
1687 // Background subtraction-------------------------------------------------------------------
1688 if (GetPlugin(kNonPhotonicElectron)) fBackgroundSubtraction->CountPoolAssociated(fInputEvent, fCentralityF);
1689 //------------------------------------------------------------------------------------------
1690
1691 fQACollection->Fill("nElectronTracksEvent", nElectronCandidates);
1692}
1693
1694//____________________________________________________________
1695Bool_t AliAnalysisTaskHFE::ProcessMCtrack(AliVParticle *track){
1696 //
1697 // Filter the Monte Carlo Track
1698 // Additionally Fill a THnSparse for Signal To Background Studies
1699 // Works for AOD and MC analysis Type
1700 //
1701 fVarManager->NewTrack(track, NULL, fCentralityF, -1, kTRUE);
1702
1703
1704 Double_t vertex[3] = {0.,0.,0.}; // Production vertex cut to mask gammas which are NOT supposed to have hits in the first ITS layer(s)
1705 if(IsESDanalysis()){
1706 AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(track);
1707 if(mctrack){
1708 vertex[0] = mctrack->Particle()->Vx();
1709 vertex[1] = mctrack->Particle()->Vy();
1710 }
1711 } else {
1712 AliAODMCParticle *aodmctrack = dynamic_cast<AliAODMCParticle *>(track);
1713 if(aodmctrack) aodmctrack->XvYvZv(vertex);
1714 }
1715
1716 //printf("MC Generated\n");
1717 if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepMCGenerated, track)) return kFALSE;
1718 //printf("MC Generated pass\n");
1719 fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCGenerated, kFALSE);
1720
1721 // Step GeneratedZOutNoPileUp
1722 if((fIdentifiedAsPileUp) || (TMath::Abs(fVz) > fCuts->GetVertexRange()) || (fCentralityF < 0)) return kFALSE;
1723 fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCGeneratedZOutNoPileUpCentralityFine, kFALSE);
1724 //printf("In ProcessMCtrack %f\n",fCentralityF);
1725
1726 // Step Generated Event Cut
1727 if(!fPassTheEventCut) return kFALSE;
1728 fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCGeneratedEventCut, kFALSE);
1729
1730 if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepMCInAcceptance, track)) return kFALSE;
1731 fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCInAcceptance, kFALSE);
1732 return kTRUE;
1733}
1734
1735//____________________________________________________________
1736Bool_t AliAnalysisTaskHFE::PreSelectTrack(AliESDtrack *track) const {
1737 //
1738 // Preselect tracks
1739 //
1740
1741
1742 Bool_t survived = kTRUE;
1743
1744 if(fCutspreselect) {
1745 //printf("test preselect\n");
1746 if(!fCutspreselect->IsSelected(track)) survived=kFALSE;
1747 }
1748 //printf("survived %d\n",(Int_t)survived);
1749
1750 if(survived && fPIDpreselect){
1751 // Apply PID
1752 AliHFEpidObject hfetrack;
1753 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1754 hfetrack.SetRecTrack(track);
1755 if(!fPIDpreselect->IsSelected(&hfetrack)) {
1756 //printf("Did not pass AliHFEcuts::kPID\n");
1757 survived = kFALSE;
1758 }
1759 //else printf("Pass AliHFEcuts::kPID\n");
1760 }
1761
1762 return survived;
1763
1764}
1765//____________________________________________________________
1766void AliAnalysisTaskHFE::MakeEventContainer(){
1767 //
1768 // Create the event container for the correction framework and link it
1769 // 1st bin: Vertex z-position
1770 // 2nd bin: V0AND decision (normalization to sigma_inel)
1771 // 3rd bin: Centrality class (for pp defined as number of contributors in vertex.)
1772 // 4th bin: Number of contributors > 0
1773 //
1774
1775 const Int_t kNvar = 4; // number of variables on the grid:
1776 Int_t nBins[kNvar] = {120, 2, 11, 2};
1777 Double_t binMin[kNvar] = {-30. , 0., 0.0, 0.};
1778 Double_t binMax[kNvar] = {30., 2., 11.0, 2.};
1779
1780 AliCFContainer *evCont = new AliCFContainer("eventContainer", "Container for events", AliHFEcuts::kNcutStepsEvent, kNvar, nBins);
1781
1782 Double_t *vertexBins = AliHFEtools::MakeLinearBinning(nBins[0], binMin[0], binMax[0]);
1783 Double_t *v0andBins = AliHFEtools::MakeLinearBinning(nBins[1], binMin[1], binMax[1]);
1784 Double_t *centralityBins = AliHFEtools::MakeLinearBinning(nBins[2], binMin[2], binMax[2]);
1785 Double_t *contributorsBins = AliHFEtools::MakeLinearBinning(nBins[3], binMin[3], binMax[3]);
1786 evCont->SetBinLimits(0, vertexBins);
1787 evCont->SetBinLimits(1, v0andBins);
1788 evCont->SetBinLimits(2, centralityBins);
1789 evCont->SetBinLimits(3, contributorsBins);
1790 delete[] vertexBins; delete[] v0andBins; delete[] centralityBins; delete[] contributorsBins;
1791
1792 fCFM->SetEventContainer(evCont);
1793}
1794
1795//____________________________________________________________
1796void AliAnalysisTaskHFE::MakeParticleContainer(){
1797 //
1798 // Create the particle container for the correction framework manager and
1799 // link it
1800 //
1801 if(!fContainer) fContainer = new AliHFEcontainer("trackContainer");
1802 fVarManager->DefineVariables(fContainer);
1803
1804 // Create Correction Framework containers
1805 fContainer->CreateContainer("MCTrackCont", "Track Container filled with MC information", AliHFEcuts::kNcutStepsMCTrack);
1806 fContainer->CreateContainer("recTrackContReco", "Track Container filled with MC information", AliHFEcuts::kNcutStepsRecTrack + fPID->GetNumberOfPIDdetectors());
1807 fContainer->CreateContainer("recTrackContMC", "Track Container filled with MC information", AliHFEcuts::kNcutStepsRecTrack + fPID->GetNumberOfPIDdetectors());
1808
1809 fContainer->CreateContainer("hadronicBackground", "Container for Hadronic Background", 3);
1810 fContainer->CreateContainer("recTrackContDEReco", "Container for displaced electron analysis with Reco information", 1);
1811 fContainer->CreateContainer("recTrackContDEMC", "Container for displaced electron analysis with MC information", 1);
1812 fContainer->CreateContainer("recTrackContSecvtxReco", "Container for secondary vertexing analysis with Reco information", 1);
1813 fContainer->CreateContainer("recTrackContSecvtxMC", "Container for secondary vertexing analysis with MC information", 1);
1814
1815 if(HasMCData()){
1816 fContainer->CreateContainer("conversionElecs", "Container for weighted conversion electrons",7);
1817 fContainer->CreateContainer("mesonElecs", "Container for weighted electrons from meson decays",7);
1818 fContainer->Sumw2("conversionElecs");
1819 fContainer->Sumw2("mesonElecs");
1820
1821 if(fisNonHFEsystematics){
1822 const Char_t *sourceName[kElecBgSpecies]={"Pion","Eta","Omega","Phi","EtaPrime","Rho"};
1823 const Char_t *levelName[kBgLevels]={"Best","Lower","Upper"};
1824 for(Int_t iSource = 0; iSource < kElecBgSpecies; iSource++){
1825 for(Int_t iLevel = 0; iLevel < kBgLevels; iLevel++){
1826 fContainer->CreateContainer(Form("conversionElecs%s%s",sourceName[iSource],levelName[iLevel]), Form("Container for weighted conversion electrons from %s grandm., %s level",sourceName[iSource],levelName[iLevel]),5);
1827 fContainer->CreateContainer(Form("mesonElecs%s%s",sourceName[iSource],levelName[iLevel]), Form("Container for weighted electrons from %s decays, %s level",sourceName[iSource],levelName[iLevel]),5);
1828 fContainer->Sumw2(Form("conversionElecs%s%s",sourceName[iSource],levelName[iLevel]));
1829 fContainer->Sumw2(Form("mesonElecs%s%s",sourceName[iSource],levelName[iLevel]));
1830 if(IsPbPb())break;
1831 }
1832 }
1833 }
1834 //fContainer->CreateContainer("charmElecs", "Container for weighted charm electrons",2);
1835 }
1836
1837 fContainer->CreateCorrelationMatrix("correlationstepafterPID","THnSparse with correlations");
1838 fContainer->CreateCorrelationMatrix("correlationstepafterDE","THnSparse with correlations");
1839 if(!fVarManager->IsVariableDefined("centrality")) {
1840 //printf("Create the two other correlation maps\n");
1841 fContainer->CreateCorrelationMatrix("correlationstepbeforePID","THnSparse with correlations");
1842 fContainer->CreateCorrelationMatrix("correlationstepafterTOF","THnSparse with correlations");
1843 }
1844
1845 // Define the step names
1846 for(UInt_t istep = 0; istep < AliHFEcuts::kNcutStepsMCTrack; istep++){
1847 fContainer->SetStepTitle("MCTrackCont", AliHFEcuts::MCCutName(istep), istep);
1848 }
1849 for(UInt_t istep = 0; istep < AliHFEcuts::kNcutStepsRecTrack; istep++){
1850 fContainer->SetStepTitle("recTrackContReco", AliHFEcuts::RecoCutName(istep), istep);
1851 fContainer->SetStepTitle("recTrackContMC", AliHFEcuts::RecoCutName(istep), istep);
1852 }
1853 for(UInt_t ipid = 0; ipid < fPID->GetNumberOfPIDdetectors(); ipid++){
1854 fContainer->SetStepTitle("recTrackContReco", fPID->SortedDetectorName(ipid), AliHFEcuts::kNcutStepsRecTrack + ipid);
1855 fContainer->SetStepTitle("recTrackContMC", fPID->SortedDetectorName(ipid), AliHFEcuts::kNcutStepsRecTrack + ipid);
1856 }
1857}
1858//____________________________________________________________
1859void AliAnalysisTaskHFE::InitContaminationQA(){
1860 //
1861 // Add QA for Impact Parameter cut
1862 //
1863
1864 TObjArray *array = fVarManager->GetVariables();
1865 Int_t nvars = array->GetEntriesFast();
1866 for(Int_t v = 0; v < nvars; v++) {
1867 AliHFEvarManager::AliHFEvariable *variable = (AliHFEvarManager::AliHFEvariable *) array->At(v);
1868 if(!variable) continue;
1869 TString name(((AliHFEvarManager::AliHFEvariable *)variable)->GetName());
1870 if(!name.CompareTo("pt")) {
1871 const Int_t nBinPt = variable->GetNumberOfBins();
1872 const Double_t *kPtRange = variable->GetBinning();
1873
1874 fQACollection->CreateTH1Farray("hadronsBeforeIPcut", "Hadrons before IP cut", nBinPt, kPtRange);
1875 fQACollection->CreateTH1Farray("hadronsAfterIPcut", "Hadrons after IP cut", nBinPt, kPtRange);
1876
1877 fQACollection->CreateTH2Farray("Ke3Kecorr", "Ke3 decay e and K correlation; Ke3K p_{t}; Ke3e p_{t}; ", nBinPt, kPtRange, 20,0.,20.);
1878 fQACollection->CreateTH2Farray("Ke3K0Lecorr", "Ke3 decay e and K0L correlation; Ke3K0L p_{t}; Ke3e p_{t}; ", nBinPt, kPtRange, 20,0.,20.);
1879 fQACollection->CreateTH1Farray("Kptspectra", "Charged Kaons: MC p_{t} ", nBinPt, kPtRange);
1880 fQACollection->CreateTH1Farray("K0Lptspectra", "K0L: MC p_{t} ", nBinPt, kPtRange);
1881
1882 const Double_t kDCAbound[2] = {-0.2, 0.2};
1883
1884 const Int_t nDimDca=6;
1885 const Int_t nBinDca[nDimDca] = { 8, nBinPt, 800, 12, 6, 2};
1886 Double_t minimaDca[nDimDca] = { -1., 0., kDCAbound[0], -1., -1, -1.1};
1887 Double_t maximaDca[nDimDca] = { 7., 20., kDCAbound[1], 11., 5, 1.1};
1888
1889 Double_t *sourceBins = AliHFEtools::MakeLinearBinning(nBinDca[0], minimaDca[0], maximaDca[0]);
1890 Double_t *dcaBins = AliHFEtools::MakeLinearBinning(nBinDca[2], minimaDca[2], maximaDca[2]);
1891 Double_t *centralityBins = AliHFEtools::MakeLinearBinning(nBinDca[3], minimaDca[3], maximaDca[3]);
1892 Double_t *v0PIDBins = AliHFEtools::MakeLinearBinning(nBinDca[4], minimaDca[4], maximaDca[4]);
1893 Double_t *chargeBins = AliHFEtools::MakeLinearBinning(nBinDca[5], minimaDca[5], maximaDca[5]);
1894
1895 fQACollection->CreateTHnSparseNoLimits("Dca", "Dca; source (0-all, 1-charm,etc); pT [GeV/c]; dca; centrality bin; v0pid; charge", nDimDca, nBinDca);
1896 ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(0, sourceBins);
1897 ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(1, kPtRange);
1898 ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(2, dcaBins);
1899 ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(3, centralityBins);
1900 ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(4, v0PIDBins);
1901 ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(5, chargeBins);
1902
1903 break;
1904 }
1905 }
1906
1907}
1908
1909//____________________________________________________________
1910void AliAnalysisTaskHFE::InitHistoITScluster(){
1911 //
1912 // Initialize a temporary histogram to monitor the chi2/ITS cluster
1913 if(IsPbPb()) {
1914 const Int_t kNDim = 7;
1915 const Int_t kNBins[kNDim] = {88, 20,90,11, 160, 2, 1000};
1916 const Double_t kMin[kNDim] = {0.1, -1,0, 0.,0., 0, 0.};
1917 const Double_t kMax[kNDim] = {20., 1, 2.*TMath::Pi(), 11.,160, 2, 100.};
1918 fQACollection->CreateTHnSparse("fChi2perITScluster", "chi2/ITS cluster; p_{T} (GeV/c);eta;phi; centrality class;nclus;sharebit; #chi^{2}/ITS cluster", kNDim, kNBins, kMin, kMax);
1919 fQACollection->BinLogAxis("fChi2perITScluster", 0);
1920 }
1921 else
1922 {
1923 const Int_t kNDim = 3;
1924 const Int_t kNBins[kNDim] = {44, 11, 1000};
1925 const Double_t kMin[kNDim] = {0.1, 0., 0.};
1926 const Double_t kMax[kNDim] = {20., 11., 100.};
1927 fQACollection->CreateTHnSparse("fChi2perITScluster", "chi2/ITS cluster; p_{T} (GeV/c); centrality class; #chi^{2}/ITS cluster", kNDim, kNBins, kMin, kMax);
1928 fQACollection->BinLogAxis("fChi2perITScluster", 0);
1929 }
1930}
1931
1932//____________________________________________________________
1933void AliAnalysisTaskHFE::SelectSpecialTrigger(const Char_t *trgclust, Int_t runMin, Int_t runMax){
1934 //
1935 // Select only events triggered by a special trigeer cluster
1936 //
1937 if(!fSpecialTrigger) fSpecialTrigger = new AliOADBContainer("SpecialTrigger");
1938 fSpecialTrigger->AppendObject(new TObjString(trgclust), runMin, runMax);
1939}
1940
1941//____________________________________________________________
1942const Char_t * AliAnalysisTaskHFE::GetSpecialTrigger(Int_t run){
1943 //
1944 // Derive selected trigger string for given run
1945 //
1946 if(!fSpecialTrigger) return NULL;
1947 TObjString *trg = dynamic_cast<TObjString *>(fSpecialTrigger->GetObject(run));
1948 if(!trg) return NULL;
1949 return trg->String().Data();
1950}
1951
1952//____________________________________________________________
1953void AliAnalysisTaskHFE::PrintStatus() const {
1954 //
1955 // Print Analysis status
1956 //
1957 printf("\n\tAnalysis Settings\n\t========================================\n\n");
1958 printf("\tSecondary Vertex finding: %s\n", GetPlugin(kSecVtx) ? "YES" : "NO");
1959 printf("\tPrimary Vertex resolution: %s\n", GetPlugin(kPriVtx) ? "YES" : "NO");
1960 printf("\tDisplaced electron analysis step: %s\n", GetPlugin(kDEstep) ? "YES" : "NO");
1961 printf("\tTagged Track Analysis: %s\n", GetPlugin(kTaggedTrackAnalysis) ? "YES" : "NO");
1962 printf("\n");
1963 printf("\tParticle Identification Detectors:\n");
1964 fPID->PrintStatus();
1965 printf("\n");
1966 printf("\tQA: \n");
1967 printf("\t\tPID: %s\n", IsQAOn(kPIDqa) ? "YES" : "NO");
1968 printf("\t\tCUTS: %s\n", (fCuts != NULL && fCuts->IsQAOn()) ? "YES" : "NO");
1969 printf("\t\tMC: %s\n", IsQAOn(kMCqa) ? "YES" : "NO");
1970 printf("\n");
1971}
1972
1973//____________________________________________________________
1974Bool_t AliAnalysisTaskHFE::FillProductionVertex(const AliVParticle * const track) const{
1975 //
1976 // Find the production vertex of the associated MC track
1977 //
1978 if(!fMCEvent) return kFALSE;
1979 const AliVParticle *mctrack = NULL;
1980 TString objectType = track->IsA()->GetName();
1981 if(objectType.CompareTo("AliESDtrack") == 0 || objectType.CompareTo("AliAODTrack") == 0){
1982 // Reconstructed track
1983 mctrack = fMCEvent->GetTrack(TMath::Abs(track->GetLabel()));
1984 } else {
1985 // MCParticle
1986 mctrack = track;
1987 }
1988
1989 if(!mctrack) return kFALSE;
1990
1991 Double_t xv = 0.0;
1992 Double_t yv = 0.0;
1993
1994 if(TString(mctrack->IsA()->GetName()).CompareTo("AliMCParticle") == 0){
1995 // case MCParticle
1996 const AliMCParticle *mcpart = dynamic_cast<const AliMCParticle *>(mctrack);
1997 if(mcpart){
1998 xv = mcpart->Xv();
1999 yv = mcpart->Yv();
2000 }
2001 } else {
2002 // case AODMCParticle
2003 const AliAODMCParticle *mcpart = dynamic_cast<const AliAODMCParticle *>(mctrack);
2004 if(mcpart){
2005 xv = mcpart->Xv();
2006 yv = mcpart->Yv();
2007 }
2008 }
2009
2010 //printf("xv %f, yv %f\n",xv,yv);
2011 fQACollection->Fill("radius", TMath::Abs(xv),TMath::Abs(yv));
2012
2013 return kTRUE;
2014
2015}
2016//__________________________________________
2017void AliAnalysisTaskHFE::SwitchOnPlugin(Int_t plug){
2018 //
2019 // Switch on Plugin
2020 // Available:
2021 // - Primary vertex studies
2022 // - Secondary vertex Studies
2023 // - Post Processing
2024 //
2025 switch(plug){
2026 case kPriVtx: SETBIT(fPlugins, plug); break;
2027 case kSecVtx: SETBIT(fPlugins, plug); break;
2028 case kIsElecBackGround: SETBIT(fPlugins, plug); break;
2029 case kPostProcess: SETBIT(fPlugins, plug); break;
2030 case kDEstep: SETBIT(fPlugins, plug); break;
2031 case kTaggedTrackAnalysis: SETBIT(fPlugins, plug); break;
2032 case kNonPhotonicElectron: SETBIT(fPlugins, plug); break;
2033 default: AliError("Unknown Plugin");
2034 };
2035}
2036//__________________________________________
2037Bool_t AliAnalysisTaskHFE::ProcessCutStep(Int_t cutStep, AliVParticle *track){
2038 //
2039 // Check single track cuts for a given cut step
2040 // Fill the particle container
2041 //
2042 const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
2043 if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
2044 if(fVarManager->IsSignalTrack()) {
2045 fVarManager->FillContainer(fContainer, "recTrackContReco", cutStep, kFALSE);
2046 fVarManager->FillContainer(fContainer, "recTrackContMC", cutStep, kTRUE);
2047 }
2048 return kTRUE;
2049}
2050//___________________________________________________
2051Bool_t AliAnalysisTaskHFE::ReadCentrality() {
2052 //
2053 // Recover the centrality of the event from ESD or AOD
2054 //
2055
2056 Float_t fCentralityLimitstemp[12];
2057 Float_t fCentralityLimitsdefault[12]= {0.,5.,10., 20., 30., 40., 50., 60.,70.,80., 90., 100.};
2058 if(!fPbPbUserCentralityBinning) memcpy(fCentralityLimitstemp,fCentralityLimitsdefault,sizeof(fCentralityLimitsdefault));
2059 else memcpy(fCentralityLimitstemp,fCentralityLimits,sizeof(fCentralityLimitsdefault));
2060
2061
2062 Int_t bin = -1;
4437a0d2 2063 if(IsPbPb()||IspPb()) {
01745870 2064 // Centrality
2065 AliCentrality *centrality = fInputEvent->GetCentrality();
2066 fCentralityPercent = centrality->GetCentralityPercentile(fCentralityEstimator.Data());
2067 //printf("centrality %f\n",fCentralityPercent);
2068
2069 for(Int_t ibin = 0; ibin < 11; ibin++){
2070 if(fCentralityPercent >= fCentralityLimitstemp[ibin] && fCentralityPercent < fCentralityLimitstemp[ibin+1]){
2071 bin = ibin;
2072 //printf("test bin %f, low %f, high %f, %d\n",fCentralityPercent,fCentralityLimitstemp[ibin],fCentralityLimitstemp[ibin+1],ibin);
2073 break;
2074 }
2075 }
2076
2077 if(bin == -1) bin = 11; // Overflow
2078 } else {
2079 // PP: Tracklet multiplicity, use common definition
2080 Int_t itsMultiplicity = GetITSMultiplicity(fInputEvent);
2081 Int_t multiplicityLimits[8] = {0, 1, 9, 17, 25, 36, 60, 500};
2082 for(Int_t ibin = 0; ibin < 7; ibin++){
2083 if(itsMultiplicity >= multiplicityLimits[ibin] && itsMultiplicity < multiplicityLimits[ibin + 1]){
2084 bin = ibin;
2085 break;
2086 }
2087 }
2088 if(bin == -1) bin = 7; // Overflow
2089 }
2090 fCentralityF = bin;
2091 AliDebug(2, Form("Centrality class %d\n", fCentralityF));
2092
2093
2094 // contributors, to be outsourced
2095 const AliVVertex *vtx;
2096 if(IsAODanalysis()){
2097 AliAODEvent *fAOD = dynamic_cast<AliAODEvent *>(fInputEvent);
2098 if(!fAOD){
2099 AliError("AOD Event required for AOD Analysis");
2100 return kFALSE;
2101 }
2102 vtx = fAOD->GetPrimaryVertex();
2103 } else {
2104 AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
2105 if(!fESD){
2106 AliError("ESD Event required for ESD Analysis");
2107 return kFALSE;
2108 }
2109 vtx = fESD->GetPrimaryVertex() ;
2110 }
2111 if(!vtx){
2112 fContributors = 0.5;
2113 return kFALSE;
2114 }
2115 else {
2116 Int_t contributorstemp = vtx->GetNContributors();
2117 if( contributorstemp <= 0) {
2118 fContributors = 0.5;
2119 //printf("Number of contributors %d and vz %f\n",contributorstemp,vtx->GetZ());
2120 }
2121 else fContributors = 1.5;
2122 //printf("Number of contributors %d\n",contributorstemp);
2123 }
2124 return kTRUE;
2125}
2126
2127//___________________________________________________
2128Int_t AliAnalysisTaskHFE::GetITSMultiplicity(AliVEvent *ev){
2129 //
2130 // Definition of the Multiplicity according to the JPSI group (F. Kramer)
2131 //
2132 Int_t nTracklets = 0;
2133 Int_t nAcc = 0;
2134 Double_t etaRange = 1.6;
2135
2136 if (ev->IsA() == AliAODEvent::Class()) {
2137 AliAODTracklets *tracklets = ((AliAODEvent*)ev)->GetTracklets();
2138 nTracklets = tracklets->GetNumberOfTracklets();
2139 for (Int_t nn = 0; nn < nTracklets; nn++) {
2140 Double_t theta = tracklets->GetTheta(nn);
2141 Double_t eta = -TMath::Log(TMath::Tan(theta/2.0));
2142 if (TMath::Abs(eta) < etaRange) nAcc++;
2143 }
2144 } else if (ev->IsA() == AliESDEvent::Class()) {
2145 nTracklets = ((AliESDEvent*)ev)->GetMultiplicity()->GetNumberOfTracklets();
2146 for (Int_t nn = 0; nn < nTracklets; nn++) {
2147 Double_t eta = ((AliESDEvent*)ev)->GetMultiplicity()->GetEta(nn);
2148 if (TMath::Abs(eta) < etaRange) nAcc++;
2149 }
2150 } else return -1;
2151
2152 return nAcc;
2153}
2154
2155//___________________________________________________
2156void AliAnalysisTaskHFE::RejectionPileUpVertexRangeEventCut() {
2157 //
2158 // Recover the centrality of the event from ESD or AOD
2159 //
2160 if(IsAODanalysis()){
2161
2162 AliAODEvent *fAOD = dynamic_cast<AliAODEvent *>(fInputEvent);
2163 if(!fAOD){
2164 AliError("AOD Event required for AOD Analysis");
2165 return;
2166 }
2167 // PileUp
2168 fIdentifiedAsPileUp = kFALSE;
2169 if(fRemovePileUp && fAOD->IsPileupFromSPD()) fIdentifiedAsPileUp = kTRUE;
2170 // Z vertex
2171 fIdentifiedAsOutInz = kFALSE;
2172 //printf("Z vertex %f and out %f\n",fAOD->GetPrimaryVertex()->GetZ(),fCuts->GetVertexRange());
2173 if(TMath::Abs(fAOD->GetPrimaryVertex()->GetZ()) > fCuts->GetVertexRange()) fIdentifiedAsOutInz = kTRUE;
2174 // Event Cut
2175 fPassTheEventCut = kTRUE;
2176 if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fAOD)) fPassTheEventCut = kFALSE;
2177
2178
2179 } else {
2180
2181 AliDebug(3, "Processing ESD Centrality");
2182 AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
2183 if(!fESD){
2184 AliError("ESD Event required for ESD Analysis");
2185 return;
2186 }
2187 // PileUp
2188 fIdentifiedAsPileUp = kFALSE;
2189 if(fRemovePileUp && fESD->IsPileupFromSPD()) fIdentifiedAsPileUp = kTRUE;
2190
2191
2192
2193 // Z vertex
2194 fIdentifiedAsOutInz = kFALSE;
2195 Bool_t findvertex = kTRUE;
2196 const AliESDVertex* vtxESD = fESD->GetPrimaryVertex();
2197 if((!vtxESD) || (vtxESD->GetNContributors() <= 0)) findvertex = kFALSE;
2198 if(findvertex) {
2199 if(TMath::Abs(vtxESD->GetZ()) > fCuts->GetVertexRange()) fIdentifiedAsOutInz = kTRUE;
2200 }
2201
2202 //Event Cut
2203 fPassTheEventCut = kTRUE;
2204 if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) fPassTheEventCut = kFALSE;
2205
2206 }
2207
2208}
8c07fb0d 2209
2210//___________________________________________________
4437a0d2 2211Bool_t AliAnalysisTaskHFE::CheckTRDTriggerESD(AliESDEvent *ev) {
8c07fb0d 2212//
2213// Check TRD trigger; pPb settings
2214//
2215 Bool_t cint8=kFALSE;
2216 Bool_t cint7=kFALSE;
2217 Bool_t cint5=kFALSE;
2218 Bool_t trdtrgevent=kFALSE;
2219
4437a0d2 2220
2221 // mb selection of WU events
8c07fb0d 2222 if(fWhichTRDTrigger==1)
2223 {
4437a0d2 2224 if(ev->IsTriggerClassFired("CINT7WU-B-NOPF-ALL"))
8c07fb0d 2225 {
2226 DrawTRDTrigger(ev);
2227 return kTRUE;
2228 }
4437a0d2 2229 else return kFALSE;
8c07fb0d 2230 }
618038b6 2231
4437a0d2 2232
2233 // HSE no cleanup
8c07fb0d 2234 if(fWhichTRDTrigger==2)
2235 {
4437a0d2 2236 cint8= ev->IsTriggerClassFired("CINT8WUHSE-B-NOPF-CENT");
2237 cint7= ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT");
2238 cint5= (ev->IsTriggerClassFired("CINT5WU-B-NOPF-ALL")) &&
2239 (ev->GetHeader()->GetL1TriggerInputs() & (1 << 10));
2240 if((cint7==kFALSE)&&(cint8==kFALSE)&&(cint5==kFALSE)) return kFALSE;
8c07fb0d 2241 else
2242 {
2243 DrawTRDTrigger(ev);
2244 return kTRUE;
2245 }
2246 }
2247
4437a0d2 2248
2249
618038b6 2250 //HQU no cleanup
8c07fb0d 2251 if(fWhichTRDTrigger==3)
2252 {
8c07fb0d 2253 cint8= ev->IsTriggerClassFired("CINT8WUHQU-B-NOPF-CENT");
2254 cint7= ev->IsTriggerClassFired("CINT7WUHQU-B-NOPF-CENT");
2255 cint5= (ev->IsTriggerClassFired("CINT5WU-B-NOPF-ALL")) &&
2256 (ev->GetHeader()->GetL1TriggerInputs() & (1 << 12));
8c07fb0d 2257 if((cint7==kFALSE)&&(cint8==kFALSE)&&(cint5==kFALSE)) return kFALSE;
2258 else
2259 {
2260 DrawTRDTrigger(ev);
2261 return kTRUE;
2262 }
8c07fb0d 2263 }
618038b6 2264
4437a0d2 2265
2266
2267 return trdtrgevent;
2268
2269}
2270
2271
2272//___________________________________________________
2273Bool_t AliAnalysisTaskHFE::CheckTRDTrigger(AliVEvent *ev) {
2274//
2275// Check TRD trigger; pPb settings
2276//
2277
2278 fTRDTriggerAnalysis->CalcTriggers(ev);
2279
2280 // HSE cleanup
2281 if(fWhichTRDTrigger==6)
2282 {
2283 if(!fTRDTriggerAnalysis->IsFired(AliTRDTriggerAnalysis::kHSE))
2284 {
2285 return kFALSE;
2286 }
2287 else
2288 {
2289 // DrawTRDTrigger(ev);
2290 return kTRUE;
2291 }
2292 }
2293
2294
2295
618038b6 2296 // HQU cleanup
4437a0d2 2297 if(fWhichTRDTrigger==7)
8c07fb0d 2298 {
618038b6 2299
2300 if(!fTRDTriggerAnalysis->IsFired(AliTRDTriggerAnalysis::kHQU))
2301 {
2302 return kFALSE;
2303 }
8c07fb0d 2304 else
2305 {
4437a0d2 2306 // DrawTRDTrigger(ev);
8c07fb0d 2307 return kTRUE;
2308 }
2309 }
4437a0d2 2310
2311 // HSE or HQU cleanup
2312 if(fWhichTRDTrigger==8)
2313 {
2314 if((fTRDTriggerAnalysis->IsFired(AliTRDTriggerAnalysis::kHSE))||(fTRDTriggerAnalysis->IsFired(AliTRDTriggerAnalysis::kHQU)))
2315 {
2316 // DrawTRDTrigger(ev);
2317 return kTRUE;
2318 }
2319 else
2320 {
2321 return kFALSE;
2322 }
2323 }
2324
2325 return kFALSE;
8c07fb0d 2326
2327}
2328
2329void AliAnalysisTaskHFE::DrawTRDTrigger(AliESDEvent *ev) {
2330
2331 Int_t ntriggerbit=0;
2332 fQACollection->Fill("nTriggerBit",ntriggerbit);
2333 if(ev->IsTriggerClassFired("CINT7-B-NOPF-ALLNOTRD"))
2334 {
2335 ntriggerbit=2;
2336 fQACollection->Fill("nTriggerBit",ntriggerbit);
2337 }
2338 if(ev->IsTriggerClassFired("CINT7WU-B-NOPF-ALL"))
2339 {
2340 ntriggerbit=3;
2341 fQACollection->Fill("nTriggerBit",ntriggerbit);
2342 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT")) {
2343 ntriggerbit=18;
2344 fQACollection->Fill("nTriggerBit",ntriggerbit);
2345 }
2346 if(ev->IsTriggerClassFired("CINT7WUHQU-B-NOPF-CENT")) {
2347 ntriggerbit=19;
2348 fQACollection->Fill("nTriggerBit",ntriggerbit);
2349 }
2350 }
2351 if(ev->IsTriggerClassFired("CINT7WUHJT-B-NOPF-CENT"))
2352 {
2353 ntriggerbit=4;
2354 fQACollection->Fill("nTriggerBit",ntriggerbit);
2355
2356 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT")) {
2357 ntriggerbit=13;
2358 fQACollection->Fill("nTriggerBit",ntriggerbit);
2359 }
2360 if(ev->IsTriggerClassFired("CINT7WUHQU-B-NOPF-CENT")) {
2361 ntriggerbit=14;
2362 fQACollection->Fill("nTriggerBit",ntriggerbit);
2363 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT")) {
2364 ntriggerbit=17;
2365 fQACollection->Fill("nTriggerBit",ntriggerbit);
2366 }
2367 }
2368 }
2369 if(ev->IsTriggerClassFired("CINT7WUHQU-B-NOPF-CENT"))
2370 {
2371 ntriggerbit=5;
2372 fQACollection->Fill("nTriggerBit",ntriggerbit);
2373 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT")) {
2374 ntriggerbit=11;
2375 fQACollection->Fill("nTriggerBit",ntriggerbit);
2376 }
2377 if((!(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT")))&&(!(ev->IsTriggerClassFired("CINT7WUHJT-B-NOPF-CENT")))) {
2378 ntriggerbit=21;
2379 fQACollection->Fill("nTriggerBit",ntriggerbit);
2380
2381 /*
2382 Int_t nTrdTracks = ev->GetNumberOfTrdTracks();
2383 for (Int_t iTrack = 0; iTrack < nTrdTracks; ++iTrack) {
2384 AliESDTrdTrack* trdTrack = ev->GetTrdTrack(iTrack);
2385 printf("GTU track %3i: pt = %5.1f, PID = %3i\n", iTrack, trdTrack->Pt(), trdTrack->GetPID());
2386 }*/
2387
2388
2389 }
2390
2391 }
2392 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT"))
2393 {
2394 ntriggerbit=6;
2395 fQACollection->Fill("nTriggerBit",ntriggerbit);
2396 if(ev->IsTriggerClassFired("CINT7WUHQU-B-NOPF-CENT")) {
2397 ntriggerbit=12;
2398 fQACollection->Fill("nTriggerBit",ntriggerbit);
2399 }
2400 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-FAST")){
2401 ntriggerbit=15;
2402 fQACollection->Fill("nTriggerBit",ntriggerbit);
2403
2404 if((!(ev->IsTriggerClassFired("CINT7WUHQU-B-NOPF-CENT")))&&(!(ev->IsTriggerClassFired("CINT7WUHJT-B-NOPF-CENT")))) {
2405 ntriggerbit=20;
2406 fQACollection->Fill("nTriggerBit",ntriggerbit);
2407 /*
2408 Int_t nTrdTracks = ev->GetNumberOfTrdTracks();
2409 for (Int_t iTrack = 0; iTrack < nTrdTracks; ++iTrack) {
2410 AliESDTrdTrack* trdTrack = ev->GetTrdTrack(iTrack);
2411 printf("HSE GTU track %3i: pt = %5.1f, PID = %3i\n", iTrack, trdTrack->Pt(), trdTrack->GetPID());
2412 } */
2413
2414 }
2415
2416 }
2417
2418 }
2419 if(ev->IsTriggerClassFired("CEMC7WUHEE-B-NOPF-CENT")) {
2420 ntriggerbit=7;
2421 fQACollection->Fill("nTriggerBit",ntriggerbit);
2422 }
2423 if(ev->IsTriggerClassFired("CINT7WUHJT-B-NOPF-FAST")){
2424 ntriggerbit=8;
2425 fQACollection->Fill("nTriggerBit",ntriggerbit);
2426 }
2427 if(ev->IsTriggerClassFired("CINT7WUHQU-B-NOPF-FAST")){
2428 ntriggerbit=9;
2429 fQACollection->Fill("nTriggerBit",ntriggerbit);
2430 }
2431 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-FAST")){
2432 ntriggerbit=10;
2433 fQACollection->Fill("nTriggerBit",ntriggerbit);
2434 if(ev->IsTriggerClassFired("CINT7WUHSE-B-NOPF-CENT")) {
2435 ntriggerbit=16;
2436 fQACollection->Fill("nTriggerBit",ntriggerbit);
2437 }
2438 }
2439 if(ntriggerbit==0) fQACollection->Fill("nTriggerBit",1);
2440
01745870 2441}