PWGDQ/dielectron/AliDielectronCF.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 Correction framework manager                         //
  18 //                                                                       //
  19 /*
  20
  21
  22
  23
  24
  25
  26
  27
  28
  29 */
  30 //                                                                       //
  31 ///////////////////////////////////////////////////////////////////////////
  32
  33 #include <TList.h>
  34 #include <TObjArray.h>
  35 #include <TVectorD.h>
  36 #include <TString.h>
  37 #include <TObjString.h>
  38
  39 #include <AliCFContainer.h>
  40 #include <AliAnalysisFilter.h>
  41 #include <AliAnalysisCuts.h>
  42 #include <AliVParticle.h>
  43 #include <AliLog.h>
  44
  45 #include "AliDielectronCF.h"
  46 #include "AliDielectronMC.h"
  47 #include "AliDielectronPair.h"
  48 #include "AliDielectronSignalMC.h"
  49
  50 ClassImp(AliDielectronCF)
  51
  52 AliDielectronCF::AliDielectronCF() :
  53   TNamed("DielectronCF","DielectronCF"),
  54   fNSteps(0),
  55   fNVars(0),
  56   fVarBinLimits(0x0),
  57   fNVarsLeg(0),
  58   fVarBinLimitsLeg(0x0),
  59   fNCuts(0),
  60   fValues(0x0),
  61   fStepForMCtruth(kFALSE),
  62   fStepForNoCutsMCmotherPid(kFALSE),
  63   fStepForAfterAllCuts(kTRUE),
  64   fStepForPreFilter(kFALSE),
  65   fStepsForEachCut(kFALSE),
  66   fStepsForCutsIncreasing(kFALSE),
  67   fStepsForSignal(kTRUE),
  68   fStepsForBackground(kFALSE),
  69   fNStepMasks(0),
  70   fPdgMother(-1),
  71   fSignalsMC(0x0),
  72   fCfContainer(0x0),
  73   fHasMC(kFALSE),
  74   fNAddSteps(0)
  75 {
  76   //
  77   // Default constructor
  78   //
  79   for (Int_t i=0; i<AliDielectronVarManager::kNMaxValues; ++i){
  80     fVariables[i]=0;
  81     fVariablesLeg[i]=0;
  82   }
  83
  84   for (Int_t i=0; i<kNmaxAddSteps; ++i){
  85     fStepMasks[i]=0xFFFFFF;
  86   }
  87 }
  88
  89 //________________________________________________________________
  90 AliDielectronCF::AliDielectronCF(const char* name, const char* title) :
  91   TNamed(name, title),
  92   fNSteps(0),
  93   fNVars(0),
  94   fVarBinLimits(0x0),
  95   fNVarsLeg(0),
  96   fVarBinLimitsLeg(0x0),
  97   fNCuts(0),
  98   fValues(0x0),
  99   fStepForMCtruth(kFALSE),
 100   fStepForNoCutsMCmotherPid(kFALSE),
 101   fStepForAfterAllCuts(kTRUE),
 102   fStepForPreFilter(kFALSE),
 103   fStepsForEachCut(kFALSE),
 104   fStepsForCutsIncreasing(kFALSE),
 105   fStepsForSignal(kTRUE),
 106   fStepsForBackground(kFALSE),
 107   fNStepMasks(0),
 108   fPdgMother(-1),
 109   fSignalsMC(0x0),
 110   fCfContainer(0x0),
 111   fHasMC(kFALSE),
 112   fNAddSteps(0)
 113 {
 114   //
 115   // Named constructor
 116   //
 117   for (Int_t i=0; i<AliDielectronVarManager::kNMaxValues; ++i){
 118     fVariables[i]=0;
 119     fVariablesLeg[i]=0;
 120   }
 121
 122   for (Int_t i=0; i<kNmaxAddSteps; ++i){
 123     fStepMasks[i]=0xFFFFFF;
 124   }
 125 }
 126
 127 //________________________________________________________________
 128 AliDielectronCF::~AliDielectronCF()
 129 {
 130   //
 131   // Destructor
 132   //
 133   if (fValues) delete [] fValues;
 134   if (fVarBinLimits) delete fVarBinLimits;
 135   if (fVarBinLimitsLeg) delete fVarBinLimitsLeg;
 136 }
 137
 138 //________________________________________________________________
 139 void AliDielectronCF::AddVariable(AliDielectronVarManager::ValueTypes type, Int_t nbins,
 140                                   Double_t min, Double_t max, Bool_t leg, Bool_t log)
 141 {
 142   //
 143   // Add a variable to the CF configuration
 144   // if leg is true it will add the variables of the leg
 145   // if log is true log binning will be created
 146   //
 147
 148   TVectorD *binLimits=0x0;
 149   if (!log) binLimits=MakeLinBinning(nbins,min,max);
 150   else binLimits=MakeLogBinning(nbins,min,max);
 151   AddVariable(type,binLimits,leg);
 152 }
 153
 154 //________________________________________________________________
 155 void AliDielectronCF::AddVariable(AliDielectronVarManager::ValueTypes type, const char* binLimitStr, Bool_t leg/*=kFALSE*/)
 156 {
 157   //
 158   // Add a variable to the CF configuration
 159   // specify arbitrary binning in a string.
 160   // Bin limits need to be separated by a ","
 161   //
 162   TString limits(binLimitStr);
 163   if (limits.IsNull()){
 164     AliError(Form("Bin Limit string is empty, cannot add the variable '%s'",AliDielectronVarManager::GetValueName(type)));
 165     return;
 166   }
 167
 168   TObjArray *arr=limits.Tokenize(",");
 169   Int_t nLimits=arr->GetEntries();
 170   if (nLimits<2){
 171     AliError(Form("Need at leas 2 bin limits, cannot add the variable '%s'",AliDielectronVarManager::GetValueName(type)));
 172     delete arr;
 173     return;
 174   }
 175
 176   TVectorD *binLimits=new TVectorD(nLimits);
 177   for (Int_t iLim=0; iLim<nLimits; ++iLim){
 178     (*binLimits)[iLim]=(static_cast<TObjString*>(arr->At(iLim)))->GetString().Atof();
 179   }
 180
 181   delete arr;
 182   AddVariable(type,binLimits,leg);
 183 }
 184
 185 //________________________________________________________________
 186 void AliDielectronCF::AddVariable(AliDielectronVarManager::ValueTypes type, TVectorD *binLimits, Bool_t leg/*=kFALSE*/)
 187 {
 188   //
 189   // Add variable with the binning given in the TVectorD
 190   //
 191   if (!leg){
 192     if (!fVarBinLimits){
 193       fVarBinLimits=new TObjArray;
 194       fVarBinLimits->SetOwner();
 195     }
 196     fVarBinLimits->Add(binLimits);
 197     fVariables[fNVars]  = (UInt_t)type;
 198     ++fNVars;
 199   } else {
 200     if (!fVarBinLimitsLeg){
 201       fVarBinLimitsLeg=new TObjArray;
 202       fVarBinLimitsLeg->SetOwner();
 203     }
 204     fVarBinLimitsLeg->Add(binLimits);
 205     fVariablesLeg[fNVarsLeg]  = (UInt_t)type;
 206     ++fNVarsLeg;
 207   }
 208 }
 209
 210 //________________________________________________________________
 211 void AliDielectronCF::InitialiseContainer(const AliAnalysisFilter& filter)
 212 {
 213   //
 214   // Initialise container based on the cuts in the analysis filter
 215   //
 216
 217   fNCuts=filter.GetCuts()->GetEntries();
 218
 219   fHasMC=AliDielectronMC::Instance()->HasMC();
 220   fNAddSteps=1;
 221   if (fHasMC){
 222     if (fStepsForSignal && fSignalsMC) fNAddSteps+=fSignalsMC->GetEntries();
 223     if (fStepsForBackground) ++fNAddSteps;
 224   } else {
 225     //if
 226     fStepForMCtruth=kFALSE;
 227     fStepForNoCutsMCmotherPid=kFALSE;
 228     fStepsForSignal=kFALSE;
 229     fStepsForBackground=kFALSE;
 230   }
 231
 232   fNSteps=0;
 233   if (fStepForMCtruth && fSignalsMC) fNSteps+=fSignalsMC->GetEntries();
 234   if (fStepForNoCutsMCmotherPid && fSignalsMC) fNSteps+=fSignalsMC->GetEntries();
 235   if (fStepForAfterAllCuts) fNSteps+=fNAddSteps;
 236
 237   if (fStepsForEachCut&&fNCuts>1)        fNSteps+=(fNAddSteps*fNCuts);     //one step for each cut + Signal (MC)
 238   if (fStepsForCutsIncreasing&&fNCuts>2) fNSteps+=(fNAddSteps*(fNCuts-2)); //one step for the increasing cuts + Signal (MC)
 239                                                       // e.g. cut2&cut3, cut2&cut3&cut4
 240   fNSteps+=(fNAddSteps*fNStepMasks);                            // cuts for the additional cut masks
 241
 242   if (fStepForPreFilter) fNSteps+=fNAddSteps; //Add at the end for Prefilter (maxcutmask+1)
 243
 244   // create the container
 245
 246   Int_t *nbins=new Int_t[fNVars+2*fNVarsLeg];
 247   for (Int_t i=0;i<fNVars;++i) {
 248     Int_t nBins=(static_cast<TVectorD*>(fVarBinLimits->At(i)))->GetNrows()-1;
 249     nbins[i]=nBins;
 250   }
 251   for (Int_t i=0;i<fNVarsLeg;++i){
 252     Int_t nBins=(static_cast<TVectorD*>(fVarBinLimitsLeg->At(i)))->GetNrows()-1;
 253     nbins[i+fNVars]=nBins;
 254     nbins[i+fNVars+fNVarsLeg]=nBins;
 255   }
 256
 257   fCfContainer = new AliCFContainer(GetName(), GetTitle(), fNSteps, fNVars+2*fNVarsLeg, nbins);
 258   delete [] nbins;
 259
 260   // initialize the variables and their bin limits
 261   for (Int_t iVar=0; iVar<fNVars; iVar++) {
 262     UInt_t type=fVariables[iVar];
 263     Double_t *binLim = (static_cast<TVectorD*>(fVarBinLimits->At(iVar)))->GetMatrixArray();
 264
 265     fCfContainer->SetBinLimits(iVar, binLim);
 266     fCfContainer->SetVarTitle(iVar, AliDielectronVarManager::GetValueName(type));
 267   }
 268
 269   // initialize the variables and their bin limits for the Legs
 270   for (Int_t iVar=0; iVar<fNVarsLeg; iVar++) {
 271     UInt_t type=fVariablesLeg[iVar];
 272     Double_t *binLim=(static_cast<TVectorD*>(fVarBinLimitsLeg->At(iVar)))->GetMatrixArray();
 273
 274     //Leg1
 275     fCfContainer->SetBinLimits(iVar+fNVars, binLim);
 276     fCfContainer->SetVarTitle(iVar+fNVars, Form("Leg1_%s",AliDielectronVarManager::GetValueName(type)));
 277
 278     //Leg2
 279     fCfContainer->SetBinLimits(iVar+fNVars+fNVarsLeg, binLim);
 280     fCfContainer->SetVarTitle(iVar+fNVars+fNVarsLeg, Form("Leg2_%s",AliDielectronVarManager::GetValueName(type)));
 281   }
 282
 283   // array for storing values
 284   fValues = new Double_t[fNVars+2*fNVarsLeg];
 285
 286   //=================//
 287   // Set step titles //
 288   //=================//
 289   Int_t step=0;
 290
 291   //Pure MC truth
 292   if(fStepForMCtruth && fSignalsMC) {
 293     for(Int_t i=0; i<fSignalsMC->GetEntries(); i++)
 294       fCfContainer->SetStepTitle(step++, Form("MC truth (Signal: %s)", fSignalsMC->At(i)->GetTitle()));
 295   }
 296
 297   //before cuts (MC truth)
 298   if (fStepForNoCutsMCmotherPid && fSignalsMC){
 299     for(Int_t i=0; i<fSignalsMC->GetEntries(); i++)
 300       fCfContainer->SetStepTitle(step++,Form("No cuts (Signal: %s)",fSignalsMC->At(i)->GetTitle()));
 301   }
 302
 303   TString cutName;
 304   //Steps for each of the cuts
 305   if (fStepsForEachCut&&fNCuts>1){
 306     for (Int_t iCut=0; iCut<fNCuts;++iCut) {
 307       cutName=filter.GetCuts()->At(iCut)->GetName(); //TODO: User GetTitle???
 308
 309       fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
 310
 311       if (fHasMC){
 312         if (fStepsForSignal && fSignalsMC) {
 313           for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
 314             fCfContainer->SetStepTitle(step++, Form("%s (Signal: %s)", cutName.Data(), fSignalsMC->At(i)->GetTitle())); //Step for the cut with MC truth
 315           }
 316         }
 317         if (fStepsForBackground)
 318           fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
 319       }
 320     }
 321   }
 322
 323   //Steps for increasing cut match
 324   if (fStepsForCutsIncreasing&&fNCuts>2){
 325     cutName=filter.GetCuts()->At(0)->GetName(); //TODO: User GetTitle???
 326     for (Int_t iCut=1; iCut<fNCuts-1;++iCut) {
 327       cutName+="&";
 328       cutName+=filter.GetCuts()->At(iCut)->GetName();
 329
 330       fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
 331
 332       if (fHasMC){
 333         if (fStepsForSignal && fSignalsMC)
 334           for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
 335             fCfContainer->SetStepTitle(step++, Form("%s (Signal: %s)", cutName.Data(), fSignalsMC->At(i)->GetTitle())); //Step for the cut with MC truth
 336           }
 337         if (fStepsForBackground)
 338           fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
 339       }
 340     }
 341   }
 342
 343   //Steps of user defined cut combinations
 344   for (UInt_t iComb=0; iComb<fNStepMasks; ++iComb){
 345     cutName="";
 346     UInt_t mask=fStepMasks[iComb];
 347     for (Int_t iCut=0; iCut<fNCuts;++iCut) {
 348       if (mask&(1<<iCut)){
 349         if (cutName.IsNull()){
 350           cutName=filter.GetCuts()->At(iCut)->GetName();
 351         }else{
 352           cutName+="&";
 353           cutName+=filter.GetCuts()->At(iCut)->GetName();
 354         }
 355       }
 356     }
 357
 358     fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
 359
 360     if (fHasMC){
 361       if (fStepsForSignal && fSignalsMC)
 362         for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
 363           fCfContainer->SetStepTitle(step++, Form("%s (Signal: %s)", cutName.Data(), fSignalsMC->At(i)->GetTitle())); //Step for the cut with MC truth
 364         }
 365       if (fStepsForBackground)
 366         fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
 367     }
 368   }
 369
 370   //After All cuts
 371   if (fStepForAfterAllCuts){
 372     cutName="No pair cuts";
 373     if (filter.GetCuts()->At(0)){
 374       cutName=filter.GetCuts()->At(0)->GetName(); //TODO: User GetTitle???
 375       for (Int_t iCut=1; iCut<fNCuts;++iCut) {
 376         cutName+="&";
 377         cutName+=filter.GetCuts()->At(iCut)->GetName();
 378       }
 379     }
 380     fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
 381     if (fHasMC){
 382       if (fStepsForSignal && fSignalsMC)
 383         for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
 384           fCfContainer->SetStepTitle(step++, Form("%s (Signal: %s)", cutName.Data(), fSignalsMC->At(i)->GetTitle())); //Step for the cut with MC truth
 385         }
 386       if (fStepsForBackground)
 387         fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
 388     }
 389   }
 390
 391   //Additional Step for result after PreFilter
 392   if (fStepForPreFilter){
 393     cutName="PreFilter";
 394     fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
 395     if (fHasMC){
 396       if (fStepsForSignal && fSignalsMC)
 397         for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
 398           fCfContainer->SetStepTitle(step++, Form("%s (Signal %s)", cutName.Data(), fSignalsMC->At(i)->GetTitle())); //Step for the cut with MC truth
 399         }
 400       if (fStepsForBackground)
 401         fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
 402     }
 403   }
 404
 405
 406
 407   if (step!=fNSteps) {
 408     AliError(Form("Something went wrong in the naming of the steps!!! (%d != %d)",step,fNSteps));
 409   }
 410 }
 411
 412 //________________________________________________________________
 413 void AliDielectronCF::Fill(UInt_t mask, const AliDielectronPair *particle)
 414 {
 415   //
 416   // Fill the containers
 417   //
 418
 419   // Check the MC truths
 420   Bool_t* isMCTruth=0x0;
 421   if(fHasMC && fSignalsMC && fSignalsMC->GetEntries()>0) {
 422     isMCTruth=new Bool_t[fSignalsMC->GetEntries()];
 423     for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) isMCTruth[i]=kFALSE;
 424   }
 425
 426   Bool_t isBackground = kFALSE;
 427   if(fHasMC && isMCTruth) {
 428     for(Int_t i=0; i<fSignalsMC->GetEntries(); ++i) {
 429       isMCTruth[i] = AliDielectronMC::Instance()->IsMCTruth(particle, (AliDielectronSignalMC*)fSignalsMC->At(i));
 430       isBackground = (isBackground || isMCTruth[i]);
 431     }
 432     // background is considered that pair which does not fulfill any of the signals
 433     isBackground = !isBackground;
 434   }
 435
 436   Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
 437   AliDielectronVarManager::Fill(particle,valuesPair);
 438
 439   for (Int_t iVar=0; iVar<fNVars; ++iVar){
 440     Int_t var=fVariables[iVar];
 441     fValues[iVar]=valuesPair[var];
 442   }
 443
 444   if (fNVarsLeg>0){
 445     Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues];
 446     AliDielectronVarManager::Fill(particle->GetFirstDaughter(),valuesLeg1);
 447     Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues];
 448     AliDielectronVarManager::Fill(particle->GetSecondDaughter(),valuesLeg2);
 449
 450     for (Int_t iVar=0; iVar<fNVarsLeg; ++iVar){
 451       Int_t var=fVariablesLeg[iVar];
 452       fValues[iVar+fNVars]=valuesLeg1[var];
 453       fValues[iVar+fNVars+fNVarsLeg]=valuesLeg2[var];
 454     }
 455   }
 456
 457   UInt_t selectedMask=(1<<fNCuts)-1;
 458
 459   //============//
 460   // Fill steps //
 461   //============//
 462   // Pure MC steps are handled in FillMC
 463   Int_t step=0;
 464   if (fStepForMCtruth && isMCTruth) step+=fSignalsMC->GetEntries();
 465
 466   //No cuts (MC truth)
 467   if (fStepForNoCutsMCmotherPid && isMCTruth){
 468     for(Int_t i=0; i<fSignalsMC->GetEntries(); ++i) {
 469       if(isMCTruth[i]) {
 470         fCfContainer->Fill(fValues,step);
 471       }
 472       ++step;
 473     }
 474   }
 475
 476   //Steps for each of the cuts
 477   if (fStepsForEachCut&&fNCuts>1){
 478     for (Int_t iCut=0; iCut<fNCuts;++iCut) {
 479       if (mask&(1<<iCut)) {
 480         fCfContainer->Fill(fValues,step);
 481         ++step;
 482
 483         if (fHasMC){
 484           if ( fStepsForSignal && isMCTruth){
 485             for(Int_t i=0; i<fSignalsMC->GetEntries(); ++i) {
 486               if(isMCTruth[i]) {
 487                 fCfContainer->Fill(fValues,step);
 488               }
 489               ++step;
 490             }
 491           }
 492           if ( fStepsForBackground ){
 493             if (isBackground) fCfContainer->Fill(fValues,step);
 494             ++step;
 495           }
 496         }
 497       } else {
 498         step+=fNAddSteps;
 499       }
 500     }
 501   }
 502
 503
 504   //Steps for increasing cut match
 505   if (fStepsForCutsIncreasing&&fNCuts>2){
 506     for (Int_t iCut=1; iCut<fNCuts-1;++iCut) {
 507       if (mask&(1<<((iCut+1)-1))) {
 508         fCfContainer->Fill(fValues,step);
 509         ++step;
 510
 511         if (fHasMC){
 512           if ( fStepsForSignal && isMCTruth){
 513             for(Int_t i=0; i<fSignalsMC->GetEntries(); ++i) {
 514               if(isMCTruth[i]) {
 515                 fCfContainer->Fill(fValues,step);
 516               }
 517               ++step;
 518             }
 519           }
 520           if ( fStepsForBackground ){
 521             if (isBackground) fCfContainer->Fill(fValues,step);
 522             ++step;
 523           }
 524         }
 525       } else {
 526         step+=fNAddSteps;
 527       }
 528     }
 529   }
 530
 531   //Steps of user defined cut combinations
 532   for (UInt_t iComb=0; iComb<fNStepMasks; ++iComb){
 533     UInt_t userMask=fStepMasks[iComb];
 534     if (mask&userMask) {
 535       fCfContainer->Fill(fValues,step);
 536       ++step;
 537
 538       if (fHasMC){
 539         if ( fStepsForSignal && isMCTruth){
 540           for(Int_t i=0; i<fSignalsMC->GetEntries(); ++i) {
 541             if(isMCTruth[i]) {
 542               fCfContainer->Fill(fValues,step);
 543             }
 544             ++step;
 545           }
 546         }
 547         if ( fStepsForBackground ){
 548           if (isBackground) fCfContainer->Fill(fValues,step);
 549           ++step;
 550         }
 551       }
 552     } else {
 553       step+=fNAddSteps;
 554     }
 555   }
 556
 557   //All cuts
 558   if (fStepForAfterAllCuts){
 559     if (mask == selectedMask){
 560       fCfContainer->Fill(fValues,step);
 561       ++step;
 562
 563       if (fHasMC){
 564         if ( fStepsForSignal && isMCTruth){
 565           for(Int_t i=0; i<fSignalsMC->GetEntries(); ++i) {
 566             if(isMCTruth[i]) {
 567               fCfContainer->Fill(fValues,step);
 568             }
 569             ++step;
 570           }
 571         }
 572         if ( fStepsForBackground ){
 573           if (isBackground) fCfContainer->Fill(fValues,step);
 574           ++step;
 575         }
 576       }
 577     } else {
 578       step+=fNAddSteps;
 579     }
 580   }
 581   if (fStepForPreFilter) {
 582         if (mask&(1<<fNCuts)) {
 583           fCfContainer->Fill(fValues,step);
 584           ++step;
 585
 586           if (fHasMC){
 587             if ( fStepsForSignal && isMCTruth){
 588               for(Int_t i=0; i<fSignalsMC->GetEntries(); ++i) {
 589                 if(isMCTruth[i]) {
 590                   fCfContainer->Fill(fValues,step);
 591                 }
 592                 ++step;
 593               }
 594             }
 595             if ( fStepsForBackground ){
 596               if (isBackground) fCfContainer->Fill(fValues,step);
 597               ++step;
 598             }
 599           }
 600         }
 601         else {
 602           step+=fNAddSteps;
 603         }
 604   }
 605   if (step!=fNSteps) {
 606     AliError("Something went wrong in the step filling!!!");
 607   }
 608   if(isMCTruth) delete [] isMCTruth;
 609 }
 610
 611 //________________________________________________________________
 612 void AliDielectronCF::FillMC(const TObject *particle)
 613 {
 614   //
 615   // fill MC part of the Container
 616   //
 617   if (!fStepForMCtruth) return;
 618
 619   Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
 620   AliDielectronVarManager::Fill(particle,valuesPair);
 621
 622   AliVParticle *d1=0x0;
 623   AliVParticle *d2=0x0;
 624   AliDielectronMC::Instance()->GetDaughters(particle,d1,d2);
 625
 626   valuesPair[AliDielectronVarManager::kThetaHE]=AliDielectronPair::ThetaPhiCM(d1,d2,kTRUE,kTRUE);
 627   valuesPair[AliDielectronVarManager::kPhiHE]=AliDielectronPair::ThetaPhiCM(d1,d2,kTRUE,kFALSE);
 628   valuesPair[AliDielectronVarManager::kThetaCS]=AliDielectronPair::ThetaPhiCM(d1,d2,kFALSE,kTRUE);
 629   valuesPair[AliDielectronVarManager::kPhiCS]=AliDielectronPair::ThetaPhiCM(d1,d2,kFALSE,kFALSE);
 630
 631   //TODO: temporary solution, set manually the pair type to 1: unlikesign SE
 632   valuesPair[AliDielectronVarManager::kPairType]=1;
 633
 634   for (Int_t iVar=0; iVar<fNVars; ++iVar){
 635     Int_t var=fVariables[iVar];
 636     fValues[iVar]=valuesPair[var];
 637   }
 638
 639   if (fNVarsLeg>0){
 640     Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues];
 641     Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues];
 642     if (d1->Pt()>d2->Pt()){
 643       AliDielectronVarManager::Fill(d1,valuesLeg1);
 644       AliDielectronVarManager::Fill(d2,valuesLeg2);
 645     } else {
 646       AliDielectronVarManager::Fill(d2,valuesLeg1);
 647       AliDielectronVarManager::Fill(d1,valuesLeg2);
 648     }
 649
 650     for (Int_t iVar=0; iVar<fNVarsLeg; ++iVar){
 651       Int_t var=fVariablesLeg[iVar];
 652       fValues[iVar+fNVars]=valuesLeg1[var];
 653       fValues[iVar+fNVars+fNVarsLeg]=valuesLeg2[var];
 654     }
 655   }
 656
 657   fCfContainer->Fill(fValues,0);
 658 }
 659
 660
 661 //________________________________________________________________
 662 void AliDielectronCF::FillMC(Int_t label1, Int_t label2, Int_t nSignal) {
 663   //
 664   // fill the pure MC part of the container starting from a pair of 2 particles (part1 and part2 are legs)
 665   //
 666   if (!fStepForMCtruth) return;
 667
 668   AliVParticle* part1 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label1);
 669   AliVParticle* part2 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label2);
 670
 671   AliDielectronMC* dieMC = AliDielectronMC::Instance();
 672
 673   Int_t mLabel1 = dieMC->GetMothersLabel(label1);    // should work for both ESD and AOD
 674   Int_t mLabel2 = dieMC->GetMothersLabel(label2);
 675   // check the same mother option
 676   AliDielectronSignalMC* sigMC = (AliDielectronSignalMC*)fSignalsMC->At(nSignal);
 677   if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kSame && mLabel1!=mLabel2) return;
 678   if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kDifferent && mLabel1==mLabel2) return;
 679
 680   // fill the leg variables
 681   if (fNVarsLeg>0){
 682     Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues];
 683     Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues];
 684     if (part1->Pt()>part2->Pt()){
 685       AliDielectronVarManager::Fill(part1,valuesLeg1);
 686       AliDielectronVarManager::Fill(part2,valuesLeg2);
 687     } else {
 688       AliDielectronVarManager::Fill(part2,valuesLeg1);
 689       AliDielectronVarManager::Fill(part1,valuesLeg2);
 690     }
 691
 692     for (Int_t iVar=0; iVar<fNVarsLeg; ++iVar){
 693       Int_t var=fVariablesLeg[iVar];
 694       fValues[iVar+fNVars]=valuesLeg1[var];
 695       fValues[iVar+fNVars+fNVarsLeg]=valuesLeg2[var];
 696     }
 697   }
 698
 699   Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
 700   AliDielectronVarManager::FillVarMCParticle2(part1,part2,valuesPair);
 701
 702   valuesPair[AliDielectronVarManager::kThetaHE]=AliDielectronPair::ThetaPhiCM(part1,part2,kTRUE,kTRUE);
 703   valuesPair[AliDielectronVarManager::kPhiHE]=AliDielectronPair::ThetaPhiCM(part1,part2,kTRUE,kFALSE);
 704   valuesPair[AliDielectronVarManager::kThetaCS]=AliDielectronPair::ThetaPhiCM(part1,part2,kFALSE,kTRUE);
 705   valuesPair[AliDielectronVarManager::kPhiCS]=AliDielectronPair::ThetaPhiCM(part1,part2,kFALSE,kFALSE);
 706
 707   if(part1->Charge()*part2->Charge()<0)
 708     valuesPair[AliDielectronVarManager::kPairType]=1;
 709   else if(part1->Charge()>0)
 710     valuesPair[AliDielectronVarManager::kPairType]=0;
 711   else
 712     valuesPair[AliDielectronVarManager::kPairType]=2;
 713
 714   for(Int_t iVar=0; iVar<fNVars; ++iVar){
 715     Int_t var=fVariables[iVar];
 716     fValues[iVar]=valuesPair[var];
 717   }
 718
 719   fCfContainer->Fill(fValues,nSignal);
 720 }
 721
 722
 723 //_____________________________________________________________________________
 724 TVectorD* AliDielectronCF::MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax) const
 725 {
 726   //
 727   // Make logarithmic binning
 728   // the user has to delete the array afterwards!!!
 729   //
 730
 731   //check limits
 732   if (xmin<1e-20 || xmax<1e-20){
 733     AliError("For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
 734     return MakeLinBinning(nbinsX, xmin, xmax);
 735   }
 736   if (xmax<xmin){
 737     Double_t tmp=xmin;
 738     xmin=xmax;
 739     xmax=tmp;
 740   }
 741   TVectorD *binLim=new TVectorD(nbinsX+1);
 742   Double_t first=xmin;
 743   Double_t last=xmax;
 744   Double_t expMax=TMath::Log(last/first);
 745   for (Int_t i=0; i<nbinsX+1; ++i){
 746     (*binLim)[i]=first*TMath::Exp(expMax/nbinsX*(Double_t)i);
 747   }
 748   return binLim;
 749 }
 750
 751 //_____________________________________________________________________________
 752 TVectorD* AliDielectronCF::MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax) const
 753 {
 754   //
 755   // Make logarithmic binning
 756   // the user has to delete the array afterwards!!!
 757   //
 758   if (xmax<xmin){
 759     Double_t tmp=xmin;
 760     xmin=xmax;
 761     xmax=tmp;
 762   }
 763   TVectorD *binLim=new TVectorD(nbinsX+1);
 764   Double_t first=xmin;
 765   Double_t last=xmax;
 766   Double_t binWidth=(last-first)/nbinsX;
 767   for (Int_t i=0; i<nbinsX+1; ++i){
 768     (*binLim)[i]=first+binWidth*(Double_t)i;
 769   }
 770   return binLim;
 771 }
 772