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