ROOT logo
/**************************************************************************
 * Copyright(c) 1998-2008, 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.                  *
 **************************************************************************/

/* $Id$ */

/////////////////////////////////////////////////////////////
//
// AliAnalysisTaskSE for the comparison of heavy flavor
// decay candidates with the MC truth.
//
// Author: A.Dainese, andrea.dainese@lnl.infn.it
/////////////////////////////////////////////////////////////

#include <TClonesArray.h>
#include <TNtuple.h>
#include <TList.h>
#include <TH1F.h>

#include "AliAnalysisManager.h"
#include "AliAODHandler.h"
#include "AliAODEvent.h"
#include "AliAODVertex.h"
#include "AliESDtrack.h"
#include "AliAODTrack.h"
#include "AliAODMCHeader.h"
#include "AliAODMCParticle.h"
#include "AliAODRecoDecayHF2Prong.h"
#include "AliAODRecoDecayHF3Prong.h"
#include "AliAODRecoDecayHF4Prong.h"
#include "AliAODRecoCascadeHF.h"
#include "AliAnalysisVertexingHF.h"
#include "AliAnalysisTaskSE.h"
#include "AliAnalysisTaskSECompareHF.h"

ClassImp(AliAnalysisTaskSECompareHF)


//________________________________________________________________________
AliAnalysisTaskSECompareHF::AliAnalysisTaskSECompareHF():
AliAnalysisTaskSE(),
fOutput(0), 
fNtupleCmp(0),
fHistMass(0),
fHistNEvents(0)
{
  // Default constructor

  // NO DefineOutput() HERE (ONLY IN STANDARD CONSTRUCTOR)
}

//________________________________________________________________________
AliAnalysisTaskSECompareHF::AliAnalysisTaskSECompareHF(const char *name):
AliAnalysisTaskSE(name),
fOutput(0), 
fNtupleCmp(0),
fHistMass(0),
fHistNEvents(0)
{
  // Standard constructor

  // Output slot #1 writes into a TList container
  DefineOutput(1,TList::Class());  //My private output
  // Output slot #2 writes into a TNtuple container
  DefineOutput(2,TNtuple::Class());  //My private output
}

//________________________________________________________________________
AliAnalysisTaskSECompareHF::~AliAnalysisTaskSECompareHF()
{
  // Destructor
  if (fOutput) {
    delete fOutput;
    fOutput = 0;
  }
}  

//________________________________________________________________________
void AliAnalysisTaskSECompareHF::Init()
{
  // Initialization

  if(fDebug > 1) printf("AnalysisTaskSECompareHF::Init() \n");
  
  return;
}

//________________________________________________________________________
void AliAnalysisTaskSECompareHF::UserCreateOutputObjects()
{
  // Create the output container
  //
  if(fDebug > 1) printf("AnalysisTaskSECompareHF::UserCreateOutputObjects() \n");

  // Several histograms are more conveniently managed in a TList
  fOutput = new TList();
  fOutput->SetOwner();

  fHistMass = new TH1F("fHistMass", "D^{0} invariant mass; M [GeV]; Entries",200,1.765,1.965);
  fHistMass->Sumw2();
  fHistMass->SetMinimum(0);
  fOutput->Add(fHistMass);

  fHistNEvents = new TH1F("fHistNEvents", "Number of processed events; ; Events",3,-1.5,1.5);
  fHistNEvents->Sumw2();
  fHistNEvents->SetMinimum(0);
  fOutput->Add(fHistNEvents);

  OpenFile(2); // 2 is the slot number of the ntuple
  fNtupleCmp = new TNtuple("fNtupleCmp","Charm comparison","pdg:nprongs:VxRec:VxTrue:ErrVx:VyRec:VyTrue:ErrVy:VzRec:VzTrue:ErrVz:Chi2toNDF:PtRec:Mrec:CPta:Prodd0");

  return;
}

//________________________________________________________________________
void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/)
{
  // Execute analysis for current event:
  // heavy flavor candidates association to MC truth

  
  AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());

  TClonesArray *inputArrayVertices = 0;
  TClonesArray *inputArrayD0toKpi = 0;
  TClonesArray *inputArrayDstar = 0;

  if(!aod && AODEvent() && IsStandardAOD()) {
    // In case there is an AOD handler writing a standard AOD, use the AOD 
    // event in memory rather than the input (ESD) event.    
    aod = dynamic_cast<AliAODEvent*> (AODEvent());
    // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
    // have to taken from the AOD event hold by the AliAODExtension
    AliAODHandler* aodHandler = (AliAODHandler*) 
      ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
    if(aodHandler->GetExtensions()) {
      AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
      AliAODEvent *aodFromExt = ext->GetAOD();
      // load HF vertices                
      inputArrayVertices = (TClonesArray*)aodFromExt->GetList()->FindObject("VerticesHF");
      // load D0->Kpi candidates
      inputArrayD0toKpi = (TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
      // load D*+ candidates                                                   
      inputArrayDstar = (TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
    }
  } else if(aod) {
    // load HF vertices                
    inputArrayVertices = (TClonesArray*)aod->GetList()->FindObject("VerticesHF");
    // load D0->Kpi candidates                                                 
    inputArrayD0toKpi = (TClonesArray*)aod->GetList()->FindObject("D0toKpi");
    // load D*+ candidates                                                   
    inputArrayDstar = (TClonesArray*)aod->GetList()->FindObject("Dstar");
  }


  if(!inputArrayVertices || !aod) {
    printf("AliAnalysisTaskSECompareHF::UserExec: Vertices branch not found!\n");
    return;
  }
  if(!inputArrayD0toKpi) {
    printf("AliAnalysisTaskSECompareHF::UserExec: D0toKpi branch not found!\n");
    return;
  }
  if(!inputArrayDstar) {
    printf("AliAnalysisTaskSECompareHF::UserExec: Dstar branch not found!\n");
    return;
  }
  

  fHistNEvents->Fill(0); // count event
  // Post the data already here
  PostData(1,fOutput);

  // AOD primary vertex
  AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
  //vtx1->Print();

  // load MC particles
  TClonesArray *mcArray = 
    (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
  if(!mcArray) {
    printf("AliAnalysisTaskSECompareHF::UserExec: MC particles branch not found!\n");
    return;
  }

  // load MC header
  AliAODMCHeader *mcHeader = 
    (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
  if(!mcHeader) {
    printf("AliAnalysisTaskSECompareHF::UserExec: MC header branch not found!\n");
    return;
  }
  
  Int_t nprongs,lab,okD0,okD0bar,pdg;
  Bool_t unsetvtx;    
  Double_t invmass,posRec[3],posTrue[3],covRec[6],errx,erry,errz;
  AliAODRecoDecayHF2Prong *d2=0;
  AliAODRecoDecayHF3Prong *d3=0;
  AliAODRecoDecayHF4Prong *d4=0;

  Int_t pdgDgD0toKpi[2]={321,211};
  Int_t pdgDgDplustoKpipi[3]={321,211,211};
  Int_t pdgDgD0toKpipipi[4]={321,211,211,211};

  // loop over vertices
  Int_t nVertices = inputArrayVertices->GetEntriesFast();
  if(fDebug>1) printf("Number of vertices: %d\n",nVertices);

  for (Int_t iVtx = 0; iVtx < nVertices; iVtx++) {
    AliAODVertex *vtx = (AliAODVertex*)inputArrayVertices->UncheckedAt(iVtx);

    vtx->GetXYZ(posRec);  
    vtx->GetCovarianceMatrix(covRec);
    errx=1.; erry=1.; errz=1.;
    if(covRec[0]>0) errx = TMath::Sqrt(covRec[0]);
    if(covRec[2]>0) erry = TMath::Sqrt(covRec[2]);
    if(covRec[5]>0) errz = TMath::Sqrt(covRec[5]);
	  

    // get parent AliAODRecoDecayHF
    nprongs= vtx->GetNDaughters();
    // check that the daughters have kITSrefit and kTPCrefit
    Bool_t allDgOK=kTRUE;
    for(Int_t idg=0; idg<nprongs; idg++) {
      AliAODTrack *track = (AliAODTrack*)vtx->GetDaughter(idg);
      if(!(track->GetStatus()&AliESDtrack::kITSrefit)) allDgOK=kFALSE;
      if(!(track->GetStatus()&AliESDtrack::kTPCrefit)) allDgOK=kFALSE;
    }
    if(!allDgOK) continue;


    switch(nprongs) {
    case 2: // look for D0->Kpi
      d2 = (AliAODRecoDecayHF2Prong*)vtx->GetParent();
      if(d2->IsLikeSign()) continue;
      if(d2->Charge() != 0) continue; // these are D* 
      lab = d2->MatchToMC(421,mcArray,2,pdgDgD0toKpi);
      if(lab>=0) {
	unsetvtx=kFALSE;
	if(!d2->GetOwnPrimaryVtx()) {
	  d2->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
	  unsetvtx=kTRUE;
	}
	okD0=1; okD0bar=1; 
	AliAODMCParticle *dMC = (AliAODMCParticle*)mcArray->At(lab);
	pdg = dMC->GetPdgCode();
	invmass = (pdg==421 ? d2->InvMassD0() : d2->InvMassD0bar());
	// get a daughter for true pos of decay vertex
	AliAODMCParticle *dg0MC = (AliAODMCParticle*)mcArray->At(dMC->GetDaughter(0));
	dg0MC->XvYvZv(posTrue);
	fHistMass->Fill(invmass);
	// Post the data already here
	PostData(1,fOutput);
	Float_t tmp[16]={(Float_t)pdg,(Float_t)nprongs,
			 (Float_t)posRec[0],(Float_t)posTrue[0],(Float_t)errx,
			 (Float_t)posRec[1],(Float_t)posTrue[1],(Float_t)erry,
			 (Float_t)posRec[2],(Float_t)posTrue[2],(Float_t)errz,
			 (Float_t)d2->GetReducedChi2(),(Float_t)d2->Pt(),(Float_t)invmass,
			 (Float_t)d2->CosPointingAngle(),(Float_t)d2->Prodd0d0()};
	fNtupleCmp->Fill(tmp);
	PostData(2,fNtupleCmp);
      
	if(unsetvtx) d2->UnsetOwnPrimaryVtx();
      }
      break;
    case 3: // look for D+
      d3 = (AliAODRecoDecayHF3Prong*)vtx->GetParent();
      if(d3->IsLikeSign()) continue;
      lab = d3->MatchToMC(411,mcArray,3,pdgDgDplustoKpipi);
      if(lab>=0) {
	unsetvtx=kFALSE;
	if(!d3->GetOwnPrimaryVtx()) {
	  d3->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
	  unsetvtx=kTRUE;
	}
	AliAODMCParticle *dMC = (AliAODMCParticle*)mcArray->At(lab);
	pdg = dMC->GetPdgCode();
	invmass = d3->InvMassDplus();
	// get a daughter for true pos of decay vertex
	AliAODMCParticle *dg0MC = (AliAODMCParticle*)mcArray->At(dMC->GetDaughter(0));
	dg0MC->XvYvZv(posTrue);
	Float_t tmp[16]={(Float_t)pdg,(Float_t)nprongs,
			 (Float_t)posRec[0],(Float_t)posTrue[0],(Float_t)errx,
			 (Float_t)posRec[1],(Float_t)posTrue[1],(Float_t)erry,
			 (Float_t)posRec[2],(Float_t)posTrue[2],(Float_t)errz,
			 (Float_t)d3->GetReducedChi2(),(Float_t)d3->Pt(),(Float_t)invmass,
			 (Float_t)d3->CosPointingAngle(),(Float_t)(d3->Getd0Prong(0)*d3->Getd0Prong(1)*d3->Getd0Prong(2))};
	fNtupleCmp->Fill(tmp);
	PostData(2,fNtupleCmp);
	if(unsetvtx) d3->UnsetOwnPrimaryVtx();
      }
      break;
    case 4: // look for D0->Kpipipi
      d4 = (AliAODRecoDecayHF4Prong*)vtx->GetParent();
      if(d4->IsLikeSign()) continue;
      lab = d4->MatchToMC(421,mcArray,4,pdgDgD0toKpipipi);
      if(lab>=0) {
	unsetvtx=kFALSE;
	if(!d4->GetOwnPrimaryVtx()) {
	  d4->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
	  unsetvtx=kTRUE;
	}
	okD0=0; okD0bar=0; 
	AliAODMCParticle *dMC = (AliAODMCParticle*)mcArray->At(lab);
	pdg = dMC->GetPdgCode();
	//invmass = (pdg==421 ? d->InvMassD0() : d->InvMassD0bar());
	invmass = 10.; 	  // dummy
	// get a daughter for true pos of decay vertex
	AliAODMCParticle *dg0MC = (AliAODMCParticle*)mcArray->At(dMC->GetDaughter(0));
	dg0MC->XvYvZv(posTrue);
	Float_t tmp[16]={(Float_t)pdg,(Float_t)nprongs,
			 (Float_t)posRec[0],(Float_t)posTrue[0],(Float_t)errx,
			 (Float_t)posRec[1],(Float_t)posTrue[1],(Float_t)erry,
			 (Float_t)posRec[2],(Float_t)posTrue[2],(Float_t)errz,
			 (Float_t)d4->GetReducedChi2(),(Float_t)d4->Pt(),(Float_t)invmass,
			 (Float_t)d4->CosPointingAngle(),(Float_t)(d4->Getd0Prong(0)*d4->Getd0Prong(1)*d4->Getd0Prong(2)*d4->Getd0Prong(3))};
	fNtupleCmp->Fill(tmp);
	PostData(2,fNtupleCmp);
	if(unsetvtx) d4->UnsetOwnPrimaryVtx();
      }
      break;
    }

  } // end loop on vertices


  // loop over D*+ candidates
  /*
  for (Int_t iDstar = 0; iDstar < inputArrayDstar->GetEntries(); iDstar++) {
    AliAODRecoCascadeHF *c = (AliAODRecoCascadeHF*)inputArrayDstar->UncheckedAt(iDstar);
    Int_t labDstar = c->MatchToMC(413,421,mcArray);
    if(labDstar>=0) printf("GOOD MATCH FOR D*+\n"); 
  }
  */

  // Post the data already here
  PostData(1,fOutput);
  PostData(2,fNtupleCmp);

  return;
}
//________________________________________________________________________
void AliAnalysisTaskSECompareHF::Terminate(Option_t */*option*/)
{
  // Terminate analysis
  //
  if(fDebug > 1) printf("AnalysisTaskSECompareHF: Terminate() \n");

  fOutput = dynamic_cast<TList*> (GetOutputData(1));
  if (!fOutput) {     
    printf("ERROR: fOutput not available\n");
    return;
  }

  fHistMass = dynamic_cast<TH1F*>(fOutput->FindObject("fHistMass"));
  fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));

  //fNtupleCmp = dynamic_cast<TNtuple*> (GetOutputData(2));

  return;
}

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