ROOT logo
/**************************************************************************
* Copyright(c) 1998-1999, 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.                  *
**************************************************************************/
//
// Task for Efficiency studies
// Used for testing classes AliHFEcontainer and AliHFEfilter
// Creates Efficiency Histograms
//
// Author:
//   Markus Fasel <M.Fasel@gsi.de>
//
#include <TAxis.h>
#include <TCanvas.h>
#include <TFile.h>
#include <TGraphErrors.h>
#include <TH1F.h>
#include <TH2.h>
#include <TIterator.h>
#include <TLegend.h>
#include <TObjArray.h>
#include <TPad.h>

#include "AliAnalysisManager.h"
#include "AliCFAcceptanceCuts.h"
#include "AliCFTrackIsPrimaryCuts.h"
#include "AliCFContainer.h"
#include "AliCFEffGrid.h"
#include "AliESDEvent.h"
#include "AliHFEcollection.h"
#include "AliHFEcontainer.h"
#include "AliHFEcutStep.h"
#include "AliHFEefficiency.h"
#include "AliHFEextraCuts.h"
#include "AliHFEtrackFilter.h"
#include "AliHFEtools.h"
#include "AliMCEvent.h"
#include "AliMCEventHandler.h"

ClassImp(AliHFEefficiency)

AliHFEefficiency::AliHFEefficiency():
  AliAnalysisTaskSE(),
  fFilter(NULL),
  fMCcut(NULL),
  fAcceptanceCuts(NULL),
  fEfficiency(NULL),
  fOutput(NULL),
  fCutTRD(kFALSE)
{
  //
  // Default constructor
  //
}

AliHFEefficiency::AliHFEefficiency(const Char_t *name):
  AliAnalysisTaskSE(name),
  fFilter(NULL),
  fMCcut(NULL),
  fAcceptanceCuts(NULL),
  fEfficiency(NULL),
  fOutput(NULL),
  fCutTRD(kFALSE)
{
  //
  // Main Constructor
  //
  DefineOutput(1, AliHFEcontainer::Class());
  DefineOutput(2, TList::Class());

  SetRunTerminate();
}

AliHFEefficiency::~AliHFEefficiency(){
  //
  // Destructor
  //
  if(fFilter) delete fFilter;
  if(fEfficiency) delete fEfficiency;
  if(fAcceptanceCuts) delete fAcceptanceCuts;
  if(fOutput) delete fOutput;
}

void AliHFEefficiency::UserCreateOutputObjects(){
  //
  // Create the output objects
  //
  fEfficiency = new AliHFEcontainer("Efficiency");
  fEfficiency->SetNumberOfVariables(4);
  // InitializeVariables
  fEfficiency->SetVariableName(0, "pt");
  fEfficiency->MakeLogarithmicBinning(0, 40, 0.1, 10);
  fEfficiency->SetVariableName(1, "Eta");
  fEfficiency->MakeLinearBinning(1, 20, -0.9, 0.9);
  fEfficiency->SetVariableName(2, "phi");
  fEfficiency->MakeLinearBinning(2 , 18, 0, 2 * TMath::Pi());
  fEfficiency->SetVariableName(3, "Charge");
  fEfficiency->MakeLinearBinning(3, 2, -1.1, 1.1);
  // New container
  fEfficiency->CreateContainer("MCFilter", "Efficiency for MC Studies", 2);
  AliCFContainer *cont = fEfficiency->GetCFContainer("MCFilter");
  cont->SetStepTitle(0, "MC Signal");
  cont->SetStepTitle(1, "MC Signal in acceptance");

  fFilter = new AliHFEtrackFilter("testfilter");
  fFilter->GenerateCutSteps();
  fFilter->InitCF(fEfficiency);
  fMCcut = fFilter->GetMCSignalCuts();
  if(fCutTRD){
    AliHFEcutStep *hfeTRD = fFilter->GetCutStep("HFETRD");
    AliHFEextraCuts *hfecuts = dynamic_cast<AliHFEextraCuts *>(hfeTRD->GetCut("HFETRDCuts"));
    if(hfecuts) hfecuts->SetMinTrackletsTRD(4);
  }
  AliHFEcutStep *hfeITS = fFilter->GetCutStep("HFEITS");
  AliHFEextraCuts *hfeitscuts = dynamic_cast<AliHFEextraCuts *>(hfeITS->GetCut("HFEPixelsCuts"));
  if(hfeitscuts)hfeitscuts->SetRequireITSpixel(AliHFEextraCuts::kFirst);

  AliHFEcutStep *primary = fFilter->GetCutStep("Primary");
  AliCFTrackIsPrimaryCuts *primcuts = dynamic_cast<AliCFTrackIsPrimaryCuts *>(primary->GetCut("PrimaryCuts"));
  if(primcuts){ 
    primcuts->SetMaxDCAToVertexXY(3);
    primcuts->SetMaxDCAToVertexZ(5);
  }

  fAcceptanceCuts = new AliCFAcceptanceCuts("Acceptance", "MC Acceptance Cuts");
  fAcceptanceCuts->SetMinNHitITS(3);
  fAcceptanceCuts->SetMinNHitTPC(2);
  fAcceptanceCuts->SetMinNHitTRD(0);

  // create additional histos for testing purpose
  fOutput = new AliHFEcollection("histos", "QA histograms");
  fOutput->CreateTH1F("hNtracks","Number of Tracks per Event", 100, 0, 100);
  fOutput->CreateTH1F("hPt","Pt of the tracks", 40, 0.1, 10, 0);
  fOutput->CreateTH1F("kinkIndex", "Kink Index", 400, -200, 200);
  fOutput->CreateTH1F("itspixel", "ITS PIXEL", 2, 0, 2);
  fOutput->CreateTH1F("mcmother", "Mother PDG", 1000, -500, 500);
  fOutput->CreateTH2F("ptres", "Pt Resolution", 40, 0.1, 10, 200, -1.5, 1.5, 0);
}

void AliHFEefficiency::UserExec(Option_t *){
  //
  // Event Loop
  // Filter track, fill Efficiency container
  //
  AliESDEvent *esdevent = dynamic_cast<AliESDEvent *>(fInputEvent);
  if(!esdevent){
    AliError("ESD Event required");
    return;
  }
  fEfficiency->NewEvent();
  fFilter->SetRecEvent(fInputEvent);
  if(fMCEvent){
    AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
   if ( ! mcH ) {
     AliError("Cannot get MC truth event handler");
     return;
    }  
    if(mcH &&(!mcH->InitOk())) return;
    if(mcH &&(!mcH->TreeK())) return;
    if(mcH &&(!mcH->TreeTR())) return;
    fFilter->SetMC(fMCEvent);
    FilterMC();
  }
  fFilter->FilterTracks(esdevent);
  TObjArray *tracks = fFilter->GetFilteredTracks();
  TIterator *iter = tracks->MakeIterator();
  fOutput->Fill("hNtracks", tracks->GetEntriesFast());
  AliESDtrack *track = NULL;
  while((track = dynamic_cast<AliESDtrack *>(iter->Next()))){
    fOutput->Fill("hPt", track->Pt());  
    fOutput->Fill("kinkIndex", track->GetKinkIndex(0));
    if(TESTBIT(track->GetITSClusterMap(), 0))
      fOutput->Fill("itspixel",1);
    else
      fOutput->Fill("itspixel",0);
    if(fMCEvent){
      AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(track->GetLabel())));
      if(mctrack){
        Int_t motherLabel = mctrack->Particle()->GetFirstMother();
        if(motherLabel){
          AliMCParticle *mother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(motherLabel));
          if(mother)fOutput->Fill("mcmother", mother->Particle()->GetPdgCode());
        }
        fOutput->Fill("ptres", mctrack->Pt(), (track->Pt() - mctrack->Pt())/mctrack->Pt());
      }
    }
  }
  delete iter;
  fFilter->Flush();
  PostData(1, fEfficiency);
  PostData(2, fOutput->GetList());
}

void AliHFEefficiency::Terminate(Option_t *){
  //
  // Evaluate Results
  //
  fEfficiency = dynamic_cast<AliHFEcontainer *>(GetOutputData(1));
  if(!fEfficiency){
    AliError("No Output data available");
    return;
  }

  if(!IsRunningTerminate()) return;
  PostProcess();

  TList *l = dynamic_cast<TList *>(GetOutputData(2));
  if(l) DrawPtResolution(l);
}

void AliHFEefficiency::FilterMC(){
  //
  // Cut MC tracks
  //
  Double_t cont[4] = {0., 0., 0., 0.};
  AliMCParticle *track = NULL;
  for(Int_t itrack = 0; itrack < fMCEvent->GetNumberOfTracks(); itrack++){
    track = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(itrack));
    if(!track) continue;
    if(!fMCcut->IsSelected(track)) continue;
    cont[0] = track->Pt();
    cont[1] = track->Eta();
    cont[2] = track->Phi();
    cont[3] = track->Charge()/3;
    //AliMCParticle *mother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(track->Particle()->GetFirstMother()));
    //if(TMath::Abs(mother->Particle()->GetPdgCode()) != 443) continue;
    fEfficiency->FillCFContainer("MCFilter", 0, cont);
    if(!fAcceptanceCuts->IsSelected(track)) continue;
    fEfficiency->FillCFContainer("MCFilter", 1, cont);
  }
}

void AliHFEefficiency::Load(const char* filename){
  //
  // Load results for post processing
  //
  TFile *input = TFile::Open(filename);
  AliHFEcontainer *ctmp = dynamic_cast<AliHFEcontainer *>(input->Get("Efficiency"));
  if(ctmp) fEfficiency = dynamic_cast<AliHFEcontainer *>(ctmp->Clone());
  input->Close();
  delete input;
}

void AliHFEefficiency::PostProcess(){
  //
  // Calculate Efficiencies
  //
  fEfficiency->Print();

  // Prepare containers
  AliCFContainer *crec  = fEfficiency->MakeMergedCFContainer("merged", "Merged Container", "MCFilter:testfilter_container_signal");
  AliCFContainer *cmc   = fEfficiency->MakeMergedCFContainer("merged", "Merged Container", "MCFilter:testfilter_container_signalMC");
  AliCFEffGrid *gridRec = new AliCFEffGrid("effrec", "Efficiency Calculation using rec vars", *crec);
  AliCFEffGrid *gridMC  = new AliCFEffGrid("effmc", "Efficiency Calculation using mc vars", *cmc);

  TCanvas *ctmp = NULL;
  for(Int_t ivar = 0; ivar < crec->GetNVar(); ivar++){
    ctmp = new TCanvas(Form("cEffSignal%s", crec->GetVarTitle(ivar)), Form("Signal Efficiency (%s)", crec->GetVarTitle(ivar)));
    ctmp->cd();
    DrawSignalEfficiency(gridRec, crec, ivar);
    
    ctmp = new TCanvas(Form("cCutEff%s", crec->GetVarTitle(ivar)), Form("Cut Step Efficiency (%s)", crec->GetVarTitle(ivar)));
    ctmp->cd();
    DrawCutEfficiency(gridRec, crec, ivar);
  
    ctmp = new TCanvas(Form("cEffSignalMC%s", crec->GetVarTitle(ivar)), Form("Signal Efficiency (%s, MC)", crec->GetVarTitle(ivar)));
    ctmp->cd();
    DrawSignalEfficiency(gridMC, cmc, ivar);

    ctmp = new TCanvas(Form("cCutEffMC%s", crec->GetVarTitle(ivar)), Form("Cut Step Efficiency (%s, MC)", crec->GetVarTitle(ivar)));
    ctmp->cd();
    DrawCutEfficiency(gridMC, cmc, ivar);
  } 

  CalculatePTsmearing();

  delete gridRec;
  delete gridMC;
  delete crec;
  delete cmc;
}

void AliHFEefficiency::DrawSignalEfficiency(AliCFEffGrid *eff, AliCFContainer *cont, Int_t var){
  //
  // Efficiency of a cut step with respect to the signal
  //
  TH1 *hEff = NULL;
  Bool_t kFirst = kTRUE;
  TLegend *leg = new TLegend(0.5, 0.7, 0.89, 0.89);
  for(Int_t istep = 1; istep < cont->GetNStep(); istep++){
    eff->CalculateEfficiency(istep, 0);
    hEff = eff->Project(var);
    hEff->SetTitle(Form("Signal Efficiency (%s)", cont->GetVarTitle(var)));
    hEff->GetYaxis()->SetRangeUser(0., 1.6);
    hEff->GetYaxis()->SetTitle("Efficiency");
    hEff->SetStats(kFALSE);
    hEff->SetMarkerStyle(20+istep);
    hEff->SetMarkerColor(kBlue - 5 + istep);
    hEff->Draw(kFirst ? "ep" : "epsame");
    kFirst = kFALSE;
    leg->AddEntry(hEff, cont->GetStepTitle(istep), "p");
  }
  leg->SetBorderSize(0);
  leg->SetFillStyle(0);
  leg->Draw();
  gPad->Update();
}

void AliHFEefficiency::DrawCutEfficiency(AliCFEffGrid *eff, AliCFContainer *cont, Int_t var){
  //
  // Efficiency of a cut step with respect to the previous one
  //
  Bool_t kFirst = kTRUE;
  TH1 *hEff = NULL;
  TLegend *leg = new TLegend(0.5, 0.7, 0.89, 0.89);
  TString effTitle;
  Int_t start = 0, length = 0;
  for(Int_t istep = 1; istep < cont->GetNStep(); istep++){
    eff->CalculateEfficiency(istep, istep -1);
    hEff = eff->Project(var);
    effTitle = hEff->GetTitle();
    start = effTitle.Index("projected"); length = effTitle.Length() - start;
    effTitle.Remove(start - 1, length + 1);
    effTitle.ReplaceAll("Efficiency: ", "");
    hEff->SetTitle(Form("Cut Efficiency (%s)", cont->GetVarTitle(var)));
    hEff->GetYaxis()->SetRangeUser(0., 1.6);
    hEff->GetYaxis()->SetTitle("Efficiency");
    hEff->SetStats(kFALSE);
    hEff->SetMarkerStyle(20+istep);
    hEff->SetMarkerColor(kBlue - 5 + istep);
    hEff->Draw(kFirst ? "ep" : "epsame");
    kFirst = kFALSE;
    leg->AddEntry(hEff, cont->GetStepTitle(istep), "p");
  }
  leg->SetBorderSize(0);
  leg->SetFillStyle(0);
  leg->Draw();
  gPad->Update();
}

void AliHFEefficiency::CalculatePTsmearing(){
  //
  // Efficiency of a cut step with respect to the same cut step filled with MC 
  // Information
  //
  AliCFContainer *cont = fEfficiency->MakeMergedCFContainer("merged", "Merged Container", "testfilter_container_signal:testfilter_container_signalMC");
  AliCFEffGrid *grid = new AliCFEffGrid("effc", "Efficiency Calculation", *cont);
  Bool_t kFirst = kTRUE;
  TH1 *hEff = NULL;
  TCanvas *c2 = new TCanvas("cSmear", "pT modification");
  c2->cd();
  TLegend *leg1 = new TLegend(0.5, 0.7, 0.89, 0.89);
  Int_t nStep = cont->GetNStep()/2;
  for(Int_t istep = 0; istep < nStep; istep++){
    grid->CalculateEfficiency(istep, istep + nStep);
    hEff = grid->Project(0);
    hEff->SetStats(kFALSE);
    hEff->SetMarkerStyle(20+istep);
    hEff->SetMarkerColor(kBlue - 5 + istep);
    hEff->Draw(kFirst ? "ep" : "epsame");
    kFirst = kFALSE;
    leg1->AddEntry(hEff, cont->GetStepTitle(istep), "p");
  }
  leg1->Draw();
  gPad->Update();

  delete cont;
  delete grid;
}

void AliHFEefficiency::DrawPtResolution(const TList * const l){
  //
  // Draw pt resolution
  //
  TH2 *h2 = dynamic_cast<TH2 *>(l->FindObject("ptres"));
  if(!h2) return;
  TGraphErrors *mean = new TGraphErrors(h2->GetNbinsX());
  TGraphErrors *rms = new TGraphErrors(h2->GetNbinsX());

  Double_t pt = 0., pterr = 0., mymean = 0., myrms = 0., mymeanerr = 0., myrmserr = 0.;
  TH1 *htmp = NULL;
  for(Int_t imom = 1; imom <= h2->GetXaxis()->GetLast(); imom++){
    pt = h2->GetXaxis()->GetBinCenter(imom);
    pterr = h2->GetXaxis()->GetBinWidth(imom)/2.;
    htmp = h2->ProjectionY("py", imom, imom);
    mymean = htmp->GetMean();
    myrms = htmp->GetRMS();
    mymeanerr = htmp->GetMeanError();
    myrmserr = htmp->GetRMSError();
    delete htmp;
    mean->SetPoint(imom -1, pt, mymean);
    rms->SetPoint(imom -1, pt, myrms);
    mean->SetPointError(imom-1,pterr,mymeanerr);
    rms->SetPointError(imom-1,pterr,myrmserr);
  }

  TCanvas *c1 = new TCanvas("cPtShift", "pT Shift");
  c1->cd();
  mean->SetMarkerStyle(22);
  mean->SetMarkerColor(kBlue);
  mean->Draw("ape");
  rms->SetMarkerStyle(22);
  rms->SetMarkerColor(kBlack);
  rms->Draw("pesame");
  TLegend *leg = new TLegend(0.5, 0.7, 0.89, 0.89);
  leg->AddEntry(mean, "Mean", "lp");
  leg->AddEntry(rms, "RMS", "lp");
  leg->Draw();
  gPad->Update();
} 

 AliHFEefficiency.cxx:1
 AliHFEefficiency.cxx:2
 AliHFEefficiency.cxx:3
 AliHFEefficiency.cxx:4
 AliHFEefficiency.cxx:5
 AliHFEefficiency.cxx:6
 AliHFEefficiency.cxx:7
 AliHFEefficiency.cxx:8
 AliHFEefficiency.cxx:9
 AliHFEefficiency.cxx:10
 AliHFEefficiency.cxx:11
 AliHFEefficiency.cxx:12
 AliHFEefficiency.cxx:13
 AliHFEefficiency.cxx:14
 AliHFEefficiency.cxx:15
 AliHFEefficiency.cxx:16
 AliHFEefficiency.cxx:17
 AliHFEefficiency.cxx:18
 AliHFEefficiency.cxx:19
 AliHFEefficiency.cxx:20
 AliHFEefficiency.cxx:21
 AliHFEefficiency.cxx:22
 AliHFEefficiency.cxx:23
 AliHFEefficiency.cxx:24
 AliHFEefficiency.cxx:25
 AliHFEefficiency.cxx:26
 AliHFEefficiency.cxx:27
 AliHFEefficiency.cxx:28
 AliHFEefficiency.cxx:29
 AliHFEefficiency.cxx:30
 AliHFEefficiency.cxx:31
 AliHFEefficiency.cxx:32
 AliHFEefficiency.cxx:33
 AliHFEefficiency.cxx:34
 AliHFEefficiency.cxx:35
 AliHFEefficiency.cxx:36
 AliHFEefficiency.cxx:37
 AliHFEefficiency.cxx:38
 AliHFEefficiency.cxx:39
 AliHFEefficiency.cxx:40
 AliHFEefficiency.cxx:41
 AliHFEefficiency.cxx:42
 AliHFEefficiency.cxx:43
 AliHFEefficiency.cxx:44
 AliHFEefficiency.cxx:45
 AliHFEefficiency.cxx:46
 AliHFEefficiency.cxx:47
 AliHFEefficiency.cxx:48
 AliHFEefficiency.cxx:49
 AliHFEefficiency.cxx:50
 AliHFEefficiency.cxx:51
 AliHFEefficiency.cxx:52
 AliHFEefficiency.cxx:53
 AliHFEefficiency.cxx:54
 AliHFEefficiency.cxx:55
 AliHFEefficiency.cxx:56
 AliHFEefficiency.cxx:57
 AliHFEefficiency.cxx:58
 AliHFEefficiency.cxx:59
 AliHFEefficiency.cxx:60
 AliHFEefficiency.cxx:61
 AliHFEefficiency.cxx:62
 AliHFEefficiency.cxx:63
 AliHFEefficiency.cxx:64
 AliHFEefficiency.cxx:65
 AliHFEefficiency.cxx:66
 AliHFEefficiency.cxx:67
 AliHFEefficiency.cxx:68
 AliHFEefficiency.cxx:69
 AliHFEefficiency.cxx:70
 AliHFEefficiency.cxx:71
 AliHFEefficiency.cxx:72
 AliHFEefficiency.cxx:73
 AliHFEefficiency.cxx:74
 AliHFEefficiency.cxx:75
 AliHFEefficiency.cxx:76
 AliHFEefficiency.cxx:77
 AliHFEefficiency.cxx:78
 AliHFEefficiency.cxx:79
 AliHFEefficiency.cxx:80
 AliHFEefficiency.cxx:81
 AliHFEefficiency.cxx:82
 AliHFEefficiency.cxx:83
 AliHFEefficiency.cxx:84
 AliHFEefficiency.cxx:85
 AliHFEefficiency.cxx:86
 AliHFEefficiency.cxx:87
 AliHFEefficiency.cxx:88
 AliHFEefficiency.cxx:89
 AliHFEefficiency.cxx:90
 AliHFEefficiency.cxx:91
 AliHFEefficiency.cxx:92
 AliHFEefficiency.cxx:93
 AliHFEefficiency.cxx:94
 AliHFEefficiency.cxx:95
 AliHFEefficiency.cxx:96
 AliHFEefficiency.cxx:97
 AliHFEefficiency.cxx:98
 AliHFEefficiency.cxx:99
 AliHFEefficiency.cxx:100
 AliHFEefficiency.cxx:101
 AliHFEefficiency.cxx:102
 AliHFEefficiency.cxx:103
 AliHFEefficiency.cxx:104
 AliHFEefficiency.cxx:105
 AliHFEefficiency.cxx:106
 AliHFEefficiency.cxx:107
 AliHFEefficiency.cxx:108
 AliHFEefficiency.cxx:109
 AliHFEefficiency.cxx:110
 AliHFEefficiency.cxx:111
 AliHFEefficiency.cxx:112
 AliHFEefficiency.cxx:113
 AliHFEefficiency.cxx:114
 AliHFEefficiency.cxx:115
 AliHFEefficiency.cxx:116
 AliHFEefficiency.cxx:117
 AliHFEefficiency.cxx:118
 AliHFEefficiency.cxx:119
 AliHFEefficiency.cxx:120
 AliHFEefficiency.cxx:121
 AliHFEefficiency.cxx:122
 AliHFEefficiency.cxx:123
 AliHFEefficiency.cxx:124
 AliHFEefficiency.cxx:125
 AliHFEefficiency.cxx:126
 AliHFEefficiency.cxx:127
 AliHFEefficiency.cxx:128
 AliHFEefficiency.cxx:129
 AliHFEefficiency.cxx:130
 AliHFEefficiency.cxx:131
 AliHFEefficiency.cxx:132
 AliHFEefficiency.cxx:133
 AliHFEefficiency.cxx:134
 AliHFEefficiency.cxx:135
 AliHFEefficiency.cxx:136
 AliHFEefficiency.cxx:137
 AliHFEefficiency.cxx:138
 AliHFEefficiency.cxx:139
 AliHFEefficiency.cxx:140
 AliHFEefficiency.cxx:141
 AliHFEefficiency.cxx:142
 AliHFEefficiency.cxx:143
 AliHFEefficiency.cxx:144
 AliHFEefficiency.cxx:145
 AliHFEefficiency.cxx:146
 AliHFEefficiency.cxx:147
 AliHFEefficiency.cxx:148
 AliHFEefficiency.cxx:149
 AliHFEefficiency.cxx:150
 AliHFEefficiency.cxx:151
 AliHFEefficiency.cxx:152
 AliHFEefficiency.cxx:153
 AliHFEefficiency.cxx:154
 AliHFEefficiency.cxx:155
 AliHFEefficiency.cxx:156
 AliHFEefficiency.cxx:157
 AliHFEefficiency.cxx:158
 AliHFEefficiency.cxx:159
 AliHFEefficiency.cxx:160
 AliHFEefficiency.cxx:161
 AliHFEefficiency.cxx:162
 AliHFEefficiency.cxx:163
 AliHFEefficiency.cxx:164
 AliHFEefficiency.cxx:165
 AliHFEefficiency.cxx:166
 AliHFEefficiency.cxx:167
 AliHFEefficiency.cxx:168
 AliHFEefficiency.cxx:169
 AliHFEefficiency.cxx:170
 AliHFEefficiency.cxx:171
 AliHFEefficiency.cxx:172
 AliHFEefficiency.cxx:173
 AliHFEefficiency.cxx:174
 AliHFEefficiency.cxx:175
 AliHFEefficiency.cxx:176
 AliHFEefficiency.cxx:177
 AliHFEefficiency.cxx:178
 AliHFEefficiency.cxx:179
 AliHFEefficiency.cxx:180
 AliHFEefficiency.cxx:181
 AliHFEefficiency.cxx:182
 AliHFEefficiency.cxx:183
 AliHFEefficiency.cxx:184
 AliHFEefficiency.cxx:185
 AliHFEefficiency.cxx:186
 AliHFEefficiency.cxx:187
 AliHFEefficiency.cxx:188
 AliHFEefficiency.cxx:189
 AliHFEefficiency.cxx:190
 AliHFEefficiency.cxx:191
 AliHFEefficiency.cxx:192
 AliHFEefficiency.cxx:193
 AliHFEefficiency.cxx:194
 AliHFEefficiency.cxx:195
 AliHFEefficiency.cxx:196
 AliHFEefficiency.cxx:197
 AliHFEefficiency.cxx:198
 AliHFEefficiency.cxx:199
 AliHFEefficiency.cxx:200
 AliHFEefficiency.cxx:201
 AliHFEefficiency.cxx:202
 AliHFEefficiency.cxx:203
 AliHFEefficiency.cxx:204
 AliHFEefficiency.cxx:205
 AliHFEefficiency.cxx:206
 AliHFEefficiency.cxx:207
 AliHFEefficiency.cxx:208
 AliHFEefficiency.cxx:209
 AliHFEefficiency.cxx:210
 AliHFEefficiency.cxx:211
 AliHFEefficiency.cxx:212
 AliHFEefficiency.cxx:213
 AliHFEefficiency.cxx:214
 AliHFEefficiency.cxx:215
 AliHFEefficiency.cxx:216
 AliHFEefficiency.cxx:217
 AliHFEefficiency.cxx:218
 AliHFEefficiency.cxx:219
 AliHFEefficiency.cxx:220
 AliHFEefficiency.cxx:221
 AliHFEefficiency.cxx:222
 AliHFEefficiency.cxx:223
 AliHFEefficiency.cxx:224
 AliHFEefficiency.cxx:225
 AliHFEefficiency.cxx:226
 AliHFEefficiency.cxx:227
 AliHFEefficiency.cxx:228
 AliHFEefficiency.cxx:229
 AliHFEefficiency.cxx:230
 AliHFEefficiency.cxx:231
 AliHFEefficiency.cxx:232
 AliHFEefficiency.cxx:233
 AliHFEefficiency.cxx:234
 AliHFEefficiency.cxx:235
 AliHFEefficiency.cxx:236
 AliHFEefficiency.cxx:237
 AliHFEefficiency.cxx:238
 AliHFEefficiency.cxx:239
 AliHFEefficiency.cxx:240
 AliHFEefficiency.cxx:241
 AliHFEefficiency.cxx:242
 AliHFEefficiency.cxx:243
 AliHFEefficiency.cxx:244
 AliHFEefficiency.cxx:245
 AliHFEefficiency.cxx:246
 AliHFEefficiency.cxx:247
 AliHFEefficiency.cxx:248
 AliHFEefficiency.cxx:249
 AliHFEefficiency.cxx:250
 AliHFEefficiency.cxx:251
 AliHFEefficiency.cxx:252
 AliHFEefficiency.cxx:253
 AliHFEefficiency.cxx:254
 AliHFEefficiency.cxx:255
 AliHFEefficiency.cxx:256
 AliHFEefficiency.cxx:257
 AliHFEefficiency.cxx:258
 AliHFEefficiency.cxx:259
 AliHFEefficiency.cxx:260
 AliHFEefficiency.cxx:261
 AliHFEefficiency.cxx:262
 AliHFEefficiency.cxx:263
 AliHFEefficiency.cxx:264
 AliHFEefficiency.cxx:265
 AliHFEefficiency.cxx:266
 AliHFEefficiency.cxx:267
 AliHFEefficiency.cxx:268
 AliHFEefficiency.cxx:269
 AliHFEefficiency.cxx:270
 AliHFEefficiency.cxx:271
 AliHFEefficiency.cxx:272
 AliHFEefficiency.cxx:273
 AliHFEefficiency.cxx:274
 AliHFEefficiency.cxx:275
 AliHFEefficiency.cxx:276
 AliHFEefficiency.cxx:277
 AliHFEefficiency.cxx:278
 AliHFEefficiency.cxx:279
 AliHFEefficiency.cxx:280
 AliHFEefficiency.cxx:281
 AliHFEefficiency.cxx:282
 AliHFEefficiency.cxx:283
 AliHFEefficiency.cxx:284
 AliHFEefficiency.cxx:285
 AliHFEefficiency.cxx:286
 AliHFEefficiency.cxx:287
 AliHFEefficiency.cxx:288
 AliHFEefficiency.cxx:289
 AliHFEefficiency.cxx:290
 AliHFEefficiency.cxx:291
 AliHFEefficiency.cxx:292
 AliHFEefficiency.cxx:293
 AliHFEefficiency.cxx:294
 AliHFEefficiency.cxx:295
 AliHFEefficiency.cxx:296
 AliHFEefficiency.cxx:297
 AliHFEefficiency.cxx:298
 AliHFEefficiency.cxx:299
 AliHFEefficiency.cxx:300
 AliHFEefficiency.cxx:301
 AliHFEefficiency.cxx:302
 AliHFEefficiency.cxx:303
 AliHFEefficiency.cxx:304
 AliHFEefficiency.cxx:305
 AliHFEefficiency.cxx:306
 AliHFEefficiency.cxx:307
 AliHFEefficiency.cxx:308
 AliHFEefficiency.cxx:309
 AliHFEefficiency.cxx:310
 AliHFEefficiency.cxx:311
 AliHFEefficiency.cxx:312
 AliHFEefficiency.cxx:313
 AliHFEefficiency.cxx:314
 AliHFEefficiency.cxx:315
 AliHFEefficiency.cxx:316
 AliHFEefficiency.cxx:317
 AliHFEefficiency.cxx:318
 AliHFEefficiency.cxx:319
 AliHFEefficiency.cxx:320
 AliHFEefficiency.cxx:321
 AliHFEefficiency.cxx:322
 AliHFEefficiency.cxx:323
 AliHFEefficiency.cxx:324
 AliHFEefficiency.cxx:325
 AliHFEefficiency.cxx:326
 AliHFEefficiency.cxx:327
 AliHFEefficiency.cxx:328
 AliHFEefficiency.cxx:329
 AliHFEefficiency.cxx:330
 AliHFEefficiency.cxx:331
 AliHFEefficiency.cxx:332
 AliHFEefficiency.cxx:333
 AliHFEefficiency.cxx:334
 AliHFEefficiency.cxx:335
 AliHFEefficiency.cxx:336
 AliHFEefficiency.cxx:337
 AliHFEefficiency.cxx:338
 AliHFEefficiency.cxx:339
 AliHFEefficiency.cxx:340
 AliHFEefficiency.cxx:341
 AliHFEefficiency.cxx:342
 AliHFEefficiency.cxx:343
 AliHFEefficiency.cxx:344
 AliHFEefficiency.cxx:345
 AliHFEefficiency.cxx:346
 AliHFEefficiency.cxx:347
 AliHFEefficiency.cxx:348
 AliHFEefficiency.cxx:349
 AliHFEefficiency.cxx:350
 AliHFEefficiency.cxx:351
 AliHFEefficiency.cxx:352
 AliHFEefficiency.cxx:353
 AliHFEefficiency.cxx:354
 AliHFEefficiency.cxx:355
 AliHFEefficiency.cxx:356
 AliHFEefficiency.cxx:357
 AliHFEefficiency.cxx:358
 AliHFEefficiency.cxx:359
 AliHFEefficiency.cxx:360
 AliHFEefficiency.cxx:361
 AliHFEefficiency.cxx:362
 AliHFEefficiency.cxx:363
 AliHFEefficiency.cxx:364
 AliHFEefficiency.cxx:365
 AliHFEefficiency.cxx:366
 AliHFEefficiency.cxx:367
 AliHFEefficiency.cxx:368
 AliHFEefficiency.cxx:369
 AliHFEefficiency.cxx:370
 AliHFEefficiency.cxx:371
 AliHFEefficiency.cxx:372
 AliHFEefficiency.cxx:373
 AliHFEefficiency.cxx:374
 AliHFEefficiency.cxx:375
 AliHFEefficiency.cxx:376
 AliHFEefficiency.cxx:377
 AliHFEefficiency.cxx:378
 AliHFEefficiency.cxx:379
 AliHFEefficiency.cxx:380
 AliHFEefficiency.cxx:381
 AliHFEefficiency.cxx:382
 AliHFEefficiency.cxx:383
 AliHFEefficiency.cxx:384
 AliHFEefficiency.cxx:385
 AliHFEefficiency.cxx:386
 AliHFEefficiency.cxx:387
 AliHFEefficiency.cxx:388
 AliHFEefficiency.cxx:389
 AliHFEefficiency.cxx:390
 AliHFEefficiency.cxx:391
 AliHFEefficiency.cxx:392
 AliHFEefficiency.cxx:393
 AliHFEefficiency.cxx:394
 AliHFEefficiency.cxx:395
 AliHFEefficiency.cxx:396
 AliHFEefficiency.cxx:397
 AliHFEefficiency.cxx:398
 AliHFEefficiency.cxx:399
 AliHFEefficiency.cxx:400
 AliHFEefficiency.cxx:401
 AliHFEefficiency.cxx:402
 AliHFEefficiency.cxx:403
 AliHFEefficiency.cxx:404
 AliHFEefficiency.cxx:405
 AliHFEefficiency.cxx:406
 AliHFEefficiency.cxx:407
 AliHFEefficiency.cxx:408
 AliHFEefficiency.cxx:409
 AliHFEefficiency.cxx:410
 AliHFEefficiency.cxx:411
 AliHFEefficiency.cxx:412
 AliHFEefficiency.cxx:413
 AliHFEefficiency.cxx:414
 AliHFEefficiency.cxx:415
 AliHFEefficiency.cxx:416
 AliHFEefficiency.cxx:417
 AliHFEefficiency.cxx:418
 AliHFEefficiency.cxx:419
 AliHFEefficiency.cxx:420
 AliHFEefficiency.cxx:421