[u/mrichter/AliRoot.git] / PWGDQ / dielectron / AliDielectronHF.cxx

/*************************************************************************
* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
*                                                                        *
* Author: The ALICE Off-line Project.                                    *
* Contributors are mentioned in the code where appropriate.              *
*                                                                        *
* Permission to use, copy, modify and distribute this software and its   *
* documentation strictly for non-commercial purposes is hereby granted   *
* without fee, provided that the above copyright notice appears in all   *
* copies and that both the copyright notice and this permission notice   *
* appear in the supporting documentation. The authors make no claims     *
* about the suitability of this software for any purpose. It is          *
* provided "as is" without express or implied warranty.                  *
**************************************************************************/

///////////////////////////////////////////////////////////////////////////
//                Dielectron HF                                  //
//                                                                       //
//                                                                       //
/*
Detailed description


*/
//                                                                       //
///////////////////////////////////////////////////////////////////////////

#include <TVectorD.h>
#include <TH1.h>
#include <TAxis.h>
#include <AliVParticle.h>

#include <AliLog.h>

#include "AliDielectron.h"
#include "AliDielectronHelper.h"
#include "AliDielectronMC.h"
#include "AliDielectronPair.h"
#include "AliDielectronSignalMC.h"

#include "AliDielectronHistos.h"
#include "AliDielectronHF.h"

ClassImp(AliDielectronHF)

AliDielectronHF::AliDielectronHF() :
  TNamed(),
  fArrPairType(0x0),
  fPairType(kOSonly),
  fSignalsMC(0x0),
  fAxes(kMaxCuts),
  fHasMC(kFALSE),
  fStepGenerated(kFALSE),
  fRefObj(0x0)
{
  //
  // Default Constructor
  //
  for (Int_t i=0; i<kMaxCuts; ++i){
    fVarCuts[i]=0;
    fVarCutType[i]=0;
    fBinType[i]=kStdBin;
  }
  fAxes.SetOwner(kTRUE);
}

//______________________________________________
AliDielectronHF::AliDielectronHF(const char* name, const char* title) :
  TNamed(name, title),
  fArrPairType(0x0),
  fPairType(kOSonly),
  fSignalsMC(0x0),
  fAxes(kMaxCuts),
  fHasMC(kFALSE),
  fStepGenerated(kFALSE),
  fRefObj(0x0)
{
  //
  // Named Constructor
  //
  for (Int_t i=0; i<kMaxCuts; ++i){
    fVarCuts[i]=0;
    fVarCutType[i]=0;
    fBinType[i]=kStdBin;
  }
  fAxes.SetOwner(kTRUE);
}

//______________________________________________
AliDielectronHF::~AliDielectronHF()
{
  //
  // Default Destructor
  //
  fAxes.Delete();
}

//________________________________________________________________
void AliDielectronHF::SetRefHist(TH1 *obj, UInt_t vars[4])
{
  //
  // store reference object and its varaibles
  //

  //  UInt_t val[2]={AliDielectronVarManager::kM,AliDielectronVarManager::kPt};
  AliDielectronHistos::StoreVariables(obj,vars);
  AliDielectronHistos::AdaptNameTitle(obj,"Pair");
  obj->SetName("");
  fRefObj     = obj;

}

//________________________________________________________________
void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
                                             Int_t nbins, Double_t min, Double_t max, Bool_t log, Bool_t leg, EBinType btype)
{
  //
  // Add a variable to the mixing handler
  //

  // limit number of variables to kMaxCuts
  if (fAxes.GetEntriesFast()>=kMaxCuts) return;
  
  TVectorD *binLimits=0x0;
  if (!log) binLimits=AliDielectronHelper::MakeLinBinning(nbins,min,max);
  else binLimits=AliDielectronHelper::MakeLogBinning(nbins,min,max);
  if (!binLimits) return;

  Int_t size=fAxes.GetEntriesFast();
  fVarCuts[size]=(UShort_t)type;
  fVarCutType[size]=leg;
  fAxes.Add(binLimits->Clone());
  fBinType[size]=btype;
}

//________________________________________________________________
void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
                                             const char* binLimitStr, Bool_t leg, EBinType btype)
{
  //
  // Add a variable to the mixing handler with arbitrary binning
  //

  // limit number of variables to kMaxCuts
  if (fAxes.GetEntriesFast()>=kMaxCuts) return;
  
  TVectorD *binLimits=AliDielectronHelper::MakeArbitraryBinning(binLimitStr);
  if (!binLimits) return;
  
  Int_t size=fAxes.GetEntriesFast();
  fVarCuts[size]=(UShort_t)type;
  fVarCutType[size]=leg;
  fAxes.Add(binLimits);
  fBinType[size]=btype;
}

//________________________________________________________________
void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
                                             TVectorD * binLimits, Bool_t leg, EBinType btype)
{
  //
  // Add a variable to the mixing handler with a vector
  // the TVectorD is assumed to be surplus after the creation and will be deleted!!!
  //

  // limit number of variables to kMaxCuts
  if (fAxes.GetEntriesFast()>=kMaxCuts) return;
  
  if (!binLimits) return;
  
  Int_t size=fAxes.GetEntriesFast();
  fVarCuts[size]=(UShort_t)type;
  fVarCutType[size]=leg;
  fAxes.Add(binLimits);
  fBinType[size]=btype;
}

//______________________________________________
void AliDielectronHF::Fill(Int_t label1, Int_t label2, Int_t nSignal) 
{
  //
  // fill the pure MC part of the container starting from a pair of 2 particles (part1 and part2 are legs)
  //
  // fill only if we have asked for these steps
  if(!fStepGenerated) return;

  AliVParticle* part1 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label1);
  AliVParticle* part2 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label2);

  AliDielectronMC* dieMC = AliDielectronMC::Instance();
  
  Int_t mLabel1 = dieMC->GetMothersLabel(label1);    // should work for both ESD and AOD
  Int_t mLabel2 = dieMC->GetMothersLabel(label2);

  // check the same mother option
  AliDielectronSignalMC* sigMC = (AliDielectronSignalMC*)fSignalsMC->At(nSignal);
  if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kSame && mLabel1!=mLabel2) return;
  if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kDifferent && mLabel1==mLabel2) return;
    
  // fill the leg variables
  Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues];
  Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues];
  AliDielectronVarManager::Fill(part1,valuesLeg1);
  AliDielectronVarManager::Fill(part2,valuesLeg2);
    
  // fill the pair and event variables
  Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
  AliDielectronVarManager::Fill(dieMC->GetMCEvent(), valuesPair);
  AliDielectronVarManager::FillVarMCParticle2(part1,part2,valuesPair);

  if(part1->Charge()*part2->Charge()<0) {
    //valuesPair[AliDielectronVarManager::kPairType]=1;
//     for(Int_t i=0; i<fAxes.GetEntriesFast(); i++) {
//       printf("[D]: %s %f %f %f \n",
// 	     AliDielectronVarManager::GetValueName(fVarCuts[i]),
// 	     valuesLeg1[fVarCuts[i]], valuesLeg2[fVarCuts[i]], valuesPair[fVarCuts[i]]); 
//     }
    Fill(nSignal, valuesPair,  valuesLeg1, valuesLeg2);
  }
  // on OS at the moment
  // else if(part1->Charge()>0)
  //   valuesPair[AliDielectronVarManager::kPairType]=0;
  // else
  //   valuesPair[AliDielectronVarManager::kPairType]=2; // if one of the two particles is neutral, the pair will go here

  return;
}
//______________________________________________
void AliDielectronHF::Fill(Int_t pairIndex, const AliDielectronPair *particle)
{
  //
  // fill histograms for event, pair and daughter cuts and pair types
  //
  
  // only OS pairs in case of MC
  if(fHasMC && pairIndex!=AliDielectron::kEv1PM) return;

  // only selected pair types in case of data
  if(!IsPairTypeSelected(pairIndex)) return;

  // get event and pair variables
  Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
  AliDielectronVarManager::Fill(particle,valuesPair);

  // get leg variables
  Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues]={0};
  AliDielectronVarManager::Fill(particle->GetFirstDaughter(),valuesLeg1);
  Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues]={0};
  AliDielectronVarManager::Fill(particle->GetSecondDaughter(),valuesLeg2);

  // fill
  if(!fHasMC) { Fill(pairIndex, valuesPair,  valuesLeg1, valuesLeg2); }
  if(fHasMC && fSignalsMC) {
    for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
      if(AliDielectronMC::Instance()->IsMCTruth(particle, (AliDielectronSignalMC*)fSignalsMC->At(i))) 
	Fill(i, valuesPair,  valuesLeg1, valuesLeg2);
    }
  }
  
  return;
  
}

//______________________________________________
void AliDielectronHF::Fill(Int_t Index, Double_t * const valuesPair, Double_t * const valuesLeg1, Double_t * const valuesLeg2)
{
  //
  // main fill function using index and values as input
  //

  TObjArray *histArr = static_cast<TObjArray*>(fArrPairType.At(Index));
  if(!histArr) return;

  Int_t size  = GetNumberOfBins();
  // loop over all histograms
  for(Int_t ihist=0; ihist<size; ihist++) {

    Int_t sizeAdd   = 1;
    Bool_t selected = kTRUE;
    
    // loop over all cut variables
    Int_t nvars = fAxes.GetEntriesFast();
    for(Int_t ivar=0; ivar<nvars; ivar++) {
      
      // get bin limits
      TVectorD *bins = static_cast<TVectorD*>(fAxes.At(ivar));
      Int_t nbins    = bins->GetNrows()-1;

      // bin limits for current ivar bin
      Int_t ibin   = (ihist/sizeAdd)%nbins;
      Double_t lowEdge = (*bins)[ibin];
      Double_t upEdge  = (*bins)[ibin+1];
      switch(fBinType[ivar]) {
      case kStdBin:     upEdge=(*bins)[ibin+1];     break;
      case kBinToMax:   upEdge=(*bins)[nbins];      break;
      case kBinFromMin: lowEdge=(*bins)[0];         break;
      case kSymBin:     upEdge=(*bins)[nbins-ibin];
	if(ibin>=((Double_t)(nbins+1))/2) upEdge=(*bins)[nbins]; // to avoid low>up
	break;
      }

      // leg variable
      if(fVarCutType[ivar]) {
	if( (valuesLeg1[fVarCuts[ivar]] < lowEdge || valuesLeg1[fVarCuts[ivar]] >= upEdge) ||
	    (valuesLeg2[fVarCuts[ivar]] < lowEdge || valuesLeg2[fVarCuts[ivar]] >= upEdge) ) {
	  selected=kFALSE;
	  break;
	}
      }
      else { // pair and event variables
	if( (valuesPair[fVarCuts[ivar]] < lowEdge || valuesPair[fVarCuts[ivar]] >= upEdge) ) {
	  selected=kFALSE;
	  break;
	}
      }

      sizeAdd*=nbins;
    } //end of var cut loop

    // do not fill the histogram
    if(!selected) continue;

    // fill the object with Pair and event values (TODO: this needs to be changed)
    TH1 *tmp=static_cast<TH1*>(histArr->At(ihist));
    TString title = tmp->GetName();
    AliDebug(10,title.Data());

    AliDielectronHistos::FillValues(tmp, valuesPair);
    //    AliDebug(10,Form("Fill var %d %s value %f in %s \n",fVar,AliDielectronVarManager::GetValueName(fVar),valuesPair[fVar],tmp->GetName()));
  } //end of hist loop

}

//______________________________________________
void AliDielectronHF::Init()
{
  //
  // initialise event buffers
  //

  // has MC signals
  fHasMC=AliDielectronMC::Instance()->HasMC();
  Int_t steps = 0;
  if(fHasMC) steps=fSignalsMC->GetEntries();
  if(fStepGenerated) steps*=2; 

  // init pair type array
  fArrPairType.SetName(Form("%s_HF",GetName()));
  if(fHasMC) fArrPairType.Expand(steps);
  else fArrPairType.Expand(AliDielectron::kEv1PMRot+1);

  Int_t size  = GetNumberOfBins();
  AliDebug(10,Form("Creating a histo array with size %d \n",size));

  Int_t sizeAdd  = 1; 

  // fill object array with the histograms
  TObjArray *histArr = new TObjArray();
  histArr->Expand(size);

  //  printf("fRefObj %p \n",fRefObj);
  for(Int_t ihist=0; ihist<size; ihist++) {
    histArr->AddAt(fRefObj->Clone(""), ihist);
    //histArr->AddAt(fRefObj->Clone(Form("h%04d",ihist)), ihist);
  }

  // loop over all cut variables
  Int_t nvars = fAxes.GetEntriesFast();
  for(Int_t ivar=0; ivar<nvars; ivar++) {
    
    // get bin limits
    TVectorD *bins = static_cast<TVectorD*>(fAxes.At(ivar));
    Int_t nbins    = bins->GetNrows()-1;


    // loop over all histograms an set unique titles
    for(Int_t ihist=0; ihist<size; ihist++) {

      // get the lower limit for current ivar bin
      Int_t ibin   = (ihist/sizeAdd)%nbins; 
      Double_t lowEdge = (*bins)[ibin];
      Double_t upEdge  = (*bins)[ibin+1];
      switch(fBinType[ivar]) {
      case kStdBin:     upEdge=(*bins)[ibin+1];     break;
      case kBinToMax:   upEdge=(*bins)[nbins];      break;
      case kBinFromMin: lowEdge=(*bins)[0];         break;
      case kSymBin:     upEdge=(*bins)[nbins-ibin];
	if(ibin>=((Double_t)(nbins+1))/2) upEdge=(*bins)[nbins]; // to avoid low>up
	break;
      }

      TH1 *tmp=static_cast<TH1*>(histArr->At(ihist));
      TString title = tmp->GetName();
      if(ivar!=0)           title+=":";
      if(fVarCutType[ivar]) title+="Leg";
      title+=AliDielectronVarManager::GetValueName(fVarCuts[ivar]);
      title+=Form("#%.2f#%.2f",lowEdge,upEdge);
      tmp->SetName(title.Data());
      AliDebug(10,title.Data());
    }
    sizeAdd*=nbins;
  }

  // copy array to the selected pair types/ MC sources
  if(fHasMC) {
    for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
	TString title = Form("(Signal: %s)",fSignalsMC->At(i)->GetTitle());
	fArrPairType[i]=(TObjArray*)histArr->Clone(title.Data());
	if(fStepGenerated)  {
	  title+=" MC truth";
	  fArrPairType[i+fSignalsMC->GetEntries()]=(TObjArray*)histArr->Clone(title.Data());
	}
      }
   }
   else {
    for(Int_t i=0; i<AliDielectron::kEv1PMRot+1; i++) {
      if(IsPairTypeSelected(i)) fArrPairType[i]=(TObjArray*)histArr->Clone(Form("%s",AliDielectron::PairClassName(i)));
      else fArrPairType[i]=0x0;
    }
  }
  
  // clean up
  if(histArr) {
    delete histArr;
    histArr=0;
  }

}

//______________________________________________
Int_t AliDielectronHF::GetNumberOfBins() const
{
  //
  // return the number of bins this mixing handler has
  //
  Int_t size=1;
  for (Int_t i=0; i<fAxes.GetEntriesFast(); ++i)
    size*=((static_cast<TVectorD*>(fAxes.At(i)))->GetNrows()-1);
  return size;
}

//______________________________________________
Bool_t AliDielectronHF::IsPairTypeSelected(Int_t itype)
{
  //
  // check whether a pair type was selected
  //
  
  Bool_t selected = kFALSE;

  // fill all
  if(fPairType==kAll) selected = kTRUE;

  switch(itype) {
  case AliDielectron::kEv1PP: 
  case AliDielectron::kEv1MM:
    if(fPairType==kOSandLS)   selected = kTRUE;
    break;
  case AliDielectron::kEv1PM: selected = kTRUE;
    break;
  case AliDielectron::kEv1PEv2P:
  case AliDielectron::kEv1MEv2P:
  case AliDielectron::kEv1PEv2M:
  case AliDielectron::kEv1MEv2M:
    if(fPairType==kOSandMIX)  selected = kTRUE;
    break;
  case AliDielectron::kEv2PP:
  case AliDielectron::kEv2PM: 
  case AliDielectron::kEv2MM:
    selected = kFALSE;
    break;
  case AliDielectron::kEv1PMRot:
    if(fPairType==kOSandROT)  selected = kTRUE;
    break;
  }
  
  return selected;

}
Commit	Line	Data
5e2cf960	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 <TAxis.h>
443a091c	31	#include <AliVParticle.h>
5e2cf960	32
	33	#include <AliLog.h>
	34
	35	#include "AliDielectron.h"
	36	#include "AliDielectronHelper.h"
443a091c	37	#include "AliDielectronMC.h"
	38	#include "AliDielectronPair.h"
	39	#include "AliDielectronSignalMC.h"
5e2cf960	40
35363fb9	41	#include "AliDielectronHistos.h"
5e2cf960	42	#include "AliDielectronHF.h"
	43
	44	ClassImp(AliDielectronHF)
	45
	46	AliDielectronHF::AliDielectronHF() :
	47	TNamed(),
	48	fArrPairType(0x0),
	49	fPairType(kOSonly),
443a091c	50	fSignalsMC(0x0),
5e2cf960	51	fAxes(kMaxCuts),
4d7704c5	52	fHasMC(kFALSE),
35363fb9	53	fStepGenerated(kFALSE),
35363fb9	54	fRefObj(0x0)
5e2cf960	55	{
	56	//
	57	// Default Constructor
	58	//
	59	for (Int_t i=0; i<kMaxCuts; ++i){
	60	fVarCuts[i]=0;
	61	fVarCutType[i]=0;
	62	fBinType[i]=kStdBin;
	63	}
	64	fAxes.SetOwner(kTRUE);
	65	}
	66
	67	//______________________________________________
	68	AliDielectronHF::AliDielectronHF(const char* name, const char* title) :
	69	TNamed(name, title),
	70	fArrPairType(0x0),
	71	fPairType(kOSonly),
443a091c	72	fSignalsMC(0x0),
5e2cf960	73	fAxes(kMaxCuts),
4d7704c5	74	fHasMC(kFALSE),
35363fb9	75	fStepGenerated(kFALSE),
35363fb9	76	fRefObj(0x0)
5e2cf960	77	{
	78	//
	79	// Named Constructor
	80	//
	81	for (Int_t i=0; i<kMaxCuts; ++i){
	82	fVarCuts[i]=0;
	83	fVarCutType[i]=0;
	84	fBinType[i]=kStdBin;
	85	}
	86	fAxes.SetOwner(kTRUE);
	87	}
	88
	89	//______________________________________________
	90	AliDielectronHF::~AliDielectronHF()
	91	{
	92	//
	93	// Default Destructor
	94	//
	95	fAxes.Delete();
	96	}
	97
	98	//________________________________________________________________
35363fb9	99	void AliDielectronHF::SetRefHist(TH1 *obj, UInt_t vars[4])
5e2cf960	100	{
5e2cf960	101	//
35363fb9	102	// store reference object and its varaibles
5e2cf960	103	//
5e2cf960	104
35363fb9	105	// UInt_t val[2]={AliDielectronVarManager::kM,AliDielectronVarManager::kPt};
	106	AliDielectronHistos::StoreVariables(obj,vars);
	107	AliDielectronHistos::AdaptNameTitle(obj,"Pair");
	108	obj->SetName("");
	109	fRefObj = obj;
5e2cf960	110
	111	}
	112
	113	//________________________________________________________________
	114	void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
	115	Int_t nbins, Double_t min, Double_t max, Bool_t log, Bool_t leg, EBinType btype)
	116	{
	117	//
	118	// Add a variable to the mixing handler
	119	//
	120
	121	// limit number of variables to kMaxCuts
	122	if (fAxes.GetEntriesFast()>=kMaxCuts) return;
	123
	124	TVectorD *binLimits=0x0;
	125	if (!log) binLimits=AliDielectronHelper::MakeLinBinning(nbins,min,max);
	126	else binLimits=AliDielectronHelper::MakeLogBinning(nbins,min,max);
	127	if (!binLimits) return;
	128
	129	Int_t size=fAxes.GetEntriesFast();
	130	fVarCuts[size]=(UShort_t)type;
	131	fVarCutType[size]=leg;
	132	fAxes.Add(binLimits->Clone());
	133	fBinType[size]=btype;
	134	}
	135
	136	//________________________________________________________________
	137	void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
	138	const char* binLimitStr, Bool_t leg, EBinType btype)
	139	{
	140	//
	141	// Add a variable to the mixing handler with arbitrary binning
	142	//
	143
	144	// limit number of variables to kMaxCuts
	145	if (fAxes.GetEntriesFast()>=kMaxCuts) return;
	146
	147	TVectorD *binLimits=AliDielectronHelper::MakeArbitraryBinning(binLimitStr);
	148	if (!binLimits) return;
	149
	150	Int_t size=fAxes.GetEntriesFast();
	151	fVarCuts[size]=(UShort_t)type;
	152	fVarCutType[size]=leg;
	153	fAxes.Add(binLimits);
	154	fBinType[size]=btype;
	155	}
	156
	157	//________________________________________________________________
	158	void AliDielectronHF::AddCutVariable(AliDielectronVarManager::ValueTypes type,
	159	TVectorD * binLimits, Bool_t leg, EBinType btype)
	160	{
	161	//
	162	// Add a variable to the mixing handler with a vector
	163	// the TVectorD is assumed to be surplus after the creation and will be deleted!!!
	164	//
	165
	166	// limit number of variables to kMaxCuts
	167	if (fAxes.GetEntriesFast()>=kMaxCuts) return;
	168
	169	if (!binLimits) return;
	170
	171	Int_t size=fAxes.GetEntriesFast();
	172	fVarCuts[size]=(UShort_t)type;
	173	fVarCutType[size]=leg;
174	fAxes.Add(binLimits);
175	fBinType[size]=btype;
176	}
177
178	//______________________________________________
443a091c	179	void AliDielectronHF::Fill(Int_t label1, Int_t label2, Int_t nSignal)
5e2cf960	180	{
5e2cf960	181	//
443a091c	182	// fill the pure MC part of the container starting from a pair of 2 particles (part1 and part2 are legs)
5e2cf960	183	//
4d7704c5	184	// fill only if we have asked for these steps
4d7704c5	185	if(!fStepGenerated) return;
5e2cf960	186
443a091c	187	AliVParticle* part1 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label1);
443a091c	188	AliVParticle* part2 = AliDielectronMC::Instance()->GetMCTrackFromMCEvent(label2);
5e2cf960	189
443a091c	190	AliDielectronMC* dieMC = AliDielectronMC::Instance();
	191
	192	Int_t mLabel1 = dieMC->GetMothersLabel(label1); // should work for both ESD and AOD
	193	Int_t mLabel2 = dieMC->GetMothersLabel(label2);
5e2cf960	194
443a091c	195	// check the same mother option
	196	AliDielectronSignalMC* sigMC = (AliDielectronSignalMC*)fSignalsMC->At(nSignal);
	197	if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kSame && mLabel1!=mLabel2) return;
	198	if(sigMC->GetMothersRelation()==AliDielectronSignalMC::kDifferent && mLabel1==mLabel2) return;
	199
	200	// fill the leg variables
	201	Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues];
	202	Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues];
	203	AliDielectronVarManager::Fill(part1,valuesLeg1);
	204	AliDielectronVarManager::Fill(part2,valuesLeg2);
	205
	206	// fill the pair and event variables
	207	Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
	208	AliDielectronVarManager::Fill(dieMC->GetMCEvent(), valuesPair);
	209	AliDielectronVarManager::FillVarMCParticle2(part1,part2,valuesPair);
	210
	211	if(part1->Charge()*part2->Charge()<0) {
88204efa	212	//valuesPair[AliDielectronVarManager::kPairType]=1;
	213	// for(Int_t i=0; i<fAxes.GetEntriesFast(); i++) {
	214	// printf("[D]: %s %f %f %f \n",
	215	// AliDielectronVarManager::GetValueName(fVarCuts[i]),
	216	// valuesLeg1[fVarCuts[i]], valuesLeg2[fVarCuts[i]], valuesPair[fVarCuts[i]]);
	217	// }
443a091c	218	Fill(nSignal, valuesPair, valuesLeg1, valuesLeg2);
	219	}
	220	// on OS at the moment
	221	// else if(part1->Charge()>0)
	222	// valuesPair[AliDielectronVarManager::kPairType]=0;
	223	// else
	224	// valuesPair[AliDielectronVarManager::kPairType]=2; // if one of the two particles is neutral, the pair will go here
	225
	226	return;
	227	}
	228	//______________________________________________
	229	void AliDielectronHF::Fill(Int_t pairIndex, const AliDielectronPair *particle)
	230	{
	231	//
	232	// fill histograms for event, pair and daughter cuts and pair types
	233	//
5e2cf960	234
443a091c	235	// only OS pairs in case of MC
	236	if(fHasMC && pairIndex!=AliDielectron::kEv1PM) return;
	237
	238	// only selected pair types in case of data
	239	if(!IsPairTypeSelected(pairIndex)) return;
5e2cf960	240
	241	// get event and pair variables
	242	Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
	243	AliDielectronVarManager::Fill(particle,valuesPair);
	244
	245	// get leg variables
	246	Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues]={0};
	247	AliDielectronVarManager::Fill(particle->GetFirstDaughter(),valuesLeg1);
	248	Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues]={0};
	249	AliDielectronVarManager::Fill(particle->GetSecondDaughter(),valuesLeg2);
	250
443a091c	251	// fill
	252	if(!fHasMC) { Fill(pairIndex, valuesPair, valuesLeg1, valuesLeg2); }
	253	if(fHasMC && fSignalsMC) {
	254	for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
	255	if(AliDielectronMC::Instance()->IsMCTruth(particle, (AliDielectronSignalMC*)fSignalsMC->At(i)))
	256	Fill(i, valuesPair, valuesLeg1, valuesLeg2);
	257	}
	258	}
	259
	260	return;
	261
	262	}
	263
	264	//______________________________________________
	265	void AliDielectronHF::Fill(Int_t Index, Double_t * const valuesPair, Double_t * const valuesLeg1, Double_t * const valuesLeg2)
	266	{
	267	//
	268	// main fill function using index and values as input
	269	//
	270
	271	TObjArray histArr = static_cast<TObjArray>(fArrPairType.At(Index));
	272	if(!histArr) return;
	273
	274	Int_t size = GetNumberOfBins();
5e2cf960	275	// loop over all histograms
	276	for(Int_t ihist=0; ihist<size; ihist++) {
	277
35363fb9	278	Int_t sizeAdd = 1;
5e2cf960	279	Bool_t selected = kTRUE;
	280
	281	// loop over all cut variables
	282	Int_t nvars = fAxes.GetEntriesFast();
	283	for(Int_t ivar=0; ivar<nvars; ivar++) {
	284
	285	// get bin limits
	286	TVectorD bins = static_cast<TVectorD>(fAxes.At(ivar));
	287	Int_t nbins = bins->GetNrows()-1;
	288
	289	// bin limits for current ivar bin
35363fb9	290	Int_t ibin = (ihist/sizeAdd)%nbins;
5e2cf960	291	Double_t lowEdge = (*bins)[ibin];
	292	Double_t upEdge = (*bins)[ibin+1];
	293	switch(fBinType[ivar]) {
	294	case kStdBin: upEdge=(*bins)[ibin+1]; break;
	295	case kBinToMax: upEdge=(*bins)[nbins]; break;
	296	case kBinFromMin: lowEdge=(*bins)[0]; break;
35363fb9	297	case kSymBin: upEdge=(*bins)[nbins-ibin];
5e2cf960	298	if(ibin>=((Double_t)(nbins+1))/2) upEdge=(*bins)[nbins]; // to avoid low>up
	299	break;
	300	}
443a091c	301
5e2cf960	302	// leg variable
443a091c	303	if(fVarCutType[ivar]) {
5e2cf960	304	if( (valuesLeg1[fVarCuts[ivar]] < lowEdge \|\| valuesLeg1[fVarCuts[ivar]] >= upEdge) \|\|
	305	(valuesLeg2[fVarCuts[ivar]] < lowEdge \|\| valuesLeg2[fVarCuts[ivar]] >= upEdge) ) {
	306	selected=kFALSE;
	307	break;
	308	}
35363fb9	309	}
5e2cf960	310	else { // pair and event variables
	311	if( (valuesPair[fVarCuts[ivar]] < lowEdge \|\| valuesPair[fVarCuts[ivar]] >= upEdge) ) {
	312	selected=kFALSE;
	313	break;
	314	}
	315	}
35363fb9	316
5e2cf960	317	sizeAdd*=nbins;
5e2cf960	318	} //end of var cut loop
35363fb9	319
5e2cf960	320	// do not fill the histogram
5e2cf960	321	if(!selected) continue;
35363fb9	322
	323	// fill the object with Pair and event values (TODO: this needs to be changed)
	324	TH1 tmp=static_cast<TH1>(histArr->At(ihist));
	325	TString title = tmp->GetName();
5e2cf960	326	AliDebug(10,title.Data());
35363fb9	327
	328	AliDielectronHistos::FillValues(tmp, valuesPair);
	329	// AliDebug(10,Form("Fill var %d %s value %f in %s \n",fVar,AliDielectronVarManager::GetValueName(fVar),valuesPair[fVar],tmp->GetName()));
5e2cf960	330	} //end of hist loop
35363fb9	331
5e2cf960	332	}
	333
	334	//______________________________________________
	335	void AliDielectronHF::Init()
	336	{
	337	//
	338	// initialise event buffers
	339	//
	340
443a091c	341	// has MC signals
443a091c	342	fHasMC=AliDielectronMC::Instance()->HasMC();
4d7704c5	343	Int_t steps = 0;
	344	if(fHasMC) steps=fSignalsMC->GetEntries();
	345	if(fStepGenerated) steps*=2;
443a091c	346
5e2cf960	347	// init pair type array
5e2cf960	348	fArrPairType.SetName(Form("%s_HF",GetName()));
4d7704c5	349	if(fHasMC) fArrPairType.Expand(steps);
443a091c	350	else fArrPairType.Expand(AliDielectron::kEv1PMRot+1);
5e2cf960	351
5e2cf960	352	Int_t size = GetNumberOfBins();
	353	AliDebug(10,Form("Creating a histo array with size %d \n",size));
	354
	355	Int_t sizeAdd = 1;
	356
	357	// fill object array with the histograms
	358	TObjArray *histArr = new TObjArray();
	359	histArr->Expand(size);
	360
35363fb9	361	// printf("fRefObj %p \n",fRefObj);
5e2cf960	362	for(Int_t ihist=0; ihist<size; ihist++) {
35363fb9	363	histArr->AddAt(fRefObj->Clone(""), ihist);
35363fb9	364	//histArr->AddAt(fRefObj->Clone(Form("h%04d",ihist)), ihist);
5e2cf960	365	}
	366
	367	// loop over all cut variables
	368	Int_t nvars = fAxes.GetEntriesFast();
	369	for(Int_t ivar=0; ivar<nvars; ivar++) {
	370
	371	// get bin limits
	372	TVectorD bins = static_cast<TVectorD>(fAxes.At(ivar));
	373	Int_t nbins = bins->GetNrows()-1;
	374
	375
	376	// loop over all histograms an set unique titles
	377	for(Int_t ihist=0; ihist<size; ihist++) {
	378
	379	// get the lower limit for current ivar bin
	380	Int_t ibin = (ihist/sizeAdd)%nbins;
	381	Double_t lowEdge = (*bins)[ibin];
	382	Double_t upEdge = (*bins)[ibin+1];
	383	switch(fBinType[ivar]) {
	384	case kStdBin: upEdge=(*bins)[ibin+1]; break;
	385	case kBinToMax: upEdge=(*bins)[nbins]; break;
	386	case kBinFromMin: lowEdge=(*bins)[0]; break;
	387	case kSymBin: upEdge=(*bins)[nbins-ibin];
	388	if(ibin>=((Double_t)(nbins+1))/2) upEdge=(*bins)[nbins]; // to avoid low>up
	389	break;
	390	}
	391
35363fb9	392	TH1 tmp=static_cast<TH1>(histArr->At(ihist));
35363fb9	393	TString title = tmp->GetName();
5e2cf960	394	if(ivar!=0) title+=":";
	395	if(fVarCutType[ivar]) title+="Leg";
	396	title+=AliDielectronVarManager::GetValueName(fVarCuts[ivar]);
	397	title+=Form("#%.2f#%.2f",lowEdge,upEdge);
35363fb9	398	tmp->SetName(title.Data());
5e2cf960	399	AliDebug(10,title.Data());
	400	}
	401	sizeAdd*=nbins;
	402	}
	403
443a091c	404	// copy array to the selected pair types/ MC sources
	405	if(fHasMC) {
	406	for(Int_t i=0; i<fSignalsMC->GetEntries(); i++) {
4d7704c5	407	TString title = Form("(Signal: %s)",fSignalsMC->At(i)->GetTitle());
	408	fArrPairType[i]=(TObjArray*)histArr->Clone(title.Data());
	409	if(fStepGenerated) {
	410	title+=" MC truth";
	411	fArrPairType[i+fSignalsMC->GetEntries()]=(TObjArray*)histArr->Clone(title.Data());
	412	}
	413	}
	414	}
	415	else {
443a091c	416	for(Int_t i=0; i<AliDielectron::kEv1PMRot+1; i++) {
	417	if(IsPairTypeSelected(i)) fArrPairType[i]=(TObjArray*)histArr->Clone(Form("%s",AliDielectron::PairClassName(i)));
	418	else fArrPairType[i]=0x0;
	419	}
5e2cf960	420	}
	421
	422	// clean up
	423	if(histArr) {
	424	delete histArr;
	425	histArr=0;
	426	}
	427
	428	}
	429
	430	//______________________________________________
	431	Int_t AliDielectronHF::GetNumberOfBins() const
	432	{
	433	//
	434	// return the number of bins this mixing handler has
	435	//
	436	Int_t size=1;
	437	for (Int_t i=0; i<fAxes.GetEntriesFast(); ++i)
	438	size=((static_cast<TVectorD>(fAxes.At(i)))->GetNrows()-1);
	439	return size;
	440	}
	441
	442	//______________________________________________
	443	Bool_t AliDielectronHF::IsPairTypeSelected(Int_t itype)
	444	{
	445	//
	446	// check whether a pair type was selected
	447	//
	448
	449	Bool_t selected = kFALSE;
	450
	451	// fill all
	452	if(fPairType==kAll) selected = kTRUE;
	453
	454	switch(itype) {
	455	case AliDielectron::kEv1PP:
	456	case AliDielectron::kEv1MM:
	457	if(fPairType==kOSandLS) selected = kTRUE;
	458	break;
	459	case AliDielectron::kEv1PM: selected = kTRUE;
	460	break;
	461	case AliDielectron::kEv1PEv2P:
	462	case AliDielectron::kEv1MEv2P:
	463	case AliDielectron::kEv1PEv2M:
	464	case AliDielectron::kEv1MEv2M:
	465	if(fPairType==kOSandMIX) selected = kTRUE;
	466	break;
	467	case AliDielectron::kEv2PP:
	468	case AliDielectron::kEv2PM:
	469	case AliDielectron::kEv2MM:
	470	selected = kFALSE;
	471	break;
	472	case AliDielectron::kEv1PMRot:
	473	if(fPairType==kOSandROT) selected = kTRUE;
	474	break;
	475	}
	476
	477	return selected;
	478
	479	}
	480
	481