From 86c32a36d08368003da0aec8a82c204a464ea17d Mon Sep 17 00:00:00 2001 From: rvernet Date: Tue, 26 Aug 2008 15:13:50 +0000 Subject: [PATCH] Support for AOD in example tasks and macros. --- CORRFW/test/AliCFRsnTask.C | 40 +++++- CORRFW/test/AliCFRsnTask.cxx | 197 +++++++++++++++++------------ CORRFW/test/AliCFSingleTrackTask.C | 2 +- CORRFW/test/AliCFV0Task.C | 64 ++++------ CORRFW/test/AliCFV0Task.cxx | 129 +++++++++++++++---- CORRFW/test/AliCFV0Task.h | 4 +- 6 files changed, 278 insertions(+), 158 deletions(-) diff --git a/CORRFW/test/AliCFRsnTask.C b/CORRFW/test/AliCFRsnTask.C index f3c19076462..55dbee97620 100644 --- a/CORRFW/test/AliCFRsnTask.C +++ b/CORRFW/test/AliCFRsnTask.C @@ -13,6 +13,7 @@ const Double32_t nsigmavtx = 3. ; //max track sigma to PVertex Bool_t AliCFRsnTask( const Bool_t useGrid = 1, + const Bool_t readAOD = 0, const char * kTagXMLFile="wn.xml", // XML file containing tags ) { @@ -44,9 +45,23 @@ Bool_t AliCFRsnTask( analysisChain = tagAna->QueryTags(runCuts,lhcCuts,detCuts,eventCuts); } else {// local data - analysisChain = new TChain("esdTree"); //here put your input data path - analysisChain->Add("AliESDs.root"); + if (readAOD) { + analysisChain = new TChain("aodTree"); + //analysisChain->Add("AliAOD.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/001/AliAOD.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/002/AliAOD.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/003/AliAOD.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/004/AliAOD.root"); + } + else { + analysisChain = new TChain("esdTree"); + //analysisChain->Add("AliESDs.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/001/AliESDs.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/002/AliESDs.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/003/AliESDs.root"); + analysisChain->Add("/home/vernet/Data/LHC08b2/300000/004/AliESDs.root"); + } } @@ -102,11 +117,15 @@ Bool_t AliCFRsnTask( recKineCuts->SetChargeRec(charge); AliCFPairQualityCuts *recQualityCuts = new AliCFPairQualityCuts("recQualityCuts","rec-level quality cuts"); - recQualityCuts->SetMinNClusterTPC(minclustersTPC,minclustersTPC); - recQualityCuts->SetRequireITSRefit(kTRUE,kTRUE); + if (!readAOD) recQualityCuts->SetMinNClusterTPC(minclustersTPC,minclustersTPC); + recQualityCuts->SetStatus(AliESDtrack::kTPCrefit & AliESDtrack::kITSrefit, + AliESDtrack::kTPCrefit & AliESDtrack::kITSrefit) ; + + AliCFPairIsPrimaryCuts *recIsPrimaryCuts = new AliCFPairIsPrimaryCuts("recIsPrimaryCuts","rec-level isPrimary cuts"); - recIsPrimaryCuts->SetMaxNSigmaToVertex(nsigmavtx,nsigmavtx); + if (readAOD) recIsPrimaryCuts->SetAODType(AliAODTrack::kPrimary,AliAODTrack::kPrimary); + else recIsPrimaryCuts->SetMaxNSigmaToVertex(nsigmavtx,nsigmavtx); AliCFPairPidCut* cutPID = new AliCFPairPidCut("cutPID","ESD_PID") ; Double_t prior_pp[AliPID::kSPECIES] = {0.0244519, @@ -125,6 +144,9 @@ Bool_t AliCFRsnTask( cutPID->SetPriors(prior_pbpb); cutPID->SetProbabilityCut(0.,0.); cutPID->SetDetectors("TPC ITS TOF TRD","TPC ITS TOF TRD"); + if (readAOD) cutPID->SetAODmode(kTRUE); + else cutPID->SetAODmode(kFALSE); + switch(PDG) { case -313 : cutPID->SetParticleType(AliPID::kKaon ,kTRUE,AliPID::kPion ,kTRUE); break; case 313 : cutPID->SetParticleType(AliPID::kPion ,kTRUE,AliPID::kKaon ,kTRUE); break; @@ -182,9 +204,13 @@ Bool_t AliCFRsnTask( else mgr->SetAnalysisType(AliAnalysisManager::kLocalAnalysis); AliMCEventHandler* mcHandler = new AliMCEventHandler(); - AliESDInputHandler* esdHandler = new AliESDInputHandler(); mgr->SetMCtruthEventHandler(mcHandler); - mgr->SetInputEventHandler(esdHandler); + + AliInputEventHandler* dataHandler ; + if (readAOD) dataHandler = new AliAODInputHandler(); + else dataHandler = new AliESDInputHandler(); + mgr->SetInputEventHandler(dataHandler); + // Create and connect containers for input/output diff --git a/CORRFW/test/AliCFRsnTask.cxx b/CORRFW/test/AliCFRsnTask.cxx index fefb018f648..2538e1d59ad 100644 --- a/CORRFW/test/AliCFRsnTask.cxx +++ b/CORRFW/test/AliCFRsnTask.cxx @@ -1,3 +1,4 @@ + /************************************************************************** * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * * * @@ -33,7 +34,9 @@ #include "AliLog.h" #include "AliRsnDaughter.h" #include "AliCFPair.h" -#include "AliRsnParticle.h" +#include "AliRsnMCInfo.h" +#include "AliRsnPairParticle.h" +#include "AliAODEvent.h" //__________________________________________________________________________ AliCFRsnTask::AliCFRsnTask() : @@ -111,9 +114,8 @@ void AliCFRsnTask::UserExec(Option_t *) // Info("UserExec","") ; - AliESDEvent* fESD = dynamic_cast(fInputEvent); - if (!fESD) { - Error("UserExec","NO ESD FOUND!"); + if (!fInputEvent) { + Error("UserExec","NO EVENT FOUND!"); return; } @@ -122,6 +124,9 @@ void AliCFRsnTask::UserExec(Option_t *) AliStack* stack = fMCEvent->Stack(); + Bool_t isESDEvent = strcmp(fInputEvent->ClassName(),"AliESDEvent") == 0 ? kTRUE : kFALSE ; + Bool_t isAODEvent = strcmp(fInputEvent->ClassName(),"AliAODEvent") == 0 ? kTRUE : kFALSE ; + // MC-event selection Double_t containerInput[2] ; @@ -145,86 +150,96 @@ void AliCFRsnTask::UserExec(Option_t *) //Now go to rec level - Info("UserExec","Looping on ESD event"); - - //SET THE ESD AS EVENT INFO IN RECONSTRUCTION CUTS - TObjArray* fCutsReco = fCFManager->GetParticleCutsList(AliCFManager::kPartRecCuts); - TObjArray* fCutsPID = fCFManager->GetParticleCutsList(AliCFManager::kPartSelCuts); - TObjArrayIter iter1(fCutsReco); - TObjArrayIter iter2(fCutsPID); - AliCFCutBase *cut = 0; - while ( (cut = (AliCFCutBase*)iter1.Next()) ) { - cut->SetEvtInfo(fESD); - } - while ( (cut = (AliCFCutBase*)iter2.Next()) ) { - cut->SetEvtInfo(fESD); - } - + Info("UserExec","Looping on %s",fInputEvent->ClassName()); + // Loop on negative tracks - for (Int_t iTrack1 = 0; iTrack1GetNumberOfTracks(); iTrack1++) { - AliESDtrack* esdTrack1 = fESD->GetTrack(iTrack1); + for (Int_t iTrack1 = 0; iTrack1GetNumberOfTracks(); iTrack1++) { + AliVParticle* track1 = fInputEvent->GetTrack(iTrack1); //track1 is negative - if (esdTrack1->Charge()>=0) continue; - Int_t esdLabel1 = esdTrack1->GetLabel(); - if (esdLabel1<0) continue; + if (track1->Charge()>=0) continue; + Int_t label1 = track1->GetLabel(); + if (label1<0) continue; //Loop on positive tracks - for (Int_t iTrack2 = 0; iTrack2GetNumberOfTracks(); iTrack2++) { - AliESDtrack* esdTrack2 = fESD->GetTrack(iTrack2); + for (Int_t iTrack2 = 0; iTrack2GetNumberOfTracks(); iTrack2++) { + AliVParticle* track2 = fInputEvent->GetTrack(iTrack2); //track2 is positive - if (esdTrack2->Charge()<=0) continue; - Int_t esdLabel2 = esdTrack2->GetLabel(); - if (esdLabel2<0) continue; - + if (track2->Charge()<=0) continue; + Int_t label2 = track2->GetLabel(); + if (label2<0) continue; + //Create Resonance daughter objects - AliRsnDaughter* tmp1 = new AliRsnDaughter(esdTrack1); - AliRsnDaughter* tmp2 = new AliRsnDaughter(esdTrack2); - AliRsnDaughter track1(*tmp1); - AliRsnDaughter track2(*tmp2); - delete tmp1; - delete tmp2; - - //Set MC information to resonance daughters - TParticle *part1 = stack->Particle(esdLabel1); - track1.InitParticle(part1); - track1.GetParticle()->SetPDG(part1->GetPdgCode()); + AliRsnDaughter* daughter1tmp = 0x0 ; + AliRsnDaughter* daughter2tmp = 0x0 ; + if (isESDEvent) { + daughter1tmp = new AliRsnDaughter((AliESDtrack*)track1) ; + daughter2tmp = new AliRsnDaughter((AliESDtrack*)track2) ; + } + else if (isAODEvent) { + daughter1tmp = new AliRsnDaughter((AliAODTrack*)track1) ; + daughter2tmp = new AliRsnDaughter((AliAODTrack*)track2) ; + } + else { + Error("UserExec","Error: input data file is not an ESD nor an AOD"); + return; + } + + AliRsnDaughter daughter1(*daughter1tmp); + AliRsnDaughter daughter2(*daughter2tmp); + delete daughter1tmp; + delete daughter2tmp; + + AliCFPair pair(track1,track2); // This object is used for cuts + // (to be replaced when AliRsnPairParticle + // inherits from AliVParticle) + + //Set MC PDG information to resonance daughters + TParticle *part1 = stack->Particle(label1); + daughter1.InitMCInfo(part1); + daughter1.GetMCInfo()->SetPDG(part1->GetPdgCode()); + daughter1.SetM(part1->GetCalcMass()); // assign true mass Int_t mother1 = part1->GetFirstMother(); - track1.GetParticle()->SetMother(mother1); + daughter1.GetMCInfo()->SetMother(mother1); if (mother1 >= 0) { TParticle *mum = stack->Particle(mother1); - track1.GetParticle()->SetMotherPDG(mum->GetPdgCode()); + daughter1.GetMCInfo()->SetMotherPDG(mum->GetPdgCode()); } - - TParticle *part2 = stack->Particle(esdLabel2); - track2.InitParticle(part2); - track2.GetParticle()->SetPDG(part2->GetPdgCode()); + + TParticle *part2 = stack->Particle(label2); + daughter2.InitMCInfo(part2); + daughter2.GetMCInfo()->SetPDG(part2->GetPdgCode()); + daughter2.SetM(part2->GetCalcMass()); // assign true mass Int_t mother2 = part2->GetFirstMother(); - track2.GetParticle()->SetMother(mother2); + daughter2.GetMCInfo()->SetMother(mother2); if (mother2 >= 0) { TParticle *mum = stack->Particle(mother2); - track2.GetParticle()->SetMotherPDG(mum->GetPdgCode()); + daughter2.GetMCInfo()->SetMotherPDG(mum->GetPdgCode()); } //make a mother resonance from the 2 candidate daughters - AliRsnDaughter rsn = AliRsnDaughter::Sum(track1,track2) ; - AliCFPair pair(esdTrack1,esdTrack2); // This object is used for cuts (to be replaced) + AliRsnPairParticle rsn; + rsn.SetPair(&daughter1,&daughter2); //check if true resonance - if (rsn.GetParticle()->PDG() != fRsnPDG) continue; + if (!rsn.IsTruePair(fRsnPDG)) continue; if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,&pair)) continue; //check if associated MC resonance passes the cuts - Int_t motherLabel=rsn.Label(); + Int_t motherLabel = rsn.GetDaughter(0)->GetMCInfo()->Mother() ; if (motherLabel<0) continue ; + AliMCParticle* mcRsn = fMCEvent->GetTrack(motherLabel); if (!mcRsn) continue; if (!fCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,mcRsn)) continue; - //fill the container - containerInput[0] = rsn.Pt() ; - containerInput[1] = GetRapidity(rsn.E(),rsn.Pz()); + // fill the container + Double_t rsnEnergy = rsn.GetDaughter(0)->E() + rsn.GetDaughter(1)->E() ; + Double_t rsnPz = rsn.GetDaughter(0)->Pz() + rsn.GetDaughter(1)->Pz() ; + + containerInput[0] = rsn.GetPt() ; + containerInput[1] = GetRapidity(rsnEnergy,rsnPz); fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed) ; if (!fCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,&pair)) continue ; @@ -251,53 +266,67 @@ void AliCFRsnTask::Terminate(Option_t*) // the results graphically or save the results to file. Info("Terminate",""); + AliAnalysisTaskSE::Terminate(); + - Double_t max1 = fCFManager->GetParticleContainer()->ShowProjection(0,0)->GetMaximum(); - Double_t max2 = fCFManager->GetParticleContainer()->ShowProjection(1,0)->GetMaximum(); + //draw some example plots.... - fCFManager->GetParticleContainer()->ShowProjection(0,0)->GetYaxis()->SetRangeUser(0,max1*1.2); - fCFManager->GetParticleContainer()->ShowProjection(0,1)->GetYaxis()->SetRangeUser(0,max1*1.2); - fCFManager->GetParticleContainer()->ShowProjection(0,2)->GetYaxis()->SetRangeUser(0,max1*1.2); - fCFManager->GetParticleContainer()->ShowProjection(0,3)->GetYaxis()->SetRangeUser(0,max1*1.2); + AliCFContainer *cont= dynamic_cast (GetOutputData(2)); - fCFManager->GetParticleContainer()->ShowProjection(1,0)->GetYaxis()->SetRangeUser(0,max2*1.2); - fCFManager->GetParticleContainer()->ShowProjection(1,1)->GetYaxis()->SetRangeUser(0,max2*1.2); - fCFManager->GetParticleContainer()->ShowProjection(1,2)->GetYaxis()->SetRangeUser(0,max2*1.2); - fCFManager->GetParticleContainer()->ShowProjection(1,3)->GetYaxis()->SetRangeUser(0,max2*1.2); + TH1D* h00 = cont->ShowProjection(0,0) ; + TH1D* h01 = cont->ShowProjection(0,1) ; + TH1D* h02 = cont->ShowProjection(0,2) ; + TH1D* h03 = cont->ShowProjection(0,3) ; - fCFManager->GetParticleContainer()->ShowProjection(0,0)->SetMarkerStyle(23) ; - fCFManager->GetParticleContainer()->ShowProjection(0,1)->SetMarkerStyle(24) ; - fCFManager->GetParticleContainer()->ShowProjection(0,2)->SetMarkerStyle(25) ; - fCFManager->GetParticleContainer()->ShowProjection(0,3)->SetMarkerStyle(26) ; + TH1D* h10 = cont->ShowProjection(1,0) ; + TH1D* h11 = cont->ShowProjection(1,1) ; + TH1D* h12 = cont->ShowProjection(1,2) ; + TH1D* h13 = cont->ShowProjection(1,3) ; - fCFManager->GetParticleContainer()->ShowProjection(1,0)->SetMarkerStyle(23) ; - fCFManager->GetParticleContainer()->ShowProjection(1,1)->SetMarkerStyle(24) ; - fCFManager->GetParticleContainer()->ShowProjection(1,2)->SetMarkerStyle(25) ; - fCFManager->GetParticleContainer()->ShowProjection(1,3)->SetMarkerStyle(26) ; + Double_t max1 = h00->GetMaximum(); + Double_t max2 = h10->GetMaximum(); + + h00->GetYaxis()->SetRangeUser(0,max1*1.2); + h01->GetYaxis()->SetRangeUser(0,max1*1.2); + h02->GetYaxis()->SetRangeUser(0,max1*1.2); + h03->GetYaxis()->SetRangeUser(0,max1*1.2); + + h10->GetYaxis()->SetRangeUser(0,max2*1.2); + h11->GetYaxis()->SetRangeUser(0,max2*1.2); + h12->GetYaxis()->SetRangeUser(0,max2*1.2); + h13->GetYaxis()->SetRangeUser(0,max2*1.2); + + h00->SetMarkerStyle(23) ; + h01->SetMarkerStyle(24) ; + h02->SetMarkerStyle(25) ; + h03->SetMarkerStyle(26) ; + + h10->SetMarkerStyle(23) ; + h11->SetMarkerStyle(24) ; + h12->SetMarkerStyle(25) ; + h13->SetMarkerStyle(26) ; TCanvas * c =new TCanvas("c","",1400,800); c->Divide(4,2); c->cd(1); - fCFManager->GetParticleContainer()->ShowProjection(0,0)->Draw("p"); + h00->Draw("p"); c->cd(2); - fCFManager->GetParticleContainer()->ShowProjection(0,1)->Draw("p"); + h01->Draw("p"); c->cd(3); - fCFManager->GetParticleContainer()->ShowProjection(0,2)->Draw("p"); + h02->Draw("p"); c->cd(4); - fCFManager->GetParticleContainer()->ShowProjection(0,3)->Draw("p"); + h03->Draw("p"); c->cd(5); - fCFManager->GetParticleContainer()->ShowProjection(1,0)->Draw("p"); + h10->Draw("p"); c->cd(6); - fCFManager->GetParticleContainer()->ShowProjection(1,1)->Draw("p"); + h11->Draw("p"); c->cd(7); - fCFManager->GetParticleContainer()->ShowProjection(1,2)->Draw("p"); + h12->Draw("p"); c->cd(8); - fCFManager->GetParticleContainer()->ShowProjection(1,3)->Draw("p"); + h13->Draw("p"); c->SaveAs("plots.eps"); - - delete fHistEventsProcessed ; } //___________________________________________________________________________ diff --git a/CORRFW/test/AliCFSingleTrackTask.C b/CORRFW/test/AliCFSingleTrackTask.C index ebcb7f026e7..d3ccf5409a3 100644 --- a/CORRFW/test/AliCFSingleTrackTask.C +++ b/CORRFW/test/AliCFSingleTrackTask.C @@ -268,6 +268,6 @@ void Load() { gSystem->Load("libCORRFW.so") ; //compile online the task class - gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ALICE_ROOT/TPC -I$ROOTSYS/include"); + gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ROOTSYS/include"); gROOT->LoadMacro("./AliCFSingleTrackTask.cxx+"); } diff --git a/CORRFW/test/AliCFV0Task.C b/CORRFW/test/AliCFV0Task.C index eac6970ba41..a34a5a88a65 100644 --- a/CORRFW/test/AliCFV0Task.C +++ b/CORRFW/test/AliCFV0Task.C @@ -14,6 +14,7 @@ const Int_t chargeV0 = 0 ; Bool_t AliCFV0Task( const Bool_t useGrid = 1, + const Bool_t readAOD = 0, const char * kTagXMLFile="wn.xml", // XML file containing tags ) { @@ -48,9 +49,15 @@ Bool_t AliCFV0Task( analysisChain = tagAna->QueryTags(runCuts,lhcCuts,detCuts,eventCuts); } else {// local data - analysisChain = new TChain("esdTree"); //here put your input data path - analysisChain->Add("AliESDs.root"); + if (readAOD) { + analysisChain = new TChain("aodTree"); + analysisChain->Add("AliAOD.root"); + } + else { + analysisChain = new TChain("esdTree"); + analysisChain->Add("AliESDs.root"); + } } @@ -105,7 +112,9 @@ Bool_t AliCFV0Task( recKineCuts->SetChargeRec(0); AliCFPairQualityCuts *recQualityCuts = new AliCFPairQualityCuts("recQualityCuts","V0 rec-level quality cuts"); - recQualityCuts->SetMinNClusterTPC(minclustersTPC,minclustersTPC); + if (!readAOD) recQualityCuts->SetMinNClusterTPC(minclustersTPC,minclustersTPC); + recQualityCuts->SetStatus(AliESDtrack::kTPCrefit & AliESDtrack::kITSrefit, + AliESDtrack::kTPCrefit & AliESDtrack::kITSrefit) ; AliCFV0TopoCuts *recTopoCuts = new AliCFV0TopoCuts("recTopoCuts","V0 Topological Cuts"); recTopoCuts->SetMaxDcaDaughters(0.1); @@ -126,8 +135,10 @@ Bool_t AliCFV0Task( 0.7152 , 0.0442, 0.0733 }; - cutPID->SetPriors(prior_pbpb); + cutPID->SetPriors(prior_pp); cutPID->SetDetectors("ITS TPC TRD","ITS TPC TRD"); + if (readAOD) cutPID->SetAODmode(kTRUE); + else cutPID->SetAODmode(kFALSE); cutPID->SetProbabilityCut(0,0); switch(TMath::Abs(PDG)) { case 310 : cutPID->SetParticleType(AliPID::kPion , kTRUE, AliPID::kPion , kTRUE); break; @@ -170,7 +181,7 @@ Bool_t AliCFV0Task( // create the task AliCFV0Task *task = new AliCFV0Task("AliCFV0Task"); task->SetCFManager(man); //here is set the CF manager - task->SetRebuildV0s(kTRUE); + if (!readAOD) task->SetRebuildV0s(kTRUE); task->SetV0PDG(PDG); //SETUP THE ANALYSIS MANAGER TO READ INPUT CHAIN AND WRITE DESIRED OUTPUTS @@ -182,10 +193,13 @@ Bool_t AliCFV0Task( else mgr->SetAnalysisType(AliAnalysisManager::kLocalAnalysis); AliMCEventHandler* mcHandler = new AliMCEventHandler(); - AliESDInputHandler* esdHandler = new AliESDInputHandler(); mgr->SetMCtruthEventHandler(mcHandler); - mgr->SetInputEventHandler(esdHandler); + AliInputEventHandler* dataHandler ; + + if (readAOD) dataHandler = new AliAODInputHandler(); + else dataHandler = new AliESDInputHandler(); + mgr->SetInputEventHandler(dataHandler); // Create and connect containers for input/output @@ -224,41 +238,13 @@ Bool_t AliCFV0Task( } void Load(Bool_t useGrid) { - //remove this file which can cause problems - if (!useGrid) - gSystem->Exec("rm $ALICE_ROOT/ANALYSIS/AliAnalysisSelector_cxx.so"); - + //load the required aliroot libraries gSystem->Load("libANALYSIS") ; gSystem->Load("libANALYSISalice") ; - - gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ROOTSYS/include -I$ALICE_ROOT/PWG2/RESONANCES -I$ALICE_ROOT/CORRFW"); + gSystem->Load("libCORRFW.so") ; - if (!useGrid) { - //load local CF library - gSystem->Load("libCORRFW"); - } - else { - //compile online the task class - gSystem->Exec("tar zxvf cf_classes.tgz"); - gROOT->LoadMacro("./AliCFFrame.cxx+"); - gROOT->LoadMacro("./AliCFVGrid.cxx+"); - gROOT->LoadMacro("./AliCFGrid.cxx+"); - gROOT->LoadMacro("./AliCFGridSparse.cxx+"); - gROOT->LoadMacro("./AliCFContainer.cxx+"); - gROOT->LoadMacro("./AliCFCutBase.cxx+"); - gROOT->LoadMacro("./AliCFTrackKineCuts.cxx+"); - gROOT->LoadMacro("./AliCFTrackQualityCuts.cxx+"); - gROOT->LoadMacro("./AliCFParticleGenCuts.cxx+"); - gROOT->LoadMacro("./AliCFAcceptanceCuts.cxx+"); - gROOT->LoadMacro("./AliCFTrackCutPid.cxx+"); - gROOT->LoadMacro("./AliCFManager.cxx+"); - gSystem->Exec("tar zxvf cf_V0.tgz"); - gROOT->LoadMacro("./AliCFPair.cxx+"); - gROOT->LoadMacro("./AliCFPairAcceptanceCuts.cxx+"); - gROOT->LoadMacro("./AliCFPairQualityCuts.cxx+"); - gROOT->LoadMacro("./AliCFPairPidCut.cxx+"); - gROOT->LoadMacro("./AliCFV0TopoCuts.cxx+"); - } + //compile online the task class + gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ROOTSYS/include"); gROOT->LoadMacro("./AliCFV0Task.cxx+"); } diff --git a/CORRFW/test/AliCFV0Task.cxx b/CORRFW/test/AliCFV0Task.cxx index 67a76eedf42..8653fd2342c 100644 --- a/CORRFW/test/AliCFV0Task.cxx +++ b/CORRFW/test/AliCFV0Task.cxx @@ -41,8 +41,10 @@ #include "AliV0vertexer.h" #include "AliCFPair.h" #include "AliESDEvent.h" +#include "AliAODEvent.h" #include "TChain.h" #include "AliCFParticleGenCuts.h" +#include "AliAODv0.h" //__________________________________________________________________________ AliCFV0Task::AliCFV0Task() : @@ -124,15 +126,21 @@ void AliCFV0Task::UserExec(Option_t *) // Info("UserExec","") ; - AliESDEvent* fESD = dynamic_cast(fInputEvent); - if (!fESD) { - Error("UserExec","NO ESD FOUND!"); + if (!fInputEvent) { + Error("UserExec","NO EVENT FOUND!"); return; } if (!fMCEvent) Error("UserExec","NO MC INFO FOUND!"); fCFManager->SetEventInfo(fMCEvent); + Bool_t isESDEvent = strcmp(fInputEvent->ClassName(),"AliESDEvent") == 0 ? kTRUE : kFALSE ; + Bool_t isAODEvent = strcmp(fInputEvent->ClassName(),"AliAODEvent") == 0 ? kTRUE : kFALSE ; + + AliESDEvent* fESD = dynamic_cast(fInputEvent); + AliAODEvent* fAOD = dynamic_cast(fInputEvent); + + // MC-event selection Double_t containerInput[2] ; @@ -156,23 +164,49 @@ void AliCFV0Task::UserExec(Option_t *) //Now go to rec level - Info("UserExec","Looping on ESD event"); - - if (fRebuildV0s) RebuildV0s(fESD) ; - Info("UserExec","There are %d V0s in event",fESD->GetNumberOfV0s()); - for (Int_t iV0 = 0; iV0GetNumberOfV0s(); iV0++) { - - AliESDv0* esdV0 = fESD->GetV0(iV0); - - //check if mother reconstructed V0 can be associated to a MC V0 - Int_t labMCV0 = IsMcV0(esdV0) ; - if (labMCV0 <0) continue; - - esdV0->ChangeMassHypothesis(fV0PDG); //important to do that before entering the cut check + Info("UserExec","Looping on %s",fInputEvent->ClassName()); + + //SET THE ESD AS EVENT INFO IN RECONSTRUCTION CUTS + TObjArray* fCutsReco = fCFManager->GetParticleCutsList(AliCFManager::kPartRecCuts); + TObjArrayIter iter1(fCutsReco); + AliCFCutBase *cut = 0; + while ( (cut = (AliCFCutBase*)iter1.Next()) ) { + cut->SetEvtInfo(fInputEvent); + } - AliCFPair pair(esdV0,fESD); - if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,&pair)) continue; + if (isESDEvent && fRebuildV0s) RebuildV0s(fESD) ; + Info("UserExec","There are %d V0s in event",fInputEvent->GetNumberOfV0s()); + + AliESDv0* esdV0 = 0x0; + AliAODv0* aodV0 = 0x0; + AliCFPair* pair = 0x0; + Int_t labMCV0 = 0; + + for (Int_t iV0 = 0; iV0GetNumberOfV0s(); iV0++) { + + if (isESDEvent) { + esdV0 = fESD->GetV0(iV0); + //check if mother reconstructed V0 can be associated to a MC V0 + labMCV0 = IsMcV0(esdV0); + if (labMCV0 <0) continue; + pair = new AliCFPair(esdV0,fESD); + } + else if (isAODEvent) { + aodV0 = fAOD->GetV0(iV0); + labMCV0 = IsMcV0(aodV0); + if (labMCV0 <0) continue; + pair = new AliCFPair(aodV0); + } + else { + Error("UserExec","Error: input data file is not an ESD nor an AOD"); + return; + } + + pair->SetV0PDG(fV0PDG); + pair->SetLabel(labMCV0); + if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pair)) continue; + //check if associated MC v0 passes the cuts AliMCParticle* mcV0 = fMCEvent->GetTrack(labMCV0); if (!mcV0) continue; @@ -180,17 +214,28 @@ void AliCFV0Task::UserExec(Option_t *) if (!fCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,mcV0)) continue; //fill the container - Double_t mom[3]; - esdV0->GetPxPyPz(mom[0],mom[1],mom[2]); - Double_t pt2 = mom[0]*mom[0]+mom[1]*mom[1] ; - Double_t pt = TMath::Sqrt(pt2); - Double_t energy = TMath::Sqrt(pt2 + mom[2]*mom[2] + TMath::Power(TDatabasePDG::Instance()->GetParticle(fV0PDG)->Mass(),2)); - + Double_t pt, rapidity; + + if (isESDEvent) { + Double_t mom[3]; + esdV0->GetPxPyPz(mom[0],mom[1],mom[2]); + Double_t pt2 = mom[0]*mom[0]+mom[1]*mom[1] ; + pt = TMath::Sqrt(pt2); + Double_t energy = TMath::Sqrt(pt2 + mom[2]*mom[2] + TMath::Power(TDatabasePDG::Instance()->GetParticle(fV0PDG)->Mass(),2)); + rapidity = GetRapidity(energy,mom[2]) ; + } + else { + pt = aodV0->Pt(); + rapidity = aodV0->Y(fV0PDG); + } + containerInput[0] = pt ; - containerInput[1] = GetRapidity(energy,mom[2]); + containerInput[1] = rapidity ; fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed) ; - if (!fCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,&pair)) continue ; + if (!fCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pair)) continue ; fCFManager->GetParticleContainer()->Fill(containerInput,kStepSelected); + + delete pair; } fHistEventsProcessed->Fill(0); @@ -325,6 +370,10 @@ Int_t AliCFV0Task::GetV0Label(UInt_t lab1, UInt_t lab2) const { //___________________________________________________________________________ Double_t AliCFV0Task::GetRapidity(Double_t energy, Double_t pz) { + // + // calculates rapidity, checking energy is larger than pz, otherwise returns 999. + // + if (energy == pz || energy == -pz) { printf("GetRapidity : ERROR : rapidity for 4-vector with E = Pz -- infinite result"); return 999; @@ -339,6 +388,11 @@ Double_t AliCFV0Task::GetRapidity(Double_t energy, Double_t pz) { //___________________________________________________________________________ Int_t AliCFV0Task::IsMcV0(AliESDv0* v0) const { + // + // check if the passed V0 is associated to a MC one, + // and returns the corresponding geant label. + // returns -1 if the V0 is fake (i.e. label<0). + // AliESDEvent* fESD = dynamic_cast(fInputEvent); @@ -357,6 +411,29 @@ Int_t AliCFV0Task::IsMcV0(AliESDv0* v0) const { return GetV0Label((UInt_t)nlab,(UInt_t)plab) ; } +//___________________________________________________________________________ +Int_t AliCFV0Task::IsMcV0(AliAODv0* v0) const { + // + // check if the passed V0 is associated to a MC one, + // and returns the corresponding geant label. + // returns -1 if the V0 is fake (i.e. label<0). + // + + AliAODVertex * vtx = v0->GetSecondaryVtx(); + + AliAODTrack *nTrack = (AliAODTrack*)vtx->GetDaughter(1); //neg is filled after pos in AliAnalysisTaskStrange + AliAODTrack *pTrack = (AliAODTrack*)vtx->GetDaughter(0); + + if (!nTrack || !pTrack) return -1 ; + + Int_t nlab = nTrack->GetLabel() ; + Int_t plab = pTrack->GetLabel() ; + + if (nlab <0 || plab <0) return -1 ; + + return GetV0Label((UInt_t)nlab,(UInt_t)plab) ; +} + //___________________________________________________________________________ void AliCFV0Task::RebuildV0s(AliESDEvent* fESD) { diff --git a/CORRFW/test/AliCFV0Task.h b/CORRFW/test/AliCFV0Task.h index c96148bf47b..005b0b8703d 100644 --- a/CORRFW/test/AliCFV0Task.h +++ b/CORRFW/test/AliCFV0Task.h @@ -29,6 +29,7 @@ class AliMCEventHandler; class AliCFManager; class AliESDv0; class AliESDEvent; +class AliAODv0; class AliCFV0Task : public AliAnalysisTaskSE { public: @@ -58,6 +59,7 @@ class AliCFV0Task : public AliAnalysisTaskSE { void SetRebuildV0s(Bool_t flag) {fRebuildV0s = flag;} // setter for V0 on-the-fly reconstruction void SetV0PDG(Int_t code) {fV0PDG = code; } // defines the PDG code of searched V0's Int_t IsMcV0(AliESDv0*) const ; // checks if the AliESDv0 can be associated, returns mother label + Int_t IsMcV0(AliAODv0*) const ; // checks if the AliAODv0 can be associated, returns mother label Int_t GetV0Label(UInt_t, UInt_t) const ; // returns label of V0 given 2 daughter labels static Double_t GetRapidity(Double_t, Double_t) ; // returns the rapidity of the V0 (assuming PDG code) @@ -67,7 +69,7 @@ class AliCFV0Task : public AliAnalysisTaskSE { Bool_t fRebuildV0s; // flag for on-the-fly V0 reconstruction Int_t fV0PDG; // PDG code of searched V0's - AliCFManager *fCFManager ; // pointer to the CF manager + AliCFManager *fCFManager ; // pointer to the CF manager // Histograms //Number of events -- 2.39.3