]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGDQ/dielectron/AliDielectron.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / PWGDQ / dielectron / AliDielectron.cxx
CommitLineData
b2a297fa 1/*************************************************************************
2* Copyright(c) 1998-2009, 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///////////////////////////////////////////////////////////////////////////
17// Dielectron Analysis Main class //
18// //
19/*
20Framework to perform event selectoin, single track selection and track pair
21selection.
22
23Convention for the signs of the pair in fPairCandidates:
24The names are available via the function PairClassName(Int_t i)
25
260: ev1+ ev1+ (same event like sign +)
271: ev1+ ev1- (same event unlike sign)
282: ev1- ev1- (same event like sign -)
29
303: ev1+ ev2+ (mixed event like sign +)
314: ev1- ev2+ (mixed event unlike sign -+)
326: ev1+ ev2- (mixed event unlike sign +-)
337: ev1- ev2- (mixed event like sign -)
34
355: ev2+ ev2+ (same event like sign +)
368: ev2+ ev2- (same event unlike sign)
379: ev2- ev2- (same event like sign -)
38
2a14a7b1 3910: ev1+ ev1- (same event track rotation)
b2a297fa 40
41*/
42// //
43///////////////////////////////////////////////////////////////////////////
44
45#include <TString.h>
46#include <TList.h>
47#include <TMath.h>
5720c765 48#include <TObject.h>
794ba2d4 49#include <TGrid.h>
b2a297fa 50
5720c765 51#include <AliKFParticle.h>
b2a297fa 52
4533e78e 53#include <AliESDInputHandler.h>
54#include <AliAnalysisManager.h>
55#include <AliEPSelectionTask.h>
5720c765 56#include <AliEventplane.h>
b2a297fa 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"
2a14a7b1 65#include "AliDielectronTrackRotator.h"
572b0139 66#include "AliDielectronDebugTree.h"
ba15fdfb 67#include "AliDielectronSignalMC.h"
5720c765 68#include "AliDielectronMixingHandler.h"
42c76bea 69#include "AliDielectronPairLegCuts.h"
061ca303 70#include "AliDielectronV0Cuts.h"
30f9a393 71#include "AliDielectronPID.h"
7238fa88 72#include "AliDielectronHistos.h"
572b0139 73
b2a297fa 74#include "AliDielectron.h"
75
76ClassImp(AliDielectron)
77
78const char* AliDielectron::fgkTrackClassNames[4] = {
79 "ev1+",
80 "ev1-",
81 "ev2+",
82 "ev2-"
83};
84
2a14a7b1 85const char* AliDielectron::fgkPairClassNames[11] = {
b2a297fa 86 "ev1+_ev1+",
87 "ev1+_ev1-",
88 "ev1-_ev1-",
89 "ev1+_ev2+",
90 "ev1-_ev2+",
91 "ev2+_ev2+",
92 "ev1+_ev2-",
93 "ev1-_ev2-",
94 "ev2+_ev2-",
2a14a7b1 95 "ev2-_ev2-",
96 "ev1+_ev1-_TR"
b2a297fa 97};
98
99//________________________________________________________________
100AliDielectron::AliDielectron() :
101 TNamed("AliDielectron","AliDielectron"),
d8bb1abe 102 fCutQA(kFALSE),
103 fQAmonitor(0x0),
30f9a393
JB		 104 fPostPIDCntrdCorr(0x0),
105 fPostPIDWdthCorr(0x0),
7965afb6 106 fLegEffMap(0x0),
0b369e55 107 fPairEffMap(0x0),
b2a297fa 108 fEventFilter("EventFilter"),
109 fTrackFilter("TrackFilter"),
61d106d3 110 fPairPreFilter("PairPreFilter"),
111 fPairPreFilterLegs("PairPreFilterLegs"),
b2a297fa 112 fPairFilter("PairFilter"),
5720c765 113 fEventPlanePreFilter("EventPlanePreFilter"),
114 fEventPlanePOIPreFilter("EventPlanePOIPreFilter"),
b2a297fa 115 fPdgMother(443),
8df8e382 116 fPdgLeg1(11),
117 fPdgLeg2(11),
ba15fdfb 118 fSignalsMC(0x0),
554e40f8 119 fNoPairing(kFALSE),
5156b880 120 fProcessLS(kTRUE),
08b801a6 121 fUseKF(kTRUE),
620b34a5 122 fHistoArray(0x0),
b2a297fa 123 fHistos(0x0),
a94c2e7e 124 fUsedVars(new TBits(AliDielectronVarManager::kNMaxValues)),
5720c765 125 fPairCandidates(new TObjArray(11)),
572b0139 126 fCfManagerPair(0x0),
2a14a7b1 127 fTrackRotator(0x0),
554e40f8 128 fDebugTree(0x0),
5720c765 129 fMixing(0x0),
130 fPreFilterEventPlane(kFALSE),
131 fLikeSignSubEvents(kFALSE),
fb7d2d99 132 fPreFilterUnlikeOnly(kFALSE),
ba15fdfb 133 fPreFilterAllSigns(kFALSE),
5720c765 134 fHasMC(kFALSE),
135 fStoreRotatedPairs(kFALSE),
2e02dba4 136 fDontClearArrays(kFALSE),
42c76bea 137 fEventProcess(kTRUE),
5720c765 138 fEstimatorFilename(""),
a823f01b 139 fTRDpidCorrectionFilename(""),
140 fVZEROCalibrationFilename(""),
794ba2d4 141 fVZERORecenteringFilename(""),
36a1e727
JW		 142 fZDCRecenteringFilename("")
143
b2a297fa 144{
145 //
146 // Default constructor
147 //
148
149}
150
151//________________________________________________________________
152AliDielectron::AliDielectron(const char* name, const char* title) :
153 TNamed(name,title),
d8bb1abe 154 fCutQA(kFALSE),
155 fQAmonitor(0x0),
30f9a393
JB		 156 fPostPIDCntrdCorr(0x0),
157 fPostPIDWdthCorr(0x0),
7f7b25bf
JB		 158 fLegEffMap(0x0),
159 fPairEffMap(0x0),
b2a297fa 160 fEventFilter("EventFilter"),
161 fTrackFilter("TrackFilter"),
61d106d3 162 fPairPreFilter("PairPreFilter"),
163 fPairPreFilterLegs("PairPreFilterLegs"),
b2a297fa 164 fPairFilter("PairFilter"),
5720c765 165 fEventPlanePreFilter("EventPlanePreFilter"),
166 fEventPlanePOIPreFilter("EventPlanePOIPreFilter"),
b2a297fa 167 fPdgMother(443),
8df8e382 168 fPdgLeg1(11),
169 fPdgLeg2(11),
ba15fdfb 170 fSignalsMC(0x0),
554e40f8 171 fNoPairing(kFALSE),
5156b880 172 fProcessLS(kTRUE),
08b801a6 173 fUseKF(kTRUE),
620b34a5 174 fHistoArray(0x0),
b2a297fa 175 fHistos(0x0),
a94c2e7e 176 fUsedVars(new TBits(AliDielectronVarManager::kNMaxValues)),
5720c765 177 fPairCandidates(new TObjArray(11)),
572b0139 178 fCfManagerPair(0x0),
2a14a7b1 179 fTrackRotator(0x0),
554e40f8 180 fDebugTree(0x0),
5720c765 181 fMixing(0x0),
182 fPreFilterEventPlane(kFALSE),
183 fLikeSignSubEvents(kFALSE),
fb7d2d99 184 fPreFilterUnlikeOnly(kFALSE),
ba15fdfb 185 fPreFilterAllSigns(kFALSE),
5720c765 186 fHasMC(kFALSE),
187 fStoreRotatedPairs(kFALSE),
2e02dba4 188 fDontClearArrays(kFALSE),
42c76bea 189 fEventProcess(kTRUE),
5720c765 190 fEstimatorFilename(""),
a823f01b 191 fTRDpidCorrectionFilename(""),
192 fVZEROCalibrationFilename(""),
794ba2d4 193 fVZERORecenteringFilename(""),
36a1e727 194 fZDCRecenteringFilename("")
b2a297fa 195{
196 //
197 // Named constructor
198 //
199
200}
201
202//________________________________________________________________
203AliDielectron::~AliDielectron()
204{
205 //
206 // Default destructor
207 //
d8bb1abe 208 if (fQAmonitor) delete fQAmonitor;
30f9a393
JB		 209 if (fPostPIDCntrdCorr) delete fPostPIDCntrdCorr;
210 if (fPostPIDWdthCorr) delete fPostPIDWdthCorr;
7965afb6 211 if (fLegEffMap) delete fLegEffMap;
0b369e55 212 if (fPairEffMap) delete fPairEffMap;
b2a297fa 213 if (fHistos) delete fHistos;
a94c2e7e 214 if (fUsedVars) delete fUsedVars;
42c76bea 215 if (fPairCandidates && fEventProcess) delete fPairCandidates;
572b0139 216 if (fDebugTree) delete fDebugTree;
5720c765 217 if (fMixing) delete fMixing;
ba15fdfb 218 if (fSignalsMC) delete fSignalsMC;
5720c765 219 if (fCfManagerPair) delete fCfManagerPair;
d4619886 220 if (fHistoArray) delete fHistoArray;
b2a297fa 221}
222
223//________________________________________________________________
224void AliDielectron::Init()
225{
226 //
227 // Initialise objects
228 //
fb7d2d99 229
230 if(GetHasMC()) AliDielectronMC::Instance()->SetHasMC(GetHasMC());
5720c765 231
c1bd0a44 232 if(fEventProcess) InitPairCandidateArrays();
5720c765 233
ba15fdfb 234 if (fCfManagerPair) {
235 fCfManagerPair->SetSignalsMC(fSignalsMC);
236 fCfManagerPair->InitialiseContainer(fPairFilter);
237 }
1201a1a9 238 if (fTrackRotator) {
239 fTrackRotator->SetTrackArrays(&fTracks[0],&fTracks[1]);
240 fTrackRotator->SetPdgLegs(fPdgLeg1,fPdgLeg2);
241 }
8df8e382 242 if (fDebugTree) fDebugTree->SetDielectron(this);
30f9a393 243
30f9a393 244 if(fEstimatorFilename.Contains(".root")) AliDielectronVarManager::InitEstimatorAvg(fEstimatorFilename.Data());
5720c765 245 if(fTRDpidCorrectionFilename.Contains(".root")) AliDielectronVarManager::InitTRDpidEffHistograms(fTRDpidCorrectionFilename.Data());
a823f01b 246 if(fVZEROCalibrationFilename.Contains(".root")) AliDielectronVarManager::SetVZEROCalibrationFile(fVZEROCalibrationFilename.Data());
247 if(fVZERORecenteringFilename.Contains(".root")) AliDielectronVarManager::SetVZERORecenteringFile(fVZERORecenteringFilename.Data());
36a1e727 248 if(fZDCRecenteringFilename.Contains(".root")) AliDielectronVarManager::SetZDCRecenteringFile(fZDCRecenteringFilename.Data());
794ba2d4 249
5720c765 250 if (fMixing) fMixing->Init(this);
443a091c 251 if (fHistoArray) {
252 fHistoArray->SetSignalsMC(fSignalsMC);
253 fHistoArray->Init();
254 }
d8bb1abe 255
30f9a393
JB		 256 if(fPostPIDCntrdCorr) AliDielectronPID::SetCentroidCorrFunction(fPostPIDCntrdCorr);
257 if(fPostPIDWdthCorr) AliDielectronPID::SetWidthCorrFunction(fPostPIDWdthCorr);
258
42c76bea 259 if(!fEventProcess) {
260 AliDielectronPairLegCuts *trk2leg = new AliDielectronPairLegCuts("trk2leg","trk2leg");
261 // move all track cuts (if any) into pair leg cuts
262 TIter listIterator(fTrackFilter.GetCuts());
263 while (AliAnalysisCuts *thisCut = (AliAnalysisCuts*) listIterator()) {
264 trk2leg->GetLeg1Filter().AddCuts((AliAnalysisCuts*)thisCut->Clone());
265 trk2leg->GetLeg2Filter().AddCuts((AliAnalysisCuts*)thisCut->Clone());
266 }
267 // add pair leg cuts to pair filter
268 fPairFilter.AddCuts(trk2leg);
269 }
270
d8bb1abe 271 if (fCutQA) {
272 fQAmonitor = new AliDielectronCutQA(Form("QAcuts_%s",GetName()),"QAcuts");
273 fQAmonitor->AddTrackFilter(&fTrackFilter);
2ef1eb25 274 if(!fNoPairing) fQAmonitor->AddPairFilter(&fPairFilter);
d8bb1abe 275 fQAmonitor->AddEventFilter(&fEventFilter);
276 fQAmonitor->Init();
277 }
a94c2e7e 278
279 if(fHistos) {
280 (*fUsedVars)|= (*fHistos->GetUsedVars());
281 }
282
d8bb1abe 283}
b2a297fa 284
285//________________________________________________________________
42c76bea 286
287void AliDielectron::Process(TObjArray *arr)
288{
289 //
290 // Process the pair array
291 //
292
293 // set pair arrays
294 fPairCandidates = arr;
295
296 //fill debug tree if a manager is attached
297 // if (fDebugTree) FillDebugTree();
298 //in case there is a histogram manager, fill the QA histograms
299 // if (fHistos && fSignalsMC) FillMCHistograms(ev1);
300
301 // apply cuts and fill output
302 if (fHistos) FillHistogramsFromPairArray();
303
304 // never clear arrays !!!!
305
306
307}
308
309//________________________________________________________________
310Bool_t AliDielectron::Process(AliVEvent *ev1, AliVEvent *ev2)
b2a297fa 311{
312 //
313 // Process the events
314 //
315
45b2b1b8 316 //at least first event is needed!
317 if (!ev1){
318 AliError("At least first event must be set!");
42c76bea 319 return 0;
45b2b1b8 320 }
5720c765 321
061ca303 322 // modify event numbers in MC so that we can identify new events
323 // in AliDielectronV0Cuts (not neeeded for collision data)
324 if(GetHasMC()) {
325 ev1->SetBunchCrossNumber(1);
326 ev1->SetOrbitNumber(1);
327 ev1->SetPeriodNumber(1);
328 }
329
c315310a 330 // set event
331 AliDielectronVarManager::SetEvent(ev1);
332 if (fMixing){
333 //set mixing bin to event data
334 Int_t bin=fMixing->FindBin(AliDielectronVarManager::GetData());
335 AliDielectronVarManager::SetValue(AliDielectronVarManager::kMixingBin,bin);
336 }
337
e4942d79 338 // set efficiency maps
339 AliDielectronVarManager::SetLegEffMap(fLegEffMap);
340 AliDielectronVarManager::SetPairEffMap(fPairEffMap);
341
b2a297fa 342 //in case we have MC load the MC event and process the MC particles
d8bb1abe 343 // why do not apply the event cuts first ????
40875e45 344 if (AliDielectronMC::Instance()->ConnectMCEvent()){
5720c765 345 ProcessMC(ev1);
572b0139 346 }
0caf5fbb 347
b2a297fa 348 //if candidate array doesn't exist, create it
349 if (!fPairCandidates->UncheckedAt(0)) {
350 InitPairCandidateArrays();
351 } else {
352 ClearArrays();
353 }
354
355 //mask used to require that all cuts are fulfilled
356 UInt_t selectedMask=(1<<fEventFilter.GetCuts()->GetEntries())-1;
357
358 //apply event cuts
d8bb1abe 359 UInt_t cutmask = fEventFilter.IsSelected(ev1);
360 if(fCutQA) fQAmonitor->FillAll(ev1);
361 if(fCutQA) fQAmonitor->Fill(cutmask,ev1);
362 if ((ev1&&cutmask!=selectedMask) ||
42c76bea 363 (ev2&&fEventFilter.IsSelected(ev2)!=selectedMask)) return 0;
5720c765 364
b2a297fa 365 //fill track arrays for the first event
554e40f8 366 if (ev1){
367 FillTrackArrays(ev1);
ba15fdfb 368 if (((fPreFilterAllSigns)||(fPreFilterUnlikeOnly)) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(0, 1, fTracks[0], fTracks[1]);
554e40f8 369 }
370
b2a297fa 371
372 //fill track arrays for the second event
554e40f8 373 if (ev2) {
374 FillTrackArrays(ev2,1);
ba15fdfb 375 if (((fPreFilterAllSigns)||(fPreFilterUnlikeOnly)) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(2, 3, fTracks[2], fTracks[3]);
554e40f8 376 }
b2a297fa 377
5720c765 378 // TPC event plane correction
4533e78e 379 if (ev1 && fPreFilterEventPlane && ( fEventPlanePreFilter.GetCuts()->GetEntries()>0 || fEventPlanePOIPreFilter.GetCuts()->GetEntries()>0))
380 EventPlanePreFilter(0, 1, fTracks[0], fTracks[1], ev1);
f0351f9c 381
554e40f8 382 if (!fNoPairing){
383 // create pairs and fill pair candidate arrays
384 for (Int_t itrackArr1=0; itrackArr1<4; ++itrackArr1){
385 for (Int_t itrackArr2=itrackArr1; itrackArr2<4; ++itrackArr2){
5156b880
JB		 386 if(!fProcessLS && GetPairIndex(itrackArr1,itrackArr2)!=kEv1PM) continue;
387 FillPairArrays(itrackArr1, itrackArr2);
554e40f8 388 }
b2a297fa 389 }
6551594b 390
554e40f8 391 //track rotation
1201a1a9 392 if (fTrackRotator) {
393 fTrackRotator->SetEvent(ev1);
394 FillPairArrayTR();
395 }
554e40f8 396 }
572b0139 397
398 //fill debug tree if a manager is attached
399 if (fDebugTree) FillDebugTree();
5720c765 400
401 //process event mixing
402 if (fMixing) {
403 fMixing->Fill(ev1,this);
4533e78e 404 // FillHistograms(0x0,kTRUE);
5720c765 405 }
406
b829fbd2 407 // fill candidate variables
408 Double_t ntracks = fTracks[0].GetEntriesFast() + fTracks[1].GetEntriesFast();
409 Double_t npairs = PairArray(AliDielectron::kEv1PM)->GetEntriesFast();
410 AliDielectronVarManager::SetValue(AliDielectronVarManager::kTracks, ntracks);
411 AliDielectronVarManager::SetValue(AliDielectronVarManager::kPairs, npairs);
412
5720c765 413 //in case there is a histogram manager, fill the QA histograms
93d4ec8a 414 if (fHistos && fSignalsMC) FillMCHistograms(ev1);
5720c765 415 if (fHistos) FillHistograms(ev1);
416
f0351f9c 417
5720c765 418 // clear arrays
2e02dba4 419 if (!fDontClearArrays) ClearArrays();
061ca303 420
421 // reset TPC EP and unique identifiers for v0 cut class
5720c765 422 AliDielectronVarManager::SetTPCEventPlane(0x0);
061ca303 423 if(GetHasMC()) { // only for MC needed
424 for (Int_t iCut=0; iCut<fTrackFilter.GetCuts()->GetEntries();++iCut) {
425 if ( fTrackFilter.GetCuts()->At(iCut)->IsA() == AliDielectronV0Cuts::Class() )
426 ((AliDielectronV0Cuts*)fTrackFilter.GetCuts()->At(iCut))->ResetUniqueEventNumbers();
427 }
428 }
429
42c76bea 430 return 1;
431
b2a297fa 432}
433
434//________________________________________________________________
5720c765 435void AliDielectron::ProcessMC(AliVEvent *ev1)
b2a297fa 436{
437 //
438 // Process the MC data
439 //
440
ba15fdfb 441 AliDielectronMC *dieMC=AliDielectronMC::Instance();
442
5720c765 443 if (fHistos) FillHistogramsMC(dieMC->GetMCEvent(), ev1);
ba15fdfb 444
b829fbd2 445 // mc tracks
446 if(!dieMC->GetNMCTracks()) return;
88204efa 447
ba15fdfb 448 // signals to be studied
b829fbd2 449 if(!fSignalsMC) return;
ba15fdfb 450 Int_t nSignals = fSignalsMC->GetEntries();
b829fbd2 451 if(!nSignals) return;
452
453 //loop over all MC data and Fill the HF, CF containers and histograms if they exist
454 if(fCfManagerPair) fCfManagerPair->SetPdgMother(fPdgMother);
ba15fdfb 455
833da63d 456 Bool_t bFillCF = (fCfManagerPair ? fCfManagerPair->GetStepForMCtruth() : kFALSE);
457 Bool_t bFillHF = (fHistoArray ? fHistoArray->GetStepForMCGenerated() : kFALSE);
458 Bool_t bFillHist = kFALSE;
b829fbd2 459 if(fHistos) {
833da63d 460 const THashList *histlist = fHistos->GetHistogramList();
b829fbd2 461 for(Int_t isig=0;isig<nSignals;isig++) {
462 TString sigName = fSignalsMC->At(isig)->GetName();
463 bFillHist |= histlist->FindObject(Form("Pair_%s_MCtruth",sigName.Data()))!=0x0;
464 bFillHist |= histlist->FindObject(Form("Track_Leg_%s_MCtruth",sigName.Data()))!=0x0;
465 bFillHist |= histlist->FindObject(Form("Track_%s_%s_MCtruth",fgkPairClassNames[1],sigName.Data()))!=0x0;
466 if(bFillHist) break;
467 }
833da63d 468 }
469 // check if there is anything to fill
470 if(!bFillCF && !bFillHF && !bFillHist) return;
471
472
ba15fdfb 473 // initialize 2D arrays of labels for particles from each MC signal
474 Int_t** labels1; // labels for particles satisfying branch 1
475 Int_t** labels2; // labels for particles satisfying branch 2
476 Int_t** labels12; // labels for particles satisfying both branches
477 labels1 = new Int_t*[nSignals];
478 labels2 = new Int_t*[nSignals];
479 labels12 = new Int_t*[nSignals];
480 Int_t* indexes1=new Int_t[nSignals];
481 Int_t* indexes2=new Int_t[nSignals];
482 Int_t* indexes12=new Int_t[nSignals];
483 for(Int_t isig=0;isig<nSignals;++isig) {
484 *(labels1+isig) = new Int_t[dieMC->GetNMCTracks()];
485 *(labels2+isig) = new Int_t[dieMC->GetNMCTracks()];
486 *(labels12+isig) = new Int_t[dieMC->GetNMCTracks()];
487 for(Int_t ip=0; ip<dieMC->GetNMCTracks();++ip) {
488 labels1[isig][ip] = -1;
489 labels2[isig][ip] = -1;
490 labels12[isig][ip] = -1;
491 }
492 indexes1[isig]=0;
493 indexes2[isig]=0;
494 indexes12[isig]=0;
b2a297fa 495 }
ba15fdfb 496
497 Bool_t truth1=kFALSE;
498 Bool_t truth2=kFALSE;
499 // loop over the MC tracks
500 for(Int_t ipart=0; ipart<dieMC->GetNMCTracks(); ++ipart) {
501 for(Int_t isig=0; isig<nSignals; ++isig) { // loop over signals
502 // Proceed only if this signal is required in the pure MC step
503 // NOTE: Some signals can be satisfied by many particles and this leads to high
504 // computation times (e.g. secondary electrons from the GEANT transport). Be aware of this!!
505 if(!((AliDielectronSignalMC*)fSignalsMC->At(isig))->GetFillPureMCStep()) continue;
506
507 truth1 = dieMC->IsMCTruth(ipart, (AliDielectronSignalMC*)fSignalsMC->At(isig), 1);
508 truth2 = dieMC->IsMCTruth(ipart, (AliDielectronSignalMC*)fSignalsMC->At(isig), 2);
509
510 // particles satisfying both branches are treated separately to avoid double counting during pairing
511 if(truth1 && truth2) {
512 labels12[isig][indexes12[isig]] = ipart;
513 ++indexes12[isig];
514 }
515 else {
516 if(truth1) {
517 labels1[isig][indexes1[isig]] = ipart;
518 ++indexes1[isig];
519 }
520 if(truth2) {
521 labels2[isig][indexes2[isig]] = ipart;
522 ++indexes2[isig];
523 }
524 }
525 }
526 } // end loop over MC particles
527
528 // Do the pairing and fill the CF container with pure MC info
529 for(Int_t isig=0; isig<nSignals; ++isig) {
bc3cacd4 530 // printf("INDEXES: %d-%d both%d\n",indexes1[isig],indexes2[isig],indexes12[isig]);
ba15fdfb 531 // mix the particles which satisfy only one of the signal branches
532 for(Int_t i1=0;i1<indexes1[isig];++i1) {
bc3cacd4 533 if(!indexes2[isig]) FillMCHistograms(labels1[isig][i1], -1, isig); // (e.g. single electrons only, no pairs)
534 for(Int_t i2=0;i2<indexes2[isig];++i2) {
88204efa 535 if(bFillCF) fCfManagerPair->FillMC(labels1[isig][i1], labels2[isig][i2], isig);
536 if(bFillHF) fHistoArray->Fill(labels1[isig][i1], labels2[isig][i2], isig);
833da63d 537 FillMCHistograms(labels1[isig][i1], labels2[isig][i2], isig);
ba15fdfb 538 }
539 }
540 // mix the particles which satisfy both branches
541 for(Int_t i1=0;i1<indexes12[isig];++i1) {
542 for(Int_t i2=0; i2<i1; ++i2) {
88204efa 543 if(bFillCF) fCfManagerPair->FillMC(labels12[isig][i1], labels12[isig][i2], isig);
544 if(bFillHF) fHistoArray->Fill(labels12[isig][i1], labels12[isig][i2], isig);
833da63d 545 FillMCHistograms(labels12[isig][i1], labels12[isig][i2], isig);
ba15fdfb 546 }
547 }
548 } // end loop over signals
549
550 // release the memory
551 for(Int_t isig=0;isig<nSignals;++isig) {
552 delete [] *(labels1+isig);
553 delete [] *(labels2+isig);
554 delete [] *(labels12+isig);
555 }
556 delete [] labels1;
557 delete [] labels2;
558 delete [] labels12;
559 delete [] indexes1;
560 delete [] indexes2;
561 delete [] indexes12;
b2a297fa 562}
563
554e40f8 564//________________________________________________________________
565void AliDielectron::FillHistogramsTracks(TObjArray **tracks)
566{
567 //
568 // Fill Histogram information for tracks after prefilter
569 // ignore mixed events - for prefilter, only single tracks +/- are relevant
570 //
571
572 TString className,className2;
573 Double_t values[AliDielectronVarManager::kNMaxValues];
a94c2e7e 574 AliDielectronVarManager::SetFillMap(fUsedVars);
575
554e40f8 576 //Fill track information, separately for the track array candidates
577 for (Int_t i=0; i<2; ++i){
578 className.Form("Pre_%s",fgkTrackClassNames[i]);
579 if (!fHistos->GetHistogramList()->FindObject(className.Data())) continue;
580 Int_t ntracks=tracks[i]->GetEntriesFast();
581 for (Int_t itrack=0; itrack<ntracks; ++itrack){
582 AliDielectronVarManager::Fill(tracks[i]->UncheckedAt(itrack), values);
583 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
584 }
585 }
586}
ba15fdfb 587
588
589//________________________________________________________________
5720c765 590void AliDielectron::FillHistogramsMC(const AliMCEvent *ev, AliVEvent *ev1)
ba15fdfb 591{
592 //
593 // Fill Histogram information for MCEvents
594 //
595
5720c765 596 Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
a94c2e7e 597 AliDielectronVarManager::SetFillMap(fUsedVars);
598
ba15fdfb 599 // Fill event information
5720c765 600 AliDielectronVarManager::Fill(ev1, values); // ESD/AOD information
601 AliDielectronVarManager::Fill(ev, values); // MC truth info
ba15fdfb 602 if (fHistos->GetHistogramList()->FindObject("MCEvent"))
603 fHistos->FillClass("MCEvent", AliDielectronVarManager::kNMaxValues, values);
604}
605
606
b2a297fa 607//________________________________________________________________
c315310a 608void AliDielectron::FillHistograms(const AliVEvent *ev, Bool_t pairInfoOnly)
b2a297fa 609{
610 //
611 // Fill Histogram information for tracks and pairs
612 //
613
61d106d3 614 TString className,className2;
5720c765 615 Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
a94c2e7e 616 AliDielectronVarManager::SetFillMap(fUsedVars);
d327d9cd 617
0caf5fbb 618 //Fill event information
619 if (ev){
606603d9 620 if (fHistos->GetHistogramList()->FindObject("Event")) {
0caf5fbb 621 fHistos->FillClass("Event", AliDielectronVarManager::kNMaxValues, AliDielectronVarManager::GetData());
606603d9 622 }
5720c765 623 }
0caf5fbb 624
b2a297fa 625 //Fill track information, separately for the track array candidates
5720c765 626 if (!pairInfoOnly){
4533e78e 627 className2.Form("Track_%s",fgkPairClassNames[1]); // unlike sign, SE only
5720c765 628 for (Int_t i=0; i<4; ++i){
629 className.Form("Track_%s",fgkTrackClassNames[i]);
4533e78e 630 Bool_t mergedtrkClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
631 Bool_t trkClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
632 if (!trkClass && !mergedtrkClass) continue;
5720c765 633 Int_t ntracks=fTracks[i].GetEntriesFast();
634 for (Int_t itrack=0; itrack<ntracks; ++itrack){
635 AliDielectronVarManager::Fill(fTracks[i].UncheckedAt(itrack), values);
4533e78e 636 if(trkClass)
637 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
638 if(mergedtrkClass && i<2)
639 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values); //only ev1
5720c765 640 }
b2a297fa 641 }
642 }
643
61d106d3 644 //Fill Pair information, separately for all pair candidate arrays and the legs
645 TObjArray arrLegs(100);
b2a297fa 646 for (Int_t i=0; i<10; ++i){
647 className.Form("Pair_%s",fgkPairClassNames[i]);
61d106d3 648 className2.Form("Track_Legs_%s",fgkPairClassNames[i]);
649 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
650 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
651 if (!pairClass&&!legClass) continue;
b2a297fa 652 Int_t ntracks=PairArray(i)->GetEntriesFast();
653 for (Int_t ipair=0; ipair<ntracks; ++ipair){
61d106d3 654 AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair));
655
656 //fill pair information
657 if (pairClass){
658 AliDielectronVarManager::Fill(pair, values);
659 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
660 }
661
662 //fill leg information, don't fill the information twice
663 if (legClass){
54ce9dc4 664 AliVParticle *d1=pair->GetFirstDaughterP();
665 AliVParticle *d2=pair->GetSecondDaughterP();
61d106d3 666 if (!arrLegs.FindObject(d1)){
667 AliDielectronVarManager::Fill(d1, values);
668 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
669 arrLegs.Add(d1);
670 }
671 if (!arrLegs.FindObject(d2)){
672 AliDielectronVarManager::Fill(d2, values);
673 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
674 arrLegs.Add(d2);
675 }
676 }
b2a297fa 677 }
61d106d3 678 if (legClass) arrLegs.Clear();
b2a297fa 679 }
680
681}
2a14a7b1 682//________________________________________________________________
ffbede40 683void AliDielectron::FillHistogramsPair(AliDielectronPair *pair,Bool_t fromPreFilter/*=kFALSE*/)
2a14a7b1 684{
685 //
686 // Fill Histogram information for pairs and the track in the pair
687 // NOTE: in this funtion the leg information may be filled multiple
688 // times. This funtion is used in the track rotation pairing
689 // and those legs are not saved!
ffbede40 690 //
2a14a7b1 691 TString className,className2;
692 Double_t values[AliDielectronVarManager::kNMaxValues];
a94c2e7e 693 AliDielectronVarManager::SetFillMap(fUsedVars);
694
2a14a7b1 695 //Fill Pair information, separately for all pair candidate arrays and the legs
696 TObjArray arrLegs(100);
697 const Int_t type=pair->GetType();
ffbede40 698 if (fromPreFilter) {
699 className.Form("RejPair_%s",fgkPairClassNames[type]);
700 className2.Form("RejTrack_%s",fgkPairClassNames[type]);
701 } else {
702 className.Form("Pair_%s",fgkPairClassNames[type]);
703 className2.Form("Track_Legs_%s",fgkPairClassNames[type]);
704 }
2a14a7b1 705
706 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
707 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
708
709 //fill pair information
710 if (pairClass){
711 AliDielectronVarManager::Fill(pair, values);
712 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
713 }
714
715 if (legClass){
54ce9dc4 716 AliVParticle *d1=pair->GetFirstDaughterP();
2a14a7b1 717 AliDielectronVarManager::Fill(d1, values);
718 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
719
54ce9dc4 720 AliVParticle *d2=pair->GetSecondDaughterP();
2a14a7b1 721 AliDielectronVarManager::Fill(d2, values);
722 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
723 }
724}
b2a297fa 725
726//________________________________________________________________
727void AliDielectron::FillTrackArrays(AliVEvent * const ev, Int_t eventNr)
728{
729 //
730 // select tracks and fill track candidate arrays
731 // eventNr = 0: First event, use track arrays 0 and 1
732 // eventNr = 1: Second event, use track arrays 2 and 3
733 //
99345a64 734
7ab4dc00 735 Int_t ntracks=ev->GetNumberOfTracks();
99345a64 736
b2a297fa 737 UInt_t selectedMask=(1<<fTrackFilter.GetCuts()->GetEntries())-1;
738 for (Int_t itrack=0; itrack<ntracks; ++itrack){
739 //get particle
740 AliVParticle *particle=ev->GetTrack(itrack);
67fd1119 741
b2a297fa 742 //apply track cuts
d8bb1abe 743 UInt_t cutmask=fTrackFilter.IsSelected(particle);
744 //fill cut QA
745 if(fCutQA) fQAmonitor->FillAll(particle);
746 if(fCutQA) fQAmonitor->Fill(cutmask,particle);
747
748 if (cutmask!=selectedMask) continue;
5720c765 749
b2a297fa 750 //fill selected particle into the corresponding track arrays
751 Short_t charge=particle->Charge();
752 if (charge>0) fTracks[eventNr*2].Add(particle);
753 else if (charge<0) fTracks[eventNr*2+1].Add(particle);
d8bb1abe 754
755
b2a297fa 756 }
757}
758
61d106d3 759//________________________________________________________________
4533e78e 760void AliDielectron::EventPlanePreFilter(Int_t arr1, Int_t arr2, TObjArray arrTracks1, TObjArray arrTracks2, const AliVEvent *ev)
2a14a7b1 761{
61d106d3 762 //
5720c765 763 // Prefilter tracks and tracks from pairs
764 // Needed for rejection in the Q-Vector of the event plane
765 // remove contribution of all tracks to the Q-vector that are in invariant mass window
61d106d3 766 //
0caf5fbb 767
5720c765 768 AliEventplane *evplane = const_cast<AliVEvent *>(ev)->GetEventplane();
4533e78e 769 if(!evplane) { // nanoAODs , here we do NOT have sub event reaction planes
770 // if(1) {
771 // get the EPselectionTask for recalculation of weighting factors
772 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
773 AliEPSelectionTask *eptask = dynamic_cast<AliEPSelectionTask *>(man->GetTask("EventplaneSelection"));
774 if(!eptask) return;
775
776 // track mapping
777 TMap mapRemovedTracks;
778
5720c765 779
4533e78e 780 Double_t cQX=0., cQY=0.;
781 // apply cuts to the tracks, e.g. etagap
782 if(fEventPlanePreFilter.GetCuts()->GetEntries()) {
783 UInt_t selectedMask=(1<<fEventPlanePreFilter.GetCuts()->GetEntries())-1;
784 Int_t ntracks=ev->GetNumberOfTracks();
785 for (Int_t itrack=0; itrack<ntracks; ++itrack){
786 AliVParticle *particle=ev->GetTrack(itrack);
787 AliVTrack *track= static_cast<AliVTrack*>(particle);
788 if (!track) continue;
789 //event plane cuts
790 UInt_t cutMask=fEventPlanePreFilter.IsSelected(track);
791 //apply cut
792 if (cutMask==selectedMask) continue;
793
794 mapRemovedTracks.Add(track,track);
795 cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
796 cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
797 }
798 }
5720c765 799
4533e78e 800 // POI (particle of interest) rejection
801 Int_t pairIndex=GetPairIndex(arr1,arr2);
802
803 Int_t ntrack1=arrTracks1.GetEntriesFast();
804 Int_t ntrack2=arrTracks2.GetEntriesFast();
805 AliDielectronPair candidate;
806 candidate.SetKFUsage(fUseKF);
807
808 UInt_t selectedMask=(1<<fEventPlanePOIPreFilter.GetCuts()->GetEntries())-1;
809 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
810 Int_t end=ntrack2;
811 if (arr1==arr2) end=itrack1;
812 Bool_t accepted=kFALSE;
813 for (Int_t itrack2=0; itrack2<end; ++itrack2){
814 TObject *track1=arrTracks1.UncheckedAt(itrack1);
815 TObject *track2=arrTracks2.UncheckedAt(itrack2);
816 if (!track1 || !track2) continue;
817 //create the pair
818 candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
819 static_cast<AliVTrack*>(track2), fPdgLeg2);
820 candidate.SetType(pairIndex);
821 candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
822
823 //event plane pair cuts
824 UInt_t cutMask=fEventPlanePOIPreFilter.IsSelected(&candidate);
825 //apply cut
826 if (cutMask==selectedMask) continue;
827
828 accepted=kTRUE;
829 //remove the tracks from the Track arrays
830 arrTracks2.AddAt(0x0,itrack2);
831 }
832 if ( accepted ) arrTracks1.AddAt(0x0,itrack1);
833 }
834 //compress the track arrays
835 arrTracks1.Compress();
836 arrTracks2.Compress();
837
838 //Modify the components: subtract the tracks
839 ntrack1=arrTracks1.GetEntriesFast();
840 ntrack2=arrTracks2.GetEntriesFast();
841 // remove leg1 contribution
842 for (Int_t itrack=0; itrack<ntrack1; ++itrack){
843 AliVTrack *track= static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack));
5720c765 844 if (!track) continue;
4533e78e 845 // track contribution was already removed
846 if (mapRemovedTracks.FindObject(track)) continue;
847 else mapRemovedTracks.Add(track,track);
5720c765 848
4533e78e 849 cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
850 cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
5720c765 851 }
4533e78e 852 // remove leg2 contribution
853 for (Int_t itrack=0; itrack<ntrack2; ++itrack){
854 AliVTrack *track= static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack));
5720c765 855 if (!track) continue;
4533e78e 856 // track contribution was already removed
857 if (mapRemovedTracks.FindObject(track)) continue;
858 else mapRemovedTracks.Add(track,track);
859
860 cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
861 cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
5720c765 862 }
5720c765 863
4533e78e 864 // build a corrected alieventplane using the values from the var manager
865 // these uncorrected values are filled using the stored magnitude and angle in the header
866 TVector2 qcorr;
867 qcorr.Set(AliDielectronVarManager::GetValue(AliDielectronVarManager::kTPCxH2uc)-cQX,
868 AliDielectronVarManager::GetValue(AliDielectronVarManager::kTPCyH2uc)-cQY);
869 // fill alieventplane
870 AliEventplane cevplane;
871 cevplane.SetQVector(&qcorr);
872 AliDielectronVarManager::SetTPCEventPlane(&cevplane);
873 cevplane.SetQVector(0);
874 return;
875 } //end: nanoAODs
876 else
877 {
878 // this is done in case of ESDs or AODs
879 Bool_t isESD=(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliESDInputHandler::Class());
880 // copy event plane object
881 AliEventplane cevplane(*evplane);
882 // Int_t nMaxID=cevplane->GetQContributionXArray()->GetSize();
883
884 TVector2 *qcorr = cevplane.GetQVector();
885 if(!qcorr) return;
886 TVector2 *qcsub1 = 0x0;
887 TVector2 *qcsub2 = 0x0;
888
889 // eta gap ?
890 Bool_t etagap = kFALSE;
891 for (Int_t iCut=0; iCut<fEventPlanePreFilter.GetCuts()->GetEntries();++iCut) {
892 TString cutName=fEventPlanePreFilter.GetCuts()->At(iCut)->GetName();
893 if(cutName.Contains("eta") || cutName.Contains("Eta")) etagap=kTRUE;
61d106d3 894 }
4533e78e 895
896 // subevent configuration for eta gap or LS (default is rndm)
897 if(fLikeSignSubEvents && etagap) {
898 // start with the full Qvector/event in both sub events
899 qcsub1 = new TVector2(*qcorr);
900 qcsub2 = new TVector2(*qcorr);
901 cevplane.SetQsub(qcsub1,qcsub2);
902
903 Int_t ntracks=ev->GetNumberOfTracks();
904 // track removals
905 for (Int_t itrack=0; itrack<ntracks; ++itrack){
906 AliVParticle *particle=ev->GetTrack(itrack);
907 AliVTrack *track= static_cast<AliVTrack*>(particle);
908 if (!track) continue;
909 if (track->GetID()>=0 && !isESD) continue;
910 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
911
912 // set contributions to zero
913 // charge sub1+ sub2-
914 if(fLikeSignSubEvents) {
915 Short_t charge=track->Charge();
916 if (charge<0) {
917 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
918 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
919 }
920 if (charge>0) {
921 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
922 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
923 }
924 }
925 // eta sub1+ sub2-
926 if(etagap) {
927 Double_t eta=track->Eta();
928 if (eta<0.0) {
929 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
930 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
931 }
932 if (eta>0.0) {
933 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
934 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
935 }
936 }
937 } // end: loop over tracks
938 } // end: sub event configuration
939
940 // apply cuts, e.g. etagap
941 if(fEventPlanePreFilter.GetCuts()->GetEntries()) {
942 UInt_t selectedMask=(1<<fEventPlanePreFilter.GetCuts()->GetEntries())-1;
943 Int_t ntracks=ev->GetNumberOfTracks();
944 for (Int_t itrack=0; itrack<ntracks; ++itrack){
945 AliVParticle *particle=ev->GetTrack(itrack);
946 AliVTrack *track= static_cast<AliVTrack*>(particle);
947 if (!track) continue;
948 if (track->GetID()>=0 && !isESD) continue;
949 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
950
951 //event plane cuts
952 UInt_t cutMask=fEventPlanePreFilter.IsSelected(track);
953 //apply cut
954 if (cutMask==selectedMask) continue;
955
956 // set contributions to zero
957 cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
958 cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
959 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
960 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
961 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
962 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
963 }
964 } // end: track cuts
965
966 // POI (particle of interest) rejection
967 Int_t pairIndex=GetPairIndex(arr1,arr2);
968 Int_t ntrack1=arrTracks1.GetEntriesFast();
969 Int_t ntrack2=arrTracks2.GetEntriesFast();
970 AliDielectronPair candidate;
971 candidate.SetKFUsage(fUseKF);
972
973 UInt_t selectedMask=(1<<fEventPlanePOIPreFilter.GetCuts()->GetEntries())-1;
974 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
975 Int_t end=ntrack2;
976 if (arr1==arr2) end=itrack1;
977 Bool_t accepted=kFALSE;
978 for (Int_t itrack2=0; itrack2<end; ++itrack2){
979 TObject *track1=arrTracks1.UncheckedAt(itrack1);
980 TObject *track2=arrTracks2.UncheckedAt(itrack2);
981 if (!track1 || !track2) continue;
982 //create the pair
983 candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
984 static_cast<AliVTrack*>(track2), fPdgLeg2);
985
986 candidate.SetType(pairIndex);
987 candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
988
989 //event plane cuts
990 UInt_t cutMask=fEventPlanePOIPreFilter.IsSelected(&candidate);
991 //apply cut
992 if (cutMask==selectedMask) continue;
993
994 accepted=kTRUE;
995 //remove the tracks from the Track arrays
996 arrTracks2.AddAt(0x0,itrack2);
997 }
5720c765 998 if ( accepted ) arrTracks1.AddAt(0x0,itrack1);
4533e78e 999 }
1000 //compress the track arrays
1001 arrTracks1.Compress();
1002 arrTracks2.Compress();
1003
1004 //Modify the components: subtract the tracks
1005 ntrack1=arrTracks1.GetEntriesFast();
1006 ntrack2=arrTracks2.GetEntriesFast();
1007 // remove leg1 contribution
1008 for (Int_t itrack=0; itrack<ntrack1; ++itrack){
1009 AliVTrack *track= static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack));
1010 if (!track) continue;
1011 if (track->GetID()>=0 && !isESD) continue;
1012 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
1013 // set contributions to zero
1014 cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
1015 cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
1016 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
1017 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
1018 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
1019 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
1020 }
1021 // remove leg2 contribution
1022 for (Int_t itrack=0; itrack<ntrack2; ++itrack){
1023 AliVTrack *track= static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack));
1024 if (!track) continue;
1025 if (track->GetID()>=0 && !isESD) continue;
1026 Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
1027 // set contributions to zero
1028 cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
1029 cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
1030 cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
1031 cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
1032 cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
1033 cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
1034 }
554e40f8 1035
4533e78e 1036 // set corrected AliEventplane and fill variables with corrected values
1037 AliDielectronVarManager::SetTPCEventPlane(&cevplane);
1038 delete qcsub1;
1039 delete qcsub2;
1040 } // end: ESD or AOD case
554e40f8 1041
4533e78e 1042}
5720c765 1043//________________________________________________________________
1044void AliDielectron::PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2)
1045{
1046 //
1047 // Prefilter tracks from pairs
1048 // Needed for datlitz rejections
1049 // remove all tracks from the Single track arrays that pass the cuts in this filter
1050 //
1051
1052 Int_t ntrack1=arrTracks1.GetEntriesFast();
1053 Int_t ntrack2=arrTracks2.GetEntriesFast();
1054 AliDielectronPair candidate;
08b801a6 1055 candidate.SetKFUsage(fUseKF);
5720c765 1056 // flag arrays for track removal
1057 Bool_t *bTracks1 = new Bool_t[ntrack1];
1058 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1) bTracks1[itrack1]=kFALSE;
1059 Bool_t *bTracks2 = new Bool_t[ntrack2];
1060 for (Int_t itrack2=0; itrack2<ntrack2; ++itrack2) bTracks2[itrack2]=kFALSE;
1061
1062 UInt_t selectedMask=(1<<fPairPreFilter.GetCuts()->GetEntries())-1;
1063 UInt_t selectedMaskPair=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1064
1065 Int_t nRejPasses = 1; //for fPreFilterUnlikeOnly and no set flag
1066 if (fPreFilterAllSigns) nRejPasses = 3;
1067
1068 for (Int_t iRP=0; iRP < nRejPasses; ++iRP) {
1069 Int_t arr1RP=arr1, arr2RP=arr2;
1070 TObjArray *arrTracks1RP=&arrTracks1;
1071 TObjArray *arrTracks2RP=&arrTracks2;
1072 Bool_t *bTracks1RP = bTracks1;
1073 Bool_t *bTracks2RP = bTracks2;
1074 switch (iRP) {
1075 case 1: arr1RP=arr1;arr2RP=arr1;
1076 arrTracks1RP=&arrTracks1;
1077 arrTracks2RP=&arrTracks1;
1078 bTracks1RP = bTracks1;
1079 bTracks2RP = bTracks1;
1080 break;
1081 case 2: arr1RP=arr2;arr2RP=arr2;
1082 arrTracks1RP=&arrTracks2;
1083 arrTracks2RP=&arrTracks2;
1084 bTracks1RP = bTracks2;
1085 bTracks2RP = bTracks2;
1086 break;
1087 default: ;//nothing to do
1088 }
1089 Int_t ntrack1RP=(*arrTracks1RP).GetEntriesFast();
1090 Int_t ntrack2RP=(*arrTracks2RP).GetEntriesFast();
1091
1092 Int_t pairIndex=GetPairIndex(arr1RP,arr2RP);
1093
1094 for (Int_t itrack1=0; itrack1<ntrack1RP; ++itrack1){
1095 Int_t end=ntrack2RP;
1096 if (arr1RP==arr2RP) end=itrack1;
1097 for (Int_t itrack2=0; itrack2<end; ++itrack2){
1098 TObject *track1=(*arrTracks1RP).UncheckedAt(itrack1);
1099 TObject *track2=(*arrTracks2RP).UncheckedAt(itrack2);
1100 if (!track1 || !track2) continue;
1101 //create the pair
1102 candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
1103 static_cast<AliVTrack*>(track2), fPdgLeg2);
1104
1105 candidate.SetType(pairIndex);
1106 candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
1107 //relate to the production vertex
1108 // if (AliDielectronVarManager::GetKFVertex()) candidate.SetProductionVertex(*AliDielectronVarManager::GetKFVertex());
1109
1110 //pair cuts
1111 UInt_t cutMask=fPairPreFilter.IsSelected(&candidate);
1112
1113 //apply cut
1114 if (cutMask!=selectedMask) continue;
1115 if (fCfManagerPair) fCfManagerPair->Fill(selectedMaskPair+1 ,&candidate);
1116 if (fHistos) FillHistogramsPair(&candidate,kTRUE);
1117 //set flags for track removal
1118 bTracks1RP[itrack1]=kTRUE;
1119 bTracks2RP[itrack2]=kTRUE;
1120 }
1121 }
1122 }
1123
1124 //remove the tracks from the Track arrays
1125 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
1126 if(bTracks1[itrack1]) arrTracks1.AddAt(0x0, itrack1);
1127 }
1128 for (Int_t itrack2=0; itrack2<ntrack2; ++itrack2){
1129 if(bTracks2[itrack2]) arrTracks2.AddAt(0x0, itrack2);
1130 }
1131
1132 // clean up
1133 delete [] bTracks1;
1134 delete [] bTracks2;
1135
1136 //compress the track arrays
61d106d3 1137 arrTracks1.Compress();
1138 arrTracks2.Compress();
1139
1140 //apply leg cuts after the pre filter
1141 if ( fPairPreFilterLegs.GetCuts()->GetEntries()>0 ) {
1142 selectedMask=(1<<fPairPreFilterLegs.GetCuts()->GetEntries())-1;
1143 //loop over tracks from array 1
1144 for (Int_t itrack=0; itrack<arrTracks1.GetEntriesFast();++itrack){
1145 //test cuts
1146 UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks1.UncheckedAt(itrack));
554e40f8 1147
61d106d3 1148 //apply cut
d4619886 1149 if (cutMask!=selectedMask) arrTracks1.AddAt(0x0,itrack);
61d106d3 1150 }
1151 arrTracks1.Compress();
1152
1153 //in case of like sign don't loop over second array
1154 if (arr1==arr2) {
1155 arrTracks2=arrTracks1;
1156 } else {
1157
1158 //loop over tracks from array 2
1159 for (Int_t itrack=0; itrack<arrTracks2.GetEntriesFast();++itrack){
1160 //test cuts
1161 UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks2.UncheckedAt(itrack));
1162 //apply cut
1163 if (cutMask!=selectedMask) arrTracks2.AddAt(0x0,itrack);
1164 }
1165 arrTracks2.Compress();
1166
1167 }
1168 }
554e40f8 1169 //For unlike-sign monitor track-cuts:
ba15fdfb 1170 if (arr1!=arr2&&fHistos) {
554e40f8 1171 TObjArray *unlikesignArray[2] = {&arrTracks1,&arrTracks2};
1172 FillHistogramsTracks(unlikesignArray);
1173 }
61d106d3 1174}
1175
b2a297fa 1176//________________________________________________________________
2a14a7b1 1177void AliDielectron::FillPairArrays(Int_t arr1, Int_t arr2)
1178{
b2a297fa 1179 //
1180 // select pairs and fill pair candidate arrays
1181 //
61d106d3 1182
1183 TObjArray arrTracks1=fTracks[arr1];
1184 TObjArray arrTracks2=fTracks[arr2];
1185
1186 //process pre filter if set
ba15fdfb 1187 if ((!fPreFilterAllSigns) && (!fPreFilterUnlikeOnly) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(arr1, arr2, arrTracks1, arrTracks2);
61d106d3 1188
b2a297fa 1189 Int_t pairIndex=GetPairIndex(arr1,arr2);
1190
61d106d3 1191 Int_t ntrack1=arrTracks1.GetEntriesFast();
1192 Int_t ntrack2=arrTracks2.GetEntriesFast();
b2a297fa 1193
1194 AliDielectronPair *candidate=new AliDielectronPair;
08b801a6 1195 candidate->SetKFUsage(fUseKF);
91a4a151 1196 // switch OFF the KF usage in case of AODs && ME (since there is no MoveToSameVertex functionality)
1197 Bool_t isESD=(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliESDInputHandler::Class());
1198 if(!isESD && pairIndex>AliDielectron::kEv1MM) candidate->SetKFUsage(kFALSE);
b2a297fa 1199
1200 UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1201
1202 for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
1203 Int_t end=ntrack2;
1204 if (arr1==arr2) end=itrack1;
1205 for (Int_t itrack2=0; itrack2<end; ++itrack2){
7070f483 1206 //create the pair (direct pointer to the memory by this daughter reference are kept also for ME)
1207 candidate->SetTracks(&(*static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1))), fPdgLeg1,
1208 &(*static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2))), fPdgLeg2);
b2a297fa 1209 candidate->SetType(pairIndex);
7070f483 1210
e4339752 1211 Int_t label=AliDielectronMC::Instance()->GetLabelMotherWithPdg(candidate,fPdgMother);
1212 candidate->SetLabel(label);
1213 if (label>-1) candidate->SetPdgCode(fPdgMother);
88db3193 1214 else candidate->SetPdgCode(0);
b2a297fa 1215
4fae8ef9 1216 // check for gamma kf particle
1217 label=AliDielectronMC::Instance()->GetLabelMotherWithPdg(candidate,22);
1218 if (label>-1) {
1219 candidate->SetGammaTracks(static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1)), fPdgLeg1,
1220 static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2)), fPdgLeg2);
1221 // should we set the pdgmothercode and the label
1222 }
1223
b2a297fa 1224 //pair cuts
1225 UInt_t cutMask=fPairFilter.IsSelected(candidate);
d4619886 1226
b2a297fa 1227 //CF manager for the pair
1228 if (fCfManagerPair) fCfManagerPair->Fill(cutMask,candidate);
d4619886 1229
d8bb1abe 1230 // cut qa
dd718b9c 1231 if(pairIndex==kEv1PM && fCutQA) {
d8bb1abe 1232 fQAmonitor->FillAll(candidate);
1233 fQAmonitor->Fill(cutMask,candidate);
1234 }
b2a297fa 1235
1236 //apply cut
1237 if (cutMask!=selectedMask) continue;
1238
d4619886 1239 //histogram array for the pair
1240 if (fHistoArray) fHistoArray->Fill(pairIndex,candidate);
1241
b2a297fa 1242 //add the candidate to the candidate array
1243 PairArray(pairIndex)->Add(candidate);
1244 //get a new candidate
1245 candidate=new AliDielectronPair;
08b801a6 1246 candidate->SetKFUsage(fUseKF);
b2a297fa 1247 }
1248 }
1249 //delete the surplus candidate
1250 delete candidate;
1251}
1252
2a14a7b1 1253//________________________________________________________________
1254void AliDielectron::FillPairArrayTR()
1255{
1256 //
1257 // select pairs and fill pair candidate arrays
1258 //
1259 UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1260
1261 while ( fTrackRotator->NextCombination() ){
1262 AliDielectronPair candidate;
08b801a6 1263 candidate.SetKFUsage(fUseKF);
1201a1a9 1264 candidate.SetTracks(&fTrackRotator->GetKFTrackP(), &fTrackRotator->GetKFTrackN(),
1265 fTrackRotator->GetVTrackP(),fTrackRotator->GetVTrackN());
ffbede40 1266 candidate.SetType(kEv1PMRot);
2a14a7b1 1267
1268 //pair cuts
1269 UInt_t cutMask=fPairFilter.IsSelected(&candidate);
1270
1271 //CF manager for the pair
1272 if (fCfManagerPair) fCfManagerPair->Fill(cutMask,&candidate);
1273
1274 //apply cut
5720c765 1275 if (cutMask==selectedMask) {
d4619886 1276
1277 //histogram array for the pair
1278 if (fHistoArray) fHistoArray->Fill((Int_t)kEv1PMRot,&candidate);
1279
1280 if(fHistos) FillHistogramsPair(&candidate);
1281 if(fStoreRotatedPairs) PairArray(kEv1PMRot)->Add(new AliDielectronPair(candidate));
1282 }
2a14a7b1 1283 }
1284}
1285
572b0139 1286//________________________________________________________________
1287void AliDielectron::FillDebugTree()
1288{
1289 //
1290 // Fill Histogram information for tracks and pairs
1291 //
1292
1293 //Fill Debug tree
1294 for (Int_t i=0; i<10; ++i){
1295 Int_t ntracks=PairArray(i)->GetEntriesFast();
1296 for (Int_t ipair=0; ipair<ntracks; ++ipair){
1297 fDebugTree->Fill(static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair)));
1298 }
1299 }
1300}
b2a297fa 1301
572b0139 1302//________________________________________________________________
1303void AliDielectron::SaveDebugTree()
1304{
1305 //
1306 // delete the debug tree, this will also write the tree
1307 //
1308 if (fDebugTree) fDebugTree->DeleteStreamer();
1309}
b2a297fa 1310
ba15fdfb 1311
1312//__________________________________________________________________
1313void AliDielectron::AddSignalMC(AliDielectronSignalMC* signal) {
1314 //
1315 // Add an MC signal to the signals list
1316 //
1317 if(!fSignalsMC) {
1318 fSignalsMC = new TObjArray();
1319 fSignalsMC->SetOwner();
1320 }
1321 fSignalsMC->Add(signal);
1322}
833da63d 1323
1324//________________________________________________________________
1325void AliDielectron::FillMCHistograms(Int_t label1, Int_t label2, Int_t nSignal) {
1326 //
1327 // fill QA MC TRUTH histograms for pairs and legs of all added mc signals
1328 //
1329
1330 TString className,className2,className3;
1331 className.Form("Pair_%s_MCtruth",fSignalsMC->At(nSignal)->GetName());
1332 className2.Form("Track_Legs_%s_MCtruth",fSignalsMC->At(nSignal)->GetName());
1333 className3.Form("Track_%s_%s_MCtruth",fgkPairClassNames[1],fSignalsMC->At(nSignal)->GetName());
1334 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
1335 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
1336 Bool_t trkClass=fHistos->GetHistogramList()->FindObject(className3.Data())!=0x0;
bc3cacd4 1337 // printf("fill signal %d: pair %d legs %d trk %d \n",nSignal,pairClass,legClass,trkClass);
833da63d 1338 if(!pairClass && !legClass && !trkClass) return;
1339
1340 // printf("leg labels: %d-%d \n",label1,label2);
1341 AliVParticle* part1 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label1);
1342 AliVParticle* part2 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label2);
1343 if(!part1 && !part2) return;
1344 if(part1&&part2) {
1345 // fill only unlike sign (and only SE)
833da63d 1346 if(part1->Charge()*part2->Charge()>=0) return;
1347 }
1348
1349
1350 AliDielectronMC* dieMC = AliDielectronMC::Instance();
1351
1352 Int_t mLabel1 = dieMC->GetMothersLabel(label1); // should work for both ESD and AOD
1353 Int_t mLabel2 = dieMC->GetMothersLabel(label2);
1354
1355 // check the same mother option
1356 AliDielectronSignalMC* sigMC = (AliDielectronSignalMC*)fSignalsMC->At(nSignal);
1357 if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kSame && mLabel1!=mLabel2) return;
1358 if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kDifferent && mLabel1==mLabel2) return;
1359
bc3cacd4 1360 // fill event values
833da63d 1361 Double_t values[AliDielectronVarManager::kNMaxValues];
a94c2e7e 1362 AliDielectronVarManager::SetFillMap(fUsedVars);
833da63d 1363 AliDielectronVarManager::Fill(dieMC->GetMCEvent(), values); // get event informations
bc3cacd4 1364
1365 // fill the leg variables
1366 // printf("leg:%d trk:%d part1:%p part2:%p \n",legClass,trkClass,part1,part2);
833da63d 1367 if (legClass || trkClass) {
1368 if(part1) AliDielectronVarManager::Fill(part1,values);
1369 if(part1 && trkClass) fHistos->FillClass(className3, AliDielectronVarManager::kNMaxValues, values);
1370 if(part1 && part2 && legClass) fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1371 if(part2) AliDielectronVarManager::Fill(part2,values);
1372 if(part2 && trkClass) fHistos->FillClass(className3, AliDielectronVarManager::kNMaxValues, values);
1373 if(part1 && part2 && legClass) fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1374 }
1375
1376 //fill pair information
1377 if (pairClass && part1 && part2) {
1378 AliDielectronVarManager::FillVarMCParticle2(part1,part2,values);
1379 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
1380 }
1381
1382}
1383
f0351f9c 1384//________________________________________________________________
1385void AliDielectron::FillMCHistograms(const AliVEvent *ev) {
1386 //
1387 // fill QA MC histograms for pairs and legs of all added mc signals
1388 //
1389
4c06df48 1390 if (!fSignalsMC) return;
833da63d 1391 TString className,className2,className3;
f0351f9c 1392 Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
a94c2e7e 1393 AliDielectronVarManager::SetFillMap(fUsedVars);
f0351f9c 1394 AliDielectronVarManager::Fill(ev, values); // get event informations
1395 //loop over all added mc signals
1396 for(Int_t isig=0; isig<fSignalsMC->GetEntries(); isig++) {
1397
bc3cacd4 1398 //check if and what to fill
f0351f9c 1399 className.Form("Pair_%s",fSignalsMC->At(isig)->GetName());
1400 className2.Form("Track_Legs_%s",fSignalsMC->At(isig)->GetName());
536a997a 1401 className3.Form("Track_%s_%s",fgkPairClassNames[1],fSignalsMC->At(isig)->GetName()); // unlike sign, SE only
f0351f9c 1402 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
1403 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
536a997a 1404 Bool_t mergedtrkClass=fHistos->GetHistogramList()->FindObject(className3.Data())!=0x0;
1405 if(!pairClass && !legClass && !mergedtrkClass) continue;
1406
1407 // fill pair and/or their leg variables
1408 if(pairClass || legClass) {
1409 Int_t npairs=PairArray(AliDielectron::kEv1PM)->GetEntriesFast(); // only SE +-
1410 for (Int_t ipair=0; ipair<npairs; ++ipair){
1411 AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(AliDielectron::kEv1PM)->UncheckedAt(ipair));
1412
1413 Bool_t isMCtruth = AliDielectronMC::Instance()->IsMCTruth(pair, (AliDielectronSignalMC*)fSignalsMC->At(isig));
1414 if(isMCtruth) {
1415 //fill pair information
1416 if (pairClass){
1417 AliDielectronVarManager::Fill(pair, values);
1418 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
1419 }
1420 //fill leg information, both + and - in the same histo
1421 if (legClass){
54ce9dc4 1422 AliDielectronVarManager::Fill(pair->GetFirstDaughterP(),values);
536a997a 1423 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
54ce9dc4 1424 AliDielectronVarManager::Fill(pair->GetSecondDaughterP(),values);
536a997a 1425 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1426 }
1427 } //is signal
1428 } //loop: pairs
1429 }
833da63d 1430
1431 // fill single tracks of signals
bc3cacd4 1432 if(!mergedtrkClass) continue;
1433 // loop over SE track arrays
1434 for (Int_t i=0; i<2; ++i){
1435 Int_t ntracks=fTracks[i].GetEntriesFast();
1436 for (Int_t itrack=0; itrack<ntracks; ++itrack){
1437 Int_t label=((AliVParticle*)fTracks[i].UncheckedAt(itrack))->GetLabel();
1438 Bool_t isMCtruth1 = AliDielectronMC::Instance()->IsMCTruth(label, (AliDielectronSignalMC*)fSignalsMC->At(isig), 1);
1439 Bool_t isMCtruth2 = AliDielectronMC::Instance()->IsMCTruth(label, (AliDielectronSignalMC*)fSignalsMC->At(isig), 2);
1440 // skip if track does not correspond to the signal
1441 if(!isMCtruth1 && !isMCtruth2) continue;
1442 AliDielectronVarManager::Fill(fTracks[i].UncheckedAt(itrack), values);
1443 fHistos->FillClass(className3, AliDielectronVarManager::kNMaxValues, values);
536a997a 1444 } //loop: tracks
bc3cacd4 1445 } //loop: arrays
833da63d 1446
f0351f9c 1447 } //loop: MCsignals
1448
1449}
30f9a393
JB		 1450
1451//______________________________________________
1452void AliDielectron::SetCentroidCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
1453{
7238fa88 1454 UInt_t valType[20] = {0};
1455 valType[0]=varx; valType[1]=vary; valType[2]=varz;
1456 AliDielectronHistos::StoreVariables(fun->GetHistogram(), valType);
1457 // clone temporare histogram, otherwise it will not be streamed to file!
9ad9e048 1458 TString key = Form("cntrd%d%d%d",varx,vary,varz);
1459 fPostPIDCntrdCorr = (TH1*)fun->GetHistogram()->Clone(key.Data());
3f15dc00 1460 if(fPostPIDCntrdCorr) {
aeedbade 1461 fPostPIDCntrdCorr->GetListOfFunctions()->AddAt(fun,0);
1462 // check for corrections and add their variables to the fill map
3f15dc00 1463 printf("POST TPC PID CORRECTION added for centroids: ");
1464 switch(fPostPIDCntrdCorr->GetDimension()) {
1465 case 3: printf(" %s, ",fPostPIDCntrdCorr->GetZaxis()->GetName());
1466 case 2: printf(" %s, ",fPostPIDCntrdCorr->GetYaxis()->GetName());
1467 case 1: printf(" %s ",fPostPIDCntrdCorr->GetXaxis()->GetName());
1468 }
1469 printf("\n");
1470 fUsedVars->SetBitNumber(varx, kTRUE);
1471 fUsedVars->SetBitNumber(vary, kTRUE);
1472 fUsedVars->SetBitNumber(varz, kTRUE);
1473 }
1474
1475}
1476//______________________________________________
1477void AliDielectron::SetCentroidCorrFunction(TH1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
1478{
1479 UInt_t valType[20] = {0};
1480 valType[0]=varx; valType[1]=vary; valType[2]=varz;
1481 AliDielectronHistos::StoreVariables(fun, valType);
1482 // clone temporare histogram, otherwise it will not be streamed to file!
1483 TString key = Form("cntrd%d%d%d",varx,vary,varz);
1484 fPostPIDCntrdCorr = (TH1*)fun->Clone(key.Data());
7238fa88 1485 // check for corrections and add their variables to the fill map
1486 if(fPostPIDCntrdCorr) {
1487 printf("POST TPC PID CORRECTION added for centroids: ");
1488 switch(fPostPIDCntrdCorr->GetDimension()) {
1489 case 3: printf(" %s, ",fPostPIDCntrdCorr->GetZaxis()->GetName());
1490 case 2: printf(" %s, ",fPostPIDCntrdCorr->GetYaxis()->GetName());
1491 case 1: printf(" %s ",fPostPIDCntrdCorr->GetXaxis()->GetName());
1492 }
1493 printf("\n");
1494 fUsedVars->SetBitNumber(varx, kTRUE);
1495 fUsedVars->SetBitNumber(vary, kTRUE);
1496 fUsedVars->SetBitNumber(varz, kTRUE);
1497 }
1498
30f9a393
JB		 1499}
1500//______________________________________________
1501void AliDielectron::SetWidthCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
1502{
7238fa88 1503 UInt_t valType[20] = {0};
1504 valType[0]=varx; valType[1]=vary; valType[2]=varz;
1505 AliDielectronHistos::StoreVariables(fun->GetHistogram(), valType);
1506 // clone temporare histogram, otherwise it will not be streamed to file!
9ad9e048 1507 TString key = Form("wdth%d%d%d",varx,vary,varz);
1508 fPostPIDWdthCorr = (TH1*)fun->GetHistogram()->Clone(key.Data());
7238fa88 1509 if(fPostPIDWdthCorr) {
aeedbade 1510 fPostPIDWdthCorr->GetListOfFunctions()->AddAt(fun,0);
1511 // check for corrections and add their variables to the fill map
7238fa88 1512 printf("POST TPC PID CORRECTION added for widths: ");
1513 switch(fPostPIDWdthCorr->GetDimension()) {
1514 case 3: printf(" %s, ",fPostPIDWdthCorr->GetZaxis()->GetName());
1515 case 2: printf(" %s, ",fPostPIDWdthCorr->GetYaxis()->GetName());
1516 case 1: printf(" %s ",fPostPIDWdthCorr->GetXaxis()->GetName());
1517 }
1518 printf("\n");
1519 fUsedVars->SetBitNumber(varx, kTRUE);
1520 fUsedVars->SetBitNumber(vary, kTRUE);
1521 fUsedVars->SetBitNumber(varz, kTRUE);
1522 }
30f9a393 1523}
3f15dc00 1524//______________________________________________
1525void AliDielectron::SetWidthCorrFunction(TH1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
1526{
1527 UInt_t valType[20] = {0};
1528 valType[0]=varx; valType[1]=vary; valType[2]=varz;
1529 AliDielectronHistos::StoreVariables(fun, valType);
1530 // clone temporare histogram, otherwise it will not be streamed to file!
1531 TString key = Form("cntrd%d%d%d",varx,vary,varz);
1532 fPostPIDWdthCorr = (TH1*)fun->Clone(key.Data());
1533 // check for corrections and add their variables to the fill map
1534 if(fPostPIDWdthCorr) {
1535 printf("POST TPC PID CORRECTION added for widths: ");
1536 switch(fPostPIDWdthCorr->GetDimension()) {
1537 case 3: printf(" %s, ",fPostPIDWdthCorr->GetZaxis()->GetName());
1538 case 2: printf(" %s, ",fPostPIDWdthCorr->GetYaxis()->GetName());
1539 case 1: printf(" %s ",fPostPIDWdthCorr->GetXaxis()->GetName());
1540 }
1541 printf("\n");
1542 fUsedVars->SetBitNumber(varx, kTRUE);
1543 fUsedVars->SetBitNumber(vary, kTRUE);
1544 fUsedVars->SetBitNumber(varz, kTRUE);
1545 }
1546
1547}
9ad9e048 1548
7965afb6 1549//______________________________________________
746e66a9 1550TObject* AliDielectron::InitEffMap(TString filename)
7965afb6 1551{
1552 // init an efficiency object for on-the-fly correction calculations
7965afb6 1553 if(filename.Contains("alien://") && !gGrid) TGrid::Connect("alien://",0,0,"t");
1554
1555 TFile* file=TFile::Open(filename.Data());
0b369e55 1556 if(!file) return 0x0;
c9cd45da 1557 else printf("[I] AliDielectron::InitEffMap efficiency maps %s loaded! \n",filename.Data());
746e66a9 1558
1559 // NOTE: the spline must have the 'variable name' stored in its fHistogram
1560 TSpline3 *hEff = (TSpline3*) file->Get("hEfficiency");
1561 //if(hEff) printf("we use a TSpline!!!!!!!!!!! \n");
1562 if(hEff) return (hEff->Clone("effMap"));
1563
1564
7965afb6 1565 THnBase *hGen = (THnBase*) file->Get("hGenerated");
1566 THnBase *hFnd = (THnBase*) file->Get("hFound");
0b369e55 1567 if(!hFnd || !hGen) return 0x0;
7965afb6 1568
0b369e55 1569 hFnd->Divide(hGen);
746e66a9 1570 return (hFnd->Clone("effMap"));
7965afb6 1571}
42c76bea 1572
1573//________________________________________________________________
1574void AliDielectron::FillHistogramsFromPairArray(Bool_t pairInfoOnly/*=kFALSE*/)
1575{
1576 //
1577 // Fill Histogram information for tracks and pairs
1578 //
1579
1580 TString className,className2;
1581 Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
e4942d79 1582 AliDielectronVarManager::SetFillMap(fUsedVars);
1583 AliDielectronVarManager::SetLegEffMap(fLegEffMap);
1584 AliDielectronVarManager::SetPairEffMap(fPairEffMap);
42c76bea 1585
1586 //Fill event information
1587 if(!pairInfoOnly) {
1588 if(fHistos->GetHistogramList()->FindObject("Event")) {
1589 fHistos->FillClass("Event", AliDielectronVarManager::kNMaxValues, AliDielectronVarManager::GetData());
1590 }
1591 }
1592
1593 UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1594
1595 //Fill Pair information, separately for all pair candidate arrays and the legs
1596 TObjArray arrLegs(100);
1597 for (Int_t i=0; i<10; ++i){ // ROT pairs??
1598 Int_t npairs=PairArray(i)->GetEntriesFast();
1599 if(npairs<1) continue;
1600
1601 className.Form("Pair_%s",fgkPairClassNames[i]);
1602 className2.Form("Track_Legs_%s",fgkPairClassNames[i]);
1603 Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
1604 Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
1605
1606 // if (!pairClass&&!legClass) continue;
1607 for (Int_t ipair=0; ipair<npairs; ++ipair){
1608 AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair));
1609
1610 // apply cuts
1611 UInt_t cutMask=fPairFilter.IsSelected(pair);
1612
1613 // cut qa
1614 if(i==kEv1PM && fCutQA) {
1615 fQAmonitor->FillAll(pair);
1616 fQAmonitor->Fill(cutMask,pair);
1617 }
1618
1619 //CF manager for the pair (TODO: check steps and if they are properly filled)
1620 // if (fCfManagerPair) fCfManagerPair->Fill(cutMask,pair);
1621
1622 //apply cut
1623 if (cutMask!=selectedMask) continue;
1624
1625 //histogram array for the pair
1626 if (fHistoArray) fHistoArray->Fill(i,pair);
1627
e17c0bd6 1628 // fill map
1629 AliDielectronVarManager::SetFillMap(fUsedVars);
1630
42c76bea 1631 //fill pair information
1632 if (pairClass){
1633 AliDielectronVarManager::Fill(pair, values);
1634 fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
1635 }
1636
1637 //fill leg information, don't fill the information twice
1638 if (legClass){
54ce9dc4 1639 AliVParticle *d1=pair->GetFirstDaughterP();
1640 AliVParticle *d2=pair->GetSecondDaughterP();
42c76bea 1641 if (!arrLegs.FindObject(d1)){
1642 AliDielectronVarManager::Fill(d1, values);
1643 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1644 arrLegs.Add(d1);
1645 }
1646 if (!arrLegs.FindObject(d2)){
1647 AliDielectronVarManager::Fill(d2, values);
1648 fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1649 arrLegs.Add(d2);
1650 }
1651 }
1652 }
1653 if (legClass) arrLegs.Clear();
1654 }
1655
1656}
1657