1 /*************************************************************************
2 * Copyright(c) 1998-2009, 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 **************************************************************************/
16 ///////////////////////////////////////////////////////////////////////////
17 // Dielectron Analysis Main class //
20 Framework to perform event selectoin, single track selection and track pair
23 Convention for the signs of the pair in fPairCandidates:
24 The names are available via the function PairClassName(Int_t i)
26 0: ev1+ ev1+ (same event like sign +)
27 1: ev1+ ev1- (same event unlike sign)
28 2: ev1- ev1- (same event like sign -)
30 3: ev1+ ev2+ (mixed event like sign +)
31 4: ev1- ev2+ (mixed event unlike sign -+)
32 6: ev1+ ev2- (mixed event unlike sign +-)
33 7: ev1- ev2- (mixed event like sign -)
35 5: ev2+ ev2+ (same event like sign +)
36 8: ev2+ ev2- (same event unlike sign)
37 9: ev2- ev2- (same event like sign -)
39 10: ev1+ ev1- (same event track rotation)
43 ///////////////////////////////////////////////////////////////////////////
51 #include <AliKFParticle.h>
53 #include <AliESDInputHandler.h>
54 #include <AliAnalysisManager.h>
55 #include <AliEPSelectionTask.h>
56 #include <AliEventplane.h>
57 #include <AliVEvent.h>
58 #include <AliVParticle.h>
59 #include <AliVTrack.h>
60 #include "AliDielectronPair.h"
61 #include "AliDielectronHistos.h"
62 #include "AliDielectronCF.h"
63 #include "AliDielectronMC.h"
64 #include "AliDielectronVarManager.h"
65 #include "AliDielectronTrackRotator.h"
66 #include "AliDielectronDebugTree.h"
67 #include "AliDielectronSignalMC.h"
68 #include "AliDielectronMixingHandler.h"
69 #include "AliDielectronV0Cuts.h"
70 #include "AliDielectronPID.h"
72 #include "AliDielectron.h"
74 ClassImp(AliDielectron)
76 const char* AliDielectron::fgkTrackClassNames[4] = {
83 const char* AliDielectron::fgkPairClassNames[11] = {
97 //________________________________________________________________
98 AliDielectron::AliDielectron() :
99 TNamed("AliDielectron","AliDielectron"),
102 fPostPIDCntrdCorr(0x0),
103 fPostPIDWdthCorr(0x0),
106 fEventFilter("EventFilter"),
107 fTrackFilter("TrackFilter"),
108 fPairPreFilter("PairPreFilter"),
109 fPairPreFilterLegs("PairPreFilterLegs"),
110 fPairFilter("PairFilter"),
111 fEventPlanePreFilter("EventPlanePreFilter"),
112 fEventPlanePOIPreFilter("EventPlanePOIPreFilter"),
122 fPairCandidates(new TObjArray(11)),
127 fPreFilterEventPlane(kFALSE),
128 fLikeSignSubEvents(kFALSE),
129 fPreFilterUnlikeOnly(kFALSE),
130 fPreFilterAllSigns(kFALSE),
132 fStoreRotatedPairs(kFALSE),
133 fDontClearArrays(kFALSE),
134 fEstimatorFilename(""),
135 fTRDpidCorrectionFilename(""),
136 fVZEROCalibrationFilename(""),
137 fVZERORecenteringFilename(""),
138 fZDCRecenteringFilename("")
142 // Default constructor
147 //________________________________________________________________
148 AliDielectron::AliDielectron(const char* name, const char* title) :
152 fPostPIDCntrdCorr(0x0),
153 fPostPIDWdthCorr(0x0),
156 fEventFilter("EventFilter"),
157 fTrackFilter("TrackFilter"),
158 fPairPreFilter("PairPreFilter"),
159 fPairPreFilterLegs("PairPreFilterLegs"),
160 fPairFilter("PairFilter"),
161 fEventPlanePreFilter("EventPlanePreFilter"),
162 fEventPlanePOIPreFilter("EventPlanePOIPreFilter"),
172 fPairCandidates(new TObjArray(11)),
177 fPreFilterEventPlane(kFALSE),
178 fLikeSignSubEvents(kFALSE),
179 fPreFilterUnlikeOnly(kFALSE),
180 fPreFilterAllSigns(kFALSE),
182 fStoreRotatedPairs(kFALSE),
183 fDontClearArrays(kFALSE),
184 fEstimatorFilename(""),
185 fTRDpidCorrectionFilename(""),
186 fVZEROCalibrationFilename(""),
187 fVZERORecenteringFilename(""),
188 fZDCRecenteringFilename("")
196 //________________________________________________________________
197 AliDielectron::~AliDielectron()
200 // Default destructor
202 if (fQAmonitor) delete fQAmonitor;
203 if (fPostPIDCntrdCorr) delete fPostPIDCntrdCorr;
204 if (fPostPIDWdthCorr) delete fPostPIDWdthCorr;
205 if (fLegEffMap) delete fLegEffMap;
206 if (fPairEffMap) delete fPairEffMap;
207 if (fHistos) delete fHistos;
208 if (fPairCandidates) delete fPairCandidates;
209 if (fDebugTree) delete fDebugTree;
210 if (fMixing) delete fMixing;
211 if (fSignalsMC) delete fSignalsMC;
212 if (fCfManagerPair) delete fCfManagerPair;
213 if (fHistoArray) delete fHistoArray;
216 //________________________________________________________________
217 void AliDielectron::Init()
220 // Initialise objects
223 if(GetHasMC()) AliDielectronMC::Instance()->SetHasMC(GetHasMC());
225 InitPairCandidateArrays();
227 if (fCfManagerPair) {
228 fCfManagerPair->SetSignalsMC(fSignalsMC);
229 fCfManagerPair->InitialiseContainer(fPairFilter);
232 fTrackRotator->SetTrackArrays(&fTracks[0],&fTracks[1]);
233 fTrackRotator->SetPdgLegs(fPdgLeg1,fPdgLeg2);
235 if (fDebugTree) fDebugTree->SetDielectron(this);
237 if(fEstimatorFilename.Contains(".root")) AliDielectronVarManager::InitEstimatorAvg(fEstimatorFilename.Data());
238 if(fTRDpidCorrectionFilename.Contains(".root")) AliDielectronVarManager::InitTRDpidEffHistograms(fTRDpidCorrectionFilename.Data());
239 if(fVZEROCalibrationFilename.Contains(".root")) AliDielectronVarManager::SetVZEROCalibrationFile(fVZEROCalibrationFilename.Data());
240 if(fVZERORecenteringFilename.Contains(".root")) AliDielectronVarManager::SetVZERORecenteringFile(fVZERORecenteringFilename.Data());
241 if(fZDCRecenteringFilename.Contains(".root")) AliDielectronVarManager::SetZDCRecenteringFile(fZDCRecenteringFilename.Data());
243 if(fLegEffMap) AliDielectronVarManager::SetLegEffMap(fLegEffMap);
244 if(fPairEffMap) AliDielectronVarManager::SetPairEffMap(fPairEffMap);
247 if (fMixing) fMixing->Init(this);
249 fHistoArray->SetSignalsMC(fSignalsMC);
253 if(fPostPIDCntrdCorr) AliDielectronPID::SetCentroidCorrFunction(fPostPIDCntrdCorr);
254 if(fPostPIDWdthCorr) AliDielectronPID::SetWidthCorrFunction(fPostPIDWdthCorr);
257 fQAmonitor = new AliDielectronCutQA(Form("QAcuts_%s",GetName()),"QAcuts");
258 fQAmonitor->AddTrackFilter(&fTrackFilter);
259 fQAmonitor->AddPairFilter(&fPairFilter);
260 fQAmonitor->AddEventFilter(&fEventFilter);
265 //________________________________________________________________
266 void AliDielectron::Process(AliVEvent *ev1, AliVEvent *ev2)
269 // Process the events
272 //at least first event is needed!
274 AliError("At least first event must be set!");
278 // modify event numbers in MC so that we can identify new events
279 // in AliDielectronV0Cuts (not neeeded for collision data)
281 ev1->SetBunchCrossNumber(1);
282 ev1->SetOrbitNumber(1);
283 ev1->SetPeriodNumber(1);
287 AliDielectronVarManager::SetEvent(ev1);
289 //set mixing bin to event data
290 Int_t bin=fMixing->FindBin(AliDielectronVarManager::GetData());
291 AliDielectronVarManager::SetValue(AliDielectronVarManager::kMixingBin,bin);
294 //in case we have MC load the MC event and process the MC particles
295 // why do not apply the event cuts first ????
296 if (AliDielectronMC::Instance()->ConnectMCEvent()){
300 //if candidate array doesn't exist, create it
301 if (!fPairCandidates->UncheckedAt(0)) {
302 InitPairCandidateArrays();
307 //mask used to require that all cuts are fulfilled
308 UInt_t selectedMask=(1<<fEventFilter.GetCuts()->GetEntries())-1;
311 UInt_t cutmask = fEventFilter.IsSelected(ev1);
312 if(fCutQA) fQAmonitor->FillAll(ev1);
313 if(fCutQA) fQAmonitor->Fill(cutmask,ev1);
314 if ((ev1&&cutmask!=selectedMask) ||
315 (ev2&&fEventFilter.IsSelected(ev2)!=selectedMask)) return;
317 //fill track arrays for the first event
319 FillTrackArrays(ev1);
320 if (((fPreFilterAllSigns)||(fPreFilterUnlikeOnly)) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(0, 1, fTracks[0], fTracks[1]);
324 //fill track arrays for the second event
326 FillTrackArrays(ev2,1);
327 if (((fPreFilterAllSigns)||(fPreFilterUnlikeOnly)) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(2, 3, fTracks[2], fTracks[3]);
330 // TPC event plane correction
331 if (ev1 && fPreFilterEventPlane && ( fEventPlanePreFilter.GetCuts()->GetEntries()>0 || fEventPlanePOIPreFilter.GetCuts()->GetEntries()>0))
332 EventPlanePreFilter(0, 1, fTracks[0], fTracks[1], ev1);
335 // create pairs and fill pair candidate arrays
336 for (Int_t itrackArr1=0; itrackArr1<4; ++itrackArr1){
337 for (Int_t itrackArr2=itrackArr1; itrackArr2<4; ++itrackArr2){
338 if(!fProcessLS && GetPairIndex(itrackArr1,itrackArr2)!=kEv1PM) continue;
339 FillPairArrays(itrackArr1, itrackArr2);
345 fTrackRotator->SetEvent(ev1);
350 //fill debug tree if a manager is attached
351 if (fDebugTree) FillDebugTree();
353 //process event mixing
355 fMixing->Fill(ev1,this);
356 // FillHistograms(0x0,kTRUE);
359 // fill candidate variables
360 Double_t ntracks = fTracks[0].GetEntriesFast() + fTracks[1].GetEntriesFast();
361 Double_t npairs = PairArray(AliDielectron::kEv1PM)->GetEntriesFast();
362 AliDielectronVarManager::SetValue(AliDielectronVarManager::kTracks, ntracks);
363 AliDielectronVarManager::SetValue(AliDielectronVarManager::kPairs, npairs);
365 //in case there is a histogram manager, fill the QA histograms
366 if (fHistos && fSignalsMC) FillMCHistograms(ev1);
367 if (fHistos) FillHistograms(ev1);
371 if (!fDontClearArrays) ClearArrays();
373 // reset TPC EP and unique identifiers for v0 cut class
374 AliDielectronVarManager::SetTPCEventPlane(0x0);
375 if(GetHasMC()) { // only for MC needed
376 for (Int_t iCut=0; iCut<fTrackFilter.GetCuts()->GetEntries();++iCut) {
377 if ( fTrackFilter.GetCuts()->At(iCut)->IsA() == AliDielectronV0Cuts::Class() )
378 ((AliDielectronV0Cuts*)fTrackFilter.GetCuts()->At(iCut))->ResetUniqueEventNumbers();
384 //________________________________________________________________
385 void AliDielectron::ProcessMC(AliVEvent *ev1)
388 // Process the MC data
391 AliDielectronMC *dieMC=AliDielectronMC::Instance();
393 if (fHistos) FillHistogramsMC(dieMC->GetMCEvent(), ev1);
396 if(!dieMC->GetNMCTracks()) return;
398 // signals to be studied
399 if(!fSignalsMC) return;
400 Int_t nSignals = fSignalsMC->GetEntries();
401 if(!nSignals) return;
403 //loop over all MC data and Fill the HF, CF containers and histograms if they exist
404 if(fCfManagerPair) fCfManagerPair->SetPdgMother(fPdgMother);
406 Bool_t bFillCF = (fCfManagerPair ? fCfManagerPair->GetStepForMCtruth() : kFALSE);
407 Bool_t bFillHF = (fHistoArray ? fHistoArray->GetStepForMCGenerated() : kFALSE);
408 Bool_t bFillHist = kFALSE;
410 const THashList *histlist = fHistos->GetHistogramList();
411 for(Int_t isig=0;isig<nSignals;isig++) {
412 TString sigName = fSignalsMC->At(isig)->GetName();
413 bFillHist |= histlist->FindObject(Form("Pair_%s_MCtruth",sigName.Data()))!=0x0;
414 bFillHist |= histlist->FindObject(Form("Track_Leg_%s_MCtruth",sigName.Data()))!=0x0;
415 bFillHist |= histlist->FindObject(Form("Track_%s_%s_MCtruth",fgkPairClassNames[1],sigName.Data()))!=0x0;
419 // check if there is anything to fill
420 if(!bFillCF && !bFillHF && !bFillHist) return;
423 // initialize 2D arrays of labels for particles from each MC signal
424 Int_t** labels1; // labels for particles satisfying branch 1
425 Int_t** labels2; // labels for particles satisfying branch 2
426 Int_t** labels12; // labels for particles satisfying both branches
427 labels1 = new Int_t*[nSignals];
428 labels2 = new Int_t*[nSignals];
429 labels12 = new Int_t*[nSignals];
430 Int_t* indexes1=new Int_t[nSignals];
431 Int_t* indexes2=new Int_t[nSignals];
432 Int_t* indexes12=new Int_t[nSignals];
433 for(Int_t isig=0;isig<nSignals;++isig) {
434 *(labels1+isig) = new Int_t[dieMC->GetNMCTracks()];
435 *(labels2+isig) = new Int_t[dieMC->GetNMCTracks()];
436 *(labels12+isig) = new Int_t[dieMC->GetNMCTracks()];
437 for(Int_t ip=0; ip<dieMC->GetNMCTracks();++ip) {
438 labels1[isig][ip] = -1;
439 labels2[isig][ip] = -1;
440 labels12[isig][ip] = -1;
447 Bool_t truth1=kFALSE;
448 Bool_t truth2=kFALSE;
449 // loop over the MC tracks
450 for(Int_t ipart=0; ipart<dieMC->GetNMCTracks(); ++ipart) {
451 for(Int_t isig=0; isig<nSignals; ++isig) { // loop over signals
452 // Proceed only if this signal is required in the pure MC step
453 // NOTE: Some signals can be satisfied by many particles and this leads to high
454 // computation times (e.g. secondary electrons from the GEANT transport). Be aware of this!!
455 if(!((AliDielectronSignalMC*)fSignalsMC->At(isig))->GetFillPureMCStep()) continue;
457 truth1 = dieMC->IsMCTruth(ipart, (AliDielectronSignalMC*)fSignalsMC->At(isig), 1);
458 truth2 = dieMC->IsMCTruth(ipart, (AliDielectronSignalMC*)fSignalsMC->At(isig), 2);
460 // particles satisfying both branches are treated separately to avoid double counting during pairing
461 if(truth1 && truth2) {
462 labels12[isig][indexes12[isig]] = ipart;
467 labels1[isig][indexes1[isig]] = ipart;
471 labels2[isig][indexes2[isig]] = ipart;
476 } // end loop over MC particles
478 // Do the pairing and fill the CF container with pure MC info
479 for(Int_t isig=0; isig<nSignals; ++isig) {
480 // printf("INDEXES: %d-%d both%d\n",indexes1[isig],indexes2[isig],indexes12[isig]);
481 // mix the particles which satisfy only one of the signal branches
482 for(Int_t i1=0;i1<indexes1[isig];++i1) {
483 if(!indexes2[isig]) FillMCHistograms(labels1[isig][i1], -1, isig); // (e.g. single electrons only, no pairs)
484 for(Int_t i2=0;i2<indexes2[isig];++i2) {
485 if(bFillCF) fCfManagerPair->FillMC(labels1[isig][i1], labels2[isig][i2], isig);
486 if(bFillHF) fHistoArray->Fill(labels1[isig][i1], labels2[isig][i2], isig);
487 FillMCHistograms(labels1[isig][i1], labels2[isig][i2], isig);
490 // mix the particles which satisfy both branches
491 for(Int_t i1=0;i1<indexes12[isig];++i1) {
492 for(Int_t i2=0; i2<i1; ++i2) {
493 if(bFillCF) fCfManagerPair->FillMC(labels12[isig][i1], labels12[isig][i2], isig);
494 if(bFillHF) fHistoArray->Fill(labels12[isig][i1], labels12[isig][i2], isig);
495 FillMCHistograms(labels12[isig][i1], labels12[isig][i2], isig);
498 } // end loop over signals
500 // release the memory
501 for(Int_t isig=0;isig<nSignals;++isig) {
502 delete [] *(labels1+isig);
503 delete [] *(labels2+isig);
504 delete [] *(labels12+isig);
514 //________________________________________________________________
515 void AliDielectron::FillHistogramsTracks(TObjArray **tracks)
518 // Fill Histogram information for tracks after prefilter
519 // ignore mixed events - for prefilter, only single tracks +/- are relevant
522 TString className,className2;
523 Double_t values[AliDielectronVarManager::kNMaxValues];
525 //Fill track information, separately for the track array candidates
526 for (Int_t i=0; i<2; ++i){
527 className.Form("Pre_%s",fgkTrackClassNames[i]);
528 if (!fHistos->GetHistogramList()->FindObject(className.Data())) continue;
529 Int_t ntracks=tracks[i]->GetEntriesFast();
530 for (Int_t itrack=0; itrack<ntracks; ++itrack){
531 AliDielectronVarManager::Fill(tracks[i]->UncheckedAt(itrack), values);
532 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
538 //________________________________________________________________
539 void AliDielectron::FillHistogramsMC(const AliMCEvent *ev, AliVEvent *ev1)
542 // Fill Histogram information for MCEvents
545 Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
546 // Fill event information
547 AliDielectronVarManager::Fill(ev1, values); // ESD/AOD information
548 AliDielectronVarManager::Fill(ev, values); // MC truth info
549 if (fHistos->GetHistogramList()->FindObject("MCEvent"))
550 fHistos->FillClass("MCEvent", AliDielectronVarManager::kNMaxValues, values);
554 //________________________________________________________________
555 void AliDielectron::FillHistograms(const AliVEvent *ev, Bool_t pairInfoOnly)
558 // Fill Histogram information for tracks and pairs
561 TString className,className2;
562 Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
564 //Fill event information
566 if (fHistos->GetHistogramList()->FindObject("Event")) {
567 fHistos->FillClass("Event", AliDielectronVarManager::kNMaxValues, AliDielectronVarManager::GetData());
571 //Fill track information, separately for the track array candidates
573 className2.Form("Track_%s",fgkPairClassNames[1]); // unlike sign, SE only
574 for (Int_t i=0; i<4; ++i){
575 className.Form("Track_%s",fgkTrackClassNames[i]);
576 Bool_t mergedtrkClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
577 Bool_t trkClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
578 if (!trkClass && !mergedtrkClass) continue;
579 Int_t ntracks=fTracks[i].GetEntriesFast();
580 for (Int_t itrack=0; itrack<ntracks; ++itrack){
581 AliDielectronVarManager::Fill(fTracks[i].UncheckedAt(itrack), values);
583 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
584 if(mergedtrkClass && i<2)
585 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values); //only ev1
590 //Fill Pair information, separately for all pair candidate arrays and the legs
591 TObjArray arrLegs(100);
592 for (Int_t i=0; i<10; ++i){
593 className.Form("Pair_%s",fgkPairClassNames[i]);
594 className2.Form("Track_Legs_%s",fgkPairClassNames[i]);
595 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
596 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
597 if (!pairClass&&!legClass) continue;
598 Int_t ntracks=PairArray(i)->GetEntriesFast();
599 for (Int_t ipair=0; ipair<ntracks; ++ipair){
600 AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair));
602 //fill pair information
604 AliDielectronVarManager::Fill(pair, values);
605 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
608 //fill leg information, don't fill the information twice
610 AliVParticle *d1=pair->GetFirstDaughter();
611 AliVParticle *d2=pair->GetSecondDaughter();
612 if (!arrLegs.FindObject(d1)){
613 AliDielectronVarManager::Fill(d1, values);
614 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
617 if (!arrLegs.FindObject(d2)){
618 AliDielectronVarManager::Fill(d2, values);
619 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
624 if (legClass) arrLegs.Clear();
628 //________________________________________________________________
629 void AliDielectron::FillHistogramsPair(AliDielectronPair *pair,Bool_t fromPreFilter/*=kFALSE*/)
632 // Fill Histogram information for pairs and the track in the pair
633 // NOTE: in this funtion the leg information may be filled multiple
634 // times. This funtion is used in the track rotation pairing
635 // and those legs are not saved!
637 TString className,className2;
638 Double_t values[AliDielectronVarManager::kNMaxValues];
640 //Fill Pair information, separately for all pair candidate arrays and the legs
641 TObjArray arrLegs(100);
642 const Int_t type=pair->GetType();
644 className.Form("RejPair_%s",fgkPairClassNames[type]);
645 className2.Form("RejTrack_%s",fgkPairClassNames[type]);
647 className.Form("Pair_%s",fgkPairClassNames[type]);
648 className2.Form("Track_Legs_%s",fgkPairClassNames[type]);
651 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
652 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
654 //fill pair information
656 AliDielectronVarManager::Fill(pair, values);
657 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
661 AliVParticle *d1=pair->GetFirstDaughter();
662 AliDielectronVarManager::Fill(d1, values);
663 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
665 AliVParticle *d2=pair->GetSecondDaughter();
666 AliDielectronVarManager::Fill(d2, values);
667 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
671 //________________________________________________________________
672 void AliDielectron::FillTrackArrays(AliVEvent * const ev, Int_t eventNr)
675 // select tracks and fill track candidate arrays
676 // eventNr = 0: First event, use track arrays 0 and 1
677 // eventNr = 1: Second event, use track arrays 2 and 3
680 Int_t ntracks=ev->GetNumberOfTracks();
682 UInt_t selectedMask=(1<<fTrackFilter.GetCuts()->GetEntries())-1;
683 for (Int_t itrack=0; itrack<ntracks; ++itrack){
685 AliVParticle *particle=ev->GetTrack(itrack);
688 UInt_t cutmask=fTrackFilter.IsSelected(particle);
690 if(fCutQA) fQAmonitor->FillAll(particle);
691 if(fCutQA) fQAmonitor->Fill(cutmask,particle);
693 if (cutmask!=selectedMask) continue;
695 //fill selected particle into the corresponding track arrays
696 Short_t charge=particle->Charge();
697 if (charge>0) fTracks[eventNr*2].Add(particle);
698 else if (charge<0) fTracks[eventNr*2+1].Add(particle);
704 //________________________________________________________________
705 void AliDielectron::EventPlanePreFilter(Int_t arr1, Int_t arr2, TObjArray arrTracks1, TObjArray arrTracks2, const AliVEvent *ev)
708 // Prefilter tracks and tracks from pairs
709 // Needed for rejection in the Q-Vector of the event plane
710 // remove contribution of all tracks to the Q-vector that are in invariant mass window
713 AliEventplane *evplane = const_cast<AliVEvent *>(ev)->GetEventplane();
714 if(!evplane) { // nanoAODs , here we do NOT have sub event reaction planes
716 // get the EPselectionTask for recalculation of weighting factors
717 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
718 AliEPSelectionTask *eptask = dynamic_cast<AliEPSelectionTask *>(man->GetTask("EventplaneSelection"));
722 TMap mapRemovedTracks;
725 Double_t cQX=0., cQY=0.;
726 // apply cuts to the tracks, e.g. etagap
727 if(fEventPlanePreFilter.GetCuts()->GetEntries()) {
728 UInt_t selectedMask=(1<<fEventPlanePreFilter.GetCuts()->GetEntries())-1;
729 Int_t ntracks=ev->GetNumberOfTracks();
730 for (Int_t itrack=0; itrack<ntracks; ++itrack){
731 AliVParticle *particle=ev->GetTrack(itrack);
732 AliVTrack *track= static_cast<AliVTrack*>(particle);
733 if (!track) continue;
735 UInt_t cutMask=fEventPlanePreFilter.IsSelected(track);
737 if (cutMask==selectedMask) continue;
739 mapRemovedTracks.Add(track,track);
740 cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
741 cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
745 // POI (particle of interest) rejection
746 Int_t pairIndex=GetPairIndex(arr1,arr2);
748 Int_t ntrack1=arrTracks1.GetEntriesFast();
749 Int_t ntrack2=arrTracks2.GetEntriesFast();
750 AliDielectronPair candidate;
751 candidate.SetKFUsage(fUseKF);
753 UInt_t selectedMask=(1<<fEventPlanePOIPreFilter.GetCuts()->GetEntries())-1;
754 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
756 if (arr1==arr2) end=itrack1;
757 Bool_t accepted=kFALSE;
758 for (Int_t itrack2=0; itrack2<end; ++itrack2){
759 TObject *track1=arrTracks1.UncheckedAt(itrack1);
760 TObject *track2=arrTracks2.UncheckedAt(itrack2);
761 if (!track1 || !track2) continue;
763 candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
764 static_cast<AliVTrack*>(track2), fPdgLeg2);
765 candidate.SetType(pairIndex);
766 candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
768 //event plane pair cuts
769 UInt_t cutMask=fEventPlanePOIPreFilter.IsSelected(&candidate);
771 if (cutMask==selectedMask) continue;
774 //remove the tracks from the Track arrays
775 arrTracks2.AddAt(0x0,itrack2);
777 if ( accepted ) arrTracks1.AddAt(0x0,itrack1);
779 //compress the track arrays
780 arrTracks1.Compress();
781 arrTracks2.Compress();
783 //Modify the components: subtract the tracks
784 ntrack1=arrTracks1.GetEntriesFast();
785 ntrack2=arrTracks2.GetEntriesFast();
786 // remove leg1 contribution
787 for (Int_t itrack=0; itrack<ntrack1; ++itrack){
788 AliVTrack *track= static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack));
789 if (!track) continue;
790 // track contribution was already removed
791 if (mapRemovedTracks.FindObject(track)) continue;
792 else mapRemovedTracks.Add(track,track);
794 cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
795 cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
797 // remove leg2 contribution
798 for (Int_t itrack=0; itrack<ntrack2; ++itrack){
799 AliVTrack *track= static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack));
800 if (!track) continue;
801 // track contribution was already removed
802 if (mapRemovedTracks.FindObject(track)) continue;
803 else mapRemovedTracks.Add(track,track);
805 cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
806 cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
809 // build a corrected alieventplane using the values from the var manager
810 // these uncorrected values are filled using the stored magnitude and angle in the header
812 qcorr.Set(AliDielectronVarManager::GetValue(AliDielectronVarManager::kTPCxH2uc)-cQX,
813 AliDielectronVarManager::GetValue(AliDielectronVarManager::kTPCyH2uc)-cQY);
814 // fill alieventplane
815 AliEventplane cevplane;
816 cevplane.SetQVector(&qcorr);
817 AliDielectronVarManager::SetTPCEventPlane(&cevplane);
818 cevplane.SetQVector(0);
823 // this is done in case of ESDs or AODs
824 Bool_t isESD=(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliESDInputHandler::Class());
825 // copy event plane object
826 AliEventplane cevplane(*evplane);
827 // Int_t nMaxID=cevplane->GetQContributionXArray()->GetSize();
829 TVector2 *qcorr = cevplane.GetQVector();
831 TVector2 *qcsub1 = 0x0;
832 TVector2 *qcsub2 = 0x0;
835 Bool_t etagap = kFALSE;
836 for (Int_t iCut=0; iCut<fEventPlanePreFilter.GetCuts()->GetEntries();++iCut) {
837 TString cutName=fEventPlanePreFilter.GetCuts()->At(iCut)->GetName();
838 if(cutName.Contains("eta") || cutName.Contains("Eta")) etagap=kTRUE;
841 // subevent configuration for eta gap or LS (default is rndm)
842 if(fLikeSignSubEvents && etagap) {
843 // start with the full Qvector/event in both sub events
844 qcsub1 = new TVector2(*qcorr);
845 qcsub2 = new TVector2(*qcorr);
846 cevplane.SetQsub(qcsub1,qcsub2);
848 Int_t ntracks=ev->GetNumberOfTracks();
850 for (Int_t itrack=0; itrack<ntracks; ++itrack){
851 AliVParticle *particle=ev->GetTrack(itrack);
852 AliVTrack *track= static_cast<AliVTrack*>(particle);
853 if (!track) continue;
854 if (track->GetID()>=0 && !isESD) continue;
855 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
857 // set contributions to zero
858 // charge sub1+ sub2-
859 if(fLikeSignSubEvents) {
860 Short_t charge=track->Charge();
862 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
863 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
866 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
867 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
872 Double_t eta=track->Eta();
874 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
875 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
878 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
879 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
882 } // end: loop over tracks
883 } // end: sub event configuration
885 // apply cuts, e.g. etagap
886 if(fEventPlanePreFilter.GetCuts()->GetEntries()) {
887 UInt_t selectedMask=(1<<fEventPlanePreFilter.GetCuts()->GetEntries())-1;
888 Int_t ntracks=ev->GetNumberOfTracks();
889 for (Int_t itrack=0; itrack<ntracks; ++itrack){
890 AliVParticle *particle=ev->GetTrack(itrack);
891 AliVTrack *track= static_cast<AliVTrack*>(particle);
892 if (!track) continue;
893 if (track->GetID()>=0 && !isESD) continue;
894 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
897 UInt_t cutMask=fEventPlanePreFilter.IsSelected(track);
899 if (cutMask==selectedMask) continue;
901 // set contributions to zero
902 cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
903 cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
904 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
905 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
906 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
907 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
911 // POI (particle of interest) rejection
912 Int_t pairIndex=GetPairIndex(arr1,arr2);
913 Int_t ntrack1=arrTracks1.GetEntriesFast();
914 Int_t ntrack2=arrTracks2.GetEntriesFast();
915 AliDielectronPair candidate;
916 candidate.SetKFUsage(fUseKF);
918 UInt_t selectedMask=(1<<fEventPlanePOIPreFilter.GetCuts()->GetEntries())-1;
919 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
921 if (arr1==arr2) end=itrack1;
922 Bool_t accepted=kFALSE;
923 for (Int_t itrack2=0; itrack2<end; ++itrack2){
924 TObject *track1=arrTracks1.UncheckedAt(itrack1);
925 TObject *track2=arrTracks2.UncheckedAt(itrack2);
926 if (!track1 || !track2) continue;
928 candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
929 static_cast<AliVTrack*>(track2), fPdgLeg2);
931 candidate.SetType(pairIndex);
932 candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
935 UInt_t cutMask=fEventPlanePOIPreFilter.IsSelected(&candidate);
937 if (cutMask==selectedMask) continue;
940 //remove the tracks from the Track arrays
941 arrTracks2.AddAt(0x0,itrack2);
943 if ( accepted ) arrTracks1.AddAt(0x0,itrack1);
945 //compress the track arrays
946 arrTracks1.Compress();
947 arrTracks2.Compress();
949 //Modify the components: subtract the tracks
950 ntrack1=arrTracks1.GetEntriesFast();
951 ntrack2=arrTracks2.GetEntriesFast();
952 // remove leg1 contribution
953 for (Int_t itrack=0; itrack<ntrack1; ++itrack){
954 AliVTrack *track= static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack));
955 if (!track) continue;
956 if (track->GetID()>=0 && !isESD) continue;
957 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
958 // set contributions to zero
959 cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
960 cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
961 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
962 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
963 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
964 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
966 // remove leg2 contribution
967 for (Int_t itrack=0; itrack<ntrack2; ++itrack){
968 AliVTrack *track= static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack));
969 if (!track) continue;
970 if (track->GetID()>=0 && !isESD) continue;
971 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
972 // set contributions to zero
973 cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
974 cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
975 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
976 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
977 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
978 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
981 // set corrected AliEventplane and fill variables with corrected values
982 AliDielectronVarManager::SetTPCEventPlane(&cevplane);
985 } // end: ESD or AOD case
988 //________________________________________________________________
989 void AliDielectron::PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2)
992 // Prefilter tracks from pairs
993 // Needed for datlitz rejections
994 // remove all tracks from the Single track arrays that pass the cuts in this filter
997 Int_t ntrack1=arrTracks1.GetEntriesFast();
998 Int_t ntrack2=arrTracks2.GetEntriesFast();
999 AliDielectronPair candidate;
1000 candidate.SetKFUsage(fUseKF);
1001 // flag arrays for track removal
1002 Bool_t *bTracks1 = new Bool_t[ntrack1];
1003 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1) bTracks1[itrack1]=kFALSE;
1004 Bool_t *bTracks2 = new Bool_t[ntrack2];
1005 for (Int_t itrack2=0; itrack2<ntrack2; ++itrack2) bTracks2[itrack2]=kFALSE;
1007 UInt_t selectedMask=(1<<fPairPreFilter.GetCuts()->GetEntries())-1;
1008 UInt_t selectedMaskPair=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1010 Int_t nRejPasses = 1; //for fPreFilterUnlikeOnly and no set flag
1011 if (fPreFilterAllSigns) nRejPasses = 3;
1013 for (Int_t iRP=0; iRP < nRejPasses; ++iRP) {
1014 Int_t arr1RP=arr1, arr2RP=arr2;
1015 TObjArray *arrTracks1RP=&arrTracks1;
1016 TObjArray *arrTracks2RP=&arrTracks2;
1017 Bool_t *bTracks1RP = bTracks1;
1018 Bool_t *bTracks2RP = bTracks2;
1020 case 1: arr1RP=arr1;arr2RP=arr1;
1021 arrTracks1RP=&arrTracks1;
1022 arrTracks2RP=&arrTracks1;
1023 bTracks1RP = bTracks1;
1024 bTracks2RP = bTracks1;
1026 case 2: arr1RP=arr2;arr2RP=arr2;
1027 arrTracks1RP=&arrTracks2;
1028 arrTracks2RP=&arrTracks2;
1029 bTracks1RP = bTracks2;
1030 bTracks2RP = bTracks2;
1032 default: ;//nothing to do
1034 Int_t ntrack1RP=(*arrTracks1RP).GetEntriesFast();
1035 Int_t ntrack2RP=(*arrTracks2RP).GetEntriesFast();
1037 Int_t pairIndex=GetPairIndex(arr1RP,arr2RP);
1039 for (Int_t itrack1=0; itrack1<ntrack1RP; ++itrack1){
1040 Int_t end=ntrack2RP;
1041 if (arr1RP==arr2RP) end=itrack1;
1042 for (Int_t itrack2=0; itrack2<end; ++itrack2){
1043 TObject *track1=(*arrTracks1RP).UncheckedAt(itrack1);
1044 TObject *track2=(*arrTracks2RP).UncheckedAt(itrack2);
1045 if (!track1 || !track2) continue;
1047 candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
1048 static_cast<AliVTrack*>(track2), fPdgLeg2);
1050 candidate.SetType(pairIndex);
1051 candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
1052 //relate to the production vertex
1053 // if (AliDielectronVarManager::GetKFVertex()) candidate.SetProductionVertex(*AliDielectronVarManager::GetKFVertex());
1056 UInt_t cutMask=fPairPreFilter.IsSelected(&candidate);
1059 if (cutMask!=selectedMask) continue;
1060 if (fCfManagerPair) fCfManagerPair->Fill(selectedMaskPair+1 ,&candidate);
1061 if (fHistos) FillHistogramsPair(&candidate,kTRUE);
1062 //set flags for track removal
1063 bTracks1RP[itrack1]=kTRUE;
1064 bTracks2RP[itrack2]=kTRUE;
1069 //remove the tracks from the Track arrays
1070 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
1071 if(bTracks1[itrack1]) arrTracks1.AddAt(0x0, itrack1);
1073 for (Int_t itrack2=0; itrack2<ntrack2; ++itrack2){
1074 if(bTracks2[itrack2]) arrTracks2.AddAt(0x0, itrack2);
1081 //compress the track arrays
1082 arrTracks1.Compress();
1083 arrTracks2.Compress();
1085 //apply leg cuts after the pre filter
1086 if ( fPairPreFilterLegs.GetCuts()->GetEntries()>0 ) {
1087 selectedMask=(1<<fPairPreFilterLegs.GetCuts()->GetEntries())-1;
1088 //loop over tracks from array 1
1089 for (Int_t itrack=0; itrack<arrTracks1.GetEntriesFast();++itrack){
1091 UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks1.UncheckedAt(itrack));
1094 if (cutMask!=selectedMask) arrTracks1.AddAt(0x0,itrack);
1096 arrTracks1.Compress();
1098 //in case of like sign don't loop over second array
1100 arrTracks2=arrTracks1;
1103 //loop over tracks from array 2
1104 for (Int_t itrack=0; itrack<arrTracks2.GetEntriesFast();++itrack){
1106 UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks2.UncheckedAt(itrack));
1108 if (cutMask!=selectedMask) arrTracks2.AddAt(0x0,itrack);
1110 arrTracks2.Compress();
1114 //For unlike-sign monitor track-cuts:
1115 if (arr1!=arr2&&fHistos) {
1116 TObjArray *unlikesignArray[2] = {&arrTracks1,&arrTracks2};
1117 FillHistogramsTracks(unlikesignArray);
1121 //________________________________________________________________
1122 void AliDielectron::FillPairArrays(Int_t arr1, Int_t arr2)
1125 // select pairs and fill pair candidate arrays
1128 TObjArray arrTracks1=fTracks[arr1];
1129 TObjArray arrTracks2=fTracks[arr2];
1131 //process pre filter if set
1132 if ((!fPreFilterAllSigns) && (!fPreFilterUnlikeOnly) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(arr1, arr2, arrTracks1, arrTracks2);
1134 Int_t pairIndex=GetPairIndex(arr1,arr2);
1136 Int_t ntrack1=arrTracks1.GetEntriesFast();
1137 Int_t ntrack2=arrTracks2.GetEntriesFast();
1139 AliDielectronPair *candidate=new AliDielectronPair;
1140 candidate->SetKFUsage(fUseKF);
1142 UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1144 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
1146 if (arr1==arr2) end=itrack1;
1147 for (Int_t itrack2=0; itrack2<end; ++itrack2){
1149 candidate->SetTracks(static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1)), fPdgLeg1,
1150 static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2)), fPdgLeg2);
1151 candidate->SetType(pairIndex);
1152 Int_t label=AliDielectronMC::Instance()->GetLabelMotherWithPdg(candidate,fPdgMother);
1153 candidate->SetLabel(label);
1154 if (label>-1) candidate->SetPdgCode(fPdgMother);
1155 else candidate->SetPdgCode(0);
1157 // check for gamma kf particle
1158 label=AliDielectronMC::Instance()->GetLabelMotherWithPdg(candidate,22);
1160 candidate->SetGammaTracks(static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1)), fPdgLeg1,
1161 static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2)), fPdgLeg2);
1162 // should we set the pdgmothercode and the label
1166 UInt_t cutMask=fPairFilter.IsSelected(candidate);
1168 //CF manager for the pair
1169 if (fCfManagerPair) fCfManagerPair->Fill(cutMask,candidate);
1172 if(pairIndex==kEv1PM && fCutQA) {
1173 fQAmonitor->FillAll(candidate);
1174 fQAmonitor->Fill(cutMask,candidate);
1178 if (cutMask!=selectedMask) continue;
1180 //histogram array for the pair
1181 if (fHistoArray) fHistoArray->Fill(pairIndex,candidate);
1183 //add the candidate to the candidate array
1184 PairArray(pairIndex)->Add(candidate);
1185 //get a new candidate
1186 candidate=new AliDielectronPair;
1187 candidate->SetKFUsage(fUseKF);
1190 //delete the surplus candidate
1194 //________________________________________________________________
1195 void AliDielectron::FillPairArrayTR()
1198 // select pairs and fill pair candidate arrays
1200 UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1202 while ( fTrackRotator->NextCombination() ){
1203 AliDielectronPair candidate;
1204 candidate.SetKFUsage(fUseKF);
1205 candidate.SetTracks(&fTrackRotator->GetKFTrackP(), &fTrackRotator->GetKFTrackN(),
1206 fTrackRotator->GetVTrackP(),fTrackRotator->GetVTrackN());
1207 candidate.SetType(kEv1PMRot);
1210 UInt_t cutMask=fPairFilter.IsSelected(&candidate);
1212 //CF manager for the pair
1213 if (fCfManagerPair) fCfManagerPair->Fill(cutMask,&candidate);
1216 if (cutMask==selectedMask) {
1218 //histogram array for the pair
1219 if (fHistoArray) fHistoArray->Fill((Int_t)kEv1PMRot,&candidate);
1221 if(fHistos) FillHistogramsPair(&candidate);
1222 if(fStoreRotatedPairs) PairArray(kEv1PMRot)->Add(new AliDielectronPair(candidate));
1227 //________________________________________________________________
1228 void AliDielectron::FillDebugTree()
1231 // Fill Histogram information for tracks and pairs
1235 for (Int_t i=0; i<10; ++i){
1236 Int_t ntracks=PairArray(i)->GetEntriesFast();
1237 for (Int_t ipair=0; ipair<ntracks; ++ipair){
1238 fDebugTree->Fill(static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair)));
1243 //________________________________________________________________
1244 void AliDielectron::SaveDebugTree()
1247 // delete the debug tree, this will also write the tree
1249 if (fDebugTree) fDebugTree->DeleteStreamer();
1253 //__________________________________________________________________
1254 void AliDielectron::AddSignalMC(AliDielectronSignalMC* signal) {
1256 // Add an MC signal to the signals list
1259 fSignalsMC = new TObjArray();
1260 fSignalsMC->SetOwner();
1262 fSignalsMC->Add(signal);
1265 //________________________________________________________________
1266 void AliDielectron::FillMCHistograms(Int_t label1, Int_t label2, Int_t nSignal) {
1268 // fill QA MC TRUTH histograms for pairs and legs of all added mc signals
1271 TString className,className2,className3;
1272 className.Form("Pair_%s_MCtruth",fSignalsMC->At(nSignal)->GetName());
1273 className2.Form("Track_Legs_%s_MCtruth",fSignalsMC->At(nSignal)->GetName());
1274 className3.Form("Track_%s_%s_MCtruth",fgkPairClassNames[1],fSignalsMC->At(nSignal)->GetName());
1275 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
1276 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
1277 Bool_t trkClass=fHistos->GetHistogramList()->FindObject(className3.Data())!=0x0;
1278 // printf("fill signal %d: pair %d legs %d trk %d \n",nSignal,pairClass,legClass,trkClass);
1279 if(!pairClass && !legClass && !trkClass) return;
1281 // printf("leg labels: %d-%d \n",label1,label2);
1282 AliVParticle* part1 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label1);
1283 AliVParticle* part2 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label2);
1284 if(!part1 && !part2) return;
1286 // fill only unlike sign (and only SE)
1287 if(part1->Charge()*part2->Charge()>=0) return;
1291 AliDielectronMC* dieMC = AliDielectronMC::Instance();
1293 Int_t mLabel1 = dieMC->GetMothersLabel(label1); // should work for both ESD and AOD
1294 Int_t mLabel2 = dieMC->GetMothersLabel(label2);
1296 // check the same mother option
1297 AliDielectronSignalMC* sigMC = (AliDielectronSignalMC*)fSignalsMC->At(nSignal);
1298 if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kSame && mLabel1!=mLabel2) return;
1299 if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kDifferent && mLabel1==mLabel2) return;
1301 // fill event values
1302 Double_t values[AliDielectronVarManager::kNMaxValues];
1303 AliDielectronVarManager::Fill(dieMC->GetMCEvent(), values); // get event informations
1305 // fill the leg variables
1306 // printf("leg:%d trk:%d part1:%p part2:%p \n",legClass,trkClass,part1,part2);
1307 if (legClass || trkClass) {
1308 if(part1) AliDielectronVarManager::Fill(part1,values);
1309 if(part1 && trkClass) fHistos->FillClass(className3, AliDielectronVarManager::kNMaxValues, values);
1310 if(part1 && part2 && legClass) fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1311 if(part2) AliDielectronVarManager::Fill(part2,values);
1312 if(part2 && trkClass) fHistos->FillClass(className3, AliDielectronVarManager::kNMaxValues, values);
1313 if(part1 && part2 && legClass) fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1316 //fill pair information
1317 if (pairClass && part1 && part2) {
1318 AliDielectronVarManager::FillVarMCParticle2(part1,part2,values);
1319 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
1324 //________________________________________________________________
1325 void AliDielectron::FillMCHistograms(const AliVEvent *ev) {
1327 // fill QA MC histograms for pairs and legs of all added mc signals
1330 if (!fSignalsMC) return;
1331 TString className,className2,className3;
1332 Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
1333 AliDielectronVarManager::Fill(ev, values); // get event informations
1334 //loop over all added mc signals
1335 for(Int_t isig=0; isig<fSignalsMC->GetEntries(); isig++) {
1337 //check if and what to fill
1338 className.Form("Pair_%s",fSignalsMC->At(isig)->GetName());
1339 className2.Form("Track_Legs_%s",fSignalsMC->At(isig)->GetName());
1340 className3.Form("Track_%s_%s",fgkPairClassNames[1],fSignalsMC->At(isig)->GetName()); // unlike sign, SE only
1341 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
1342 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
1343 Bool_t mergedtrkClass=fHistos->GetHistogramList()->FindObject(className3.Data())!=0x0;
1344 if(!pairClass && !legClass && !mergedtrkClass) continue;
1346 // fill pair and/or their leg variables
1347 if(pairClass || legClass) {
1348 Int_t npairs=PairArray(AliDielectron::kEv1PM)->GetEntriesFast(); // only SE +-
1349 for (Int_t ipair=0; ipair<npairs; ++ipair){
1350 AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(AliDielectron::kEv1PM)->UncheckedAt(ipair));
1352 Bool_t isMCtruth = AliDielectronMC::Instance()->IsMCTruth(pair, (AliDielectronSignalMC*)fSignalsMC->At(isig));
1354 //fill pair information
1356 AliDielectronVarManager::Fill(pair, values);
1357 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
1359 //fill leg information, both + and - in the same histo
1361 AliDielectronVarManager::Fill(pair->GetFirstDaughter(),values);
1362 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1363 AliDielectronVarManager::Fill(pair->GetSecondDaughter(),values);
1364 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1370 // fill single tracks of signals
1371 if(!mergedtrkClass) continue;
1372 // loop over SE track arrays
1373 for (Int_t i=0; i<2; ++i){
1374 Int_t ntracks=fTracks[i].GetEntriesFast();
1375 for (Int_t itrack=0; itrack<ntracks; ++itrack){
1376 Int_t label=((AliVParticle*)fTracks[i].UncheckedAt(itrack))->GetLabel();
1377 Bool_t isMCtruth1 = AliDielectronMC::Instance()->IsMCTruth(label, (AliDielectronSignalMC*)fSignalsMC->At(isig), 1);
1378 Bool_t isMCtruth2 = AliDielectronMC::Instance()->IsMCTruth(label, (AliDielectronSignalMC*)fSignalsMC->At(isig), 2);
1379 // skip if track does not correspond to the signal
1380 if(!isMCtruth1 && !isMCtruth2) continue;
1381 AliDielectronVarManager::Fill(fTracks[i].UncheckedAt(itrack), values);
1382 fHistos->FillClass(className3, AliDielectronVarManager::kNMaxValues, values);
1390 //______________________________________________
1391 void AliDielectron::SetCentroidCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
1393 fun->GetHistogram()->GetXaxis()->SetUniqueID(varx);
1394 fun->GetHistogram()->GetYaxis()->SetUniqueID(vary);
1395 fun->GetHistogram()->GetZaxis()->SetUniqueID(varz);
1396 fPostPIDCntrdCorr=fun;
1398 //______________________________________________
1399 void AliDielectron::SetWidthCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
1401 fun->GetHistogram()->GetXaxis()->SetUniqueID(varx);
1402 fun->GetHistogram()->GetYaxis()->SetUniqueID(vary);
1403 fun->GetHistogram()->GetZaxis()->SetUniqueID(varz);
1404 fPostPIDWdthCorr=fun;
1406 //______________________________________________
1407 THnBase* AliDielectron::InitEffMap(TString filename)
1409 // init an efficiency object for on-the-fly correction calculations
1410 if(filename.Contains("alien://") && !gGrid) TGrid::Connect("alien://",0,0,"t");
1412 TFile* file=TFile::Open(filename.Data());
1413 if(!file) return 0x0;
1414 THnBase *hGen = (THnBase*) file->Get("hGenerated");
1415 THnBase *hFnd = (THnBase*) file->Get("hFound");
1416 if(!hFnd || !hGen) return 0x0;
1419 printf("[I] AliDielectron::InitLegEffMap efficiency maps %s loaded! \n",filename.Data());
1420 return ((THnBase*) hFnd->Clone("effMap"));