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