1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 // study displaced electrons from beauty and charm
18 // with cut on impact parameters in various pT bins
22 // Hongyan Yang <hongyan@physi.uni-heidelberg.de>
23 // Carlo Bombonati <Carlo.Bombonati@cern.ch>
34 #include <TObjArray.h>
35 #include <TParticle.h>
39 #include "AliCFManager.h"
40 #include "AliMCEvent.h"
41 #include "AliMCEventHandler.h"
42 #include "AliMCParticle.h"
45 #include "AliESDEvent.h"
46 #include "AliESDInputHandler.h"
47 #include "AliESDtrack.h"
48 #include "AliESDpid.h"
50 #include "AliAnalysisManager.h"
52 #include "AliHFEpid.h"
53 #include "AliHFEcuts.h"
54 #include "AliHFEtools.h"
55 #include "AliHFEdisplacedElectrons.h"
56 #include "AliAnalysisTaskDisplacedElectrons.h"
59 //____________________________________________________________
60 AliAnalysisTaskDisplacedElectrons::AliAnalysisTaskDisplacedElectrons():
61 AliAnalysisTaskSE("Task for displaced electron study")
64 , fNminPrimVtxContrib(0)
72 , fDisplacedElectrons(0x0)
75 , fElectronsEsdPt(0x0)
76 , fElectronsDataPt(0x0)
79 , fHistDisplacedElectrons(0x0)
84 DefineInput(0, TChain::Class());
85 DefineOutput(1, TList::Class()); // displacedElectron
86 DefineOutput(2, TList::Class()); // correction framework
87 DefineOutput(3, TList::Class()); // QA
91 fDeDefaultPID = new AliESDpid;
92 fDePID = new AliHFEpid("DEPID");
97 //____________________________________________________________
98 AliAnalysisTaskDisplacedElectrons::AliAnalysisTaskDisplacedElectrons(const char * name):
99 AliAnalysisTaskSE(name)
102 , fNminPrimVtxContrib(0)
103 , fDePIDdetectors("")
110 , fDisplacedElectrons(0x0)
112 , fElectronsMcPt(0x0)
113 , fElectronsEsdPt(0x0)
114 , fElectronsDataPt(0x0)
117 , fHistDisplacedElectrons(0x0)
120 // Default constructor
122 DefineInput(0, TChain::Class());
123 DefineOutput(1, TList::Class());
124 DefineOutput(2, TList::Class());
125 DefineOutput(3, TList::Class());
128 fDeDefaultPID = new AliESDpid;
129 fDePID = new AliHFEpid("DEPID");
133 //____________________________________________________________
134 AliAnalysisTaskDisplacedElectrons::AliAnalysisTaskDisplacedElectrons(const AliAnalysisTaskDisplacedElectrons &ref):
135 AliAnalysisTaskSE(ref)
136 , fDeDebugLevel(ref.fDeDebugLevel)
137 , fNminITSCluster(ref.fNminITSCluster)
138 , fNminPrimVtxContrib(ref.fNminPrimVtxContrib)
139 , fDePIDdetectors(ref.fDePIDdetectors)
140 , fDePIDstrategy(ref.fDePIDstrategy)
141 , fDePlugins(ref.fDePlugins)
142 , fDeCuts(ref.fDeCuts)
143 , fDeDefaultPID(ref.fDeDefaultPID)
146 , fDisplacedElectrons(ref.fDisplacedElectrons)
147 , fDeNEvents(ref.fDeNEvents)
148 , fElectronsMcPt(ref.fElectronsMcPt)
149 , fElectronsEsdPt(ref.fElectronsEsdPt)
150 , fElectronsDataPt(ref.fElectronsDataPt)
151 , fDeCorrection(ref.fDeCorrection)
153 , fHistDisplacedElectrons(ref.fHistDisplacedElectrons)
160 //____________________________________________________________
161 AliAnalysisTaskDisplacedElectrons &AliAnalysisTaskDisplacedElectrons::operator=(const AliAnalysisTaskDisplacedElectrons &ref){
163 // Assignment operator
165 if(this == &ref) return *this;
166 AliAnalysisTask::operator=(ref);
167 fDeDebugLevel = ref.fDeDebugLevel;
168 fNminITSCluster = ref.fNminITSCluster;
169 fNminPrimVtxContrib = ref.fNminPrimVtxContrib;
170 fDePIDdetectors = ref.fDePIDdetectors;
171 fDePIDstrategy = ref.fDePIDstrategy;
172 fDePlugins = ref.fDePlugins;
173 fDeDefaultPID = ref.fDeDefaultPID;
175 fDeCuts = ref.fDeCuts;
177 fDisplacedElectrons = ref.fDisplacedElectrons;
178 fDeNEvents = ref.fDeNEvents;
179 fElectronsMcPt = ref.fElectronsMcPt;
180 fElectronsEsdPt = ref.fElectronsEsdPt;
181 fElectronsDataPt = ref.fElectronsDataPt;
182 fDeCorrection = ref.fDeCorrection;
184 fHistDisplacedElectrons = ref.fHistDisplacedElectrons;
189 //____________________________________________________________
190 AliAnalysisTaskDisplacedElectrons::~AliAnalysisTaskDisplacedElectrons(){
195 if(fDeDefaultPID) delete fDeDefaultPID;
196 if(fDePID) delete fDePID;
197 if(fDeCFM) delete fDeCFM;
198 if(fDisplacedElectrons) delete fDisplacedElectrons;
200 if(fDeNEvents) delete fDeNEvents;
201 if(fElectronsMcPt) delete fElectronsMcPt;
202 if(fElectronsEsdPt) delete fElectronsEsdPt;
203 if(fElectronsDataPt) delete fElectronsDataPt;
205 fDeCorrection->Clear();
206 delete fDeCorrection;
212 if(fHistDisplacedElectrons){
213 fHistDisplacedElectrons->Clear();
214 delete fHistDisplacedElectrons;
218 //____________________________________________________________
219 void AliAnalysisTaskDisplacedElectrons::UserCreateOutputObjects(){
220 // create output objects
222 // MC and Data containers
225 if(!fDeQA) fDeQA = new TList;
226 fDeQA->SetName("variousQAhistograms");
228 fDeNEvents = new TH1I("nDeEvents", "Number of Events in the DE Analysis", 2, 0, 2);
229 const Int_t nBins = 14;
230 const Float_t ptBins[nBins] = {0.0,0.5,1.0,1.5,2.0,2.5,3.0,4.0,5.0,7.0,9.0,12.0,16.0,20.0};
231 fElectronsMcPt = new TH1F("mcElectronPt", "MC: p_{T} distribution of identified electrons (mcpid);p_{T} (GeV/c);Counts;", nBins-1, ptBins);
232 fElectronsEsdPt = new TH1F("esdElectronPt", "ESD: p_{T} distribution of identified electrons (hfepid);p_{T} (GeV/c);Counts;", nBins-1, ptBins);
233 fElectronsDataPt = new TH1F("dataElectronPt", "DATA: p_{T} distribution of identified electrons (hfepid);p_{T} (GeV/c);Counts;", nBins-1, ptBins);
235 fDeQA->AddAt(fDeNEvents,0);
237 fDeQA->AddAt(fElectronsMcPt, 1);
238 fDeQA->AddAt(fElectronsEsdPt, 2);
241 fDeQA->AddAt(fElectronsDataPt, 1);
243 // Initialize correction Framework and Cuts
244 fDeCFM = new AliCFManager;
245 MakeEventContainer();
246 MakeParticleContainer();
248 if(!fDeCorrection) fDeCorrection = new TList();
249 fDeCorrection->SetName("deCorrections");
250 fDeCorrection->AddAt(fDeCFM->GetEventContainer(), 0);
251 fDeCorrection->AddAt(fDeCFM->GetParticleContainer(), 1);
252 fDeCorrection->Print();
254 for(Int_t istep = 0; istep < fDeCFM->GetEventContainer()->GetNStep(); istep++)
255 fDeCFM->SetEventCutsList(istep, 0x0);
256 for(Int_t istep = 0; istep < fDeCFM->GetParticleContainer()->GetNStep(); istep++)
257 fDeCFM->SetParticleCutsList(istep, 0x0);
260 AliWarning("Cuts not available. Default cuts will be used");
261 fDeCuts = new AliHFEcuts;
262 fDeCuts->CreateStandardCuts();
265 fDeCuts->Initialize(fDeCFM);
267 fDePID->SetHasMCData(HasMCData());
268 fDePID->AddDetector("TPC", 0);
269 fDePID->AddDetector("TOF", 1);
270 fDePID->ConfigureTPCrejection();
271 fDePID->InitializePID(); // Only restrictions to TPC allowed
274 // displaced electron study----------------------------------
275 if(GetPlugin(kDisplacedElectrons)){
277 fDisplacedElectrons = new AliHFEdisplacedElectrons;
278 fDisplacedElectrons->SetDebugLevel(fDeDebugLevel);
279 fDisplacedElectrons->SetHasMCData(HasMCData());
280 fDisplacedElectrons->SetMinPrimVtxContrib(fNminPrimVtxContrib);
281 fDisplacedElectrons->SetNitsCluster(fNminITSCluster);
283 if(!fHistDisplacedElectrons) fHistDisplacedElectrons = new TList();
284 fDisplacedElectrons->CreateOutputs(fHistDisplacedElectrons);
291 //____________________________________________________________
292 void AliAnalysisTaskDisplacedElectrons::UserExec(Option_t *){
297 if(fDeDebugLevel>=10) AliInfo("analyse single event");
300 AliError("Reconstructed Event not available");
305 AliESDInputHandler *inH = dynamic_cast<AliESDInputHandler *>(fInputHandler);
308 AliError("AliESDInputHandler dynamic cast failed!");
312 // from now on, only ESD are analyzed
313 // using HFE pid, using HFE cuts
315 AliESDpid *workingPID = inH->GetESDpid();
317 AliDebug(1, "Using ESD PID from the input handler");
318 //printf("\n ESD PID\n");
319 fDePID->SetESDpid(workingPID);
321 AliDebug(1, "Using default ESD PID");
322 //printf(" DEFAULT PID!\n\n");
323 fDePID->SetESDpid(AliHFEtools::GetDefaultPID(HasMCData()));
327 AliError("HFE cuts not available");
331 // ---- CHOICE OF ANALYSIS ----
333 // Protect against missing MC trees
334 AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
336 AliError("MC handler cuts not available");
340 if(!mcH->InitOk()) return;
341 if(!mcH->TreeK()) return;
342 if(!mcH->TreeTR()) return;
344 AliDebug(4, Form("MC Event: %p", fMCEvent));
346 AliError("No MC Event, but MC Data required");
350 ProcessMC(); // PURE MC
352 if(IsESDanalysis()) ProcessESD(); // ESD WITH MC
354 }else if(IsESDanalysis()) ProcessData(); // PURE ESD
358 PostData(1, fHistDisplacedElectrons);
359 PostData(2, fDeCorrection);
364 //____________________________________________________________
365 void AliAnalysisTaskDisplacedElectrons::ProcessMC(){
367 // handle pure MC analysis
370 Int_t nMCelectrons = 0;
371 AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
373 AliError("ESD event not available");
377 Double_t mcContainer[4]; // container for the output in THnSparse
378 memset(mcContainer, 0, sizeof(Double_t) * 4);
380 fDeCFM->SetMCEventInfo(fMCEvent);
382 Double_t nContributor[1] = {0};
383 const AliVVertex *mcPrimVtx = fMCEvent->GetPrimaryVertex();
384 if(mcPrimVtx) nContributor[0] = mcPrimVtx->GetNContributors();
387 // cut at MC event level
390 if(!fDeCFM->CheckEventCuts(AliHFEcuts::kEventStepGenerated, fMCEvent)) return;
391 if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(nContributor,AliHFEcuts::kEventStepGenerated);
393 AliStack *stack = 0x0;
395 if(!fMCEvent->Stack())return;
396 stack = fMCEvent->Stack();
397 Int_t nTracks = stack->GetNtrack();
399 AliMCParticle *mcTrack = 0x0;
401 for(Int_t itrack = 0; itrack<nTracks; itrack++){
402 if(!(stack->Particle(itrack))) continue;
403 if(mcTrack) mcTrack = 0x0;
404 mcTrack = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(itrack));
405 if(!mcTrack) continue;
406 //TParticle *mcPart = stack->Particle(itrack);
408 mcContainer[0] = mcTrack->Pt();
409 mcContainer[1] = mcTrack->Eta();
410 mcContainer[2] = mcTrack->Phi();
411 mcContainer[3] = mcTrack->Charge();
415 if (!stack->IsPhysicalPrimary(mcTrack->GetLabel())) continue;
416 // no cut but require primary
417 if(GetPlugin(kCorrection))fDeCFM->GetParticleContainer()->Fill(mcContainer, 0);
419 // all pions for reference
420 if(TMath::Abs(mcTrack->Particle()->GetPdgCode())==AliHFEdisplacedElectrons::kPDGpion && GetPlugin(kCorrection))
421 fDeCFM->GetParticleContainer()->Fill(mcContainer, 1);
423 // cut for signal: all MC electrons
424 if(TMath::Abs(mcTrack->Particle()->GetPdgCode())==AliHFEdisplacedElectrons::kPDGelectron && GetPlugin(kCorrection))
425 fDeCFM->GetParticleContainer()->Fill(mcContainer, 2);
427 // cut at track level kinematics: pt and eta
428 if(TMath::Abs(mcContainer[1])>=0.8 || mcContainer[0]>20 || mcContainer[0]<0.1) continue;
430 if(TMath::Abs(mcTrack->Particle()->GetPdgCode())==AliHFEdisplacedElectrons::kPDGelectron){
432 fElectronsMcPt->Fill(mcContainer[0]);
435 if(GetPlugin(kCorrection))fDeCFM->GetParticleContainer()->Fill(mcContainer, 3);
439 fDisplacedElectrons->FillMcOutput(fESD, fMCEvent, mcTrack);
443 if(fDeDebugLevel>=10) printf("there are %d electrons in this MC event", nMCelectrons);
447 //____________________________________________________________
448 void AliAnalysisTaskDisplacedElectrons::ProcessESD(){
450 // this is to handle ESD tracks with MC information
451 // MC pid is only used when HFE pid is implemented, for comparison
452 // corrections are taken into account
454 // process data: ESD tracks with MC information
456 const Int_t kStepPID = AliHFEcuts::kStepHFEcutsTRD + 1;
458 Double_t esdContainer[4]; // container for the output in THnSparse
459 memset(esdContainer, 0, sizeof(Double_t) * 4);
460 AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
462 AliError("No ESD event available");
466 fDeCFM->SetRecEventInfo(fESD);
467 Double_t nContrib[1] = {fESD->GetPrimaryVertex()->GetNContributors()};
469 Bool_t alreadyseen = kFALSE;
470 AliLabelContainer cont(fESD->GetNumberOfTracks());
472 Int_t nHFEelectrons = 0;
473 AliESDtrack *track = 0x0;
474 AliStack *stack = 0x0;
476 if(!(stack = fMCEvent->Stack())) return;
479 // cut at ESD event level
481 if(!fDeCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
482 if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(nContrib, AliHFEcuts::kEventStepReconstructed);
484 for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
485 track = fESD->GetTrack(itrack);
487 if(GetPlugin(kDisplacedElectrons)) {
489 esdContainer[0] = track->Pt();
490 esdContainer[1] = track->Eta();
491 esdContainer[2] = track->Phi();
492 esdContainer[3] = track->Charge();
495 alreadyseen = cont.Find(TMath::Abs(track->GetLabel()));
496 cont.Append(TMath::Abs(track->GetLabel())); // check double counting
497 if(alreadyseen) continue; // avoid double counting
498 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(esdContainer, 1+AliHFEcuts::kStepRecNoCut + AliHFEcuts::kNcutStepsMCTrack);
501 // RecKine: ITSTPC cuts : ITS & TPC refit, covmatrix: (2, 2, 0.5, 0.5, 2); min_tpccls: 50, chi2_tpccls: 3.5
502 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
503 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(esdContainer, 1+AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack);
506 // RecPrim: cut on track quality : DCA to vertex max: 3cm and 10cm; reject kink daughters
507 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
508 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(esdContainer, 1+AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack);
511 // HFEcuts: ITS layers cuts: ITS pixel layer: kFirst, kSecond, kBoth, kNone or kAny
512 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
513 if(GetPlugin(kCorrection))fDeCFM->GetParticleContainer()->Fill(esdContainer, 1+AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack);
517 // TRD: number of tracklets in TRD
518 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
519 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(esdContainer, 1+AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack);
523 // track accepted, do PID
524 // --> only electron candidate will be processed
525 AliHFEpidObject hfetrack;
526 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
527 hfetrack.SetRecTrack(track);
528 //if(HasMCData())hfetrack.SetMCTrack(mctrack);
530 if(!fDePID->IsSelected(&hfetrack)) continue;
531 else if(fDeDebugLevel>=10)
532 AliInfo("ESD info: this particle is identified as electron by HFEpid method \n");
533 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(esdContainer, 1+kStepPID + AliHFEcuts::kNcutStepsMCTrack);
537 fElectronsEsdPt->Fill(esdContainer[0]);
538 fDisplacedElectrons->FillEsdOutput(fESD, track, stack);
540 } // displaced electron analysis on ESD with MC plugin
543 if(fDeDebugLevel>=10) printf("there are %d HFE electrons in this ESD event", nHFEelectrons);
549 //____________________________________________________________
550 void AliAnalysisTaskDisplacedElectrons::ProcessData(){
552 // this is a track loop over real data
553 // no MC information at all
556 const Int_t kNcutStepsESDtrack = AliHFEcuts::kNcutStepsRecTrack + 1;
557 //const Int_t kNcutStepsTrack = AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack;
558 const Int_t kStepPID = AliHFEcuts::kStepHFEcutsTRD + 1;
560 AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
562 AliError("No ESD event available");
566 Double_t dataContainer[4]; // container for the output in THnSparse
567 memset(dataContainer, 0, sizeof(Double_t) * 4);
569 Bool_t alreadyseen = kFALSE;
570 AliLabelContainer cont(fESD->GetNumberOfTracks());
573 AliESDtrack *track = 0x0;
574 Int_t nHFEelectrons= 0;
576 fDeCFM->SetRecEventInfo(fESD);
577 Double_t nContrib[1] = {fESD->GetPrimaryVertex()->GetNContributors()};
578 if(!fDeCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
579 if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(nContrib, AliHFEcuts::kEventStepReconstructed);
582 for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
583 track = fESD->GetTrack(itrack);
585 if(GetPlugin(kDisplacedElectrons)) {
587 dataContainer[0] = track->Pt();
588 dataContainer[1] = track->Eta();
589 dataContainer[2] = track->Phi();
590 dataContainer[3] = track->Charge();
592 alreadyseen = cont.Find(TMath::Abs(track->GetLabel())); // double counted track
593 cont.Append(TMath::Abs(track->GetLabel()));
594 if(alreadyseen) continue; // avoid double counting
595 if(GetPlugin(kCorrection))fDeCFM->GetParticleContainer()->Fill(&dataContainer[4],
596 1+AliHFEcuts::kStepRecNoCut + AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack);
599 // RecKine: ITSTPC cuts : ITS & TPC refit, covmatrix: (2, 2, 0.5, 0.5, 2); min_tpccls: 50, chi2_tpccls: 3.5
600 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
601 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(&dataContainer[4],
602 1+AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack);
605 // RecPrim: cut on track quality : DCA to vertex max: 3cm and 10cm; reject kink daughters
606 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
607 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(&dataContainer[4],
608 1+AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack);
611 // HFEcuts: ITS layers cuts: ITS pixel layer: kFirst, kSecond, kBoth, kNone or kAny
612 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS, track)) continue;
613 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(&dataContainer[4],
614 1+AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack);
618 // TRD: number of tracklets in TRD0
619 if(!fDeCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
620 if(GetPlugin(kCorrection)) if(HasMCData())fDeCFM->GetParticleContainer()->Fill(&dataContainer[4],
621 1+AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack);
626 // track accepted, do PID --> only electron candidate will be processed
628 AliHFEpidObject hfetrack;
629 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
630 hfetrack.SetRecTrack(track);
632 if(!fDePID->IsSelected(&hfetrack)) continue;
633 else if(fDeDebugLevel>=10)
634 AliInfo("ESD info: this particle is identified as electron by HFEpid method \n");
635 if(GetPlugin(kCorrection)) fDeCFM->GetParticleContainer()->Fill(dataContainer, 1+kStepPID + AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack);
638 fElectronsDataPt->Fill(dataContainer[0]);
639 fDisplacedElectrons->FillDataOutput(fESD, track);
640 } // analyze displaced electrons plugin switched on
643 if(fDeDebugLevel>=10) printf("there are %d HFE electrons in this DATA event", nHFEelectrons);
647 //____________________________________________________________
648 void AliAnalysisTaskDisplacedElectrons::Terminate(Option_t *){
650 // Terminate not implemented at the moment
653 fHistDisplacedElectrons = dynamic_cast<TList *>(GetOutputData(1));
654 fDeCorrection = dynamic_cast<TList *>(GetOutputData(2));
655 fDeQA = dynamic_cast<TList *>(GetOutputData(3));
656 if(!fDeCorrection) AliError("correction list not available\n");
657 if(!fHistDisplacedElectrons) AliError("de list not available\n");
658 if(!fDeQA) AliError("qa list is not available\n");
660 fHistDisplacedElectrons->Print();
661 fDeCorrection->Print();
664 AliInfo("analysis done!\n");
668 //____________________________________________________________
669 void AliAnalysisTaskDisplacedElectrons::PrintStatus() const {
672 // Print Analysis status
675 printf("\t Analysis Settings\n\t========================================\n");
676 printf("\t running over %s\n", HasMCData()?"MC data":"pp collision data");
677 printf("\t displaced electrons' analysis is %s\n", GetPlugin(kDisplacedElectrons)?"ON":"OFF");
678 printf("\t correction container is %s\n", GetPlugin(kCorrection)?"ON":"OFF");
679 printf("\t hfe pid qa is %s\n", GetPlugin(kDePidQA)?"ON":"OFF");
680 printf("\t post processing is %s\n", GetPlugin(kPostProcess)?"ON":"OFF");
681 printf("\t cuts: %s\n", (fDeCuts != NULL) ? "YES" : "NO");
686 //__________________________________________
687 void AliAnalysisTaskDisplacedElectrons::SwitchOnPlugin(Int_t plug){
691 // - analyze impact parameter
696 case kDisplacedElectrons:
697 SETBIT(fDePlugins, plug);
700 SETBIT(fDePlugins, plug);
703 SETBIT(fDePlugins, plug);
706 SETBIT(fDePlugins, plug);
709 AliError("Unknown Plugin");
713 //____________________________________________________________
714 void AliAnalysisTaskDisplacedElectrons::MakeParticleContainer(){
716 // Create the particle container for the correction framework manager and
719 const Int_t kNcutStepsESDtrack = AliHFEcuts::kNcutStepsRecTrack + 1;
720 const Int_t kNcutStepsTrack = AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack;
722 const Int_t kNvar = 4;
723 //number of variables on the grid:pt,eta, phi, charge
724 const Double_t kPtbound[2] = {0.1, 10.};
725 const Double_t kEtabound[2] = {-0.8, 0.8};
726 const Double_t kPhibound[2] = {0., 2. * TMath::Pi()};
728 //arrays for the number of bins in each dimension
730 iBin[0] = 40; // bins in pt
731 iBin[1] = 8; // bins in eta
732 iBin[2] = 18; // bins in phi
733 iBin[3] = 2; // bins in charge
735 //arrays for lower bounds :
736 Double_t* binEdges[kNvar];
737 binEdges[0] = AliHFEtools::MakeLogarithmicBinning(iBin[0], kPtbound[0], kPtbound[1]);
738 binEdges[1] = AliHFEtools::MakeLinearBinning(iBin[1], kEtabound[0], kEtabound[1]);
739 binEdges[2] = AliHFEtools::MakeLinearBinning(iBin[2], kPhibound[0], kPhibound[1]);
740 binEdges[3] = AliHFEtools::MakeLinearBinning(iBin[3], -1.1, 1.1); // Numeric precision
742 //------------------------------------------------
743 // one "container" for MC+ESD+Data
744 //----------pure MC track-------------------------
746 // 1: MC pion total ---- be careful!!!!
747 // 2: MC electrons total
748 // 3: MC electrons in acceptance
749 //-------ESD track with MC info-------------------
750 // 4: ESD track with MC: no cut
751 // 5: ESD track with MC: cut on kine its tpc
752 // 6: ESD track with MC: rec prim
753 // 7: ESD track with MC: hfe cuts its
754 // 8: ESD track with MC: hfe cuts trd
755 // 9: ESD track with MC: hfe pid
756 //-----------data track---------------------------
757 // 10: DATA track wo MC: no cut
758 // 11: DATA track wo MC: cut on kine its tpc
759 // 12: DATA track wo MC: rec prim
760 // 13: DATA track wo MC: hfe cuts its
761 // 14: DATA track wo MC: hfe cuts trd
762 // 15: DATA track wo MC: hfe pid
763 //------------------------------------------------
765 AliCFContainer* container = new AliCFContainer("deTrackContainer", "Container for tracks",
766 (1 + kNcutStepsTrack + kNcutStepsESDtrack), kNvar, iBin);
768 //setting the bin limits
769 for(Int_t ivar = 0; ivar < kNvar; ivar++){
770 container -> SetBinLimits(ivar, binEdges[ivar]);
772 fDeCFM->SetParticleContainer(container);
776 //____________________________________________________________
777 void AliAnalysisTaskDisplacedElectrons::MakeEventContainer(){
779 // Create the event container for the correction framework and link it
786 const Int_t kNvar = 1; // number of variables on the grid: number of tracks per event
787 const Double_t kNTrackBound[2] = {-0.5, 200.5};
788 const Int_t kNBins = 201;
790 AliCFContainer *evCont = new AliCFContainer("deEventContainer", "Container for DE events", AliHFEcuts::kNcutStepsEvent, kNvar, &kNBins);
792 Double_t trackBins[kNBins];
793 for(Int_t ibin = 0; ibin < kNBins; ibin++) trackBins[ibin] = kNTrackBound[0] + static_cast<Double_t>(ibin);
794 evCont->SetBinLimits(0,trackBins);
796 fDeCFM->SetEventContainer(evCont);
802 //__________________________________________________________
803 void AliAnalysisTaskDisplacedElectrons::AddPIDdetector(TString detector){
805 // Adding PID detector to the task
807 if(!fDePIDdetectors.Length())
808 fDePIDdetectors = detector;
810 fDePIDdetectors += ":" + detector;
815 //____________________________________________________________
816 AliAnalysisTaskDisplacedElectrons::AliLabelContainer::AliLabelContainer(Int_t capacity):
824 // Default constructor
826 fContainer = new Int_t[capacity];
827 fBegin = &fContainer[0];
828 fEnd = &fContainer[capacity - 1];
829 fLast = fCurrent = fBegin;
832 //____________________________________________________________
833 Bool_t AliAnalysisTaskDisplacedElectrons::AliLabelContainer::Append(Int_t label){
835 // Add Label to the container
837 if(fLast > fEnd) return kFALSE;
842 //____________________________________________________________
843 Bool_t AliAnalysisTaskDisplacedElectrons::AliLabelContainer::Find(Int_t label) const {
845 // Find track in the list of labels
847 for(Int_t *entry = fBegin; entry <= fLast; entry++)
848 if(*entry == label) return kTRUE;
852 //____________________________________________________________
853 Int_t AliAnalysisTaskDisplacedElectrons::AliLabelContainer::Next() {
857 if(fCurrent > fLast) return -1;