ROOT logo
/*************************************************************************
* 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 Histogram framework helper                           //
//                                                                       //
/*

A helper class to extract objects(histograms and/or profiles) from
a AliDielctronHF array of objects.


How to use it:

  AliDielectronHFhelper *hf = new AliDielectronHFhelper("path/to/the/output/file.root", "ConfigName");
  // print the structure
  hf->Print();

  //apply some cuts and print them
  hf->SetRangeUser("cut1name",cutmin1,cutmax1);
  hf->SetRangeUser(AliDielectronVarManager::kPt,ptmin,ptmax);
  hf->PrintCuts();

  // collect 1-,2- or 3-dim histograms or profiles with error option (default:"")
  TObjArray *arrHists = hf->CollectHistos(AliDielectronVarManager::kM);
  TObjArray *arrProfs = hf->CollectProfiles("",AliDielectronVarManager::kM,AliDielectronVarManager::kPt);

  // then you are left with an array of histograms for all pair types or MC signals

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

#include <TObjArray.h>
#include <TKey.h>
#include <TList.h>
#include <TClass.h>
#include <TObject.h>
#include <TFile.h>
#include <TString.h>
#include <TObjString.h>
#include <TVectorD.h>
#include <TMath.h>
#include <TH1.h>

#include <AliLog.h>

#include "AliDielectron.h"
#include "AliDielectronHFhelper.h"
#include "AliDielectronHF.h"

ClassImp(AliDielectronHFhelper)

//const char* AliDielectronHFhelper::fCutVars[AliDielectronHFhelper::kMaxCuts] = {"","","","","","","","","",""};

//________________________________________________________________
AliDielectronHFhelper::AliDielectronHFhelper(const char* filename, const char* container) :
  TNamed(),
  fMainArr(0x0),
  fCutVars(0x0),
  fCutLowLimits(0),
  fCutUpLimits(0)
{
  //
  // get HF container(s) from file 'filename'
  //
  SetHFArray(filename, container);

}

//________________________________________________________________
AliDielectronHFhelper::~AliDielectronHFhelper()
{
  //
  // dtor
  //
  if(fMainArr) delete fMainArr;
  if(fCutVars)     delete fCutVars;

}

//________________________________________________________________
void AliDielectronHFhelper::SetHFArray(const char* filename, const char* container)
{
  //
  // get HF container from file
  //

  TFile *f = TFile::Open(filename);

  TList *l=f->GetListOfKeys();
  TIter nextKey(l);
  TKey *k=0x0;
  while ( (k=static_cast<TKey*>(nextKey())) ){

    TObject *o=k->ReadObj();
    if (o->IsA()==TList::Class()){

      TList *tlist=(TList*)o;

      TIter next(tlist);
      TObject *obj=0x0;
      while ((obj = next())) {
	TString objname(obj->GetName());

	if( objname.Contains(Form("%s_HF",container)) && obj->IsA()==TObjArray::Class()) {
	  fMainArr = new TObjArray( *(dynamic_cast<TObjArray*>(obj)) );
	  fMainArr->SetOwner();
	  //	  fMainArr->Print();
	  return;
	}
      }
    }
  }

}
//________________________________________________________________
void AliDielectronHFhelper::SetRangeUser(const char *varname, Double_t min, Double_t max, Bool_t leg)
{
  //
  // Set range from variable name
  //

  Int_t size=fCutLowLimits.GetNrows();

  // check if cut is already set
  for(Int_t icut=0; icut<size; icut++) {
    TString cutName = fCutVars->At(icut)->GetName();
    if(!cutName.CompareTo(Form("%s%s",(leg?"Leg":""),varname))) {
      UnsetRangeUser(varname,leg);
      SetRangeUser(varname, min, max, leg);
      return;
    }
  }

  if(size>=kMaxCuts) return;

  // arrays
  if(!fCutVars) {
    fCutVars = new TObjArray();
    fCutVars->SetOwner();
  }
  fCutLowLimits.ResizeTo(size+1);
  fCutUpLimits.ResizeTo(size+1);

  // fill
  TObjString *str = new TObjString(Form("%s%s",(leg?"Leg":""),varname));
  fCutVars->Add(str);
  fCutLowLimits(size) = min;
  fCutUpLimits(size)  = max;
  AliWarning(Form(" %s [%.2f,%.2f]",fCutVars->At(size)->GetName(),fCutLowLimits(size),fCutUpLimits(size)));
}

//________________________________________________________________
void AliDielectronHFhelper::SetRangeUser(AliDielectronVarManager::ValueTypes type, Double_t min, Double_t max, Bool_t leg)
{
  //
  // Set range from AliDielectronVarManager
  //
  SetRangeUser(AliDielectronVarManager::GetValueName(type), min, max, leg);
}

//________________________________________________________________
void AliDielectronHFhelper::UnsetRangeUser(const char *varname, Bool_t leg)
{
  //
  // unset range from variable name
  //

  Int_t size=fCutLowLimits.GetNrows();
  TVectorD newlow;
  TVectorD newup;

  // find cut and build new vectors w/o it
  Int_t ientries = 0;
  for(Int_t icut=0; icut<size; icut++) {

    TString cutName = fCutVars->At(icut)->GetName();
    if(!cutName.CompareTo(Form("%s%s",(leg?"Leg":""),varname))) {
      fCutVars->AddAt(0x0,icut);
      continue;
    }

    // fill new vectors
    newlow.ResizeTo(ientries+1);
    newup.ResizeTo(ientries+1);
    newlow(ientries) = fCutLowLimits(icut);
    newup(ientries)  = fCutUpLimits(icut);

    ientries++;
  }

  // adapt new arrays/vectors
  fCutVars->Compress();

  fCutLowLimits.ResizeTo(ientries);
  fCutUpLimits.ResizeTo(ientries);
  for(Int_t icut=0; icut<ientries; icut++) {
    fCutLowLimits(icut) = newlow(icut);
    fCutUpLimits(icut)  = newup(icut);
  }
  // PrintCuts();

}

//________________________________________________________________
void AliDielectronHFhelper::UnsetRangeUser(AliDielectronVarManager::ValueTypes type, Bool_t leg)
{
  //
  // Unset range from AliDielectronVarManager
  //
  UnsetRangeUser(AliDielectronVarManager::GetValueName(type), leg);
}

//________________________________________________________________
TObjArray* AliDielectronHFhelper::CollectProfiles(TString option,
						  AliDielectronVarManager::ValueTypes varx,
						  AliDielectronVarManager::ValueTypes vary,
						  AliDielectronVarManager::ValueTypes varz,
						  AliDielectronVarManager::ValueTypes vart)
{
  //
  // collect 1-3 dimensional TProfiles for all kind of pair types or sources
  //

  // reconstruct the histogram/profile name resp. key
  Int_t dim    = 0;
  if(varx < AliDielectronVarManager::kNMaxValues) dim++;
  if(vary < AliDielectronVarManager::kNMaxValues) dim++;
  if(varz < AliDielectronVarManager::kNMaxValues) dim++;
  if(vart < AliDielectronVarManager::kNMaxValues) dim++;
  Bool_t bPairClass=0;
  if( varx < AliDielectronVarManager::kPairMax ||
      vary < AliDielectronVarManager::kPairMax ||
      varz < AliDielectronVarManager::kPairMax ||
      vart < AliDielectronVarManager::kPairMax  ) bPairClass=kTRUE;
  Bool_t bwght= (option.Contains("weight",TString::kIgnoreCase));   // weighted histogram
  Bool_t bprf = !(option.Contains("hist",TString::kIgnoreCase));     // tprofile
  Bool_t bStdOpt=kTRUE;
  if(option.Contains("S",TString::kIgnoreCase) || option.Contains("rms",TString::kIgnoreCase)) bStdOpt=kFALSE;

  TString key = "";
  if(bwght) dim--;
  if(bprf)  dim--;
  switch(dim) {
  case 3:
    key+=Form("%s_",AliDielectronVarManager::GetValueName(varx));
    key+=Form("%s_",AliDielectronVarManager::GetValueName(vary));
    key+=Form("%s",AliDielectronVarManager::GetValueName(varz));
    if(bprf) key+=Form("-%s%s",AliDielectronVarManager::GetValueName(vart),(bStdOpt ? "avg" : "rms"));
    break;
  case 2:
    key+=Form("%s_",AliDielectronVarManager::GetValueName(varx));
    key+=Form("%s",AliDielectronVarManager::GetValueName(vary));
    if(bprf) key+=Form("-%s%s",AliDielectronVarManager::GetValueName(varz),(bStdOpt ? "avg" : "rms"));
    break;
  case 1:
    key+=Form("%s",AliDielectronVarManager::GetValueName(varx));
    if(bprf) key+=Form("-%s%s",AliDielectronVarManager::GetValueName(vary),(bStdOpt ? "avg" : "rms"));
    break;
  }
  // to differentiate btw. leg and pair histos
  if(bPairClass) key.Prepend("p");
  // prepend HF
  key.Prepend("HF_");

  // add the weighted part to the key
  if(bwght) {
    if(bprf) dim++;
    switch(dim) {
    case 3:   key+=Form("-wght%s",AliDielectronVarManager::GetValueName(vart)); break;
    case 2:   key+=Form("-wght%s",AliDielectronVarManager::GetValueName(varz)); break;
    case 1:   key+=Form("-wght%s",AliDielectronVarManager::GetValueName(vary)); break;
    case 4:   AliError(Form(" NO weighted 3D profiles supported by the framework"));
    }
  }

  //printf("--------> KEY: %s \n",key.Data());
  // init array of histograms of interest
  TObjArray *collection = new TObjArray(fMainArr->GetEntriesFast());
  collection->SetOwner(kTRUE);

  TObjArray *cloneArr = new TObjArray( *(dynamic_cast<TObjArray*>(fMainArr)) );
  if(!cloneArr) return 0x0;
  cloneArr->SetOwner(kTRUE);

  // loop over all pair types or sources
  for(Int_t i=0; i<cloneArr->GetEntriesFast(); i++) {
    if(!cloneArr->At(i))                             continue;
    if(!((TObjArray*)cloneArr->At(i))->GetEntries()) continue;

    // Get signal/source array with its histograms of interest
    AliDebug(1,Form(" Looking into step %s selected",cloneArr->At(i)->GetName()));
    TObjArray *arr = FindObjects((TObjArray*)cloneArr->At(i));
    if(arr) {

      // find requested histogram
      collection->AddAt(arr->FindObject(key.Data()), i);
      if(!collection->At(i)) { AliError(Form("Object %s does not exist",key.Data())); return 0x0; }

      // modify the signal/source name if needed (MConly)
      TString stepName = cloneArr->At(i)->GetName();
      stepName.ReplaceAll("(","");
      stepName.ReplaceAll(")","");
      stepName.ReplaceAll(": ","");
      stepName.ReplaceAll(" ","_");
      stepName.ReplaceAll("Signal","");
      ((TH1*)collection->At(i))->SetName(Form("%s_%s",key.Data(),stepName.Data()));
    }
    else
      AliError(Form("Step %d not found",i));

  }

  //clean up
  //delete cloneArr;
  //cloneArr=0;

  return collection;
}

//________________________________________________________________
TObjArray* AliDielectronHFhelper::FindObjects(TObjArray *stepArr)
{
  // rename DoCuts, return values is a tobjarray
  // apply cuts and exclude objects from the array for merging (CUT selection)
  //

  // debug
  // TString title    = stepArr->At(0)->GetTitle();
  // TObjArray* vars  = title.Tokenize(":");
  // AliDebug(1,Form(" number of cuts/vars: %d/%d",fCutLowLimits.GetNrows(),vars->GetEntriesFast()));
  if(!stepArr) { AliError("step is empty"); return 0x0;}

  // check for missing cuts
  CheckCuts(stepArr);

  // loop over all cuts
  for(Int_t icut=0; icut<fCutLowLimits.GetNrows(); icut++) {

    Bool_t bFndBin      = kFALSE; // bin with exact limits found
    const char *cutvar  = fCutVars->At(icut)->GetName(); // cut variable
    Double_t min        = fCutLowLimits(icut);           // lower limit
    Double_t max        = fCutUpLimits(icut);            // upper limit
    AliDebug(5,Form(" Cut %d: %s [%.2f,%.2f]",icut,cutvar,min,max));

    // loop over the full grid of the signalArray (pair,source)
    for(Int_t i=0; i<stepArr->GetEntriesFast(); i++) {

      // continue if already empty
      if(!stepArr->At(i)) continue;

      // collect bins from the name
      TString title    = stepArr->At(i)->GetName();
      if(title.IsNull()) continue;
      AliDebug(5,Form(" %03d object name: %s",i,title.Data()));

      // loop over all variables
      TObjArray *vars  = title.Tokenize(":");
      for(Int_t ivar=0; ivar<vars->GetEntriesFast(); ivar++) {
	TString binvar = vars->At(ivar)->GetName();
	AliDebug(10,Form(" --> %d check bin var %s",ivar,binvar.Data()));

	// check cuts set by the user, and compare to the bin limits
	if(binvar.Contains(cutvar)) {
	  // bin limits
	  TObjArray *limits = binvar.Tokenize("#");
	  if(!limits) continue;
	  Double_t binmin = atof(limits->At(1)->GetName()); // lower bin limit
	  Double_t binmax = atof(limits->At(2)->GetName()); // upper bin limit
	  AliDebug(10,Form(" bin %s var %s [%.2f,%.2f]",binvar.Data(),limits->At(0)->GetName(),binmin,binmax));
	  delete limits;

	  // cut and remove objects from the array
	  if(binmin < min || binmax < min || binmin > max || binmax > max ) {
	    AliDebug(10,Form(" removed, out of range! lower bin,cut: %.2f,%.2f or upper bin,cut: %.2f>%.2f",binmin,min,binmax,max));
	    stepArr->AddAt(0x0,i);
	  }
	  if(bFndBin && !(binmin == min && binmax == max)) {
	    stepArr->AddAt(0x0,i);
	    AliDebug(10,Form(" removed, within range! lower bin,cut: %.2f,%.2f or upper bin,cut: %.2f,%.2f",binmin,min,binmax,max));
	  }

	  // did we found a bin with exact the cut ranges
	  // this can happen only once per variable
	  if(binmin==min && binmax==max) bFndBin=kTRUE;

	}

      }
      // clean up
      if(vars)   delete vars;

    } //end loop: pair step

  } //end: loop cuts

  // compress the array by removing all empty entries
  stepArr->Compress();
  AliDebug(1,Form(" Compression: %d arrays left",stepArr->GetEntriesFast()));

  // merge left objects
  TObjArray* hist = Merge(stepArr);
  if(hist) AliDebug(1,Form(" final array has %d objects",hist->GetEntries()));
  return hist;
}

//________________________________________________________________
TObjArray* AliDielectronHFhelper::Merge(TObjArray *arr)
{
  //
  // merge left objects into a single one (LAST step)
  //

  if(arr->GetEntriesFast()<1) { AliError(" No more objects left!"); return 0x0; }

  TObjArray *final=(TObjArray*) arr->At(0)->Clone();
  if(!final) return 0x0;
  final->SetOwner();

  TList listH;
  TString listHargs;
  listHargs.Form("((TCollection*)0x%lx)", (ULong_t)&listH);
  Int_t error = 0;

  //  final->Reset("CE");
  //  final->SetTitle(""); //TODO: change in future
  for(Int_t i=1; i<arr->GetEntriesFast(); i++) {
    listH.Add(arr->At(i));
    //   printf("%d: ent %.0f \n",i,((TH1*)((TObjArray*)arr->At(i))->At(0))->GetEntries());
    //    final->Add((TH1*)arr->At(i));
  }
  //  arr->Clear();

  final->Execute("Merge", listHargs.Data(), &error);
  // returns a tobjarray of histograms/profiles
  return final;
}

//________________________________________________________________
void AliDielectronHFhelper::CheckCuts(TObjArray *arr)
{
  //
  // Compare binning and cut variables. Add necessary cuts (full range, no exclusion)
  //

  // build array with bin variable, minimum and maximum bin values
  TString titleFIRST       = arr->First()->GetName();
  TString titleLAST        = arr->Last()->GetName();
  TObjArray* binvarsF  = titleFIRST.Tokenize(":#");
  TObjArray* binvarsL  = titleLAST.Tokenize(":#");
  Double_t binmin[kMaxCuts]= {0.0};
  Double_t binmax[kMaxCuts]= {0.0};
  if(!binvarsF) { delete binvarsL; return; }
  if(!binvarsL) { delete binvarsF; return; }
  for(Int_t ivar=0; ivar<binvarsF->GetEntriesFast(); ivar++) {

    TString elementF=binvarsF->At(ivar)->GetName();
    TString elementL=binvarsL->At(ivar)->GetName();
    AliDebug(1,Form(" binvar %d: %s,%s",ivar,elementF.Data(),elementL.Data()));

    switch(ivar%3) {
    case 0: continue; break;
    case 1: binmin[(int)ivar/3]=atof(elementF.Data()); break;
    case 2: binmax[(int)ivar/3]=atof(elementL.Data()); break;
    }

    binvarsF->AddAt(0x0,ivar);
  }
  binvarsF->Compress();

  // loop over all vars and cuts, check for missing stuff
  for(Int_t ivar=0; ivar<binvarsF->GetEntriesFast(); ivar++) {

    TString binvar=binvarsF->At(ivar)->GetName();
    Bool_t selected=kFALSE;

    AliDebug(1,Form(" check cuts %d %s [%.2f,%.2f]",ivar,binvar.Data(),binmin[ivar],binmax[ivar]));
    // loop over all cuts and check for missing stuff
    for(Int_t icut=0; icut<fCutLowLimits.GetNrows(); icut++) {
      if(binvar.Contains(fCutVars->At(icut)->GetName())) { selected=kTRUE; break; }
      //      else break;
    }

    // add missing cut with max limits
    if(!selected) {
      AliWarning(Form(" Bin variable %s not covered. Add cut!",binvar.Data()));
      Bool_t leg = binvar.BeginsWith("Leg");
      if(leg) binvar.Remove(0,3);
      SetRangeUser(binvar.Data(),binmin[ivar],binmax[ivar], leg);
    }

  }

  // clean up
  delete binvarsF;
  delete binvarsL;
}

//________________________________________________________________
void AliDielectronHFhelper::Print(const Option_t* /*option*/) const
{

  //
  // Print out object contents
  //
  AliInfo(Form(" Container:               %s",fMainArr->GetName()));

  // pairtypes, steps and sources
  Int_t stepLast=0;
  AliInfo(Form(" Number of filled steps:  %d",fMainArr->GetEntries()));
  for(Int_t istep=0; istep<fMainArr->GetEntriesFast(); istep++) {
    if(!fMainArr->At(istep))                             continue;
    if(!((TObjArray*)fMainArr->At(istep))->GetEntries()) continue;
    AliInfo(Form(" step %d: %s",istep,fMainArr->At(istep)->GetName()));
    stepLast=istep;
  }

  AliInfo(Form(" Number of objects:    %d",
	       ((TObjArray*) ((TObjArray*)fMainArr->At(stepLast)) ->First())->GetEntriesFast()));

  TString title       = ((TObjArray*)fMainArr->At(stepLast))->First()->GetName();
  TObjArray* binvars  = title.Tokenize(":");
  AliInfo(Form(" Number of variables:     %d",binvars->GetEntriesFast()));
  delete binvars;

  /* REACTIVATE
  // check for missing cuts
  CheckCuts(((TObjArray*)fMainArr->At(stepLast)));
  // loop over all cuts
  for(Int_t icut=0; icut<fCutLowLimits.GetNrows(); icut++) {
    const char *cutvar  = fCutVars->At(icut)->GetName(); // cut variable
    Double_t min        = fCutLowLimits(icut);           // lower limit
    Double_t max        = fCutUpLimits(icut);            // upper limit
    AliInfo(Form(" variable %d: %s [%.2f,%.2f]",icut,cutvar,min,max));
  }
  */
  TObjArray* binvars2  = title.Tokenize(":#");
  for(Int_t ivar=0; ivar<binvars2->GetEntriesFast(); ivar++) {
    if(ivar%3) continue;
    AliInfo(Form(" variable %.0f: %s",((Double_t)ivar)/3+1,binvars2->At(ivar)->GetName()));
  }
  delete binvars2;

}

//________________________________________________________________
void AliDielectronHFhelper::PrintCuts()
{

  //
  // Print cuts
  //

  // loop over all cuts
  AliInfo(" Selected cuts:");
  for(Int_t icut=0; icut<fCutLowLimits.GetNrows(); icut++)
    AliInfo(Form(" %d: %s [%.2f,%.2f]",icut,fCutVars->At(icut)->GetName(),fCutLowLimits(icut),fCutUpLimits(icut)));

}

 AliDielectronHFhelper.cxx:1
 AliDielectronHFhelper.cxx:2
 AliDielectronHFhelper.cxx:3
 AliDielectronHFhelper.cxx:4
 AliDielectronHFhelper.cxx:5
 AliDielectronHFhelper.cxx:6
 AliDielectronHFhelper.cxx:7
 AliDielectronHFhelper.cxx:8
 AliDielectronHFhelper.cxx:9
 AliDielectronHFhelper.cxx:10
 AliDielectronHFhelper.cxx:11
 AliDielectronHFhelper.cxx:12
 AliDielectronHFhelper.cxx:13
 AliDielectronHFhelper.cxx:14
 AliDielectronHFhelper.cxx:15
 AliDielectronHFhelper.cxx:16
 AliDielectronHFhelper.cxx:17
 AliDielectronHFhelper.cxx:18
 AliDielectronHFhelper.cxx:19
 AliDielectronHFhelper.cxx:20
 AliDielectronHFhelper.cxx:21
 AliDielectronHFhelper.cxx:22
 AliDielectronHFhelper.cxx:23
 AliDielectronHFhelper.cxx:24
 AliDielectronHFhelper.cxx:25
 AliDielectronHFhelper.cxx:26
 AliDielectronHFhelper.cxx:27
 AliDielectronHFhelper.cxx:28
 AliDielectronHFhelper.cxx:29
 AliDielectronHFhelper.cxx:30
 AliDielectronHFhelper.cxx:31
 AliDielectronHFhelper.cxx:32
 AliDielectronHFhelper.cxx:33
 AliDielectronHFhelper.cxx:34
 AliDielectronHFhelper.cxx:35
 AliDielectronHFhelper.cxx:36
 AliDielectronHFhelper.cxx:37
 AliDielectronHFhelper.cxx:38
 AliDielectronHFhelper.cxx:39
 AliDielectronHFhelper.cxx:40
 AliDielectronHFhelper.cxx:41
 AliDielectronHFhelper.cxx:42
 AliDielectronHFhelper.cxx:43
 AliDielectronHFhelper.cxx:44
 AliDielectronHFhelper.cxx:45
 AliDielectronHFhelper.cxx:46
 AliDielectronHFhelper.cxx:47
 AliDielectronHFhelper.cxx:48
 AliDielectronHFhelper.cxx:49
 AliDielectronHFhelper.cxx:50
 AliDielectronHFhelper.cxx:51
 AliDielectronHFhelper.cxx:52
 AliDielectronHFhelper.cxx:53
 AliDielectronHFhelper.cxx:54
 AliDielectronHFhelper.cxx:55
 AliDielectronHFhelper.cxx:56
 AliDielectronHFhelper.cxx:57
 AliDielectronHFhelper.cxx:58
 AliDielectronHFhelper.cxx:59
 AliDielectronHFhelper.cxx:60
 AliDielectronHFhelper.cxx:61
 AliDielectronHFhelper.cxx:62
 AliDielectronHFhelper.cxx:63
 AliDielectronHFhelper.cxx:64
 AliDielectronHFhelper.cxx:65
 AliDielectronHFhelper.cxx:66
 AliDielectronHFhelper.cxx:67
 AliDielectronHFhelper.cxx:68
 AliDielectronHFhelper.cxx:69
 AliDielectronHFhelper.cxx:70
 AliDielectronHFhelper.cxx:71
 AliDielectronHFhelper.cxx:72
 AliDielectronHFhelper.cxx:73
 AliDielectronHFhelper.cxx:74
 AliDielectronHFhelper.cxx:75
 AliDielectronHFhelper.cxx:76
 AliDielectronHFhelper.cxx:77
 AliDielectronHFhelper.cxx:78
 AliDielectronHFhelper.cxx:79
 AliDielectronHFhelper.cxx:80
 AliDielectronHFhelper.cxx:81
 AliDielectronHFhelper.cxx:82
 AliDielectronHFhelper.cxx:83
 AliDielectronHFhelper.cxx:84
 AliDielectronHFhelper.cxx:85
 AliDielectronHFhelper.cxx:86
 AliDielectronHFhelper.cxx:87
 AliDielectronHFhelper.cxx:88
 AliDielectronHFhelper.cxx:89
 AliDielectronHFhelper.cxx:90
 AliDielectronHFhelper.cxx:91
 AliDielectronHFhelper.cxx:92
 AliDielectronHFhelper.cxx:93
 AliDielectronHFhelper.cxx:94
 AliDielectronHFhelper.cxx:95
 AliDielectronHFhelper.cxx:96
 AliDielectronHFhelper.cxx:97
 AliDielectronHFhelper.cxx:98
 AliDielectronHFhelper.cxx:99
 AliDielectronHFhelper.cxx:100
 AliDielectronHFhelper.cxx:101
 AliDielectronHFhelper.cxx:102
 AliDielectronHFhelper.cxx:103
 AliDielectronHFhelper.cxx:104
 AliDielectronHFhelper.cxx:105
 AliDielectronHFhelper.cxx:106
 AliDielectronHFhelper.cxx:107
 AliDielectronHFhelper.cxx:108
 AliDielectronHFhelper.cxx:109
 AliDielectronHFhelper.cxx:110
 AliDielectronHFhelper.cxx:111
 AliDielectronHFhelper.cxx:112
 AliDielectronHFhelper.cxx:113
 AliDielectronHFhelper.cxx:114
 AliDielectronHFhelper.cxx:115
 AliDielectronHFhelper.cxx:116
 AliDielectronHFhelper.cxx:117
 AliDielectronHFhelper.cxx:118
 AliDielectronHFhelper.cxx:119
 AliDielectronHFhelper.cxx:120
 AliDielectronHFhelper.cxx:121
 AliDielectronHFhelper.cxx:122
 AliDielectronHFhelper.cxx:123
 AliDielectronHFhelper.cxx:124
 AliDielectronHFhelper.cxx:125
 AliDielectronHFhelper.cxx:126
 AliDielectronHFhelper.cxx:127
 AliDielectronHFhelper.cxx:128
 AliDielectronHFhelper.cxx:129
 AliDielectronHFhelper.cxx:130
 AliDielectronHFhelper.cxx:131
 AliDielectronHFhelper.cxx:132
 AliDielectronHFhelper.cxx:133
 AliDielectronHFhelper.cxx:134
 AliDielectronHFhelper.cxx:135
 AliDielectronHFhelper.cxx:136
 AliDielectronHFhelper.cxx:137
 AliDielectronHFhelper.cxx:138
 AliDielectronHFhelper.cxx:139
 AliDielectronHFhelper.cxx:140
 AliDielectronHFhelper.cxx:141
 AliDielectronHFhelper.cxx:142
 AliDielectronHFhelper.cxx:143
 AliDielectronHFhelper.cxx:144
 AliDielectronHFhelper.cxx:145
 AliDielectronHFhelper.cxx:146
 AliDielectronHFhelper.cxx:147
 AliDielectronHFhelper.cxx:148
 AliDielectronHFhelper.cxx:149
 AliDielectronHFhelper.cxx:150
 AliDielectronHFhelper.cxx:151
 AliDielectronHFhelper.cxx:152
 AliDielectronHFhelper.cxx:153
 AliDielectronHFhelper.cxx:154
 AliDielectronHFhelper.cxx:155
 AliDielectronHFhelper.cxx:156
 AliDielectronHFhelper.cxx:157
 AliDielectronHFhelper.cxx:158
 AliDielectronHFhelper.cxx:159
 AliDielectronHFhelper.cxx:160
 AliDielectronHFhelper.cxx:161
 AliDielectronHFhelper.cxx:162
 AliDielectronHFhelper.cxx:163
 AliDielectronHFhelper.cxx:164
 AliDielectronHFhelper.cxx:165
 AliDielectronHFhelper.cxx:166
 AliDielectronHFhelper.cxx:167
 AliDielectronHFhelper.cxx:168
 AliDielectronHFhelper.cxx:169
 AliDielectronHFhelper.cxx:170
 AliDielectronHFhelper.cxx:171
 AliDielectronHFhelper.cxx:172
 AliDielectronHFhelper.cxx:173
 AliDielectronHFhelper.cxx:174
 AliDielectronHFhelper.cxx:175
 AliDielectronHFhelper.cxx:176
 AliDielectronHFhelper.cxx:177
 AliDielectronHFhelper.cxx:178
 AliDielectronHFhelper.cxx:179
 AliDielectronHFhelper.cxx:180
 AliDielectronHFhelper.cxx:181
 AliDielectronHFhelper.cxx:182
 AliDielectronHFhelper.cxx:183
 AliDielectronHFhelper.cxx:184
 AliDielectronHFhelper.cxx:185
 AliDielectronHFhelper.cxx:186
 AliDielectronHFhelper.cxx:187
 AliDielectronHFhelper.cxx:188
 AliDielectronHFhelper.cxx:189
 AliDielectronHFhelper.cxx:190
 AliDielectronHFhelper.cxx:191
 AliDielectronHFhelper.cxx:192
 AliDielectronHFhelper.cxx:193
 AliDielectronHFhelper.cxx:194
 AliDielectronHFhelper.cxx:195
 AliDielectronHFhelper.cxx:196
 AliDielectronHFhelper.cxx:197
 AliDielectronHFhelper.cxx:198
 AliDielectronHFhelper.cxx:199
 AliDielectronHFhelper.cxx:200
 AliDielectronHFhelper.cxx:201
 AliDielectronHFhelper.cxx:202
 AliDielectronHFhelper.cxx:203
 AliDielectronHFhelper.cxx:204
 AliDielectronHFhelper.cxx:205
 AliDielectronHFhelper.cxx:206
 AliDielectronHFhelper.cxx:207
 AliDielectronHFhelper.cxx:208
 AliDielectronHFhelper.cxx:209
 AliDielectronHFhelper.cxx:210
 AliDielectronHFhelper.cxx:211
 AliDielectronHFhelper.cxx:212
 AliDielectronHFhelper.cxx:213
 AliDielectronHFhelper.cxx:214
 AliDielectronHFhelper.cxx:215
 AliDielectronHFhelper.cxx:216
 AliDielectronHFhelper.cxx:217
 AliDielectronHFhelper.cxx:218
 AliDielectronHFhelper.cxx:219
 AliDielectronHFhelper.cxx:220
 AliDielectronHFhelper.cxx:221
 AliDielectronHFhelper.cxx:222
 AliDielectronHFhelper.cxx:223
 AliDielectronHFhelper.cxx:224
 AliDielectronHFhelper.cxx:225
 AliDielectronHFhelper.cxx:226
 AliDielectronHFhelper.cxx:227
 AliDielectronHFhelper.cxx:228
 AliDielectronHFhelper.cxx:229
 AliDielectronHFhelper.cxx:230
 AliDielectronHFhelper.cxx:231
 AliDielectronHFhelper.cxx:232
 AliDielectronHFhelper.cxx:233
 AliDielectronHFhelper.cxx:234
 AliDielectronHFhelper.cxx:235
 AliDielectronHFhelper.cxx:236
 AliDielectronHFhelper.cxx:237
 AliDielectronHFhelper.cxx:238
 AliDielectronHFhelper.cxx:239
 AliDielectronHFhelper.cxx:240
 AliDielectronHFhelper.cxx:241
 AliDielectronHFhelper.cxx:242
 AliDielectronHFhelper.cxx:243
 AliDielectronHFhelper.cxx:244
 AliDielectronHFhelper.cxx:245
 AliDielectronHFhelper.cxx:246
 AliDielectronHFhelper.cxx:247
 AliDielectronHFhelper.cxx:248
 AliDielectronHFhelper.cxx:249
 AliDielectronHFhelper.cxx:250
 AliDielectronHFhelper.cxx:251
 AliDielectronHFhelper.cxx:252
 AliDielectronHFhelper.cxx:253
 AliDielectronHFhelper.cxx:254
 AliDielectronHFhelper.cxx:255
 AliDielectronHFhelper.cxx:256
 AliDielectronHFhelper.cxx:257
 AliDielectronHFhelper.cxx:258
 AliDielectronHFhelper.cxx:259
 AliDielectronHFhelper.cxx:260
 AliDielectronHFhelper.cxx:261
 AliDielectronHFhelper.cxx:262
 AliDielectronHFhelper.cxx:263
 AliDielectronHFhelper.cxx:264
 AliDielectronHFhelper.cxx:265
 AliDielectronHFhelper.cxx:266
 AliDielectronHFhelper.cxx:267
 AliDielectronHFhelper.cxx:268
 AliDielectronHFhelper.cxx:269
 AliDielectronHFhelper.cxx:270
 AliDielectronHFhelper.cxx:271
 AliDielectronHFhelper.cxx:272
 AliDielectronHFhelper.cxx:273
 AliDielectronHFhelper.cxx:274
 AliDielectronHFhelper.cxx:275
 AliDielectronHFhelper.cxx:276
 AliDielectronHFhelper.cxx:277
 AliDielectronHFhelper.cxx:278
 AliDielectronHFhelper.cxx:279
 AliDielectronHFhelper.cxx:280
 AliDielectronHFhelper.cxx:281
 AliDielectronHFhelper.cxx:282
 AliDielectronHFhelper.cxx:283
 AliDielectronHFhelper.cxx:284
 AliDielectronHFhelper.cxx:285
 AliDielectronHFhelper.cxx:286
 AliDielectronHFhelper.cxx:287
 AliDielectronHFhelper.cxx:288
 AliDielectronHFhelper.cxx:289
 AliDielectronHFhelper.cxx:290
 AliDielectronHFhelper.cxx:291
 AliDielectronHFhelper.cxx:292
 AliDielectronHFhelper.cxx:293
 AliDielectronHFhelper.cxx:294
 AliDielectronHFhelper.cxx:295
 AliDielectronHFhelper.cxx:296
 AliDielectronHFhelper.cxx:297
 AliDielectronHFhelper.cxx:298
 AliDielectronHFhelper.cxx:299
 AliDielectronHFhelper.cxx:300
 AliDielectronHFhelper.cxx:301
 AliDielectronHFhelper.cxx:302
 AliDielectronHFhelper.cxx:303
 AliDielectronHFhelper.cxx:304
 AliDielectronHFhelper.cxx:305
 AliDielectronHFhelper.cxx:306
 AliDielectronHFhelper.cxx:307
 AliDielectronHFhelper.cxx:308
 AliDielectronHFhelper.cxx:309
 AliDielectronHFhelper.cxx:310
 AliDielectronHFhelper.cxx:311
 AliDielectronHFhelper.cxx:312
 AliDielectronHFhelper.cxx:313
 AliDielectronHFhelper.cxx:314
 AliDielectronHFhelper.cxx:315
 AliDielectronHFhelper.cxx:316
 AliDielectronHFhelper.cxx:317
 AliDielectronHFhelper.cxx:318
 AliDielectronHFhelper.cxx:319
 AliDielectronHFhelper.cxx:320
 AliDielectronHFhelper.cxx:321
 AliDielectronHFhelper.cxx:322
 AliDielectronHFhelper.cxx:323
 AliDielectronHFhelper.cxx:324
 AliDielectronHFhelper.cxx:325
 AliDielectronHFhelper.cxx:326
 AliDielectronHFhelper.cxx:327
 AliDielectronHFhelper.cxx:328
 AliDielectronHFhelper.cxx:329
 AliDielectronHFhelper.cxx:330
 AliDielectronHFhelper.cxx:331
 AliDielectronHFhelper.cxx:332
 AliDielectronHFhelper.cxx:333
 AliDielectronHFhelper.cxx:334
 AliDielectronHFhelper.cxx:335
 AliDielectronHFhelper.cxx:336
 AliDielectronHFhelper.cxx:337
 AliDielectronHFhelper.cxx:338
 AliDielectronHFhelper.cxx:339
 AliDielectronHFhelper.cxx:340
 AliDielectronHFhelper.cxx:341
 AliDielectronHFhelper.cxx:342
 AliDielectronHFhelper.cxx:343
 AliDielectronHFhelper.cxx:344
 AliDielectronHFhelper.cxx:345
 AliDielectronHFhelper.cxx:346
 AliDielectronHFhelper.cxx:347
 AliDielectronHFhelper.cxx:348
 AliDielectronHFhelper.cxx:349
 AliDielectronHFhelper.cxx:350
 AliDielectronHFhelper.cxx:351
 AliDielectronHFhelper.cxx:352
 AliDielectronHFhelper.cxx:353
 AliDielectronHFhelper.cxx:354
 AliDielectronHFhelper.cxx:355
 AliDielectronHFhelper.cxx:356
 AliDielectronHFhelper.cxx:357
 AliDielectronHFhelper.cxx:358
 AliDielectronHFhelper.cxx:359
 AliDielectronHFhelper.cxx:360
 AliDielectronHFhelper.cxx:361
 AliDielectronHFhelper.cxx:362
 AliDielectronHFhelper.cxx:363
 AliDielectronHFhelper.cxx:364
 AliDielectronHFhelper.cxx:365
 AliDielectronHFhelper.cxx:366
 AliDielectronHFhelper.cxx:367
 AliDielectronHFhelper.cxx:368
 AliDielectronHFhelper.cxx:369
 AliDielectronHFhelper.cxx:370
 AliDielectronHFhelper.cxx:371
 AliDielectronHFhelper.cxx:372
 AliDielectronHFhelper.cxx:373
 AliDielectronHFhelper.cxx:374
 AliDielectronHFhelper.cxx:375
 AliDielectronHFhelper.cxx:376
 AliDielectronHFhelper.cxx:377
 AliDielectronHFhelper.cxx:378
 AliDielectronHFhelper.cxx:379
 AliDielectronHFhelper.cxx:380
 AliDielectronHFhelper.cxx:381
 AliDielectronHFhelper.cxx:382
 AliDielectronHFhelper.cxx:383
 AliDielectronHFhelper.cxx:384
 AliDielectronHFhelper.cxx:385
 AliDielectronHFhelper.cxx:386
 AliDielectronHFhelper.cxx:387
 AliDielectronHFhelper.cxx:388
 AliDielectronHFhelper.cxx:389
 AliDielectronHFhelper.cxx:390
 AliDielectronHFhelper.cxx:391
 AliDielectronHFhelper.cxx:392
 AliDielectronHFhelper.cxx:393
 AliDielectronHFhelper.cxx:394
 AliDielectronHFhelper.cxx:395
 AliDielectronHFhelper.cxx:396
 AliDielectronHFhelper.cxx:397
 AliDielectronHFhelper.cxx:398
 AliDielectronHFhelper.cxx:399
 AliDielectronHFhelper.cxx:400
 AliDielectronHFhelper.cxx:401
 AliDielectronHFhelper.cxx:402
 AliDielectronHFhelper.cxx:403
 AliDielectronHFhelper.cxx:404
 AliDielectronHFhelper.cxx:405
 AliDielectronHFhelper.cxx:406
 AliDielectronHFhelper.cxx:407
 AliDielectronHFhelper.cxx:408
 AliDielectronHFhelper.cxx:409
 AliDielectronHFhelper.cxx:410
 AliDielectronHFhelper.cxx:411
 AliDielectronHFhelper.cxx:412
 AliDielectronHFhelper.cxx:413
 AliDielectronHFhelper.cxx:414
 AliDielectronHFhelper.cxx:415
 AliDielectronHFhelper.cxx:416
 AliDielectronHFhelper.cxx:417
 AliDielectronHFhelper.cxx:418
 AliDielectronHFhelper.cxx:419
 AliDielectronHFhelper.cxx:420
 AliDielectronHFhelper.cxx:421
 AliDielectronHFhelper.cxx:422
 AliDielectronHFhelper.cxx:423
 AliDielectronHFhelper.cxx:424
 AliDielectronHFhelper.cxx:425
 AliDielectronHFhelper.cxx:426
 AliDielectronHFhelper.cxx:427
 AliDielectronHFhelper.cxx:428
 AliDielectronHFhelper.cxx:429
 AliDielectronHFhelper.cxx:430
 AliDielectronHFhelper.cxx:431
 AliDielectronHFhelper.cxx:432
 AliDielectronHFhelper.cxx:433
 AliDielectronHFhelper.cxx:434
 AliDielectronHFhelper.cxx:435
 AliDielectronHFhelper.cxx:436
 AliDielectronHFhelper.cxx:437
 AliDielectronHFhelper.cxx:438
 AliDielectronHFhelper.cxx:439
 AliDielectronHFhelper.cxx:440
 AliDielectronHFhelper.cxx:441
 AliDielectronHFhelper.cxx:442
 AliDielectronHFhelper.cxx:443
 AliDielectronHFhelper.cxx:444
 AliDielectronHFhelper.cxx:445
 AliDielectronHFhelper.cxx:446
 AliDielectronHFhelper.cxx:447
 AliDielectronHFhelper.cxx:448
 AliDielectronHFhelper.cxx:449
 AliDielectronHFhelper.cxx:450
 AliDielectronHFhelper.cxx:451
 AliDielectronHFhelper.cxx:452
 AliDielectronHFhelper.cxx:453
 AliDielectronHFhelper.cxx:454
 AliDielectronHFhelper.cxx:455
 AliDielectronHFhelper.cxx:456
 AliDielectronHFhelper.cxx:457
 AliDielectronHFhelper.cxx:458
 AliDielectronHFhelper.cxx:459
 AliDielectronHFhelper.cxx:460
 AliDielectronHFhelper.cxx:461
 AliDielectronHFhelper.cxx:462
 AliDielectronHFhelper.cxx:463
 AliDielectronHFhelper.cxx:464
 AliDielectronHFhelper.cxx:465
 AliDielectronHFhelper.cxx:466
 AliDielectronHFhelper.cxx:467
 AliDielectronHFhelper.cxx:468
 AliDielectronHFhelper.cxx:469
 AliDielectronHFhelper.cxx:470
 AliDielectronHFhelper.cxx:471
 AliDielectronHFhelper.cxx:472
 AliDielectronHFhelper.cxx:473
 AliDielectronHFhelper.cxx:474
 AliDielectronHFhelper.cxx:475
 AliDielectronHFhelper.cxx:476
 AliDielectronHFhelper.cxx:477
 AliDielectronHFhelper.cxx:478
 AliDielectronHFhelper.cxx:479
 AliDielectronHFhelper.cxx:480
 AliDielectronHFhelper.cxx:481
 AliDielectronHFhelper.cxx:482
 AliDielectronHFhelper.cxx:483
 AliDielectronHFhelper.cxx:484
 AliDielectronHFhelper.cxx:485
 AliDielectronHFhelper.cxx:486
 AliDielectronHFhelper.cxx:487
 AliDielectronHFhelper.cxx:488
 AliDielectronHFhelper.cxx:489
 AliDielectronHFhelper.cxx:490
 AliDielectronHFhelper.cxx:491
 AliDielectronHFhelper.cxx:492
 AliDielectronHFhelper.cxx:493
 AliDielectronHFhelper.cxx:494
 AliDielectronHFhelper.cxx:495
 AliDielectronHFhelper.cxx:496
 AliDielectronHFhelper.cxx:497
 AliDielectronHFhelper.cxx:498
 AliDielectronHFhelper.cxx:499
 AliDielectronHFhelper.cxx:500
 AliDielectronHFhelper.cxx:501
 AliDielectronHFhelper.cxx:502
 AliDielectronHFhelper.cxx:503
 AliDielectronHFhelper.cxx:504
 AliDielectronHFhelper.cxx:505
 AliDielectronHFhelper.cxx:506
 AliDielectronHFhelper.cxx:507
 AliDielectronHFhelper.cxx:508
 AliDielectronHFhelper.cxx:509
 AliDielectronHFhelper.cxx:510
 AliDielectronHFhelper.cxx:511
 AliDielectronHFhelper.cxx:512
 AliDielectronHFhelper.cxx:513
 AliDielectronHFhelper.cxx:514
 AliDielectronHFhelper.cxx:515
 AliDielectronHFhelper.cxx:516
 AliDielectronHFhelper.cxx:517
 AliDielectronHFhelper.cxx:518
 AliDielectronHFhelper.cxx:519
 AliDielectronHFhelper.cxx:520
 AliDielectronHFhelper.cxx:521
 AliDielectronHFhelper.cxx:522
 AliDielectronHFhelper.cxx:523
 AliDielectronHFhelper.cxx:524
 AliDielectronHFhelper.cxx:525
 AliDielectronHFhelper.cxx:526
 AliDielectronHFhelper.cxx:527
 AliDielectronHFhelper.cxx:528
 AliDielectronHFhelper.cxx:529
 AliDielectronHFhelper.cxx:530
 AliDielectronHFhelper.cxx:531
 AliDielectronHFhelper.cxx:532
 AliDielectronHFhelper.cxx:533
 AliDielectronHFhelper.cxx:534
 AliDielectronHFhelper.cxx:535
 AliDielectronHFhelper.cxx:536
 AliDielectronHFhelper.cxx:537
 AliDielectronHFhelper.cxx:538
 AliDielectronHFhelper.cxx:539
 AliDielectronHFhelper.cxx:540
 AliDielectronHFhelper.cxx:541
 AliDielectronHFhelper.cxx:542
 AliDielectronHFhelper.cxx:543
 AliDielectronHFhelper.cxx:544
 AliDielectronHFhelper.cxx:545
 AliDielectronHFhelper.cxx:546
 AliDielectronHFhelper.cxx:547
 AliDielectronHFhelper.cxx:548
 AliDielectronHFhelper.cxx:549
 AliDielectronHFhelper.cxx:550
 AliDielectronHFhelper.cxx:551
 AliDielectronHFhelper.cxx:552
 AliDielectronHFhelper.cxx:553
 AliDielectronHFhelper.cxx:554
 AliDielectronHFhelper.cxx:555
 AliDielectronHFhelper.cxx:556
 AliDielectronHFhelper.cxx:557
 AliDielectronHFhelper.cxx:558
 AliDielectronHFhelper.cxx:559
 AliDielectronHFhelper.cxx:560
 AliDielectronHFhelper.cxx:561
 AliDielectronHFhelper.cxx:562
 AliDielectronHFhelper.cxx:563
 AliDielectronHFhelper.cxx:564
 AliDielectronHFhelper.cxx:565
 AliDielectronHFhelper.cxx:566
 AliDielectronHFhelper.cxx:567
 AliDielectronHFhelper.cxx:568
 AliDielectronHFhelper.cxx:569
 AliDielectronHFhelper.cxx:570
 AliDielectronHFhelper.cxx:571
 AliDielectronHFhelper.cxx:572
 AliDielectronHFhelper.cxx:573