86946203d73873bfbab0f742c14fe8e068dd88c9
[u/mrichter/AliRoot.git] / PWGDQ / dielectron / AliDielectron.cxx
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 /*
20 Framework to perform event selectoin, single track selection and track pair
21 selection.
22
23 Convention for the signs of the pair in fPairCandidates:
24 The names are available via the function PairClassName(Int_t i)
25
26 0: ev1+ ev1+  (same event like sign +)
27 1: ev1+ ev1-  (same event unlike sign)
28 2: ev1- ev1-  (same event like sign -)
29
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 -)
34
35 5: ev2+ ev2+  (same event like sign +)
36 8: ev2+ ev2-  (same event unlike sign)
37 9: ev2- ev2-  (same event like sign -)
38
39 10: ev1+ ev1- (same event track rotation)
40
41 */
42 //                                                                       //
43 ///////////////////////////////////////////////////////////////////////////
44
45 #include <TString.h>
46 #include <TList.h>
47 #include <TMath.h>
48 #include <TObject.h>
49
50 #include <AliKFParticle.h>
51
52 #include <AliESDInputHandler.h>
53 #include <AliAnalysisManager.h>
54 #include <AliEPSelectionTask.h>
55 #include <AliEventplane.h>
56 #include <AliVEvent.h>
57 #include <AliVParticle.h>
58 #include <AliVTrack.h>
59 #include "AliDielectronPair.h"
60 #include "AliDielectronHistos.h"
61 #include "AliDielectronCF.h"
62 #include "AliDielectronMC.h"
63 #include "AliDielectronVarManager.h"
64 #include "AliDielectronTrackRotator.h"
65 #include "AliDielectronDebugTree.h"
66 #include "AliDielectronSignalMC.h"
67 #include "AliDielectronMixingHandler.h"
68 #include "AliDielectronV0Cuts.h"
69
70 #include "AliDielectron.h"
71
72 ClassImp(AliDielectron)
73
74 const char* AliDielectron::fgkTrackClassNames[4] = {
75   "ev1+",
76   "ev1-",
77   "ev2+",
78   "ev2-"
79 };
80
81 const char* AliDielectron::fgkPairClassNames[11] = {
82   "ev1+_ev1+",
83   "ev1+_ev1-",
84   "ev1-_ev1-",
85   "ev1+_ev2+",
86   "ev1-_ev2+",
87   "ev2+_ev2+",
88   "ev1+_ev2-",
89   "ev1-_ev2-",
90   "ev2+_ev2-",
91   "ev2-_ev2-",
92   "ev1+_ev1-_TR"
93 };
94
95 //________________________________________________________________
96 AliDielectron::AliDielectron() :
97   TNamed("AliDielectron","AliDielectron"),
98   fEventFilter("EventFilter"),
99   fTrackFilter("TrackFilter"),
100   fPairPreFilter("PairPreFilter"),
101   fPairPreFilterLegs("PairPreFilterLegs"),
102   fPairFilter("PairFilter"),
103   fEventPlanePreFilter("EventPlanePreFilter"),
104   fEventPlanePOIPreFilter("EventPlanePOIPreFilter"),
105   fPdgMother(443),
106   fPdgLeg1(11),
107   fPdgLeg2(11),
108   fSignalsMC(0x0),
109   fNoPairing(kFALSE),
110   fUseKF(kTRUE),
111   fHistoArray(0x0),
112   fHistos(0x0),
113   fPairCandidates(new TObjArray(11)),
114   fCfManagerPair(0x0),
115   fTrackRotator(0x0),
116   fDebugTree(0x0),
117   fMixing(0x0),
118   fPreFilterEventPlane(kFALSE),
119   fLikeSignSubEvents(kFALSE),
120   fPreFilterUnlikeOnly(kFALSE),
121   fPreFilterAllSigns(kFALSE),
122   fHasMC(kFALSE),
123   fStoreRotatedPairs(kFALSE),
124   fDontClearArrays(kFALSE),
125   fEstimatorFilename(""),
126   fTRDpidCorrectionFilename(""),
127   fVZEROCalibrationFilename(""),
128   fVZERORecenteringFilename("")
129 {
130   //
131   // Default constructor
132   //
133
134 }
135
136 //________________________________________________________________
137 AliDielectron::AliDielectron(const char* name, const char* title) :
138   TNamed(name,title),
139   fEventFilter("EventFilter"),
140   fTrackFilter("TrackFilter"),
141   fPairPreFilter("PairPreFilter"),
142   fPairPreFilterLegs("PairPreFilterLegs"),
143   fPairFilter("PairFilter"),
144   fEventPlanePreFilter("EventPlanePreFilter"),
145   fEventPlanePOIPreFilter("EventPlanePOIPreFilter"),
146   fPdgMother(443),
147   fPdgLeg1(11),
148   fPdgLeg2(11),
149   fSignalsMC(0x0),
150   fNoPairing(kFALSE),
151   fUseKF(kTRUE),
152   fHistoArray(0x0),
153   fHistos(0x0),
154   fPairCandidates(new TObjArray(11)),
155   fCfManagerPair(0x0),
156   fTrackRotator(0x0),
157   fDebugTree(0x0),
158   fMixing(0x0),
159   fPreFilterEventPlane(kFALSE),
160   fLikeSignSubEvents(kFALSE),
161   fPreFilterUnlikeOnly(kFALSE),
162   fPreFilterAllSigns(kFALSE),
163   fHasMC(kFALSE),
164   fStoreRotatedPairs(kFALSE),
165   fDontClearArrays(kFALSE),
166   fEstimatorFilename(""),
167   fTRDpidCorrectionFilename(""),
168   fVZEROCalibrationFilename(""),
169   fVZERORecenteringFilename("")
170 {
171   //
172   // Named constructor
173   //
174   
175 }
176
177 //________________________________________________________________
178 AliDielectron::~AliDielectron()
179 {
180   //
181   // Default destructor
182   //
183   if (fHistoArray) delete fHistoArray;
184   if (fHistos) delete fHistos;
185   if (fPairCandidates) delete fPairCandidates;
186   if (fDebugTree) delete fDebugTree;
187   if (fMixing) delete fMixing;
188   if (fSignalsMC) delete fSignalsMC;
189   if (fCfManagerPair) delete fCfManagerPair;
190 }
191
192 //________________________________________________________________
193 void AliDielectron::Init()
194 {
195   //
196   // Initialise objects
197   //
198
199   if(GetHasMC()) AliDielectronMC::Instance()->SetHasMC(GetHasMC());
200    
201   InitPairCandidateArrays();
202    
203   if (fCfManagerPair) {
204     fCfManagerPair->SetSignalsMC(fSignalsMC);
205     fCfManagerPair->InitialiseContainer(fPairFilter);
206   }
207   if (fTrackRotator)  {
208     fTrackRotator->SetTrackArrays(&fTracks[0],&fTracks[1]);
209     fTrackRotator->SetPdgLegs(fPdgLeg1,fPdgLeg2);
210   }
211   if (fDebugTree) fDebugTree->SetDielectron(this);
212   if(fEstimatorFilename.Contains(".root")) AliDielectronVarManager::InitEstimatorAvg(fEstimatorFilename.Data());
213   if(fTRDpidCorrectionFilename.Contains(".root")) AliDielectronVarManager::InitTRDpidEffHistograms(fTRDpidCorrectionFilename.Data());
214   if(fVZEROCalibrationFilename.Contains(".root")) AliDielectronVarManager::SetVZEROCalibrationFile(fVZEROCalibrationFilename.Data());
215   if(fVZERORecenteringFilename.Contains(".root")) AliDielectronVarManager::SetVZERORecenteringFile(fVZERORecenteringFilename.Data());
216   
217   if (fMixing) fMixing->Init(this);
218   if (fHistoArray) {
219     fHistoArray->SetSignalsMC(fSignalsMC);
220     fHistoArray->Init();
221   }
222
223
224 //________________________________________________________________
225 void AliDielectron::Process(AliVEvent *ev1, AliVEvent *ev2)
226 {
227   //
228   // Process the events
229   //
230
231   //at least first event is needed!
232   if (!ev1){
233     AliError("At least first event must be set!");
234     return;
235   }
236
237   // modify event numbers in MC so that we can identify new events 
238   // in AliDielectronV0Cuts (not neeeded for collision data)
239   if(GetHasMC()) {
240     ev1->SetBunchCrossNumber(1);
241     ev1->SetOrbitNumber(1);
242     ev1->SetPeriodNumber(1);
243   }
244
245   // set event
246   AliDielectronVarManager::SetEvent(ev1);
247   if (fMixing){
248     //set mixing bin to event data
249     Int_t bin=fMixing->FindBin(AliDielectronVarManager::GetData());
250     AliDielectronVarManager::SetValue(AliDielectronVarManager::kMixingBin,bin);
251   }
252
253   //in case we have MC load the MC event and process the MC particles
254   if (AliDielectronMC::Instance()->ConnectMCEvent()){
255     ProcessMC(ev1);
256   }
257
258   //if candidate array doesn't exist, create it
259   if (!fPairCandidates->UncheckedAt(0)) {
260     InitPairCandidateArrays();
261   } else {
262     ClearArrays();
263   }
264
265   //mask used to require that all cuts are fulfilled
266   UInt_t selectedMask=(1<<fEventFilter.GetCuts()->GetEntries())-1;
267
268   //apply event cuts
269     if ((ev1&&fEventFilter.IsSelected(ev1)!=selectedMask) ||
270         (ev2&&fEventFilter.IsSelected(ev2)!=selectedMask)) return;
271
272   //fill track arrays for the first event
273   if (ev1){
274     FillTrackArrays(ev1);
275     if (((fPreFilterAllSigns)||(fPreFilterUnlikeOnly)) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(0, 1, fTracks[0], fTracks[1]);
276   }
277
278
279   //fill track arrays for the second event
280   if (ev2) {
281     FillTrackArrays(ev2,1);
282     if (((fPreFilterAllSigns)||(fPreFilterUnlikeOnly)) && ( fPairPreFilter.GetCuts()->GetEntries()>0 )) PairPreFilter(2, 3, fTracks[2], fTracks[3]);
283   }
284
285   // TPC event plane correction
286   if (ev1 && fPreFilterEventPlane && ( fEventPlanePreFilter.GetCuts()->GetEntries()>0 || fEventPlanePOIPreFilter.GetCuts()->GetEntries()>0)) 
287     EventPlanePreFilter(0, 1, fTracks[0], fTracks[1], ev1);
288
289   if (!fNoPairing){
290     // create pairs and fill pair candidate arrays
291     for (Int_t itrackArr1=0; itrackArr1<4; ++itrackArr1){
292       for (Int_t itrackArr2=itrackArr1; itrackArr2<4; ++itrackArr2){
293         FillPairArrays(itrackArr1, itrackArr2);
294       }
295     }
296
297     //track rotation
298     if (fTrackRotator) {
299       fTrackRotator->SetEvent(ev1);
300       FillPairArrayTR();
301     }
302   }
303
304   //fill debug tree if a manager is attached
305   if (fDebugTree) FillDebugTree();
306
307   //process event mixing
308   if (fMixing) {
309     fMixing->Fill(ev1,this);
310     //     FillHistograms(0x0,kTRUE);
311   }
312
313   //in case there is a histogram manager, fill the QA histograms
314   if (fHistos && fSignalsMC) FillMCHistograms(ev1);
315   if (fHistos) FillHistograms(ev1);
316
317
318   // clear arrays
319   if (!fDontClearArrays) ClearArrays();
320
321   // reset TPC EP and unique identifiers for v0 cut class 
322   AliDielectronVarManager::SetTPCEventPlane(0x0);
323   if(GetHasMC()) { // only for MC needed
324     for (Int_t iCut=0; iCut<fTrackFilter.GetCuts()->GetEntries();++iCut) {
325       if ( fTrackFilter.GetCuts()->At(iCut)->IsA() == AliDielectronV0Cuts::Class() )
326         ((AliDielectronV0Cuts*)fTrackFilter.GetCuts()->At(iCut))->ResetUniqueEventNumbers();
327     }
328   }
329
330 }
331
332 //________________________________________________________________
333 void AliDielectron::ProcessMC(AliVEvent *ev1)
334 {
335   //
336   // Process the MC data
337   //
338
339   AliDielectronMC *dieMC=AliDielectronMC::Instance();
340
341   if (fHistos) FillHistogramsMC(dieMC->GetMCEvent(), ev1);
342
343   if(!fSignalsMC) return;
344   //loop over all MC data and Fill the CF container if it exist
345   if (!fCfManagerPair && !fHistoArray) return;
346   if(fCfManagerPair) fCfManagerPair->SetPdgMother(fPdgMother);
347
348   Bool_t bFillCF = (fCfManagerPair ? fCfManagerPair->GetStepForMCtruth()  : kFALSE);
349   Bool_t bFillHF = (fHistoArray    ? fHistoArray->GetStepForMCGenerated() : kFALSE);
350   if(!bFillCF && !bFillHF) return;
351
352   // signals to be studied
353   Int_t nSignals = fSignalsMC->GetEntries();
354
355   // initialize 2D arrays of labels for particles from each MC signal
356   Int_t** labels1;      // labels for particles satisfying branch 1
357   Int_t** labels2;      // labels for particles satisfying branch 2
358   Int_t** labels12;     // labels for particles satisfying both branches
359   labels1 = new Int_t*[nSignals];
360   labels2 = new Int_t*[nSignals];
361   labels12 = new Int_t*[nSignals];
362   Int_t* indexes1=new Int_t[nSignals];
363   Int_t* indexes2=new Int_t[nSignals];
364   Int_t* indexes12=new Int_t[nSignals];
365   for(Int_t isig=0;isig<nSignals;++isig) {
366     *(labels1+isig) = new Int_t[dieMC->GetNMCTracks()];
367     *(labels2+isig) = new Int_t[dieMC->GetNMCTracks()];
368     *(labels12+isig) = new Int_t[dieMC->GetNMCTracks()];
369     for(Int_t ip=0; ip<dieMC->GetNMCTracks();++ip) {
370       labels1[isig][ip] = -1;
371       labels2[isig][ip] = -1;
372       labels12[isig][ip] = -1;
373     }
374     indexes1[isig]=0;
375     indexes2[isig]=0;
376     indexes12[isig]=0;
377   }
378
379   Bool_t truth1=kFALSE;
380   Bool_t truth2=kFALSE;
381   // loop over the MC tracks
382   for(Int_t ipart=0; ipart<dieMC->GetNMCTracks(); ++ipart) {
383     for(Int_t isig=0; isig<nSignals; ++isig) {       // loop over signals
384       // Proceed only if this signal is required in the pure MC step
385       // NOTE: Some signals can be satisfied by many particles and this leads to high
386       //       computation times (e.g. secondary electrons from the GEANT transport). Be aware of this!!
387       if(!((AliDielectronSignalMC*)fSignalsMC->At(isig))->GetFillPureMCStep()) continue;
388
389       truth1 = dieMC->IsMCTruth(ipart, (AliDielectronSignalMC*)fSignalsMC->At(isig), 1);
390       truth2 = dieMC->IsMCTruth(ipart, (AliDielectronSignalMC*)fSignalsMC->At(isig), 2);
391
392       // particles satisfying both branches are treated separately to avoid double counting during pairing
393       if(truth1 && truth2) {
394         labels12[isig][indexes12[isig]] = ipart;
395         ++indexes12[isig];
396       }
397       else {
398         if(truth1) {
399           labels1[isig][indexes1[isig]] = ipart;
400           ++indexes1[isig];
401         }
402         if(truth2) {
403           labels2[isig][indexes2[isig]] = ipart;
404           ++indexes2[isig];
405         }
406       }
407     }
408   }  // end loop over MC particles
409
410   // Do the pairing and fill the CF container with pure MC info
411   for(Int_t isig=0; isig<nSignals; ++isig) {
412     // mix the particles which satisfy only one of the signal branches
413     for(Int_t i1=0;i1<indexes1[isig];++i1) {
414       for(Int_t i2=0;i2<indexes2[isig];++i2) {
415         if(bFillCF) fCfManagerPair->FillMC(labels1[isig][i1], labels2[isig][i2], isig);
416         if(bFillHF) fHistoArray->Fill(labels1[isig][i1], labels2[isig][i2], isig);
417       }
418     }
419     // mix the particles which satisfy both branches
420     for(Int_t i1=0;i1<indexes12[isig];++i1) {
421       for(Int_t i2=0; i2<i1; ++i2) {
422         if(bFillCF) fCfManagerPair->FillMC(labels12[isig][i1], labels12[isig][i2], isig);
423         if(bFillHF) fHistoArray->Fill(labels12[isig][i1], labels12[isig][i2], isig);
424       }
425     }
426   }    // end loop over signals
427
428   // release the memory
429   for(Int_t isig=0;isig<nSignals;++isig) {
430     delete [] *(labels1+isig);
431     delete [] *(labels2+isig);
432     delete [] *(labels12+isig);
433   }
434   delete [] labels1;
435   delete [] labels2;
436   delete [] labels12;
437   delete [] indexes1;
438   delete [] indexes2;
439   delete [] indexes12;
440 }
441
442 //________________________________________________________________
443 void AliDielectron::FillHistogramsTracks(TObjArray **tracks)
444 {
445   //
446   // Fill Histogram information for tracks after prefilter
447   // ignore mixed events - for prefilter, only single tracks +/- are relevant 
448   //
449   
450   TString  className,className2;
451   Double_t values[AliDielectronVarManager::kNMaxValues];
452   
453   //Fill track information, separately for the track array candidates
454   for (Int_t i=0; i<2; ++i){
455     className.Form("Pre_%s",fgkTrackClassNames[i]);
456     if (!fHistos->GetHistogramList()->FindObject(className.Data())) continue;
457     Int_t ntracks=tracks[i]->GetEntriesFast();
458     for (Int_t itrack=0; itrack<ntracks; ++itrack){
459       AliDielectronVarManager::Fill(tracks[i]->UncheckedAt(itrack), values);
460       fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
461     }
462   }
463 }
464
465
466 //________________________________________________________________
467 void AliDielectron::FillHistogramsMC(const AliMCEvent *ev, AliVEvent *ev1)
468 {
469   //
470   // Fill Histogram information for MCEvents
471   //
472
473   Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
474   // Fill event information
475   AliDielectronVarManager::Fill(ev1, values);    // ESD/AOD information
476   AliDielectronVarManager::Fill(ev, values);     // MC truth info
477   if (fHistos->GetHistogramList()->FindObject("MCEvent"))
478     fHistos->FillClass("MCEvent", AliDielectronVarManager::kNMaxValues, values);
479 }
480
481
482 //________________________________________________________________
483 void AliDielectron::FillHistograms(const AliVEvent *ev, Bool_t pairInfoOnly)
484 {
485   //
486   // Fill Histogram information for tracks and pairs
487   //
488   
489   TString  className,className2;
490   Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
491
492   //Fill event information
493   if (ev){
494     if (fHistos->GetHistogramList()->FindObject("Event")) {
495       fHistos->FillClass("Event", AliDielectronVarManager::kNMaxValues, AliDielectronVarManager::GetData());
496     }
497   }
498
499   //Fill track information, separately for the track array candidates
500   if (!pairInfoOnly){
501     className2.Form("Track_%s",fgkPairClassNames[1]);  // unlike sign, SE only
502     for (Int_t i=0; i<4; ++i){
503       className.Form("Track_%s",fgkTrackClassNames[i]);
504       Bool_t mergedtrkClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
505       Bool_t trkClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
506       if (!trkClass && !mergedtrkClass) continue;
507       Int_t ntracks=fTracks[i].GetEntriesFast();
508       for (Int_t itrack=0; itrack<ntracks; ++itrack){
509         AliDielectronVarManager::Fill(fTracks[i].UncheckedAt(itrack), values);
510         if(trkClass)
511           fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
512         if(mergedtrkClass && i<2)
513           fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values); //only ev1
514       }
515     }
516   }
517
518   //Fill Pair information, separately for all pair candidate arrays and the legs
519   TObjArray arrLegs(100);
520   for (Int_t i=0; i<10; ++i){
521     className.Form("Pair_%s",fgkPairClassNames[i]);
522     className2.Form("Track_Legs_%s",fgkPairClassNames[i]);
523     Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
524     Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
525     if (!pairClass&&!legClass) continue;
526     Int_t ntracks=PairArray(i)->GetEntriesFast();
527     for (Int_t ipair=0; ipair<ntracks; ++ipair){
528       AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair));
529       
530       //fill pair information
531       if (pairClass){
532         AliDielectronVarManager::Fill(pair, values);
533         fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
534       }
535
536       //fill leg information, don't fill the information twice
537       if (legClass){
538         AliVParticle *d1=pair->GetFirstDaughter();
539         AliVParticle *d2=pair->GetSecondDaughter();
540         if (!arrLegs.FindObject(d1)){
541           AliDielectronVarManager::Fill(d1, values);
542           fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
543           arrLegs.Add(d1);
544         }
545         if (!arrLegs.FindObject(d2)){
546           AliDielectronVarManager::Fill(d2, values);
547           fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
548           arrLegs.Add(d2);
549         }
550       }
551     }
552     if (legClass) arrLegs.Clear();
553   }
554   
555 }
556 //________________________________________________________________
557 void AliDielectron::FillHistogramsPair(AliDielectronPair *pair,Bool_t fromPreFilter/*=kFALSE*/)
558 {
559   //
560   // Fill Histogram information for pairs and the track in the pair
561   // NOTE: in this funtion the leg information may be filled multiple
562   //       times. This funtion is used in the track rotation pairing
563   //       and those legs are not saved!
564   //
565   TString  className,className2;
566   Double_t values[AliDielectronVarManager::kNMaxValues];
567   
568   //Fill Pair information, separately for all pair candidate arrays and the legs
569   TObjArray arrLegs(100);
570   const Int_t type=pair->GetType();
571   if (fromPreFilter) {
572     className.Form("RejPair_%s",fgkPairClassNames[type]);
573     className2.Form("RejTrack_%s",fgkPairClassNames[type]);
574   } else {
575     className.Form("Pair_%s",fgkPairClassNames[type]);
576     className2.Form("Track_Legs_%s",fgkPairClassNames[type]);
577   }
578   
579   Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
580   Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
581   
582   //fill pair information
583   if (pairClass){
584     AliDielectronVarManager::Fill(pair, values);
585     fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
586   }
587
588   if (legClass){
589     AliVParticle *d1=pair->GetFirstDaughter();
590     AliDielectronVarManager::Fill(d1, values);
591     fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
592     
593     AliVParticle *d2=pair->GetSecondDaughter();
594     AliDielectronVarManager::Fill(d2, values);
595     fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
596   }
597 }
598
599 //________________________________________________________________
600 void AliDielectron::FillTrackArrays(AliVEvent * const ev, Int_t eventNr)
601 {
602   //
603   // select tracks and fill track candidate arrays
604   // eventNr = 0: First  event, use track arrays 0 and 1
605   // eventNr = 1: Second event, use track arrays 2 and 3
606   //
607
608   Int_t ntracks=ev->GetNumberOfTracks();
609
610   UInt_t selectedMask=(1<<fTrackFilter.GetCuts()->GetEntries())-1;
611   for (Int_t itrack=0; itrack<ntracks; ++itrack){
612     //get particle
613     AliVParticle *particle=ev->GetTrack(itrack);
614
615     //apply track cuts
616     if (fTrackFilter.IsSelected(particle)!=selectedMask) continue;
617
618     //fill selected particle into the corresponding track arrays
619     Short_t charge=particle->Charge();
620     if (charge>0)      fTracks[eventNr*2].Add(particle);
621     else if (charge<0) fTracks[eventNr*2+1].Add(particle);
622   }
623 }
624
625 //________________________________________________________________
626 void AliDielectron::EventPlanePreFilter(Int_t arr1, Int_t arr2, TObjArray arrTracks1, TObjArray arrTracks2, const AliVEvent *ev)
627 {
628   //
629   // Prefilter tracks and tracks from pairs
630   // Needed for rejection in the Q-Vector of the event plane
631   // remove contribution of all tracks to the Q-vector that are in invariant mass window 
632   //
633
634   AliEventplane *evplane = const_cast<AliVEvent *>(ev)->GetEventplane();
635   if(!evplane) { // nanoAODs , here we do NOT have sub event reaction planes
636     //  if(1) {
637     // get the EPselectionTask for recalculation of weighting factors
638     AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
639     AliEPSelectionTask *eptask = dynamic_cast<AliEPSelectionTask *>(man->GetTask("EventplaneSelection"));
640     if(!eptask) return;
641
642     // track mapping
643     TMap mapRemovedTracks;
644
645     // eta gap ?
646     Bool_t etagap = kFALSE;
647     for (Int_t iCut=0; iCut<fEventPlanePreFilter.GetCuts()->GetEntries();++iCut) {
648       TString cutName=fEventPlanePreFilter.GetCuts()->At(iCut)->GetName();
649       if(cutName.Contains("eta") || cutName.Contains("Eta"))  etagap=kTRUE;
650     }
651
652     Double_t cQX=0., cQY=0.;
653     // apply cuts to the tracks, e.g. etagap
654     if(fEventPlanePreFilter.GetCuts()->GetEntries()) {
655       UInt_t selectedMask=(1<<fEventPlanePreFilter.GetCuts()->GetEntries())-1;
656       Int_t ntracks=ev->GetNumberOfTracks();
657       for (Int_t itrack=0; itrack<ntracks; ++itrack){
658         AliVParticle *particle=ev->GetTrack(itrack);
659         AliVTrack *track= static_cast<AliVTrack*>(particle);
660         if (!track) continue;
661         //event plane cuts
662         UInt_t cutMask=fEventPlanePreFilter.IsSelected(track);
663         //apply cut
664         if (cutMask==selectedMask) continue;
665
666         mapRemovedTracks.Add(track,track);
667         cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
668         cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
669       }
670     }
671
672     // POI (particle of interest) rejection
673     Int_t pairIndex=GetPairIndex(arr1,arr2);
674
675     Int_t ntrack1=arrTracks1.GetEntriesFast();
676     Int_t ntrack2=arrTracks2.GetEntriesFast();
677     AliDielectronPair candidate;
678     candidate.SetKFUsage(fUseKF);
679
680     UInt_t selectedMask=(1<<fEventPlanePOIPreFilter.GetCuts()->GetEntries())-1;
681     for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
682       Int_t end=ntrack2;
683       if (arr1==arr2) end=itrack1;
684       Bool_t accepted=kFALSE;
685       for (Int_t itrack2=0; itrack2<end; ++itrack2){
686         TObject *track1=arrTracks1.UncheckedAt(itrack1);
687         TObject *track2=arrTracks2.UncheckedAt(itrack2);
688         if (!track1 || !track2) continue;
689         //create the pair
690         candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
691                             static_cast<AliVTrack*>(track2), fPdgLeg2);
692         candidate.SetType(pairIndex);
693         candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
694
695         //event plane pair cuts
696         UInt_t cutMask=fEventPlanePOIPreFilter.IsSelected(&candidate);
697         //apply cut
698         if (cutMask==selectedMask) continue;
699
700         accepted=kTRUE;
701         //remove the tracks from the Track arrays
702         arrTracks2.AddAt(0x0,itrack2);
703       }
704       if ( accepted ) arrTracks1.AddAt(0x0,itrack1);
705     }
706     //compress the track arrays
707     arrTracks1.Compress();
708     arrTracks2.Compress();
709
710     //Modify the components: subtract the tracks
711     ntrack1=arrTracks1.GetEntriesFast();
712     ntrack2=arrTracks2.GetEntriesFast();
713     // remove leg1 contribution
714     for (Int_t itrack=0; itrack<ntrack1; ++itrack){
715       AliVTrack *track= static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack));
716       if (!track) continue;
717       // track contribution was already removed
718       if (mapRemovedTracks.FindObject(track)) continue;
719       else mapRemovedTracks.Add(track,track);
720
721       cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
722       cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
723     }
724     // remove leg2 contribution
725     for (Int_t itrack=0; itrack<ntrack2; ++itrack){
726       AliVTrack *track= static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack));
727       if (!track) continue;
728       // track contribution was already removed
729       if (mapRemovedTracks.FindObject(track)) continue;
730       else mapRemovedTracks.Add(track,track);
731
732       cQX += (eptask->GetWeight(track) * TMath::Cos(2*track->Phi()));
733       cQY += (eptask->GetWeight(track) * TMath::Sin(2*track->Phi()));
734     }
735
736     // build a corrected alieventplane using the values from the var manager
737     // these uncorrected values are filled using the stored magnitude and angle  in the header
738     TVector2 qcorr;
739     qcorr.Set(AliDielectronVarManager::GetValue(AliDielectronVarManager::kTPCxH2uc)-cQX,
740               AliDielectronVarManager::GetValue(AliDielectronVarManager::kTPCyH2uc)-cQY);
741     // fill alieventplane
742     AliEventplane cevplane;
743     cevplane.SetQVector(&qcorr);
744     AliDielectronVarManager::SetTPCEventPlane(&cevplane);
745     cevplane.SetQVector(0);
746     return;
747   } //end: nanoAODs
748   else
749     {
750     // this is done in case of ESDs or AODs
751     Bool_t isESD=(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliESDInputHandler::Class());
752     // copy event plane object
753     AliEventplane cevplane(*evplane);
754     //    Int_t nMaxID=cevplane->GetQContributionXArray()->GetSize();
755
756     TVector2 *qcorr  = cevplane.GetQVector();
757     if(!qcorr) return;
758     TVector2 *qcsub1 = 0x0;
759     TVector2 *qcsub2 = 0x0;
760
761     // eta gap ?
762     Bool_t etagap = kFALSE;
763     for (Int_t iCut=0; iCut<fEventPlanePreFilter.GetCuts()->GetEntries();++iCut) {
764       TString cutName=fEventPlanePreFilter.GetCuts()->At(iCut)->GetName();
765       if(cutName.Contains("eta") || cutName.Contains("Eta"))  etagap=kTRUE;
766     }
767
768     // subevent configuration for eta gap or LS (default is rndm)
769     if(fLikeSignSubEvents && etagap) {
770       // start with the full Qvector/event in both sub events
771       qcsub1 = new TVector2(*qcorr);
772       qcsub2 = new TVector2(*qcorr);
773       cevplane.SetQsub(qcsub1,qcsub2);
774
775       Int_t ntracks=ev->GetNumberOfTracks();
776       // track removals
777       for (Int_t itrack=0; itrack<ntracks; ++itrack){
778         AliVParticle *particle=ev->GetTrack(itrack);
779         AliVTrack *track= static_cast<AliVTrack*>(particle);
780         if (!track) continue;
781         if (track->GetID()>=0 && !isESD) continue;
782         Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
783
784         // set contributions to zero
785         // charge sub1+ sub2-
786         if(fLikeSignSubEvents) {
787           Short_t charge=track->Charge();
788           if (charge<0) {
789             cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
790             cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
791           }
792           if (charge>0) {
793             cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
794             cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
795           }
796         }
797         // eta sub1+ sub2-
798         if(etagap) {
799           Double_t eta=track->Eta();
800           if (eta<0.0) {
801             cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
802             cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
803           }
804           if (eta>0.0) {
805             cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
806             cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
807           }
808         }
809       } // end: loop over tracks
810     } // end: sub event configuration
811
812     // apply cuts, e.g. etagap
813     if(fEventPlanePreFilter.GetCuts()->GetEntries()) {
814       UInt_t selectedMask=(1<<fEventPlanePreFilter.GetCuts()->GetEntries())-1;
815       Int_t ntracks=ev->GetNumberOfTracks();
816       for (Int_t itrack=0; itrack<ntracks; ++itrack){
817         AliVParticle *particle=ev->GetTrack(itrack);
818         AliVTrack *track= static_cast<AliVTrack*>(particle);
819         if (!track) continue;
820         if (track->GetID()>=0 && !isESD) continue;
821         Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
822
823         //event plane cuts
824         UInt_t cutMask=fEventPlanePreFilter.IsSelected(track);
825         //apply cut
826         if (cutMask==selectedMask) continue;
827
828         // set contributions to zero
829         cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
830         cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
831         cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
832         cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
833         cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
834         cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
835       }
836     } // end: track cuts
837
838     // POI (particle of interest) rejection
839     Int_t pairIndex=GetPairIndex(arr1,arr2);
840     Int_t ntrack1=arrTracks1.GetEntriesFast();
841     Int_t ntrack2=arrTracks2.GetEntriesFast();
842     AliDielectronPair candidate;
843     candidate.SetKFUsage(fUseKF);
844
845     UInt_t selectedMask=(1<<fEventPlanePOIPreFilter.GetCuts()->GetEntries())-1;
846     for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
847       Int_t end=ntrack2;
848       if (arr1==arr2) end=itrack1;
849       Bool_t accepted=kFALSE;
850       for (Int_t itrack2=0; itrack2<end; ++itrack2){
851         TObject *track1=arrTracks1.UncheckedAt(itrack1);
852         TObject *track2=arrTracks2.UncheckedAt(itrack2);
853         if (!track1 || !track2) continue;
854         //create the pair
855         candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
856                             static_cast<AliVTrack*>(track2), fPdgLeg2);
857
858         candidate.SetType(pairIndex);
859         candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
860
861         //event plane cuts
862         UInt_t cutMask=fEventPlanePOIPreFilter.IsSelected(&candidate);
863         //apply cut
864         if (cutMask==selectedMask) continue;
865
866         accepted=kTRUE;
867         //remove the tracks from the Track arrays
868         arrTracks2.AddAt(0x0,itrack2);
869       }
870       if ( accepted ) arrTracks1.AddAt(0x0,itrack1);
871     }
872     //compress the track arrays
873     arrTracks1.Compress();
874     arrTracks2.Compress();
875
876     //Modify the components: subtract the tracks
877     ntrack1=arrTracks1.GetEntriesFast();
878     ntrack2=arrTracks2.GetEntriesFast();
879     // remove leg1 contribution
880     for (Int_t itrack=0; itrack<ntrack1; ++itrack){
881       AliVTrack *track= static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack));
882       if (!track) continue;
883       if (track->GetID()>=0 && !isESD) continue;
884       Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
885       // set contributions to zero
886       cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
887       cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
888       cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
889       cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
890       cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
891       cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
892     }
893     // remove leg2 contribution
894     for (Int_t itrack=0; itrack<ntrack2; ++itrack){
895       AliVTrack *track= static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack));
896       if (!track) continue;
897       if (track->GetID()>=0 && !isESD) continue;
898       Int_t tmpID = isESD ? track->GetID() : track->GetID()*-1 - 1;
899       // set contributions to zero
900       cevplane.GetQContributionXArray()->SetAt(0.0, tmpID);
901       cevplane.GetQContributionYArray()->SetAt(0.0, tmpID);
902       cevplane.GetQContributionXArraysub1()->SetAt(0.0, tmpID);
903       cevplane.GetQContributionYArraysub1()->SetAt(0.0, tmpID);
904       cevplane.GetQContributionXArraysub2()->SetAt(0.0, tmpID);
905       cevplane.GetQContributionYArraysub2()->SetAt(0.0, tmpID);
906     }
907
908     // set corrected AliEventplane and fill variables with corrected values
909     AliDielectronVarManager::SetTPCEventPlane(&cevplane);
910     delete qcsub1;
911     delete qcsub2;
912   } // end: ESD or AOD case
913
914 }
915 //________________________________________________________________
916 void AliDielectron::PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2)
917 {
918   //
919   // Prefilter tracks from pairs
920   // Needed for datlitz rejections
921   // remove all tracks from the Single track arrays that pass the cuts in this filter
922   //
923
924   Int_t ntrack1=arrTracks1.GetEntriesFast();
925   Int_t ntrack2=arrTracks2.GetEntriesFast();
926   AliDielectronPair candidate;
927   candidate.SetKFUsage(fUseKF);
928   // flag arrays for track removal
929   Bool_t *bTracks1 = new Bool_t[ntrack1];
930   for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1) bTracks1[itrack1]=kFALSE;
931   Bool_t *bTracks2 = new Bool_t[ntrack2];
932   for (Int_t itrack2=0; itrack2<ntrack2; ++itrack2) bTracks2[itrack2]=kFALSE;
933
934   UInt_t selectedMask=(1<<fPairPreFilter.GetCuts()->GetEntries())-1;
935   UInt_t selectedMaskPair=(1<<fPairFilter.GetCuts()->GetEntries())-1;
936
937   Int_t nRejPasses = 1; //for fPreFilterUnlikeOnly and no set flag 
938   if (fPreFilterAllSigns) nRejPasses = 3;
939
940   for (Int_t iRP=0; iRP < nRejPasses; ++iRP) {
941         Int_t arr1RP=arr1, arr2RP=arr2;
942         TObjArray *arrTracks1RP=&arrTracks1;
943         TObjArray *arrTracks2RP=&arrTracks2;
944         Bool_t *bTracks1RP = bTracks1;
945         Bool_t *bTracks2RP = bTracks2;
946         switch (iRP) {
947                 case 1: arr1RP=arr1;arr2RP=arr1;
948                                 arrTracks1RP=&arrTracks1;
949                                 arrTracks2RP=&arrTracks1;
950                                 bTracks1RP = bTracks1;
951                                 bTracks2RP = bTracks1;
952                                 break;
953                 case 2: arr1RP=arr2;arr2RP=arr2;
954                                 arrTracks1RP=&arrTracks2;
955                                 arrTracks2RP=&arrTracks2;
956                                 bTracks1RP = bTracks2;
957                                 bTracks2RP = bTracks2;
958                                 break;
959                 default: ;//nothing to do
960         }
961         Int_t ntrack1RP=(*arrTracks1RP).GetEntriesFast();
962         Int_t ntrack2RP=(*arrTracks2RP).GetEntriesFast();
963
964         Int_t pairIndex=GetPairIndex(arr1RP,arr2RP);
965
966         for (Int_t itrack1=0; itrack1<ntrack1RP; ++itrack1){
967           Int_t end=ntrack2RP;
968           if (arr1RP==arr2RP) end=itrack1;
969           for (Int_t itrack2=0; itrack2<end; ++itrack2){
970                 TObject *track1=(*arrTracks1RP).UncheckedAt(itrack1);
971                 TObject *track2=(*arrTracks2RP).UncheckedAt(itrack2);
972                 if (!track1 || !track2) continue;
973                 //create the pair
974                 candidate.SetTracks(static_cast<AliVTrack*>(track1), fPdgLeg1,
975                         static_cast<AliVTrack*>(track2), fPdgLeg2);
976
977                 candidate.SetType(pairIndex);
978                 candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
979                 //relate to the production vertex
980                 //       if (AliDielectronVarManager::GetKFVertex()) candidate.SetProductionVertex(*AliDielectronVarManager::GetKFVertex());
981
982                 //pair cuts
983                 UInt_t cutMask=fPairPreFilter.IsSelected(&candidate);
984
985                 //apply cut
986                 if (cutMask!=selectedMask) continue;
987                 if (fCfManagerPair) fCfManagerPair->Fill(selectedMaskPair+1 ,&candidate);
988                 if (fHistos) FillHistogramsPair(&candidate,kTRUE);
989                 //set flags for track removal
990                 bTracks1RP[itrack1]=kTRUE;
991                 bTracks2RP[itrack2]=kTRUE;
992           }
993         }
994   }
995
996   //remove the tracks from the Track arrays
997   for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
998     if(bTracks1[itrack1]) arrTracks1.AddAt(0x0, itrack1);
999   }
1000   for (Int_t itrack2=0; itrack2<ntrack2; ++itrack2){
1001     if(bTracks2[itrack2]) arrTracks2.AddAt(0x0, itrack2);
1002   }
1003
1004   // clean up
1005   delete [] bTracks1;
1006   delete [] bTracks2;
1007
1008   //compress the track arrays
1009   arrTracks1.Compress();
1010   arrTracks2.Compress();
1011   
1012   //apply leg cuts after the pre filter
1013   if ( fPairPreFilterLegs.GetCuts()->GetEntries()>0 ) {
1014     selectedMask=(1<<fPairPreFilterLegs.GetCuts()->GetEntries())-1;
1015     //loop over tracks from array 1
1016     for (Int_t itrack=0; itrack<arrTracks1.GetEntriesFast();++itrack){
1017       //test cuts
1018       UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks1.UncheckedAt(itrack));
1019       
1020       //apply cut
1021       if (cutMask!=selectedMask) arrTracks1.AddAt(0x0,itrack);;
1022     }
1023     arrTracks1.Compress();
1024     
1025     //in case of like sign don't loop over second array
1026     if (arr1==arr2) {
1027       arrTracks2=arrTracks1;
1028     } else {
1029       
1030       //loop over tracks from array 2
1031       for (Int_t itrack=0; itrack<arrTracks2.GetEntriesFast();++itrack){
1032       //test cuts
1033         UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks2.UncheckedAt(itrack));
1034       //apply cut
1035         if (cutMask!=selectedMask) arrTracks2.AddAt(0x0,itrack);
1036       }
1037       arrTracks2.Compress();
1038       
1039     }
1040   }
1041   //For unlike-sign monitor track-cuts:
1042   if (arr1!=arr2&&fHistos) {
1043     TObjArray *unlikesignArray[2] = {&arrTracks1,&arrTracks2};
1044     FillHistogramsTracks(unlikesignArray);
1045   }
1046 }
1047
1048 //________________________________________________________________
1049 void AliDielectron::FillPairArrays(Int_t arr1, Int_t arr2)
1050 {
1051   //
1052   // select pairs and fill pair candidate arrays
1053   //
1054
1055   TObjArray arrTracks1=fTracks[arr1];
1056   TObjArray arrTracks2=fTracks[arr2];
1057
1058   //process pre filter if set
1059   if ((!fPreFilterAllSigns) && (!fPreFilterUnlikeOnly) && ( fPairPreFilter.GetCuts()->GetEntries()>0 ))  PairPreFilter(arr1, arr2, arrTracks1, arrTracks2);
1060   
1061   Int_t pairIndex=GetPairIndex(arr1,arr2);
1062
1063   Int_t ntrack1=arrTracks1.GetEntriesFast();
1064   Int_t ntrack2=arrTracks2.GetEntriesFast();
1065
1066   AliDielectronPair *candidate=new AliDielectronPair;
1067   candidate->SetKFUsage(fUseKF);
1068
1069   UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1070   
1071   for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
1072     Int_t end=ntrack2;
1073     if (arr1==arr2) end=itrack1;
1074     for (Int_t itrack2=0; itrack2<end; ++itrack2){
1075       //create the pair
1076       candidate->SetTracks(static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1)), fPdgLeg1,
1077                              static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2)), fPdgLeg2);
1078       candidate->SetType(pairIndex);
1079       Int_t label=AliDielectronMC::Instance()->GetLabelMotherWithPdg(candidate,fPdgMother);
1080       candidate->SetLabel(label);
1081       if (label>-1) candidate->SetPdgCode(fPdgMother);
1082
1083       // check for gamma kf particle
1084       label=AliDielectronMC::Instance()->GetLabelMotherWithPdg(candidate,22);
1085       if (label>-1) {
1086         candidate->SetGammaTracks(static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1)), fPdgLeg1,
1087                                   static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2)), fPdgLeg2);
1088       // should we set the pdgmothercode and the label
1089       }
1090
1091       //pair cuts
1092       UInt_t cutMask=fPairFilter.IsSelected(candidate);
1093       
1094       //CF manager for the pair
1095       if (fCfManagerPair) fCfManagerPair->Fill(cutMask,candidate);
1096       //histogram array for the pair
1097       if (fHistoArray) fHistoArray->Fill(pairIndex,candidate);
1098
1099       //apply cut
1100       if (cutMask!=selectedMask) continue;
1101
1102       //add the candidate to the candidate array 
1103       PairArray(pairIndex)->Add(candidate);
1104       //get a new candidate
1105       candidate=new AliDielectronPair;
1106           candidate->SetKFUsage(fUseKF);
1107     }
1108   }
1109   //delete the surplus candidate
1110   delete candidate;
1111 }
1112
1113 //________________________________________________________________
1114 void AliDielectron::FillPairArrayTR()
1115 {
1116   //
1117   // select pairs and fill pair candidate arrays
1118   //
1119   UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
1120   
1121   while ( fTrackRotator->NextCombination() ){
1122     AliDielectronPair candidate;
1123     candidate.SetKFUsage(fUseKF);
1124     candidate.SetTracks(&fTrackRotator->GetKFTrackP(), &fTrackRotator->GetKFTrackN(),
1125                         fTrackRotator->GetVTrackP(),fTrackRotator->GetVTrackN());
1126     candidate.SetType(kEv1PMRot);
1127     
1128     //pair cuts
1129     UInt_t cutMask=fPairFilter.IsSelected(&candidate);
1130     
1131     //CF manager for the pair
1132     if (fCfManagerPair) fCfManagerPair->Fill(cutMask,&candidate);
1133     //histogram array for the pair
1134     if (fHistoArray) fHistoArray->Fill((Int_t)kEv1PMRot,&candidate);
1135     
1136     //apply cut
1137     if (cutMask==selectedMask) {
1138      if(fHistos) FillHistogramsPair(&candidate);
1139      if(fStoreRotatedPairs) PairArray(kEv1PMRot)->Add(new AliDielectronPair(candidate));
1140     } 
1141   }
1142 }
1143
1144 //________________________________________________________________
1145 void AliDielectron::FillDebugTree()
1146 {
1147   //
1148   // Fill Histogram information for tracks and pairs
1149   //
1150   
1151   //Fill Debug tree
1152   for (Int_t i=0; i<10; ++i){
1153     Int_t ntracks=PairArray(i)->GetEntriesFast();
1154     for (Int_t ipair=0; ipair<ntracks; ++ipair){
1155       fDebugTree->Fill(static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair)));
1156     }
1157   }
1158 }
1159
1160 //________________________________________________________________
1161 void AliDielectron::SaveDebugTree()
1162 {
1163   //
1164   // delete the debug tree, this will also write the tree
1165   //
1166   if (fDebugTree) fDebugTree->DeleteStreamer();
1167 }
1168
1169
1170 //__________________________________________________________________
1171 void AliDielectron::AddSignalMC(AliDielectronSignalMC* signal) {
1172   //
1173   //  Add an MC signal to the signals list
1174   //
1175   if(!fSignalsMC) {
1176     fSignalsMC = new TObjArray();
1177     fSignalsMC->SetOwner();
1178   }
1179   fSignalsMC->Add(signal);
1180 }
1181 //________________________________________________________________
1182 void AliDielectron::FillMCHistograms(const AliVEvent *ev) {
1183   //
1184   // fill QA MC histograms for pairs and legs of all added mc signals
1185   //
1186
1187   if (!fSignalsMC) return;
1188   TString className,className2;
1189   Double_t values[AliDielectronVarManager::kNMaxValues]={0.};
1190   AliDielectronVarManager::Fill(ev, values); // get event informations
1191   //loop over all added mc signals
1192   for(Int_t isig=0; isig<fSignalsMC->GetEntries(); isig++) {
1193
1194     className.Form("Pair_%s",fSignalsMC->At(isig)->GetName());
1195     className2.Form("Track_Legs_%s",fSignalsMC->At(isig)->GetName());
1196     Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
1197     Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
1198     if(!pairClass && !legClass) return;
1199
1200     Int_t ntracks=PairArray(AliDielectron::kEv1PM)->GetEntriesFast(); // only SE +-
1201     for (Int_t ipair=0; ipair<ntracks; ++ipair){
1202       AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(AliDielectron::kEv1PM)->UncheckedAt(ipair));
1203
1204       Bool_t isMCtruth = AliDielectronMC::Instance()->IsMCTruth(pair, (AliDielectronSignalMC*)fSignalsMC->At(isig));
1205       if(isMCtruth) {
1206         //fill pair information
1207         if (pairClass){
1208           AliDielectronVarManager::Fill(pair, values);
1209           fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
1210         }
1211         //fill leg information, both + and - in the same histo
1212         if (legClass){
1213           AliDielectronVarManager::Fill(pair->GetFirstDaughter(),values);
1214           fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1215           AliDielectronVarManager::Fill(pair->GetSecondDaughter(),values);
1216           fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
1217         }
1218       } //is signal
1219     } //loop: pairs
1220   } //loop: MCsignals
1221
1222 }