ROOT logo
// $Id$
//
// Hadronic correction task.
//
// Author: R.Reed, C.Loizides, S.Aiola

#include <TChain.h>
#include <TClonesArray.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TList.h>

#include "AliAnalysisManager.h"
#include "AliAODEvent.h"
#include "AliAODCaloCluster.h"
#include "AliESDCaloCluster.h"
#include "AliVTrack.h"
#include "AliPicoTrack.h"
#include "AliVEventHandler.h"
#include "AliEmcalParticle.h"
#include "AliEMCALGeometry.h"
#include "AliParticleContainer.h"
#include "AliClusterContainer.h"

#include "AliHadCorrTask.h"

ClassImp(AliHadCorrTask)

//________________________________________________________________________
AliHadCorrTask::AliHadCorrTask() : 
  AliAnalysisTaskEmcal("AliHadCorrTask", kFALSE),
  fOutCaloName(),
  fPhiMatch(0.05),
  fEtaMatch(0.025),
  fDoTrackClus(kTRUE),
  fHadCorr(0),
  fEexclCell(0),
  fDoExact(kFALSE),
  fEsdMode(kTRUE),
  fOutClusters(0),
  fHistMatchEtaPhiAll(0),
  fHistMatchEtaPhiAllTr(0),
  fHistMatchEtaPhiAllCl(0),
  fHistNclusvsCent(0),
  fHistNclusMatchvsCent(0),
  fHistEbefore(0),
  fHistEafter(0),
  fHistEoPCent(0),
  fHistNMatchCent(0),
  fHistNClusMatchCent(0)
{
  // Default constructor.

  for(Int_t i=0; i<8; i++) {
    fHistEsubPch[i]    = 0;
    fHistEsubPchRat[i] = 0;
    fHistEsubPchRatAll[i] = 0;

    if (i<4) {
      fHistMatchEvsP[i]    = 0;
      fHistMatchdRvsEP[i]  = 0;
      fHistNMatchEnergy[i] = 0;

      fHistEmbTrackMatchesOversub[i] = 0;
      fHistNonEmbTrackMatchesOversub[i] = 0;
      fHistOversubMCClusters[i] = 0;
      fHistOversubNonMCClusters[i] = 0;
      fHistOversub[i] = 0;

      for(Int_t j=0; j<4; j++)
	fHistNCellsEnergy[i][j] = 0;
    }
    
    for(Int_t j=0; j<9; j++) {
      for(Int_t k=0; k<2; k++) 
	fHistMatchEtaPhi[i][j][k] = 0;
    }
  } 
}

//________________________________________________________________________
AliHadCorrTask::AliHadCorrTask(const char *name, Bool_t histo) : 
  AliAnalysisTaskEmcal(name, histo),
  fOutCaloName("CaloClustersCorr"),
  fPhiMatch(0.05),
  fEtaMatch(0.025),
  fDoTrackClus(kTRUE),
  fHadCorr(0),
  fEexclCell(0),
  fDoExact(kFALSE),
  fEsdMode(kTRUE),
  fOutClusters(0),
  fHistMatchEtaPhiAll(0),
  fHistMatchEtaPhiAllTr(0),
  fHistMatchEtaPhiAllCl(0),
  fHistNclusvsCent(0),
  fHistNclusMatchvsCent(0),
  fHistEbefore(0),
  fHistEafter(0),
  fHistEoPCent(0),
  fHistNMatchCent(0),
  fHistNClusMatchCent(0)
{
  // Standard constructor.

  for(Int_t i=0; i<8; i++) {
    fHistEsubPch[i]    = 0;
    fHistEsubPchRat[i] = 0;
    fHistEsubPchRatAll[i] = 0;

    if (i<4) {
      fHistMatchEvsP[i]    = 0;
      fHistMatchdRvsEP[i]  = 0;
      fHistNMatchEnergy[i] = 0;

      fHistEmbTrackMatchesOversub[i] = 0;
      fHistNonEmbTrackMatchesOversub[i] = 0;
      fHistOversubMCClusters[i] = 0;
      fHistOversubNonMCClusters[i] = 0;
      fHistOversub[i] = 0;

      for(Int_t j=0; j<4; j++)
	fHistNCellsEnergy[i][j] = 0;
    }
    
    for(Int_t j=0; j<9; j++) {
      for(Int_t k=0; k<2; k++) 
	fHistMatchEtaPhi[i][j][k] = 0;
    }
  } 
  
  SetMakeGeneralHistograms(histo);

  fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
}

//________________________________________________________________________
AliHadCorrTask::~AliHadCorrTask()
{
  // Destructor
}

//________________________________________________________________________
UInt_t AliHadCorrTask::GetMomBin(Double_t p) const
{
  // Get momenum bin.

  UInt_t pbin=0;
  if (p<0.5) 
    pbin=0;
  else if (p>=0.5 && p<1.0) 
    pbin=1;
  else if (p>=1.0 && p<1.5) 
    pbin=2;
  else if (p>=1.5 && p<2.) 
    pbin=3;
  else if (p>=2. && p<3.) 
    pbin=4;
  else if (p>=3. && p<4.) 
    pbin=5;
  else if (p>=4. && p<5.) 
    pbin=6;
  else if (p>=5. && p<8.) 
    pbin=7;
  else 
    pbin=8;

  return pbin;
}

//________________________________________________________________________
Double_t AliHadCorrTask::GetEtaSigma(Int_t pbin) const
{
  // Get sigma in eta.

  Double_t EtaSigma[9]={0.0097,0.0075,0.0059,0.0055,0.0053,0.005,0.005,0.0045,0.0042};
  return 2.0*EtaSigma[pbin];
}

//________________________________________________________________________
Double_t AliHadCorrTask::GetPhiMean(Int_t pbin, Int_t centbin) const
{
  // Get phi mean.

  if (centbin==0) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0121,
			 0.0084,
			 0.0060,
			 0.0041,
			 0.0031,
			 0.0022,
			 0.001};
    return PhiMean[pbin];
  } else if (centbin==1) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0121,
			 0.0084,
			 0.0060,
			 0.0041,
			 0.0031,
			 0.0022,
			 0.001};
    return PhiMean[pbin];
  } else if (centbin==2) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0121,
			 0.0084,
			 0.0060,
			 0.0041,
			 0.0031,
			 0.0022,
			 0.001};
    return PhiMean[pbin];
  } else if (centbin==3) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0121,
			 0.0084,
			 0.0060,
			 0.0041,
			 0.0031,
			 0.0022,
			 0.001};  
    return PhiMean[pbin];
  } else if (centbin==4) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0127,
			 0.0089,
			 0.0068,
			 0.0049,
			 0.0038,
			 0.0028,
			 0.0018};
    return PhiMean[pbin]*(-1.);
  } else if (centbin==5) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0127,
			 0.0089,
			 0.0068,
			 0.0048,
			 0.0038,
			 0.0028,
			 0.0018};
    return PhiMean[pbin]*(-1.);
  } else if (centbin==6) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0127,
			 0.0089,
			 0.0068,
			 0.0045,
			 0.0035,
			 0.0028,
			 0.0018};
    return PhiMean[pbin]*(-1.);
  } else if (centbin==7) { 
    Double_t PhiMean[9]={0.036,
			 0.021,
			 0.0127,
			 0.0089,
			 0.0068,
			 0.0043,
			 0.0035,
			 0.0028,
			 0.0018};
    return PhiMean[pbin]*(-1.);
  }

  return 0;
}

//________________________________________________________________________
Double_t AliHadCorrTask::GetPhiSigma(Int_t pbin, Int_t centbin) const
{
  // Get phi sigma.

  if (centbin==0) { 
    Double_t PhiSigma[9]={0.0221,
			  0.0128,
			  0.0074,
			  0.0064,
			  0.0059,
			  0.0055,
			  0.0052,
			  0.0049,
			  0.0045};
    return 2.*PhiSigma[pbin];
  } else if (centbin==1) { 
    Double_t PhiSigma[9]={0.0217,
			  0.0120,
			  0.0076,
			  0.0066,
			  0.0062,
			  0.0058,
			  0.0054,
			  0.0054,
0.0045};
    return 2.*PhiSigma[pbin];
  } else if (centbin==2) { 
    Double_t PhiSigma[9]={0.0211,
			  0.0124,
			  0.0080,
			  0.0070,
			  0.0067,
			  0.0061,
			  0.0059,
			  0.0054,
			  0.0047};
    return 2.*PhiSigma[pbin];
  } else if (centbin==3) { 
    Double_t PhiSigma[9]={0.0215,
			  0.0124,
			  0.0082,
			  0.0073,
			  0.0069,
			  0.0064,
			  0.0060,
			  0.0055,
			  0.0047};  
    return 2.*PhiSigma[pbin];
  } else if (centbin==4) { 
    Double_t PhiSigma[9]={0.0199,
			  0.0108,
			  0.0072,
			  0.0071,
			  0.0060,
			  0.0055,
			  0.0052,
			  0.0049,
			  0.0045};
    return 2.*PhiSigma[pbin];
  } else if (centbin==5) { 
    Double_t PhiSigma[9]={0.0200,
			  0.0110,
			  0.0074,
			  0.0071,
			  0.0064,
			  0.0059,
			  0.0055,
			  0.0052,
			  0.0045};
    return 2.*PhiSigma[pbin];
  } else if (centbin==6) { 
    Double_t PhiSigma[9]={0.0202,
			  0.0113,
			  0.0077,
			  0.0071,
			  0.0069,
			  0.0064,
			  0.0060,
			  0.0055,
			  0.0050};
    return 2.*PhiSigma[pbin];
  } else if (centbin==7) { 
    Double_t PhiSigma[9]={0.0205,
			  0.0113,
			  0.0080,
			  0.0074,
			  0.0078,
			  0.0067,
			  0.0062,
			  0.0055,
			  0.0050};
    return 2.*PhiSigma[pbin];
  }

  return 0;
}

//________________________________________________________________________
void AliHadCorrTask::UserCreateOutputObjects()
{
  // Create my user objects.

  AliAnalysisTaskEmcal::UserCreateOutputObjects();

  if (!fCreateHisto) return;

  TString name;
  TString temp;

  const Int_t nCentChBins = fNcentBins * 2;

  fHistMatchEtaPhiAll = new TH2F("fHistMatchEtaPhiAll", "fHistMatchEtaPhiAll;#Delta#eta;#Delta#phi", fNbins, -0.1, 0.1, fNbins, -0.1, 0.1);
  fOutput->Add(fHistMatchEtaPhiAll);

  fHistMatchEtaPhiAllTr = new TH2F("fHistMatchEtaPhiAllTr", "fHistMatchEtaPhiAllTr;#Delta#eta;#Delta#phi", fNbins, -0.1, 0.1, fNbins, -0.1, 0.1);
  fOutput->Add(fHistMatchEtaPhiAllTr);

  fHistMatchEtaPhiAllCl = new TH2F("fHistMatchEtaPhiAllCl", "fHistMatchEtaPhiAllCl;#Delta#eta;#Delta#phi", fNbins, -0.1, 0.1, fNbins, -0.1, 0.1);
  fOutput->Add(fHistMatchEtaPhiAllCl);

  for(Int_t icent=0; icent<nCentChBins; ++icent) {
    for(Int_t ipt=0; ipt<9; ++ipt) {
      for(Int_t ieta=0; ieta<2; ++ieta) {
	name = Form("fHistMatchEtaPhi_%i_%i_%i",icent,ipt,ieta);
	fHistMatchEtaPhi[icent][ipt][ieta] = new TH2F(name, name, fNbins, -0.1, 0.1, fNbins, -0.1, 0.1);
	fHistMatchEtaPhi[icent][ipt][ieta]->SetXTitle("#Delta#eta");
	fHistMatchEtaPhi[icent][ipt][ieta]->SetYTitle("#Delta#phi");
	fOutput->Add(fHistMatchEtaPhi[icent][ipt][ieta]);
      }
    }

    name = Form("fHistEsubPch_%i",icent);
    temp = Form("%s (Nmatches==1)",name.Data());
    fHistEsubPch[icent]=new TH1F(name, temp, fNbins, fMinBinPt, fMaxBinPt);
    fHistEsubPch[icent]->SetXTitle("#sum p (GeV) weighted with E_{sub}");
    fOutput->Add(fHistEsubPch[icent]);
    
    name = Form("fHistEsubPchRat_%i",icent);
    temp = Form("%s (Nmatches==1)",name.Data());
    fHistEsubPchRat[icent]=new TH2F(name, temp, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.);
    fHistEsubPchRat[icent]->SetXTitle("#Sigma p (GeV)");
    fHistEsubPchRat[icent]->SetYTitle("E_{sub} / #sum p");
    fOutput->Add(fHistEsubPchRat[icent]);

    name = Form("fHistEsubPchRatAll_%i",icent);
    temp = Form("%s (all Nmatches)",name.Data());
    fHistEsubPchRatAll[icent]=new TH2F(name, temp, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.);
    fHistEsubPchRatAll[icent]->SetXTitle("#Sigma p (GeV)");
    fHistEsubPchRatAll[icent]->SetYTitle("E_{sub} / #sum p");
    fOutput->Add(fHistEsubPchRatAll[icent]);
    
    if (icent<fNcentBins) {
      for(Int_t itrk=0; itrk<4; ++itrk) {
	name = Form("fHistNCellsEnergy_%i_%i",icent,itrk);
	temp = Form("%s (Nmatches==%d);N_{cells};E_{clus} (GeV)",name.Data(),itrk);
	fHistNCellsEnergy[icent][itrk]  = new TH2F(name, temp, fNbins, fMinBinPt, fMaxBinPt, 101, -0.5, 100.5);
	fOutput->Add(fHistNCellsEnergy[icent][itrk]);
      }

      name = Form("fHistMatchEvsP_%i",icent);
      temp = Form("%s (all Nmatches)",name.Data());
      fHistMatchEvsP[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.);
      fHistMatchEvsP[icent]->SetXTitle("E_{clus} (GeV)");
      fHistMatchEvsP[icent]->SetYTitle("E_{clus} / #sum p");
      fOutput->Add(fHistMatchEvsP[icent]);

      name = Form("fHistMatchdRvsEP_%i",icent);
      temp = Form("%s (all Nmatches)",name.Data());
      fHistMatchdRvsEP[icent] = new TH2F(name, temp, fNbins, 0., 0.2, fNbins*2, 0., 10.);
      fHistMatchdRvsEP[icent]->SetXTitle("#Delta R between track and cluster");
      fHistMatchdRvsEP[icent]->SetYTitle("E_{clus} / p");
      fOutput->Add(fHistMatchdRvsEP[icent]);

      name = Form("fHistNMatchEnergy_%i",icent);
      fHistNMatchEnergy[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, 101, -0.5, 100.5);
      fHistNMatchEnergy[icent]->SetXTitle("E_{clus} (GeV)");
      fHistNMatchEnergy[icent]->SetYTitle("N_{matches}");
      fOutput->Add(fHistNMatchEnergy[icent]);
    }
  }

  fHistNclusvsCent      = new TH1F("Nclusvscent",      "NclusVsCent; Cent (%)",                     100, 0, 100);
  fHistNclusMatchvsCent = new TH1F("NclusMatchvscent", "NclusMatchVsCent (all Nmatches); Cent (%)", 100, 0, 100);
  fHistEbefore          = new TH1F("Ebefore",          "Ebefore; Cent (%); E_{clus} (GeV)",         100, 0, 100);
  fHistEafter           = new TH1F("Eafter",           "Eafter;  Cent (%); E_{clus} (GeV)",         100, 0, 100);
  fHistEoPCent          = new TH2F("EoPCent",          "EoPCent; Cent (%); E_{clus} / #sum p",      100, 0, 100, fNbins*2, 0, 10);
  fHistNMatchCent       = new TH2F("NMatchesCent",     "NMatchesCent; Cent (%); Nmatches",          100, 0, 100, 11, -0.5,  10.5);
  fHistNClusMatchCent   = new TH2F("NClusMatchesCent", "NClusMatchesCent; Cent (%); Nmatches",      100, 0, 100, 11, -0.5,  10.5);

  fOutput->Add(fHistNclusMatchvsCent);
  fOutput->Add(fHistNclusvsCent);
  fOutput->Add(fHistEbefore);
  fOutput->Add(fHistEafter);
  fOutput->Add(fHistEoPCent);
  fOutput->Add(fHistNMatchCent);
  fOutput->Add(fHistNClusMatchCent);

  if (fIsEmbedded) {
    for(Int_t icent=0; icent<fNcentBins; ++icent) {
      name = Form("fHistEmbTrackMatchesOversub_%d",icent);
      fHistEmbTrackMatchesOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
      fHistEmbTrackMatchesOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
      fHistEmbTrackMatchesOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
      fOutput->Add(fHistEmbTrackMatchesOversub[icent]);

      name = Form("fHistNonEmbTrackMatchesOversub_%d",icent);
      fHistNonEmbTrackMatchesOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
      fHistNonEmbTrackMatchesOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
      fHistNonEmbTrackMatchesOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
      fOutput->Add(fHistNonEmbTrackMatchesOversub[icent]);

      name = Form("fHistOversubMCClusters_%d",icent);
      fHistOversubMCClusters[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
      fHistOversubMCClusters[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
      fHistOversubMCClusters[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
      fOutput->Add(fHistOversubMCClusters[icent]);

      name = Form("fHistOversubNonMCClusters_%d",icent);
      fHistOversubNonMCClusters[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
      fHistOversubNonMCClusters[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
      fHistOversubNonMCClusters[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
      fOutput->Add(fHistOversubNonMCClusters[icent]);

      name = Form("fHistOversub_%d",icent);
      fHistOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
      fHistOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
      fHistOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
      fOutput->Add(fHistOversub[icent]);
    }
  }

  PostData(1, fOutput);
}

//________________________________________________________________________
void AliHadCorrTask::ExecOnce() 
{
  // Initialize the analysis.

  if (fParticleCollArray.GetEntriesFast()<2) {
    AliError(Form("Wrong number of particle collections (%d), required 2",fParticleCollArray.GetEntriesFast()));
    return;
  }

  for (Int_t i = 0; i < 2; i++) {
    AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.At(i));
    cont->SetClassName("AliEmcalParticle");
  }

  // Do base class initializations and if it fails -> bail out
  AliAnalysisTaskEmcal::ExecOnce();
  if (!fInitialized) return;

  if (dynamic_cast<AliAODEvent*>(InputEvent())) fEsdMode = kFALSE;

  if (fEsdMode) 
    fOutClusters = new TClonesArray("AliESDCaloCluster");
  else 
    fOutClusters = new TClonesArray("AliAODCaloCluster");

  fOutClusters->SetName(fOutCaloName);

  // post output in event if not yet present
  if (!(InputEvent()->FindListObject(fOutCaloName))) {
    InputEvent()->AddObject(fOutClusters);
  }
  else {
    fInitialized = kFALSE;
    AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fOutCaloName.Data()));
    return;
  }
}

//________________________________________________________________________
void AliHadCorrTask::DoTrackLoop() 
{
  // Loop over tracks to provide some QA
 
  AliParticleContainer *tracks = static_cast<AliParticleContainer*>(fParticleCollArray.At(0));
  AliParticleContainer *clusters = static_cast<AliParticleContainer*>(fParticleCollArray.At(1));

  AliEmcalParticle *emctrack = 0;

  tracks->ResetCurrentID();
  while ((emctrack = static_cast<AliEmcalParticle*>(tracks->GetNextAcceptParticle()))) {
    Int_t NmatchClus = 0;

    if (!emctrack->IsEMCAL()) continue;

    AliVTrack *track = emctrack->GetTrack();
    if (!track) continue;
    
    if (track->GetTrackEtaOnEMCal() < fGeom->GetArm1EtaMin() - fEtaMatch ||
	track->GetTrackEtaOnEMCal() > fGeom->GetArm1EtaMax() + fEtaMatch || 
	track->GetTrackPhiOnEMCal() < fGeom->GetArm1PhiMin() * TMath::DegToRad() - fPhiMatch || 
	track->GetTrackPhiOnEMCal() > fGeom->GetArm1PhiMax() * TMath::DegToRad() + fPhiMatch)
      continue;
    
    Int_t Nclus = emctrack->GetNumberOfMatchedObj();

    for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
      AliEmcalParticle *emccluster = 
        static_cast<AliEmcalParticle*>(clusters->GetAcceptParticle(emctrack->GetMatchedObjId(iClus)));
      if (!emccluster) continue;

      AliVCluster *cluster = emccluster->GetCluster();
      if (!cluster) continue;
      if (!cluster->IsEMCAL()) continue;

      Double_t etadiff = 999;
      Double_t phidiff = 999;
      AliPicoTrack::GetEtaPhiDiff(track, cluster, phidiff, etadiff);
      fHistMatchEtaPhiAllTr->Fill(etadiff,phidiff);

      if (TMath::Abs(phidiff) < fPhiMatch && TMath::Abs(etadiff) < fEtaMatch) NmatchClus++;
    }
    fHistNClusMatchCent->Fill(fCent, NmatchClus);
  }
}

//________________________________________________________________________
void AliHadCorrTask::DoMatchedTracksLoop(AliEmcalParticle *emccluster, 
                                         Double_t &totalTrkP, Int_t &Nmatches, Double_t &trkPMCfrac, Int_t &NMCmatches) 
{
  // Do the loop over matched tracks for cluster emccluster.

  AliParticleContainer *tracks = static_cast<AliParticleContainer*>(fParticleCollArray.At(0));
  AliVCluster *cluster = emccluster->GetCluster();
  Int_t iClus = emccluster->IdInCollection();

  // loop over matched tracks
  const Int_t Ntrks = emccluster->GetNumberOfMatchedObj();
  for (Int_t i = 0; i < Ntrks; ++i) {
    Int_t    iTrack = emccluster->GetMatchedObjId(i);
    
    AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(tracks->GetAcceptParticle(iTrack));
    if (!emctrack) continue;

    AliVTrack *track = emctrack->GetTrack();
    if (!track) continue;

    Double_t etadiff = 999;
    Double_t phidiff = 999;
    AliPicoTrack::GetEtaPhiDiff(track, cluster, phidiff, etadiff);
    if (fCreateHisto)
      fHistMatchEtaPhiAllCl->Fill(etadiff, phidiff);

    // check if track also points to cluster
    if (fDoTrackClus && (emctrack->GetMatchedObjId(0)) != iClus) continue;

    Double_t mom       = track->P();
    UInt_t   mombin    = GetMomBin(mom); 
    Int_t    centbinch = fCentBin;
    if (track->Charge()<0) 
      centbinch += fNcentBins;

    if (fCreateHisto) {
      Int_t etabin = 0;
      if(track->Eta() > 0) etabin=1;
      fHistMatchEtaPhi[centbinch][mombin][etabin]->Fill(etadiff, phidiff);
      fHistMatchEtaPhiAll->Fill(etadiff, phidiff);
    }
    
    Double_t etaCut   = 0.0;
    Double_t phiCutlo = 0.0;
    Double_t phiCuthi = 0.0;

    if (fPhiMatch > 0) {
      phiCutlo = -fPhiMatch;
      phiCuthi = +fPhiMatch;
    } else {
      phiCutlo = GetPhiMean(mombin, centbinch) - GetPhiSigma(mombin, fCentBin);
      phiCuthi = GetPhiMean(mombin, centbinch) + GetPhiSigma(mombin, fCentBin);
    }

    if (fEtaMatch > 0) {
      etaCut = fEtaMatch;
    } else {
      etaCut = GetEtaSigma(mombin);
    }

    if ((phidiff < phiCuthi && phidiff > phiCutlo) && TMath::Abs(etadiff) < etaCut) {
      if (track->GetLabel() > fMinMCLabel) {
	++NMCmatches;
	trkPMCfrac += mom;
      }
      ++Nmatches;
      totalTrkP += mom;

      if (fCreateHisto) {
        if (fHadCorr > 1) {
          Double_t dR         = emccluster->GetMatchedObjDistance(i);
          Double_t energyclus = cluster->E();
          fHistMatchdRvsEP[fCentBin]->Fill(dR, energyclus / mom);
        }
      }
    }
  }

  if (totalTrkP > 0)
    trkPMCfrac /= totalTrkP;
}

//________________________________________________________________________
Bool_t AliHadCorrTask::Run() 
{
  // Run the hadronic correction

  AliParticleContainer *clusters = static_cast<AliParticleContainer*>(fParticleCollArray.At(1));
  AliEmcalParticle *emccluster = 0;
  
  // provide some additional histograms
  if (fCreateHisto)
    DoTrackLoop();

  // delete output
  fOutClusters->Delete();

  Int_t clusCount = 0;

   // loop over all clusters
  clusters->ResetCurrentID();
  while ((emccluster = static_cast<AliEmcalParticle*>(clusters->GetNextAcceptParticle()))) {
    if (!emccluster->IsEMCAL()) continue;

    AliVCluster *cluster = emccluster->GetCluster();
    if (!cluster) continue;

    Double_t energyclus = 0;
    if (fCreateHisto) {
      fHistEbefore->Fill(fCent, cluster->E());
      fHistNclusvsCent->Fill(fCent);
    }
  
    // apply correction / subtraction
    // to subtract only the closest track set fHadCor to a %
    // to subtract all tracks within the cut set fHadCor to %+1
    if (fHadCorr > 1)
      energyclus = ApplyHadCorrAllTracks(emccluster, fHadCorr - 1);	
    else if (fHadCorr > 0)
      energyclus = ApplyHadCorrOneTrack(emccluster, fHadCorr);	
    else 
      energyclus = cluster->E();

    if (energyclus < 0) 
      energyclus = 0;

    if (energyclus > 0) { // create corrected cluster
      AliVCluster *oc;
      if (fEsdMode) {
        AliESDCaloCluster *ec = dynamic_cast<AliESDCaloCluster*>(cluster);
        if (!ec) continue;
	oc = new ((*fOutClusters)[clusCount]) AliESDCaloCluster(*ec);
      } else { 
        AliAODCaloCluster *ac = dynamic_cast<AliAODCaloCluster*>(cluster);
        if (!ac) continue;
	oc = new ((*fOutClusters)[clusCount]) AliAODCaloCluster(*ac);
      }
      oc->SetE(energyclus);

      ++clusCount;

      if (fCreateHisto) fHistEafter->Fill(fCent, energyclus);
    }
  }
  
  return kTRUE;
}

//________________________________________________________________________
Double_t AliHadCorrTask::ApplyHadCorrOneTrack(AliEmcalParticle *emccluster, Double_t hadCorr) 
{
  // Apply the hadronic correction with one track only.

  AliVCluster *cluster = emccluster->GetCluster();
  if (!cluster) return 0;

  Double_t energyclus  = cluster->E();
  Int_t    iMin        = emccluster->GetMatchedObjId();
  if (iMin < 0)
    return energyclus;

  AliParticleContainer *tracks = static_cast<AliParticleContainer*>(fParticleCollArray.At(0));
  AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(tracks->GetParticle(iMin));
  if (!emctrack) return energyclus;

  AliVTrack *track = emctrack->GetTrack();
  if (!track) return energyclus;

  Double_t mom = track->P();
  if (mom < 1e-6) return energyclus;

  Double_t dEtaMin    = 1e9;
  Double_t dPhiMin    = 1e9;
  AliPicoTrack::GetEtaPhiDiff(track, cluster, dPhiMin, dEtaMin);
  if (fCreateHisto)
    fHistMatchEtaPhiAllCl->Fill(dEtaMin, dPhiMin);

  // check if track also points to cluster
  Int_t cid = emctrack->GetMatchedObjId();
  if (fDoTrackClus && (cid!=emccluster->IdInCollection())) return energyclus;

  UInt_t mombin = GetMomBin(mom);
  Int_t centbinch = fCentBin;
  if (track->Charge()<0) centbinch += fNcentBins;

  // plot some histograms if switched on
  if (fCreateHisto) {
    Int_t etabin = 0;
    if(track->Eta() > 0) 
      etabin = 1;
	    
    fHistMatchEtaPhi[centbinch][mombin][etabin]->Fill(dEtaMin, dPhiMin);
    fHistMatchEtaPhiAll->Fill(dEtaMin, dPhiMin);
    
    if (mom > 0) {
      Double_t dRmin      = emccluster->GetMatchedObjDistance();
      fHistMatchEvsP[fCentBin]->Fill(energyclus, energyclus / mom);
      fHistEoPCent->Fill(fCent, energyclus / mom);
      fHistMatchdRvsEP[fCentBin]->Fill(dRmin, energyclus / mom);
    }
  }
	   
  // define eta/phi cuts
  Double_t etaCut   = 0.0;
  Double_t phiCutlo = 0.0;
  Double_t phiCuthi = 0.0;
  if (fPhiMatch > 0) {
    phiCutlo = -fPhiMatch;
    phiCuthi = +fPhiMatch;
  }
  else {
    phiCutlo = GetPhiMean(mombin, centbinch) - GetPhiSigma(mombin, fCentBin);
    phiCuthi = GetPhiMean(mombin, centbinch) + GetPhiSigma(mombin, fCentBin);
  }
  if(fEtaMatch > 0) {
    etaCut = fEtaMatch;
  }
  else {
    etaCut = GetEtaSigma(mombin);
  }
  
  // apply the correction if the track is in the eta/phi window
  if ((dPhiMin < phiCuthi && dPhiMin > phiCutlo) && TMath::Abs(dEtaMin) < etaCut) {
    energyclus -= hadCorr * mom;
  }

  return energyclus;
}

//________________________________________________________________________
Double_t AliHadCorrTask::ApplyHadCorrAllTracks(AliEmcalParticle *emccluster, Double_t hadCorr) 
{
  // Apply the hadronic correction with all tracks.

  AliParticleContainer *tracks = static_cast<AliParticleContainer*>(fParticleCollArray.At(0));

  AliVCluster *cluster = emccluster->GetCluster();
  
  Double_t energyclus = cluster->E();
  Double_t cNcells = cluster->GetNCells();
  
  Double_t totalTrkP  = 0.0; // count total track momentum
  Int_t    Nmatches   = 0;   // count total number of matches

  Double_t trkPMCfrac   = 0.0; // count total track momentum
  Int_t    NMCmatches   = 0;   // count total number of matches
  
  // do the loop over the matched tracks and get the number of matches and the total momentum
  DoMatchedTracksLoop(emccluster, totalTrkP, Nmatches, trkPMCfrac, NMCmatches);

  Double_t Esub = hadCorr * totalTrkP;
	
  if (Esub > energyclus) Esub = energyclus;
	
  // applying Peter's proposed algorithm
  // never subtract the full energy of the cluster 
  Double_t clusEexcl = fEexclCell * cNcells;
  if (energyclus < clusEexcl) clusEexcl = energyclus;
  if ((energyclus - Esub) < clusEexcl) Esub = (energyclus - clusEexcl);

  // embedding
  Double_t EsubMC       = 0;
  Double_t EsubBkg      = 0;
  Double_t EclusMC      = 0;
  Double_t EclusBkg     = 0;
  Double_t EclusCorr    = 0;
  Double_t EclusMCcorr  = 0;
  Double_t EclusBkgcorr = 0;
  Double_t overSub = 0;
  if (fIsEmbedded) {
    EsubMC   = hadCorr * totalTrkP * trkPMCfrac;
    EsubBkg  = hadCorr * totalTrkP - EsubMC;
    EclusMC  = energyclus * cluster->GetMCEnergyFraction();
    EclusBkg = energyclus - EclusMC;
 
    if (energyclus > Esub)
      EclusCorr = energyclus - Esub;

    if (EclusMC > EsubMC)
      EclusMCcorr = EclusMC - EsubMC;
  
    if (EclusBkg > EsubBkg)
      EclusBkgcorr = EclusBkg - EsubBkg;
  
    overSub = EclusMCcorr + EclusBkgcorr - EclusCorr;
  }

  // plot some histograms if switched on
  if (fCreateHisto) {
    fHistNMatchCent->Fill(fCent, Nmatches);
    fHistNMatchEnergy[fCentBin]->Fill(energyclus, Nmatches);
    
    if(Nmatches > 0) 
      fHistNclusMatchvsCent->Fill(fCent);

    if(Nmatches < 3)  
      fHistNCellsEnergy[fCentBin][Nmatches]->Fill(energyclus, cNcells);
    else
      fHistNCellsEnergy[fCentBin][3]->Fill(energyclus, cNcells);
      
    if (totalTrkP>0) {
      Double_t EoP = energyclus / totalTrkP;
      fHistEoPCent->Fill(fCent, EoP);
      fHistMatchEvsP[fCentBin]->Fill(energyclus, EoP);

      fHistEsubPchRatAll[fCentBin]->Fill(totalTrkP, Esub / totalTrkP);

      if (Nmatches == 1) {
	Int_t iMin = emccluster->GetMatchedObjId();
	AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(tracks->GetParticle(iMin));
	if (emctrack) {
	  AliVTrack *track = emctrack->GetTrack();
	  if (track) {
	    Int_t centbinchm = fCentBin;
	    if (track->Charge()<0) centbinchm += fNcentBins;
	    fHistEsubPchRat[centbinchm]->Fill(totalTrkP, Esub / totalTrkP);
	    fHistEsubPch[centbinchm]->Fill(totalTrkP, Esub);
	  }
	}
      }

      if (fIsEmbedded) {
	fHistOversub[fCentBin]->Fill(energyclus, overSub / energyclus);

	if (cluster->GetMCEnergyFraction() > 0.95)
	  fHistOversubMCClusters[fCentBin]->Fill(energyclus, overSub / energyclus);
	else if (cluster->GetMCEnergyFraction() < 0.05)
	  fHistOversubNonMCClusters[fCentBin]->Fill(energyclus, overSub / energyclus);

	if (trkPMCfrac < 0.05)
	  fHistNonEmbTrackMatchesOversub[fCentBin]->Fill(energyclus, overSub / energyclus);
	else if (trkPMCfrac > 0.95)
	  fHistEmbTrackMatchesOversub[fCentBin]->Fill(energyclus, overSub / energyclus);
      }
    }
  }

  if (fIsEmbedded && fDoExact) {
    Esub -= overSub;
    if (EclusBkgcorr + EclusMCcorr > 0) {
      Double_t newfrac = EclusMCcorr / (EclusBkgcorr + EclusMCcorr);
      cluster->SetMCEnergyFraction(newfrac);
    }
  }

  // apply the correction
  energyclus -= Esub;

  return energyclus;
}
 AliHadCorrTask.cxx:1
 AliHadCorrTask.cxx:2
 AliHadCorrTask.cxx:3
 AliHadCorrTask.cxx:4
 AliHadCorrTask.cxx:5
 AliHadCorrTask.cxx:6
 AliHadCorrTask.cxx:7
 AliHadCorrTask.cxx:8
 AliHadCorrTask.cxx:9
 AliHadCorrTask.cxx:10
 AliHadCorrTask.cxx:11
 AliHadCorrTask.cxx:12
 AliHadCorrTask.cxx:13
 AliHadCorrTask.cxx:14
 AliHadCorrTask.cxx:15
 AliHadCorrTask.cxx:16
 AliHadCorrTask.cxx:17
 AliHadCorrTask.cxx:18
 AliHadCorrTask.cxx:19
 AliHadCorrTask.cxx:20
 AliHadCorrTask.cxx:21
 AliHadCorrTask.cxx:22
 AliHadCorrTask.cxx:23
 AliHadCorrTask.cxx:24
 AliHadCorrTask.cxx:25
 AliHadCorrTask.cxx:26
 AliHadCorrTask.cxx:27
 AliHadCorrTask.cxx:28
 AliHadCorrTask.cxx:29
 AliHadCorrTask.cxx:30
 AliHadCorrTask.cxx:31
 AliHadCorrTask.cxx:32
 AliHadCorrTask.cxx:33
 AliHadCorrTask.cxx:34
 AliHadCorrTask.cxx:35
 AliHadCorrTask.cxx:36
 AliHadCorrTask.cxx:37
 AliHadCorrTask.cxx:38
 AliHadCorrTask.cxx:39
 AliHadCorrTask.cxx:40
 AliHadCorrTask.cxx:41
 AliHadCorrTask.cxx:42
 AliHadCorrTask.cxx:43
 AliHadCorrTask.cxx:44
 AliHadCorrTask.cxx:45
 AliHadCorrTask.cxx:46
 AliHadCorrTask.cxx:47
 AliHadCorrTask.cxx:48
 AliHadCorrTask.cxx:49
 AliHadCorrTask.cxx:50
 AliHadCorrTask.cxx:51
 AliHadCorrTask.cxx:52
 AliHadCorrTask.cxx:53
 AliHadCorrTask.cxx:54
 AliHadCorrTask.cxx:55
 AliHadCorrTask.cxx:56
 AliHadCorrTask.cxx:57
 AliHadCorrTask.cxx:58
 AliHadCorrTask.cxx:59
 AliHadCorrTask.cxx:60
 AliHadCorrTask.cxx:61
 AliHadCorrTask.cxx:62
 AliHadCorrTask.cxx:63
 AliHadCorrTask.cxx:64
 AliHadCorrTask.cxx:65
 AliHadCorrTask.cxx:66
 AliHadCorrTask.cxx:67
 AliHadCorrTask.cxx:68
 AliHadCorrTask.cxx:69
 AliHadCorrTask.cxx:70
 AliHadCorrTask.cxx:71
 AliHadCorrTask.cxx:72
 AliHadCorrTask.cxx:73
 AliHadCorrTask.cxx:74
 AliHadCorrTask.cxx:75
 AliHadCorrTask.cxx:76
 AliHadCorrTask.cxx:77
 AliHadCorrTask.cxx:78
 AliHadCorrTask.cxx:79
 AliHadCorrTask.cxx:80
 AliHadCorrTask.cxx:81
 AliHadCorrTask.cxx:82
 AliHadCorrTask.cxx:83
 AliHadCorrTask.cxx:84
 AliHadCorrTask.cxx:85
 AliHadCorrTask.cxx:86
 AliHadCorrTask.cxx:87
 AliHadCorrTask.cxx:88
 AliHadCorrTask.cxx:89
 AliHadCorrTask.cxx:90
 AliHadCorrTask.cxx:91
 AliHadCorrTask.cxx:92
 AliHadCorrTask.cxx:93
 AliHadCorrTask.cxx:94
 AliHadCorrTask.cxx:95
 AliHadCorrTask.cxx:96
 AliHadCorrTask.cxx:97
 AliHadCorrTask.cxx:98
 AliHadCorrTask.cxx:99
 AliHadCorrTask.cxx:100
 AliHadCorrTask.cxx:101
 AliHadCorrTask.cxx:102
 AliHadCorrTask.cxx:103
 AliHadCorrTask.cxx:104
 AliHadCorrTask.cxx:105
 AliHadCorrTask.cxx:106
 AliHadCorrTask.cxx:107
 AliHadCorrTask.cxx:108
 AliHadCorrTask.cxx:109
 AliHadCorrTask.cxx:110
 AliHadCorrTask.cxx:111
 AliHadCorrTask.cxx:112
 AliHadCorrTask.cxx:113
 AliHadCorrTask.cxx:114
 AliHadCorrTask.cxx:115
 AliHadCorrTask.cxx:116
 AliHadCorrTask.cxx:117
 AliHadCorrTask.cxx:118
 AliHadCorrTask.cxx:119
 AliHadCorrTask.cxx:120
 AliHadCorrTask.cxx:121
 AliHadCorrTask.cxx:122
 AliHadCorrTask.cxx:123
 AliHadCorrTask.cxx:124
 AliHadCorrTask.cxx:125
 AliHadCorrTask.cxx:126
 AliHadCorrTask.cxx:127
 AliHadCorrTask.cxx:128
 AliHadCorrTask.cxx:129
 AliHadCorrTask.cxx:130
 AliHadCorrTask.cxx:131
 AliHadCorrTask.cxx:132
 AliHadCorrTask.cxx:133
 AliHadCorrTask.cxx:134
 AliHadCorrTask.cxx:135
 AliHadCorrTask.cxx:136
 AliHadCorrTask.cxx:137
 AliHadCorrTask.cxx:138
 AliHadCorrTask.cxx:139
 AliHadCorrTask.cxx:140
 AliHadCorrTask.cxx:141
 AliHadCorrTask.cxx:142
 AliHadCorrTask.cxx:143
 AliHadCorrTask.cxx:144
 AliHadCorrTask.cxx:145
 AliHadCorrTask.cxx:146
 AliHadCorrTask.cxx:147
 AliHadCorrTask.cxx:148
 AliHadCorrTask.cxx:149
 AliHadCorrTask.cxx:150
 AliHadCorrTask.cxx:151
 AliHadCorrTask.cxx:152
 AliHadCorrTask.cxx:153
 AliHadCorrTask.cxx:154
 AliHadCorrTask.cxx:155
 AliHadCorrTask.cxx:156
 AliHadCorrTask.cxx:157
 AliHadCorrTask.cxx:158
 AliHadCorrTask.cxx:159
 AliHadCorrTask.cxx:160
 AliHadCorrTask.cxx:161
 AliHadCorrTask.cxx:162
 AliHadCorrTask.cxx:163
 AliHadCorrTask.cxx:164
 AliHadCorrTask.cxx:165
 AliHadCorrTask.cxx:166
 AliHadCorrTask.cxx:167
 AliHadCorrTask.cxx:168
 AliHadCorrTask.cxx:169
 AliHadCorrTask.cxx:170
 AliHadCorrTask.cxx:171
 AliHadCorrTask.cxx:172
 AliHadCorrTask.cxx:173
 AliHadCorrTask.cxx:174
 AliHadCorrTask.cxx:175
 AliHadCorrTask.cxx:176
 AliHadCorrTask.cxx:177
 AliHadCorrTask.cxx:178
 AliHadCorrTask.cxx:179
 AliHadCorrTask.cxx:180
 AliHadCorrTask.cxx:181
 AliHadCorrTask.cxx:182
 AliHadCorrTask.cxx:183
 AliHadCorrTask.cxx:184
 AliHadCorrTask.cxx:185
 AliHadCorrTask.cxx:186
 AliHadCorrTask.cxx:187
 AliHadCorrTask.cxx:188
 AliHadCorrTask.cxx:189
 AliHadCorrTask.cxx:190
 AliHadCorrTask.cxx:191
 AliHadCorrTask.cxx:192
 AliHadCorrTask.cxx:193
 AliHadCorrTask.cxx:194
 AliHadCorrTask.cxx:195
 AliHadCorrTask.cxx:196
 AliHadCorrTask.cxx:197
 AliHadCorrTask.cxx:198
 AliHadCorrTask.cxx:199
 AliHadCorrTask.cxx:200
 AliHadCorrTask.cxx:201
 AliHadCorrTask.cxx:202
 AliHadCorrTask.cxx:203
 AliHadCorrTask.cxx:204
 AliHadCorrTask.cxx:205
 AliHadCorrTask.cxx:206
 AliHadCorrTask.cxx:207
 AliHadCorrTask.cxx:208
 AliHadCorrTask.cxx:209
 AliHadCorrTask.cxx:210
 AliHadCorrTask.cxx:211
 AliHadCorrTask.cxx:212
 AliHadCorrTask.cxx:213
 AliHadCorrTask.cxx:214
 AliHadCorrTask.cxx:215
 AliHadCorrTask.cxx:216
 AliHadCorrTask.cxx:217
 AliHadCorrTask.cxx:218
 AliHadCorrTask.cxx:219
 AliHadCorrTask.cxx:220
 AliHadCorrTask.cxx:221
 AliHadCorrTask.cxx:222
 AliHadCorrTask.cxx:223
 AliHadCorrTask.cxx:224
 AliHadCorrTask.cxx:225
 AliHadCorrTask.cxx:226
 AliHadCorrTask.cxx:227
 AliHadCorrTask.cxx:228
 AliHadCorrTask.cxx:229
 AliHadCorrTask.cxx:230
 AliHadCorrTask.cxx:231
 AliHadCorrTask.cxx:232
 AliHadCorrTask.cxx:233
 AliHadCorrTask.cxx:234
 AliHadCorrTask.cxx:235
 AliHadCorrTask.cxx:236
 AliHadCorrTask.cxx:237
 AliHadCorrTask.cxx:238
 AliHadCorrTask.cxx:239
 AliHadCorrTask.cxx:240
 AliHadCorrTask.cxx:241
 AliHadCorrTask.cxx:242
 AliHadCorrTask.cxx:243
 AliHadCorrTask.cxx:244
 AliHadCorrTask.cxx:245
 AliHadCorrTask.cxx:246
 AliHadCorrTask.cxx:247
 AliHadCorrTask.cxx:248
 AliHadCorrTask.cxx:249
 AliHadCorrTask.cxx:250
 AliHadCorrTask.cxx:251
 AliHadCorrTask.cxx:252
 AliHadCorrTask.cxx:253
 AliHadCorrTask.cxx:254
 AliHadCorrTask.cxx:255
 AliHadCorrTask.cxx:256
 AliHadCorrTask.cxx:257
 AliHadCorrTask.cxx:258
 AliHadCorrTask.cxx:259
 AliHadCorrTask.cxx:260
 AliHadCorrTask.cxx:261
 AliHadCorrTask.cxx:262
 AliHadCorrTask.cxx:263
 AliHadCorrTask.cxx:264
 AliHadCorrTask.cxx:265
 AliHadCorrTask.cxx:266
 AliHadCorrTask.cxx:267
 AliHadCorrTask.cxx:268
 AliHadCorrTask.cxx:269
 AliHadCorrTask.cxx:270
 AliHadCorrTask.cxx:271
 AliHadCorrTask.cxx:272
 AliHadCorrTask.cxx:273
 AliHadCorrTask.cxx:274
 AliHadCorrTask.cxx:275
 AliHadCorrTask.cxx:276
 AliHadCorrTask.cxx:277
 AliHadCorrTask.cxx:278
 AliHadCorrTask.cxx:279
 AliHadCorrTask.cxx:280
 AliHadCorrTask.cxx:281
 AliHadCorrTask.cxx:282
 AliHadCorrTask.cxx:283
 AliHadCorrTask.cxx:284
 AliHadCorrTask.cxx:285
 AliHadCorrTask.cxx:286
 AliHadCorrTask.cxx:287
 AliHadCorrTask.cxx:288
 AliHadCorrTask.cxx:289
 AliHadCorrTask.cxx:290
 AliHadCorrTask.cxx:291
 AliHadCorrTask.cxx:292
 AliHadCorrTask.cxx:293
 AliHadCorrTask.cxx:294
 AliHadCorrTask.cxx:295
 AliHadCorrTask.cxx:296
 AliHadCorrTask.cxx:297
 AliHadCorrTask.cxx:298
 AliHadCorrTask.cxx:299
 AliHadCorrTask.cxx:300
 AliHadCorrTask.cxx:301
 AliHadCorrTask.cxx:302
 AliHadCorrTask.cxx:303
 AliHadCorrTask.cxx:304
 AliHadCorrTask.cxx:305
 AliHadCorrTask.cxx:306
 AliHadCorrTask.cxx:307
 AliHadCorrTask.cxx:308
 AliHadCorrTask.cxx:309
 AliHadCorrTask.cxx:310
 AliHadCorrTask.cxx:311
 AliHadCorrTask.cxx:312
 AliHadCorrTask.cxx:313
 AliHadCorrTask.cxx:314
 AliHadCorrTask.cxx:315
 AliHadCorrTask.cxx:316
 AliHadCorrTask.cxx:317
 AliHadCorrTask.cxx:318
 AliHadCorrTask.cxx:319
 AliHadCorrTask.cxx:320
 AliHadCorrTask.cxx:321
 AliHadCorrTask.cxx:322
 AliHadCorrTask.cxx:323
 AliHadCorrTask.cxx:324
 AliHadCorrTask.cxx:325
 AliHadCorrTask.cxx:326
 AliHadCorrTask.cxx:327
 AliHadCorrTask.cxx:328
 AliHadCorrTask.cxx:329
 AliHadCorrTask.cxx:330
 AliHadCorrTask.cxx:331
 AliHadCorrTask.cxx:332
 AliHadCorrTask.cxx:333
 AliHadCorrTask.cxx:334
 AliHadCorrTask.cxx:335
 AliHadCorrTask.cxx:336
 AliHadCorrTask.cxx:337
 AliHadCorrTask.cxx:338
 AliHadCorrTask.cxx:339
 AliHadCorrTask.cxx:340
 AliHadCorrTask.cxx:341
 AliHadCorrTask.cxx:342
 AliHadCorrTask.cxx:343
 AliHadCorrTask.cxx:344
 AliHadCorrTask.cxx:345
 AliHadCorrTask.cxx:346
 AliHadCorrTask.cxx:347
 AliHadCorrTask.cxx:348
 AliHadCorrTask.cxx:349
 AliHadCorrTask.cxx:350
 AliHadCorrTask.cxx:351
 AliHadCorrTask.cxx:352
 AliHadCorrTask.cxx:353
 AliHadCorrTask.cxx:354
 AliHadCorrTask.cxx:355
 AliHadCorrTask.cxx:356
 AliHadCorrTask.cxx:357
 AliHadCorrTask.cxx:358
 AliHadCorrTask.cxx:359
 AliHadCorrTask.cxx:360
 AliHadCorrTask.cxx:361
 AliHadCorrTask.cxx:362
 AliHadCorrTask.cxx:363
 AliHadCorrTask.cxx:364
 AliHadCorrTask.cxx:365
 AliHadCorrTask.cxx:366
 AliHadCorrTask.cxx:367
 AliHadCorrTask.cxx:368
 AliHadCorrTask.cxx:369
 AliHadCorrTask.cxx:370
 AliHadCorrTask.cxx:371
 AliHadCorrTask.cxx:372
 AliHadCorrTask.cxx:373
 AliHadCorrTask.cxx:374
 AliHadCorrTask.cxx:375
 AliHadCorrTask.cxx:376
 AliHadCorrTask.cxx:377
 AliHadCorrTask.cxx:378
 AliHadCorrTask.cxx:379
 AliHadCorrTask.cxx:380
 AliHadCorrTask.cxx:381
 AliHadCorrTask.cxx:382
 AliHadCorrTask.cxx:383
 AliHadCorrTask.cxx:384
 AliHadCorrTask.cxx:385
 AliHadCorrTask.cxx:386
 AliHadCorrTask.cxx:387
 AliHadCorrTask.cxx:388
 AliHadCorrTask.cxx:389
 AliHadCorrTask.cxx:390
 AliHadCorrTask.cxx:391
 AliHadCorrTask.cxx:392
 AliHadCorrTask.cxx:393
 AliHadCorrTask.cxx:394
 AliHadCorrTask.cxx:395
 AliHadCorrTask.cxx:396
 AliHadCorrTask.cxx:397
 AliHadCorrTask.cxx:398
 AliHadCorrTask.cxx:399
 AliHadCorrTask.cxx:400
 AliHadCorrTask.cxx:401
 AliHadCorrTask.cxx:402
 AliHadCorrTask.cxx:403
 AliHadCorrTask.cxx:404
 AliHadCorrTask.cxx:405
 AliHadCorrTask.cxx:406
 AliHadCorrTask.cxx:407
 AliHadCorrTask.cxx:408
 AliHadCorrTask.cxx:409
 AliHadCorrTask.cxx:410
 AliHadCorrTask.cxx:411
 AliHadCorrTask.cxx:412
 AliHadCorrTask.cxx:413
 AliHadCorrTask.cxx:414
 AliHadCorrTask.cxx:415
 AliHadCorrTask.cxx:416
 AliHadCorrTask.cxx:417
 AliHadCorrTask.cxx:418
 AliHadCorrTask.cxx:419
 AliHadCorrTask.cxx:420
 AliHadCorrTask.cxx:421
 AliHadCorrTask.cxx:422
 AliHadCorrTask.cxx:423
 AliHadCorrTask.cxx:424
 AliHadCorrTask.cxx:425
 AliHadCorrTask.cxx:426
 AliHadCorrTask.cxx:427
 AliHadCorrTask.cxx:428
 AliHadCorrTask.cxx:429
 AliHadCorrTask.cxx:430
 AliHadCorrTask.cxx:431
 AliHadCorrTask.cxx:432
 AliHadCorrTask.cxx:433
 AliHadCorrTask.cxx:434
 AliHadCorrTask.cxx:435
 AliHadCorrTask.cxx:436
 AliHadCorrTask.cxx:437
 AliHadCorrTask.cxx:438
 AliHadCorrTask.cxx:439
 AliHadCorrTask.cxx:440
 AliHadCorrTask.cxx:441
 AliHadCorrTask.cxx:442
 AliHadCorrTask.cxx:443
 AliHadCorrTask.cxx:444
 AliHadCorrTask.cxx:445
 AliHadCorrTask.cxx:446
 AliHadCorrTask.cxx:447
 AliHadCorrTask.cxx:448
 AliHadCorrTask.cxx:449
 AliHadCorrTask.cxx:450
 AliHadCorrTask.cxx:451
 AliHadCorrTask.cxx:452
 AliHadCorrTask.cxx:453
 AliHadCorrTask.cxx:454
 AliHadCorrTask.cxx:455
 AliHadCorrTask.cxx:456
 AliHadCorrTask.cxx:457
 AliHadCorrTask.cxx:458
 AliHadCorrTask.cxx:459
 AliHadCorrTask.cxx:460
 AliHadCorrTask.cxx:461
 AliHadCorrTask.cxx:462
 AliHadCorrTask.cxx:463
 AliHadCorrTask.cxx:464
 AliHadCorrTask.cxx:465
 AliHadCorrTask.cxx:466
 AliHadCorrTask.cxx:467
 AliHadCorrTask.cxx:468
 AliHadCorrTask.cxx:469
 AliHadCorrTask.cxx:470
 AliHadCorrTask.cxx:471
 AliHadCorrTask.cxx:472
 AliHadCorrTask.cxx:473
 AliHadCorrTask.cxx:474
 AliHadCorrTask.cxx:475
 AliHadCorrTask.cxx:476
 AliHadCorrTask.cxx:477
 AliHadCorrTask.cxx:478
 AliHadCorrTask.cxx:479
 AliHadCorrTask.cxx:480
 AliHadCorrTask.cxx:481
 AliHadCorrTask.cxx:482
 AliHadCorrTask.cxx:483
 AliHadCorrTask.cxx:484
 AliHadCorrTask.cxx:485
 AliHadCorrTask.cxx:486
 AliHadCorrTask.cxx:487
 AliHadCorrTask.cxx:488
 AliHadCorrTask.cxx:489
 AliHadCorrTask.cxx:490
 AliHadCorrTask.cxx:491
 AliHadCorrTask.cxx:492
 AliHadCorrTask.cxx:493
 AliHadCorrTask.cxx:494
 AliHadCorrTask.cxx:495
 AliHadCorrTask.cxx:496
 AliHadCorrTask.cxx:497
 AliHadCorrTask.cxx:498
 AliHadCorrTask.cxx:499
 AliHadCorrTask.cxx:500
 AliHadCorrTask.cxx:501
 AliHadCorrTask.cxx:502
 AliHadCorrTask.cxx:503
 AliHadCorrTask.cxx:504
 AliHadCorrTask.cxx:505
 AliHadCorrTask.cxx:506
 AliHadCorrTask.cxx:507
 AliHadCorrTask.cxx:508
 AliHadCorrTask.cxx:509
 AliHadCorrTask.cxx:510
 AliHadCorrTask.cxx:511
 AliHadCorrTask.cxx:512
 AliHadCorrTask.cxx:513
 AliHadCorrTask.cxx:514
 AliHadCorrTask.cxx:515
 AliHadCorrTask.cxx:516
 AliHadCorrTask.cxx:517
 AliHadCorrTask.cxx:518
 AliHadCorrTask.cxx:519
 AliHadCorrTask.cxx:520
 AliHadCorrTask.cxx:521
 AliHadCorrTask.cxx:522
 AliHadCorrTask.cxx:523
 AliHadCorrTask.cxx:524
 AliHadCorrTask.cxx:525
 AliHadCorrTask.cxx:526
 AliHadCorrTask.cxx:527
 AliHadCorrTask.cxx:528
 AliHadCorrTask.cxx:529
 AliHadCorrTask.cxx:530
 AliHadCorrTask.cxx:531
 AliHadCorrTask.cxx:532
 AliHadCorrTask.cxx:533
 AliHadCorrTask.cxx:534
 AliHadCorrTask.cxx:535
 AliHadCorrTask.cxx:536
 AliHadCorrTask.cxx:537
 AliHadCorrTask.cxx:538
 AliHadCorrTask.cxx:539
 AliHadCorrTask.cxx:540
 AliHadCorrTask.cxx:541
 AliHadCorrTask.cxx:542
 AliHadCorrTask.cxx:543
 AliHadCorrTask.cxx:544
 AliHadCorrTask.cxx:545
 AliHadCorrTask.cxx:546
 AliHadCorrTask.cxx:547
 AliHadCorrTask.cxx:548
 AliHadCorrTask.cxx:549
 AliHadCorrTask.cxx:550
 AliHadCorrTask.cxx:551
 AliHadCorrTask.cxx:552
 AliHadCorrTask.cxx:553
 AliHadCorrTask.cxx:554
 AliHadCorrTask.cxx:555
 AliHadCorrTask.cxx:556
 AliHadCorrTask.cxx:557
 AliHadCorrTask.cxx:558
 AliHadCorrTask.cxx:559
 AliHadCorrTask.cxx:560
 AliHadCorrTask.cxx:561
 AliHadCorrTask.cxx:562
 AliHadCorrTask.cxx:563
 AliHadCorrTask.cxx:564
 AliHadCorrTask.cxx:565
 AliHadCorrTask.cxx:566
 AliHadCorrTask.cxx:567
 AliHadCorrTask.cxx:568
 AliHadCorrTask.cxx:569
 AliHadCorrTask.cxx:570
 AliHadCorrTask.cxx:571
 AliHadCorrTask.cxx:572
 AliHadCorrTask.cxx:573
 AliHadCorrTask.cxx:574
 AliHadCorrTask.cxx:575
 AliHadCorrTask.cxx:576
 AliHadCorrTask.cxx:577
 AliHadCorrTask.cxx:578
 AliHadCorrTask.cxx:579
 AliHadCorrTask.cxx:580
 AliHadCorrTask.cxx:581
 AliHadCorrTask.cxx:582
 AliHadCorrTask.cxx:583
 AliHadCorrTask.cxx:584
 AliHadCorrTask.cxx:585
 AliHadCorrTask.cxx:586
 AliHadCorrTask.cxx:587
 AliHadCorrTask.cxx:588
 AliHadCorrTask.cxx:589
 AliHadCorrTask.cxx:590
 AliHadCorrTask.cxx:591
 AliHadCorrTask.cxx:592
 AliHadCorrTask.cxx:593
 AliHadCorrTask.cxx:594
 AliHadCorrTask.cxx:595
 AliHadCorrTask.cxx:596
 AliHadCorrTask.cxx:597
 AliHadCorrTask.cxx:598
 AliHadCorrTask.cxx:599
 AliHadCorrTask.cxx:600
 AliHadCorrTask.cxx:601
 AliHadCorrTask.cxx:602
 AliHadCorrTask.cxx:603
 AliHadCorrTask.cxx:604
 AliHadCorrTask.cxx:605
 AliHadCorrTask.cxx:606
 AliHadCorrTask.cxx:607
 AliHadCorrTask.cxx:608
 AliHadCorrTask.cxx:609
 AliHadCorrTask.cxx:610
 AliHadCorrTask.cxx:611
 AliHadCorrTask.cxx:612
 AliHadCorrTask.cxx:613
 AliHadCorrTask.cxx:614
 AliHadCorrTask.cxx:615
 AliHadCorrTask.cxx:616
 AliHadCorrTask.cxx:617
 AliHadCorrTask.cxx:618
 AliHadCorrTask.cxx:619
 AliHadCorrTask.cxx:620
 AliHadCorrTask.cxx:621
 AliHadCorrTask.cxx:622
 AliHadCorrTask.cxx:623
 AliHadCorrTask.cxx:624
 AliHadCorrTask.cxx:625
 AliHadCorrTask.cxx:626
 AliHadCorrTask.cxx:627
 AliHadCorrTask.cxx:628
 AliHadCorrTask.cxx:629
 AliHadCorrTask.cxx:630
 AliHadCorrTask.cxx:631
 AliHadCorrTask.cxx:632
 AliHadCorrTask.cxx:633
 AliHadCorrTask.cxx:634
 AliHadCorrTask.cxx:635
 AliHadCorrTask.cxx:636
 AliHadCorrTask.cxx:637
 AliHadCorrTask.cxx:638
 AliHadCorrTask.cxx:639
 AliHadCorrTask.cxx:640
 AliHadCorrTask.cxx:641
 AliHadCorrTask.cxx:642
 AliHadCorrTask.cxx:643
 AliHadCorrTask.cxx:644
 AliHadCorrTask.cxx:645
 AliHadCorrTask.cxx:646
 AliHadCorrTask.cxx:647
 AliHadCorrTask.cxx:648
 AliHadCorrTask.cxx:649
 AliHadCorrTask.cxx:650
 AliHadCorrTask.cxx:651
 AliHadCorrTask.cxx:652
 AliHadCorrTask.cxx:653
 AliHadCorrTask.cxx:654
 AliHadCorrTask.cxx:655
 AliHadCorrTask.cxx:656
 AliHadCorrTask.cxx:657
 AliHadCorrTask.cxx:658
 AliHadCorrTask.cxx:659
 AliHadCorrTask.cxx:660
 AliHadCorrTask.cxx:661
 AliHadCorrTask.cxx:662
 AliHadCorrTask.cxx:663
 AliHadCorrTask.cxx:664
 AliHadCorrTask.cxx:665
 AliHadCorrTask.cxx:666
 AliHadCorrTask.cxx:667
 AliHadCorrTask.cxx:668
 AliHadCorrTask.cxx:669
 AliHadCorrTask.cxx:670
 AliHadCorrTask.cxx:671
 AliHadCorrTask.cxx:672
 AliHadCorrTask.cxx:673
 AliHadCorrTask.cxx:674
 AliHadCorrTask.cxx:675
 AliHadCorrTask.cxx:676
 AliHadCorrTask.cxx:677
 AliHadCorrTask.cxx:678
 AliHadCorrTask.cxx:679
 AliHadCorrTask.cxx:680
 AliHadCorrTask.cxx:681
 AliHadCorrTask.cxx:682
 AliHadCorrTask.cxx:683
 AliHadCorrTask.cxx:684
 AliHadCorrTask.cxx:685
 AliHadCorrTask.cxx:686
 AliHadCorrTask.cxx:687
 AliHadCorrTask.cxx:688
 AliHadCorrTask.cxx:689
 AliHadCorrTask.cxx:690
 AliHadCorrTask.cxx:691
 AliHadCorrTask.cxx:692
 AliHadCorrTask.cxx:693
 AliHadCorrTask.cxx:694
 AliHadCorrTask.cxx:695
 AliHadCorrTask.cxx:696
 AliHadCorrTask.cxx:697
 AliHadCorrTask.cxx:698
 AliHadCorrTask.cxx:699
 AliHadCorrTask.cxx:700
 AliHadCorrTask.cxx:701
 AliHadCorrTask.cxx:702
 AliHadCorrTask.cxx:703
 AliHadCorrTask.cxx:704
 AliHadCorrTask.cxx:705
 AliHadCorrTask.cxx:706
 AliHadCorrTask.cxx:707
 AliHadCorrTask.cxx:708
 AliHadCorrTask.cxx:709
 AliHadCorrTask.cxx:710
 AliHadCorrTask.cxx:711
 AliHadCorrTask.cxx:712
 AliHadCorrTask.cxx:713
 AliHadCorrTask.cxx:714
 AliHadCorrTask.cxx:715
 AliHadCorrTask.cxx:716
 AliHadCorrTask.cxx:717
 AliHadCorrTask.cxx:718
 AliHadCorrTask.cxx:719
 AliHadCorrTask.cxx:720
 AliHadCorrTask.cxx:721
 AliHadCorrTask.cxx:722
 AliHadCorrTask.cxx:723
 AliHadCorrTask.cxx:724
 AliHadCorrTask.cxx:725
 AliHadCorrTask.cxx:726
 AliHadCorrTask.cxx:727
 AliHadCorrTask.cxx:728
 AliHadCorrTask.cxx:729
 AliHadCorrTask.cxx:730
 AliHadCorrTask.cxx:731
 AliHadCorrTask.cxx:732
 AliHadCorrTask.cxx:733
 AliHadCorrTask.cxx:734
 AliHadCorrTask.cxx:735
 AliHadCorrTask.cxx:736
 AliHadCorrTask.cxx:737
 AliHadCorrTask.cxx:738
 AliHadCorrTask.cxx:739
 AliHadCorrTask.cxx:740
 AliHadCorrTask.cxx:741
 AliHadCorrTask.cxx:742
 AliHadCorrTask.cxx:743
 AliHadCorrTask.cxx:744
 AliHadCorrTask.cxx:745
 AliHadCorrTask.cxx:746
 AliHadCorrTask.cxx:747
 AliHadCorrTask.cxx:748
 AliHadCorrTask.cxx:749
 AliHadCorrTask.cxx:750
 AliHadCorrTask.cxx:751
 AliHadCorrTask.cxx:752
 AliHadCorrTask.cxx:753
 AliHadCorrTask.cxx:754
 AliHadCorrTask.cxx:755
 AliHadCorrTask.cxx:756
 AliHadCorrTask.cxx:757
 AliHadCorrTask.cxx:758
 AliHadCorrTask.cxx:759
 AliHadCorrTask.cxx:760
 AliHadCorrTask.cxx:761
 AliHadCorrTask.cxx:762
 AliHadCorrTask.cxx:763
 AliHadCorrTask.cxx:764
 AliHadCorrTask.cxx:765
 AliHadCorrTask.cxx:766
 AliHadCorrTask.cxx:767
 AliHadCorrTask.cxx:768
 AliHadCorrTask.cxx:769
 AliHadCorrTask.cxx:770
 AliHadCorrTask.cxx:771
 AliHadCorrTask.cxx:772
 AliHadCorrTask.cxx:773
 AliHadCorrTask.cxx:774
 AliHadCorrTask.cxx:775
 AliHadCorrTask.cxx:776
 AliHadCorrTask.cxx:777
 AliHadCorrTask.cxx:778
 AliHadCorrTask.cxx:779
 AliHadCorrTask.cxx:780
 AliHadCorrTask.cxx:781
 AliHadCorrTask.cxx:782
 AliHadCorrTask.cxx:783
 AliHadCorrTask.cxx:784
 AliHadCorrTask.cxx:785
 AliHadCorrTask.cxx:786
 AliHadCorrTask.cxx:787
 AliHadCorrTask.cxx:788
 AliHadCorrTask.cxx:789
 AliHadCorrTask.cxx:790
 AliHadCorrTask.cxx:791
 AliHadCorrTask.cxx:792
 AliHadCorrTask.cxx:793
 AliHadCorrTask.cxx:794
 AliHadCorrTask.cxx:795
 AliHadCorrTask.cxx:796
 AliHadCorrTask.cxx:797
 AliHadCorrTask.cxx:798
 AliHadCorrTask.cxx:799
 AliHadCorrTask.cxx:800
 AliHadCorrTask.cxx:801
 AliHadCorrTask.cxx:802
 AliHadCorrTask.cxx:803
 AliHadCorrTask.cxx:804
 AliHadCorrTask.cxx:805
 AliHadCorrTask.cxx:806
 AliHadCorrTask.cxx:807
 AliHadCorrTask.cxx:808
 AliHadCorrTask.cxx:809
 AliHadCorrTask.cxx:810
 AliHadCorrTask.cxx:811
 AliHadCorrTask.cxx:812
 AliHadCorrTask.cxx:813
 AliHadCorrTask.cxx:814
 AliHadCorrTask.cxx:815
 AliHadCorrTask.cxx:816
 AliHadCorrTask.cxx:817
 AliHadCorrTask.cxx:818
 AliHadCorrTask.cxx:819
 AliHadCorrTask.cxx:820
 AliHadCorrTask.cxx:821
 AliHadCorrTask.cxx:822
 AliHadCorrTask.cxx:823
 AliHadCorrTask.cxx:824
 AliHadCorrTask.cxx:825
 AliHadCorrTask.cxx:826
 AliHadCorrTask.cxx:827
 AliHadCorrTask.cxx:828
 AliHadCorrTask.cxx:829
 AliHadCorrTask.cxx:830
 AliHadCorrTask.cxx:831
 AliHadCorrTask.cxx:832
 AliHadCorrTask.cxx:833
 AliHadCorrTask.cxx:834
 AliHadCorrTask.cxx:835
 AliHadCorrTask.cxx:836
 AliHadCorrTask.cxx:837
 AliHadCorrTask.cxx:838
 AliHadCorrTask.cxx:839
 AliHadCorrTask.cxx:840
 AliHadCorrTask.cxx:841
 AliHadCorrTask.cxx:842
 AliHadCorrTask.cxx:843
 AliHadCorrTask.cxx:844
 AliHadCorrTask.cxx:845
 AliHadCorrTask.cxx:846
 AliHadCorrTask.cxx:847
 AliHadCorrTask.cxx:848
 AliHadCorrTask.cxx:849
 AliHadCorrTask.cxx:850
 AliHadCorrTask.cxx:851
 AliHadCorrTask.cxx:852
 AliHadCorrTask.cxx:853
 AliHadCorrTask.cxx:854
 AliHadCorrTask.cxx:855
 AliHadCorrTask.cxx:856
 AliHadCorrTask.cxx:857
 AliHadCorrTask.cxx:858
 AliHadCorrTask.cxx:859
 AliHadCorrTask.cxx:860
 AliHadCorrTask.cxx:861
 AliHadCorrTask.cxx:862
 AliHadCorrTask.cxx:863
 AliHadCorrTask.cxx:864
 AliHadCorrTask.cxx:865
 AliHadCorrTask.cxx:866
 AliHadCorrTask.cxx:867
 AliHadCorrTask.cxx:868
 AliHadCorrTask.cxx:869
 AliHadCorrTask.cxx:870
 AliHadCorrTask.cxx:871
 AliHadCorrTask.cxx:872
 AliHadCorrTask.cxx:873
 AliHadCorrTask.cxx:874
 AliHadCorrTask.cxx:875
 AliHadCorrTask.cxx:876
 AliHadCorrTask.cxx:877
 AliHadCorrTask.cxx:878
 AliHadCorrTask.cxx:879
 AliHadCorrTask.cxx:880
 AliHadCorrTask.cxx:881
 AliHadCorrTask.cxx:882
 AliHadCorrTask.cxx:883
 AliHadCorrTask.cxx:884
 AliHadCorrTask.cxx:885
 AliHadCorrTask.cxx:886
 AliHadCorrTask.cxx:887
 AliHadCorrTask.cxx:888
 AliHadCorrTask.cxx:889
 AliHadCorrTask.cxx:890
 AliHadCorrTask.cxx:891
 AliHadCorrTask.cxx:892
 AliHadCorrTask.cxx:893
 AliHadCorrTask.cxx:894
 AliHadCorrTask.cxx:895
 AliHadCorrTask.cxx:896
 AliHadCorrTask.cxx:897
 AliHadCorrTask.cxx:898
 AliHadCorrTask.cxx:899
 AliHadCorrTask.cxx:900
 AliHadCorrTask.cxx:901
 AliHadCorrTask.cxx:902
 AliHadCorrTask.cxx:903
 AliHadCorrTask.cxx:904
 AliHadCorrTask.cxx:905
 AliHadCorrTask.cxx:906
 AliHadCorrTask.cxx:907
 AliHadCorrTask.cxx:908
 AliHadCorrTask.cxx:909
 AliHadCorrTask.cxx:910
 AliHadCorrTask.cxx:911
 AliHadCorrTask.cxx:912
 AliHadCorrTask.cxx:913
 AliHadCorrTask.cxx:914
 AliHadCorrTask.cxx:915
 AliHadCorrTask.cxx:916
 AliHadCorrTask.cxx:917
 AliHadCorrTask.cxx:918
 AliHadCorrTask.cxx:919
 AliHadCorrTask.cxx:920
 AliHadCorrTask.cxx:921
 AliHadCorrTask.cxx:922
 AliHadCorrTask.cxx:923
 AliHadCorrTask.cxx:924
 AliHadCorrTask.cxx:925
 AliHadCorrTask.cxx:926
 AliHadCorrTask.cxx:927
 AliHadCorrTask.cxx:928
 AliHadCorrTask.cxx:929
 AliHadCorrTask.cxx:930
 AliHadCorrTask.cxx:931
 AliHadCorrTask.cxx:932
 AliHadCorrTask.cxx:933
 AliHadCorrTask.cxx:934
 AliHadCorrTask.cxx:935
 AliHadCorrTask.cxx:936
 AliHadCorrTask.cxx:937
 AliHadCorrTask.cxx:938
 AliHadCorrTask.cxx:939
 AliHadCorrTask.cxx:940
 AliHadCorrTask.cxx:941
 AliHadCorrTask.cxx:942
 AliHadCorrTask.cxx:943
 AliHadCorrTask.cxx:944
 AliHadCorrTask.cxx:945
 AliHadCorrTask.cxx:946
 AliHadCorrTask.cxx:947
 AliHadCorrTask.cxx:948
 AliHadCorrTask.cxx:949