PWGDQ/dielectron/AliDielectronHF.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 HF                                  //
  18 //                                                                       //
  19 //                                                                       //
  20 /*
  21 Detailed description
  22
  23
  24 */
  25 //                                                                       //
  26 ///////////////////////////////////////////////////////////////////////////
  27
  28 #include <TVectorD.h>
  29 #include <TH1.h>
  30 #include <TH1F.h>
  31 #include <TH2.h>
  32 #include <TH3.h>
  33 #include <TProfile.h>
  34 #include <TProfile2D.h>
  35 #include <TProfile3D.h>
  36 #include <THnSparse.h>
  37 #include <TAxis.h>
  38 #include <TString.h>
  39 #include <TObjString.h>
  40 #include <TObjArray.h>
  41
  42 #include <AliVParticle.h>
  43 #include <AliLog.h>
  44
  45 #include "AliDielectron.h"
  46 #include "AliDielectronHelper.h"
  47 #include "AliDielectronMC.h"
  48 #include "AliDielectronPair.h"
  49 #include "AliDielectronSignalMC.h"
  50
  51 #include "AliDielectronHistos.h"
  52 #include "AliDielectronHF.h"
  53
  54 ClassImp(AliDielectronHF)
  55
  56 AliDielectronHF::AliDielectronHF() :
  57   TNamed(),
  58   fUsedVars(new TBits(AliDielectronVarManager::kNMaxValues)),
  59   fArrPairType(),
  60   fPairType(kSeOnlyOS),
  61   fSignalsMC(0x0),
  62   fVarCutType(new TBits(kMaxCuts)),
  63   fAxes(kMaxCuts),
  64   fHasMC(kFALSE),
  65   fStepGenerated(kFALSE),
  66   fRefObj(1)
  67 {
  68   //
  69   // Default Constructor
  70   //
  71   for (Int_t i=0; i<kMaxCuts; ++i){
  72     fVarCuts[i]=0;
  73     //    fVarCutType[i]=0;
  74     fBinType[i]=kStdBin;
  75   }
  76   fAxes.SetOwner(kTRUE);
  77   fRefObj.SetOwner(kTRUE);
  78   fArrPairType.SetOwner(kTRUE);
  79 }
  80
  81 //______________________________________________
  82 AliDielectronHF::AliDielectronHF(const char* name, const char* title) :
  83   TNamed(name, title),
  84   fUsedVars(new TBits(AliDielectronVarManager::kNMaxValues)),
  85   fArrPairType(),
  86   fPairType(kSeOnlyOS),
  87   fSignalsMC(0x0),
  88   fVarCutType(new TBits(kMaxCuts)),
  89   fAxes(kMaxCuts),
  90   fHasMC(kFALSE),
  91   fStepGenerated(kFALSE),
  92   fRefObj(1)
  93 {
  94   //
  95   // Named Constructor
  96   //
  97   for (Int_t i=0; i<kMaxCuts; ++i){
  98     fVarCuts[i]=0;
  99     //    fVarCutType[i]=0;
 100     fBinType[i]=kStdBin;
 101   }
 102   fAxes.SetOwner(kTRUE);
 103   fRefObj.SetOwner(kTRUE);
 104   fArrPairType.SetOwner(kTRUE);
 105 }
 106
 107 //______________________________________________
 108 AliDielectronHF::~AliDielectronHF()
 109 {
 110   //
 111   // Default Destructor
 112   //
 113   if(fUsedVars)   delete fUsedVars;
 114   if(fVarCutType) delete fVarCutType;
 115   fAxes.Delete();
 116   fRefObj.Delete();
 117   fArrPairType.Delete();
 118 }
 119
 120 //_____________________________________________________________________________
 121 void AliDielectronHF::UserProfile(const char* histClass, UInt_t valTypeP,
 122                                       const TVectorD * const binsX,
 123                                       UInt_t valTypeX, TString option, UInt_t valTypeW)
 124 {
 125   //
 126   // Histogram creation 1D case with arbitraty binning X
 127   // the TVectorD is assumed to be surplus after the creation and will be deleted!!!
 128   //
 129
 130   TH1 *hist=0x0;
 131   if(valTypeP==AliDielectronHistos::kNoProfile)
 132     hist=new TH1F("","",binsX->GetNrows()-1,binsX->GetMatrixArray());
 133   else {
 134     TString opt=""; Double_t pmin=0., pmax=0.;
 135     if(!option.IsNull()) {
 136       TObjArray *arr=option.Tokenize(";");
 137       arr->SetOwner();
 138       opt=((TObjString*)arr->At(0))->GetString();
 139       if(arr->GetEntriesFast()>1) pmin=(((TObjString*)arr->At(1))->GetString()).Atof();
 140       if(arr->GetEntriesFast()>2) pmax=(((TObjString*)arr->At(2))->GetString()).Atof();
 141       delete arr;
 142     }
 143     hist=new TProfile("","",binsX->GetNrows()-1,binsX->GetMatrixArray(),pmin,pmax,opt.Data());
 144     //      printf(" name %s PROFILE options: pmin %.1f pmax %.1f err %s \n",name,((TProfile*)hist)->GetYmin(),((TProfile*)hist)->GetYmax(),((TProfile*)hist)->GetErrorOption() );
 145   }
 146
 147   // store variales in axes
 148   UInt_t valType[4] = {0};
 149   valType[0]=valTypeX;     valType[1]=valTypeP;
 150   AliDielectronHistos::StoreVariables(hist, valType);
 151   hist->SetUniqueID(valTypeW); // store weighting variable
 152
 153   for(Int_t i=0; i<4; i++)   fUsedVars->SetBitNumber(valType[i],kTRUE);
 154   fUsedVars->SetBitNumber(valTypeW,kTRUE);
 155
 156   // adapt the name and title of the histogram in case they are empty
 157   AliDielectronHistos::AdaptNameTitle(hist, histClass);
 158   hist->SetName(Form("HF_%s",hist->GetName()));
 159
 160   fRefObj.AddLast(hist);
 161   delete binsX;
 162 }
 163
 164 //_____________________________________________________________________________
 165 void AliDielectronHF::UserProfile(const char* histClass, UInt_t valTypeP,
 166                                       const TVectorD * const binsX, const TVectorD * const binsY,
 167                                       UInt_t valTypeX, UInt_t valTypeY, TString option, UInt_t valTypeW)
 168 {
 169   //
 170   // Histogram creation 2D case with arbitraty binning X and Y
 171   // the TVectorD is assumed to be surplus after the creation and will be deleted!!!
 172   //
 173
 174   TH1 *hist=0x0;
 175   if(valTypeP==AliDielectronHistos::kNoProfile) {
 176     hist=new TH2F("","",
 177                   binsX->GetNrows()-1,binsX->GetMatrixArray(),
 178                   binsY->GetNrows()-1,binsY->GetMatrixArray());
 179   }
 180   else  {
 181     TString opt=""; Double_t pmin=0., pmax=0.;
 182     if(!option.IsNull()) {
 183       TObjArray *arr=option.Tokenize(";");
 184       arr->SetOwner();
 185       opt=((TObjString*)arr->At(0))->GetString();
 186       if(arr->GetEntriesFast()>1) pmin=(((TObjString*)arr->At(1))->GetString()).Atof();
 187       if(arr->GetEntriesFast()>2) pmax=(((TObjString*)arr->At(2))->GetString()).Atof();
 188       delete arr;
 189     }
 190     hist=new TProfile2D("","",
 191                         binsX->GetNrows()-1,binsX->GetMatrixArray(),
 192                         binsY->GetNrows()-1,binsY->GetMatrixArray());
 193     ((TProfile2D*)hist)->BuildOptions(pmin,pmax,opt.Data());
 194   }
 195
 196   // store variales in axes
 197   UInt_t valType[4] = {0};
 198   valType[0]=valTypeX;     valType[1]=valTypeY; valType[2]=valTypeP;
 199   AliDielectronHistos::StoreVariables(hist, valType);
 200   hist->SetUniqueID(valTypeW); // store weighting variable
 201
 202   for(Int_t i=0; i<4; i++)   fUsedVars->SetBitNumber(valType[i],kTRUE);
 203   fUsedVars->SetBitNumber(valTypeW,kTRUE);
 204
 205   // adapt the name and title of the histogram in case they are empty
 206   AliDielectronHistos::AdaptNameTitle(hist, histClass);
 207   hist->SetName(Form("HF_%s",hist->GetName()));
 208
 209   fRefObj.AddLast(hist);
 210   delete binsX;
 211   delete binsY;
 212 }
 213
 214 //_____________________________________________________________________________
 215 void AliDielectronHF::UserProfile(const char* histClass, UInt_t valTypeP,
 216                                       const TVectorD * const binsX, const TVectorD * const binsY, const TVectorD * const binsZ,
 217                                       UInt_t valTypeX, UInt_t valTypeY, UInt_t valTypeZ, TString option, UInt_t valTypeW)
 218 {
 219   //
 220   // Histogram creation 3D case with arbitraty binning X, Y, Z
 221   // the TVectorD is assumed to be surplus after the creation and will be deleted!!!
 222   //
 223   TH1 *hist=0x0;
 224   if(valTypeP==AliDielectronHistos::kNoProfile) {
 225     hist=new TH3F("","",
 226                   binsX->GetNrows()-1,binsX->GetMatrixArray(),
 227                   binsY->GetNrows()-1,binsY->GetMatrixArray(),
 228                   binsZ->GetNrows()-1,binsZ->GetMatrixArray());
 229   }
 230   else {
 231     TString opt=""; Double_t pmin=0., pmax=0.;
 232     if(!option.IsNull()) {
 233       TObjArray *arr=option.Tokenize(";");
 234       arr->SetOwner();
 235       opt=((TObjString*)arr->At(0))->GetString();
 236       if(arr->GetEntriesFast()>1) pmin=(((TObjString*)arr->At(1))->GetString()).Atof();
 237       if(arr->GetEntriesFast()>2) pmax=(((TObjString*)arr->At(2))->GetString()).Atof();
 238       delete arr;
 239     }
 240     hist=new TProfile3D("","",
 241                         binsX->GetNrows()-1,binsX->GetMatrixArray(),
 242                         binsY->GetNrows()-1,binsY->GetMatrixArray(),
 243                         binsZ->GetNrows()-1,binsZ->GetMatrixArray());
 244     ((TProfile3D*)hist)->BuildOptions(pmin,pmax,opt.Data());
 245   }
 246
 247   // store variales in axes
 248   UInt_t valType[4] = {0};
 249   valType[0]=valTypeX;     valType[1]=valTypeY;     valType[2]=valTypeZ;     valType[3]=valTypeP;
 250   AliDielectronHistos::StoreVariables(hist, valType);
 251   hist->SetUniqueID(valTypeW); // store weighting variable
 252
 253   for(Int_t i=0; i<4; i++)   fUsedVars->SetBitNumber(valType[i],kTRUE);
 254   fUsedVars->SetBitNumber(valTypeW,kTRUE);
 255
 256   // adapt the name and title of the histogram in case they are empty
 257   AliDielectronHistos::AdaptNameTitle(hist, histClass);
 258   hist->SetName(Form("HF_%s",hist->GetName()));
 259
 260   fRefObj.AddLast(hist);
 261   delete binsX;
 262   delete binsY;
 263   delete binsZ;
 264 }
 265
 266 //_____________________________________________________________________________
 267 void AliDielectronHF::UserSparse(const char* histClass, Int_t ndim, TObjArray *limits, UInt_t *vars, UInt_t valTypeW)
 268 {
 269   //
 270   // THnSparse creation with non-linear binning
 271   //
 272
 273   THnSparseF *hist=0;
 274   Int_t bins[ndim];
 275   // get number of bins
 276   for(Int_t idim=0 ;idim<ndim; idim++) {
 277     TVectorD *vec = (TVectorD*) limits->At(idim);
 278     bins[idim]=vec->GetNrows()-1;
 279   }
 280
 281   hist=new THnSparseF("",histClass, ndim, bins, 0x0, 0x0);
 282
 283   // set binning
 284   for(Int_t idim=0 ;idim<ndim; idim++) {
 285     TVectorD *vec = (TVectorD*) limits->At(idim);
 286     hist->SetBinEdges(idim,vec->GetMatrixArray());
 287   }
 288
 289   // store variales in axes
 290   AliDielectronHistos::StoreVariables(hist, vars);
 291   hist->SetUniqueID(valTypeW); // store weighting variable
 292
 293   // store which variables are used
 294   for(Int_t i=0; i<20; i++)   fUsedVars->SetBitNumber(vars[i],kTRUE);
 295   fUsedVars->SetBitNumber(valTypeW,kTRUE);
 296
 297   // adapt the name and title of the histogram in case they are empty
 298   TString name;
 299   for(Int_t iv=0; iv < ndim; iv++) name+=Form("%s_",AliDielectronVarManager::GetValueName(vars[iv]));
 300   name.Resize(name.Length()-1);
 301   hist->SetName(Form("HF_%s",name.Data()));
 302
 303   fRefObj.AddLast(hist);
 304   delete limits;
 305
 306 }
 307
 308 //________________________________________________________________
 309 void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
 310                                      Int_t nbins, Double_t min, Double_t max, Bool_t log, Bool_t leg, EBinType btype)
 311 {
 312   //
 313   // Add a variable to the mixing handler
 314   //
 315
 316   // limit number of variables to kMaxCuts
 317   if (fAxes.GetEntriesFast()>=kMaxCuts) return;
 318
 319   TVectorD *binLimits=0x0;
 320   if (!log) binLimits=AliDielectronHelper::MakeLinBinning(nbins,min,max);
 321   else binLimits=AliDielectronHelper::MakeLogBinning(nbins,min,max);
 322   if (!binLimits) return;
 323
 324   Int_t size=fAxes.GetEntriesFast();
 325   fVarCuts[size]=(UShort_t)type;
 326   //  fVarCutType[size]=leg;
 327   fVarCutType->SetBitNumber(size,leg);
 328   fAxes.Add(binLimits->Clone());
 329   fBinType[size]=btype;
 330   fUsedVars->SetBitNumber(type,kTRUE);
 331 }
 332
 333 //________________________________________________________________
 334 void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
 335                                              const char* binLimitStr, Bool_t leg, EBinType btype)
 336 {
 337   //
 338   // Add a variable to the mixing handler with arbitrary binning
 339   //
 340
 341   // limit number of variables to kMaxCuts
 342   if (fAxes.GetEntriesFast()>=kMaxCuts) return;
 343
 344   TVectorD *binLimits=AliDielectronHelper::MakeArbitraryBinning(binLimitStr);
 345   if (!binLimits) return;
 346
 347   Int_t size=fAxes.GetEntriesFast();
 348   fVarCuts[size]=(UShort_t)type;
 349   //  fVarCutType[size]=leg;
 350   fVarCutType->SetBitNumber(size,leg);
 351   fAxes.Add(binLimits);
 352   fBinType[size]=btype;
 353   fUsedVars->SetBitNumber(type,kTRUE);
 354 }
 355
 356 //________________________________________________________________
 357 void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
 358                                              TVectorD * binLimits, Bool_t leg, EBinType btype)
 359 {
 360   //
 361   // Add a variable to the mixing handler with a vector
 362   // the TVectorD is assumed to be surplus after the creation and will be deleted!!!
 363   //
 364
 365   // limit number of variables to kMaxCuts
 366   if (fAxes.GetEntriesFast()>=kMaxCuts) return;
 367
 368   if (!binLimits) return;
 369
 370   Int_t size=fAxes.GetEntriesFast();
 371   fVarCuts[size]=(UShort_t)type;
 372   //  fVarCutType[size]=leg;
 373   fVarCutType->SetBitNumber(size,leg);
 374   fAxes.Add(binLimits);
 375   fBinType[size]=btype;
 376   fUsedVars->SetBitNumber(type,kTRUE);
 377 }
 378
 379 //______________________________________________
 380 void AliDielectronHF::Fill(Int_t label1, Int_t label2, Int_t nSignal)
 381 {
 382   //
 383   // fill the pure MC part of the container starting from a pair of 2 particles (part1 and part2 are legs)
 384   //
 385   // fill only if we have asked for these steps
 386   if(!fStepGenerated) return;
 387
 388   AliVParticle* part1 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label1);
 389   AliVParticle* part2 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label2);
 390   if(!part1 || !part2) return;
 391
 392   AliDielectronMC* dieMC = AliDielectronMC::Instance();
 393
 394   Int_t mLabel1 = dieMC->GetMothersLabel(label1);    // should work for both ESD and AOD
 395   Int_t mLabel2 = dieMC->GetMothersLabel(label2);
 396
 397   // check the same mother option
 398   AliDielectronSignalMC* sigMC = (AliDielectronSignalMC*)fSignalsMC->At(nSignal);
 399   if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kSame && mLabel1!=mLabel2) return;
 400   if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kDifferent && mLabel1==mLabel2) return;
 401
 402   AliDielectronVarManager::SetFillMap(fUsedVars);
 403   // fill the leg variables
 404   Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues];
 405   Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues];
 406   AliDielectronVarManager::Fill(part1,valuesLeg1);
 407   AliDielectronVarManager::Fill(part2,valuesLeg2);
 408
 409   // fill the pair and event variables
 410   Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
 411   AliDielectronVarManager::Fill(dieMC->GetMCEvent(), valuesPair);
 412   AliDielectronVarManager::FillVarMCParticle2(part1,part2,valuesPair);
 413
 414   // if pair types are filled, fill mc sources at the end
 415   Int_t istep=0;
 416   if(fPairType!=kMConly) istep=AliDielectron::kEv1PMRot+1;
 417
 418   // only OS at the moment
 419   if(part1->Charge()*part2->Charge()<0) {
 420     Fill(istep+nSignal+fSignalsMC->GetEntries(), valuesPair,  valuesLeg1, valuesLeg2);
 421   }
 422
 423   return;
 424 }
 425 //______________________________________________
 426 void AliDielectronHF::Fill(Int_t pairIndex, const AliDielectronPair *particle)
 427 {
 428   //
 429   // fill histograms for event, pair and daughter cuts and pair types
 430   //
 431
 432   // only OS pairs in case of MC
 433   //////////////////////////////  if(fHasMC && pairIndex!=AliDielectron::kEv1PM) return;
 434
 435   // only selected pair types in case of data
 436   if(!IsPairTypeSelected(pairIndex)) return;
 437
 438   // get event and pair variables
 439   Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
 440   AliDielectronVarManager::SetFillMap(fUsedVars);
 441   AliDielectronVarManager::Fill(particle,valuesPair);
 442
 443   // get leg variables (TODO: do not fill for the moment since leg cuts are not opened)
 444   Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues]={0};
 445   if(fVarCutType->CountBits())  AliDielectronVarManager::Fill(particle->GetFirstDaughterP(),valuesLeg1);
 446   Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues]={0};
 447   if(fVarCutType->CountBits())  AliDielectronVarManager::Fill(particle->GetSecondDaughterP(),valuesLeg2);
 448
 449   // fill
 450
 451   // if pair types are filled, fill mc sources at the end
 452   Int_t istep = 0;
 453   if(fPairType!=kMConly) istep=AliDielectron::kEv1PMRot+1;
 454
 455   // mc source steps (only OS SE pairs)
 456   if(fHasMC && fSignalsMC && pairIndex==AliDielectron::kEv1PM) {
 457     for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
 458       if(AliDielectronMC::Instance()->IsMCTruth(particle, (AliDielectronSignalMC*)fSignalsMC->At(i)))
 459         Fill(istep+i, valuesPair,  valuesLeg1, valuesLeg2);
 460     }
 461   }
 462
 463   // all pair types w/o use of mc information
 464   if(fPairType==kMConly) return;
 465
 466   // remove comments
 467   //// select correct step if we are looking at signals too
 468   ////  if(fHasMC && fSignalsMC) pairIndex += ( fSignalsMC->GetEntries() * (fStepGenerated ? 2 : 1) );
 469   Fill(pairIndex, valuesPair,  valuesLeg1, valuesLeg2);
 470
 471   return;
 472 }
 473
 474 //______________________________________________
 475 void AliDielectronHF::Fill(Int_t index, Double_t * const valuesPair, Double_t * const valuesLeg1, Double_t * const valuesLeg2)
 476 {
 477   //
 478   // main fill function using index and values as input
 479   //
 480
 481   TObjArray *histArr = static_cast<TObjArray*>(fArrPairType.At(index));
 482   if(!histArr) return;
 483
 484   Int_t size  = GetNumberOfBins();
 485   // loop over all histograms
 486   for(Int_t ihist=0; ihist<size; ihist++) {
 487
 488     Int_t sizeAdd   = 1;
 489     Bool_t selected = kTRUE;
 490
 491     // loop over all cut variables
 492     Int_t nvars = fAxes.GetEntriesFast();
 493     for(Int_t ivar=0; ivar<nvars; ivar++) {
 494
 495       // get bin limits
 496       TVectorD *bins = static_cast<TVectorD*>(fAxes.At(ivar));
 497       Int_t nbins    = bins->GetNrows()-1;
 498
 499       // bin limits for current ivar bin
 500       Int_t ibin   = (ihist/sizeAdd)%nbins;
 501       Double_t lowEdge = (*bins)[ibin];
 502       Double_t upEdge  = (*bins)[ibin+1];
 503       switch(fBinType[ivar]) {
 504       case kStdBin:     upEdge=(*bins)[ibin+1];     break;
 505       case kBinToMax:   upEdge=(*bins)[nbins];      break;
 506       case kBinFromMin: lowEdge=(*bins)[0];         break;
 507       case kSymBin:     upEdge=(*bins)[nbins-ibin];
 508         if(ibin>=((Double_t)(nbins+1))/2) upEdge=(*bins)[nbins]; // to avoid low>up
 509         break;
 510       }
 511
 512       // leg variable
 513       if(fVarCutType->TestBitNumber(ivar)) {
 514         if( (valuesLeg1[fVarCuts[ivar]] < lowEdge || valuesLeg1[fVarCuts[ivar]] >= upEdge) ||
 515             (valuesLeg2[fVarCuts[ivar]] < lowEdge || valuesLeg2[fVarCuts[ivar]] >= upEdge) ) {
 516           selected=kFALSE;
 517           break;
 518         }
 519       }
 520       else { // pair and event variables
 521         if( (valuesPair[fVarCuts[ivar]] < lowEdge || valuesPair[fVarCuts[ivar]] >= upEdge) ) {
 522           selected=kFALSE;
 523           break;
 524         }
 525       }
 526
 527       sizeAdd*=nbins;
 528     } //end of var cut loop
 529
 530     // do not fill the histogram
 531     if(!selected) continue;
 532
 533     // fill the object with Pair and event values
 534     TObjArray *tmp = (TObjArray*) histArr->At(ihist);
 535     TString title = tmp->GetName();
 536     AliDebug(10,title.Data());
 537     for(Int_t i=0; i<tmp->GetEntriesFast(); i++) {
 538       AliDielectronHistos::FillValues(tmp->At(i), valuesPair);
 539     }
 540     //    AliDebug(10,Form("Fill var %d %s value %f in %s \n",fVar,AliDielectronVarManager::GetValueName(fVar),valuesPair[fVar],tmp->GetName()));
 541   } //end of hist loop
 542
 543 }
 544
 545 //______________________________________________
 546 void AliDielectronHF::Init()
 547 {
 548   //
 549   // initialise event buffers
 550   //
 551
 552   // has MC signals
 553   fHasMC=AliDielectronMC::Instance()->HasMC();
 554   Int_t steps = 0;
 555   if(fHasMC) steps=fSignalsMC->GetEntries();
 556   if(fStepGenerated) steps*=2;
 557
 558   // init pair type array
 559   fArrPairType.SetName(Form("%s_HF",GetName()));
 560   if(fHasMC && fPairType==kMConly) fArrPairType.Expand(steps);
 561   else fArrPairType.Expand(AliDielectron::kEv1PMRot+1+steps);
 562
 563   Int_t size  = GetNumberOfBins();
 564   AliDebug(10,Form("Creating a histo array with size %d \n",size));
 565
 566   Int_t sizeAdd  = 1;
 567
 568   // fill object array with the array of bin cells
 569   TObjArray *histArr = new TObjArray(0);
 570   if(!histArr) return;
 571   histArr->SetOwner(kTRUE);
 572   histArr->Expand(size);
 573
 574   //  printf("fRefObj %p \n",fRefObj);
 575   // array of histograms to each bin cell
 576   for(Int_t ihist=0; ihist<size; ihist++) {
 577     histArr->AddAt(fRefObj.Clone(""), ihist);
 578     //histArr->AddAt(fRefObj.Clone(Form("h%04d",ihist)), ihist);
 579   }
 580
 581   // loop over all cut variables and do the naming according to its bin cell
 582   Int_t nvars = fAxes.GetEntriesFast();
 583   for(Int_t ivar=0; ivar<nvars; ivar++) {
 584
 585     // get bin limits
 586     TVectorD *bins = static_cast<TVectorD*>(fAxes.At(ivar));
 587     Int_t nbins    = bins->GetNrows()-1;
 588
 589
 590     // loop over all bin cells an set unique titles
 591     for(Int_t ihist=0; ihist<size; ihist++) {
 592
 593       // get the lower limit for current ivar bin
 594       Int_t ibin   = (ihist/sizeAdd)%nbins;
 595       Double_t lowEdge = (*bins)[ibin];
 596       Double_t upEdge  = (*bins)[ibin+1];
 597       switch(fBinType[ivar]) {
 598       case kStdBin:     upEdge=(*bins)[ibin+1];     break;
 599       case kBinToMax:   upEdge=(*bins)[nbins];      break;
 600       case kBinFromMin: lowEdge=(*bins)[0];         break;
 601       case kSymBin:     upEdge=(*bins)[nbins-ibin];
 602         if(ibin>=((Double_t)(nbins+1))/2) upEdge=(*bins)[nbins]; // to avoid low>up
 603         break;
 604       }
 605
 606       TObjArray *tmp= (TObjArray*) histArr->At(ihist);
 607       TString title = tmp->GetName();
 608       if(!ivar)             title ="";
 609       if( ivar)             title+=":";
 610       if(fVarCutType->TestBitNumber(ivar)) title+="Leg";
 611       title+=AliDielectronVarManager::GetValueName(fVarCuts[ivar]);
 612       title+=Form("#%.2f#%.2f",lowEdge,upEdge);
 613       tmp->SetName(title.Data());
 614       AliDebug(10,title.Data());
 615     } // end: array of bin cell
 616     sizeAdd*=nbins;
 617   } //end: cut loop
 618
 619   // copy array to the selected pair types/ MC sources
 620
 621   // pair types
 622   Int_t istep=0;
 623   if(fPairType != kMConly) {
 624     for(istep=0; istep<AliDielectron::kEv1PMRot+1; istep++) {
 625
 626       if(IsPairTypeSelected(istep)) {
 627         // add a deep copy of the array
 628         fArrPairType[istep]=(TObjArray*)histArr->Clone(AliDielectron::PairClassName(istep));
 629         ((TObjArray*)fArrPairType[istep])->SetOwner();
 630       }
 631       else {
 632         fArrPairType[istep]=new TObjArray(0);
 633         ((TObjArray*)fArrPairType[istep])->SetOwner();
 634         ((TObjArray*)fArrPairType[istep])->SetName(AliDielectron::PairClassName(istep));
 635       }
 636     } //end: loop over pair types
 637   }
 638
 639   // mc sources
 640   if(fHasMC) {
 641     for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
 642       TString title = Form("(Signal: %s)",fSignalsMC->At(i)->GetTitle());
 643       fArrPairType[istep+i]=(TObjArray*)histArr->Clone(title.Data());
 644       if(fStepGenerated)  {
 645         title+=" MC truth";
 646         fArrPairType[istep+i+fSignalsMC->GetEntries()]=(TObjArray*)histArr->Clone(title.Data());
 647       }
 648     } // end: loop over sources
 649   }
 650
 651
 652   // clean up
 653   if(histArr) {
 654     delete histArr;
 655     histArr=0;
 656   }
 657
 658 }
 659
 660 //______________________________________________
 661 Int_t AliDielectronHF::GetNumberOfBins() const
 662 {
 663   //
 664   // return the number of bins this histogram grid has
 665   //
 666   Int_t size=1;
 667   for (Int_t i=0; i<fAxes.GetEntriesFast(); ++i)
 668     size*=((static_cast<TVectorD*>(fAxes.At(i)))->GetNrows()-1);
 669   return size;
 670 }
 671
 672 //______________________________________________
 673 Bool_t AliDielectronHF::IsPairTypeSelected(Int_t itype)
 674 {
 675   //
 676   // check whether a pair type was selected
 677   // TODO: cross check or replace by mixinghandlers processsing
 678
 679   Bool_t selected = kFALSE;
 680
 681   // fill all
 682   if(fPairType==kAll) return kTRUE;
 683
 684   switch(itype) {
 685   case AliDielectron::kEv1PP:
 686   case AliDielectron::kEv1MM:
 687     if(fPairType==kSeAll || fPairType==kSeMeAll || fPairType==kSeReAll )   selected = kTRUE;
 688     break;
 689   case AliDielectron::kEv1PM:
 690     if(fPairType!=kMeOnlyOS)  selected = kTRUE;
 691     break;
 692   case AliDielectron::kEv1PEv2P:
 693   case AliDielectron::kEv1MEv2M:
 694     if(fPairType==kMeAll || fPairType==kSeMeAll)   selected = kTRUE;
 695     break;
 696   case AliDielectron::kEv1PEv2M:
 697     if(fPairType==kMeAll || fPairType==kSeMeAll) selected = kTRUE;
 698     break;
 699   case AliDielectron::kEv1MEv2P:
 700     if(fPairType==kMeAll || fPairType==kSeMeAll || fPairType==kMeOnlyOS || fPairType==kSeMeOnlyOS)   selected = kTRUE;
 701     break;
 702   case AliDielectron::kEv2PP:
 703   case AliDielectron::kEv2PM:
 704   case AliDielectron::kEv2MM:
 705     selected = kFALSE;
 706     break;
 707   case AliDielectron::kEv1PMRot:
 708     if(fPairType==kSeReAll || fPairType==kSeReOnlyOS)   selected = kTRUE;
 709     break;
 710   }
 711
 712   return selected;
 713
 714 }
 715
 716