--- /dev/null
+AliAnalysisTaskSpectraAllChNanoAOD* AddTaskSpectraAllChAOD(
+ Bool_t mc=kFALSE,
+ Double_t CentCutMin=0,
+ Double_t CentCutMax=100,
+ Double_t QvecCutMin=0,
+ Double_t QvecCutMax=100,
+ Double_t EtaMin=-0.8,
+ Double_t EtaMax=0.8,
+ Double_t pt=50.,
+ Double_t ptTofMatch=.6,
+ UInt_t trkbit=1,
+ Double_t DCA=100000,
+ UInt_t minNclsTPC=70,
+ Double_t nsigmacut=5.,
+ Int_t PIDtype=3,
+ TString opt=""){
+
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+ if (!mgr)
+ {
+ ::Error("AddAliAnalysisTaskSpectraAllChAOD", "No analysis manager to connect to.");
+ return NULL;
+ }
+
+ // Check the analysis type using the event handlers connected to the analysis manager.
+ //==============================================================================
+ if (!mgr->GetInputEventHandler())
+ {
+ ::Error("AliAnalysisTaskSpectraAllChAOD", "This task requires an input event handler");
+ return NULL;
+ }
+
+ TString type = mgr->GetInputEventHandler()->GetDataType(); // can be "ESD" or "AOD"
+ if(type.Contains("ESD"))
+ {
+ ::Error("AliAnalysisTaskSpectraAllChAOD", "This task requires to run on AOD");
+ return NULL;
+ }
+
+ AliSpectraAODTrackCuts * trcuts = new AliSpectraAODTrackCuts(Form("TrackCuts%s",opt.Data()));
+ trcuts->SetDCA(DCA);
+ trcuts->SetTrackBits(trkbit);
+ trcuts->SetPt(pt);
+ trcuts->SetPtTOFMatching(ptTofMatch);
+ trcuts->SetEta(EtaMin,EtaMax);
+ trcuts->SetMinTPCcls(minNclsTPC);
+ trcuts->PrintCuts();
+
+ AliSpectraAODEventCuts * evcuts = new AliSpectraAODEventCuts(Form("EventCuts%s",opt.Data()));
+ evcuts->SetQVectorCut(QvecCutMin,QvecCutMax);
+ evcuts->SetCentralityCutMax(CentCutMax);
+ evcuts->SetCentralityCutMin(CentCutMin);
+ if(mc==1)evcuts->SetIsMC(kTRUE);
+ evcuts->PrintCuts();
+
+ AliHelperPID *pid=new AliHelperPID();
+ pid->SetName(Form("HelperPID%s",opt.Data()));
+ pid->SetNSigmaCut(nsigmacut);
+ pid->SetPIDType(PIDtype);
+ if(PIDtype==3){
+ AliPIDCombined *pidc=new AliPIDCombined();
+ pidc->SetDefaultTPCPriors();
+ pid->SetPIDCombined(pidc);
+ }
+
+ AliAnalysisTaskSpectraAllChNanoAOD *task = new AliAnalysisTaskSpectraAllChNanoAOD(Form("TaskAODSpectraCent%.0fto%.0f_QVec%.1fto%.1f_Eta%.1fto%.1f_TrBit%d%s",
+ CentCutMin,
+ CentCutMax,
+ QvecCutMin,
+ QvecCutMax,
+ EtaMin,
+ EtaMax,
+ trkbit,
+ opt.Data()));
+ task->SetEventCuts(evcuts);
+ task->SetTrackCuts(trcuts);
+ task->SetHelperPID(pid);
+ if(mc==1)task->SetIsMC(kTRUE);
+
+ TString outputFileName = AliAnalysisManager::GetCommonFileName();
+
+ TString typeofdata=mc?"MC":"Data";
+
+ outputFileName += Form(":SpectraESE_%s%s",typeofdata.Data(),opt.Data());
+
+ cout<<"outputFileName: "<<outputFileName<<endl;
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ AliAnalysisDataContainer *coutputpt1 = mgr->CreateContainer(Form("chist%s",opt.Data()), TList::Class(), AliAnalysisManager::kOutputContainer,outputFileName);
+ AliAnalysisDataContainer *coutputpt2 = mgr->CreateContainer(Form("cvcut%s",opt.Data()), AliSpectraAODEventCuts::Class(), AliAnalysisManager::kOutputContainer,outputFileName);
+ AliAnalysisDataContainer *coutputpt3 = mgr->CreateContainer(Form("ctcut%s",opt.Data()), AliSpectraAODTrackCuts::Class(), AliAnalysisManager::kOutputContainer, outputFileName);
+ AliAnalysisDataContainer *coutputpt4 = mgr->CreateContainer(Form("cpid%s",opt.Data()), AliHelperPID::Class(), AliAnalysisManager::kOutputContainer, outputFileName);
+ mgr->AddTask(task);
+
+ mgr->ConnectInput(task, 0, cinput);
+ mgr->ConnectOutput(task, 1, coutputpt1);
+ mgr->ConnectOutput(task, 2, coutputpt2);
+ mgr->ConnectOutput(task, 3, coutputpt3);
+ mgr->ConnectOutput(task, 4, coutputpt4);
+
+ return task;
+}
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+//-----------------------------------------------------------------
+// AliAnalysisTaskSpectraAllChNanoAOD class
+// NanoAOD edition: this is only meant to test the developement
+// version of NanoAODs
+//
+//-----------------------------------------------------------------
+
+#include "TChain.h"
+#include "TTree.h"
+#include "TLegend.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "THnSparse.h"
+#include "TCanvas.h"
+#include "AliAnalysisTask.h"
+#include "AliAODTrack.h"
+#include "AliAODMCParticle.h"
+#include "AliVParticle.h"
+#include "AliAODEvent.h"
+#include "AliAODInputHandler.h"
+#include "AliAnalysisTaskSpectraAllChNanoAOD.h"
+#include "AliAnalysisTaskESDfilter.h"
+#include "AliAnalysisDataContainer.h"
+#include "AliSpectraAODTrackCuts.h"
+#include "AliSpectraAODEventCuts.h"
+#include "AliHelperPID.h"
+#include "AliPIDCombined.h"
+#include "AliCentrality.h"
+#include "TProof.h"
+#include "AliVEvent.h"
+#include "AliStack.h"
+#include <TMCProcess.h>
+
+#include <iostream>
+#include "AliNanoAODHeader.h"
+#include "AliNanoAODTrack.h"
+
+using namespace AliHelperPIDNameSpace;
+using namespace std;
+
+ClassImp(AliAnalysisTaskSpectraAllChNanoAOD)
+
+//________________________________________________________________________
+AliAnalysisTaskSpectraAllChNanoAOD::AliAnalysisTaskSpectraAllChNanoAOD(const char *name) : AliAnalysisTaskSE(name),
+ fAOD(0x0),
+ fTrackCuts(0x0),
+ fEventCuts(0x0),
+ fHelperPID(0x0),
+ fIsMC(0),
+ fDoDoubleCounting(0),
+ fFillOnlyEvents(0),
+ fCharge(0),
+ fVZEROside(0),
+ fOutput(0x0),
+ fnCentBins(20),
+ fnQvecBins(40),
+ fnNchBins(200)
+{
+ // Default constructor
+ DefineInput(0, TChain::Class());
+ DefineOutput(1, TList::Class());
+ DefineOutput(2, AliSpectraAODEventCuts::Class());
+ DefineOutput(3, AliSpectraAODTrackCuts::Class());
+ DefineOutput(4, AliHelperPID::Class());
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskSpectraAllChNanoAOD::UserCreateOutputObjects()
+{
+ // create output objects
+ fOutput=new TList();
+ fOutput->SetOwner();
+ fOutput->SetName("fOutput");
+
+ if (!fTrackCuts) AliFatal("Track Cuts should be set in the steering macro");
+ if (!fEventCuts) AliFatal("Event Cuts should be set in the steering macro");
+ if (!fHelperPID) AliFatal("HelperPID should be set in the steering macro");
+
+ // binning common to all the THn
+ const Double_t ptBins[] = {0.20,0.30,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.4,2.8,3.2,3.6,4.0,5.0,6.0,7.0,8.0,9.0,10.,12.,15.};
+ const Int_t nptBins=26;
+
+ //dimensions of THnSparse for tracks
+ const Int_t nvartrk=8;
+ // pt cent Q vec IDrec IDgen isph iswd y
+ Int_t binsHistRealTrk[nvartrk] = { nptBins, fnCentBins, fnQvecBins, 4, 3, 2, 2, 2};
+ Double_t xminHistRealTrk[nvartrk] = { 0., 0., 0., -.5, -0.5, -0.5, -0.5, -0.5};
+ Double_t xmaxHistRealTrk[nvartrk] = { 10., 100., 8., 3.5, 2.5, 1.5, 1.5, 0.5};
+ THnSparseF* NSparseHistTrk = new THnSparseF("NSparseHistTrk","NSparseHistTrk",nvartrk,binsHistRealTrk,xminHistRealTrk,xmaxHistRealTrk);
+ NSparseHistTrk->GetAxis(0)->SetTitle("#it{p}_{T,rec}");
+ NSparseHistTrk->GetAxis(0)->SetName("pT_rec");
+ NSparseHistTrk->SetBinEdges(0,ptBins);
+ NSparseHistTrk->GetAxis(1)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
+ NSparseHistTrk->GetAxis(1)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
+ NSparseHistTrk->GetAxis(2)->SetTitle("Q vec");
+ NSparseHistTrk->GetAxis(2)->SetName("Q_vec");
+ NSparseHistTrk->GetAxis(3)->SetTitle("ID rec");
+ NSparseHistTrk->GetAxis(3)->SetName("ID_rec");
+ NSparseHistTrk->GetAxis(4)->SetTitle("ID gen");
+ NSparseHistTrk->GetAxis(4)->SetName("ID_gen");
+ NSparseHistTrk->GetAxis(5)->SetTitle("isph");
+ NSparseHistTrk->GetAxis(5)->SetName("isph");
+ NSparseHistTrk->GetAxis(6)->SetTitle("iswd");
+ NSparseHistTrk->GetAxis(6)->SetName("iswd");
+ NSparseHistTrk->GetAxis(7)->SetTitle("y");
+ NSparseHistTrk->GetAxis(7)->SetName("y");
+ fOutput->Add(NSparseHistTrk);
+
+ //dimensions of THnSparse for stack
+ const Int_t nvarst=5;
+ // pt cent IDgen isph y
+ Int_t binsHistRealSt[nvarst] = { nptBins, fnCentBins, 3, 2, 2};
+ Double_t xminHistRealSt[nvarst] = { 0., 0., -0.5, -0.5, -0.5};
+ Double_t xmaxHistRealSt[nvarst] = { 10., 100., 2.5, 1.5, 0.5};
+ THnSparseF* NSparseHistSt = new THnSparseF("NSparseHistSt","NSparseHistSt",nvarst,binsHistRealSt,xminHistRealSt,xmaxHistRealSt);
+ NSparseHistSt->GetAxis(0)->SetTitle("#it{p}_{T,gen}");
+ NSparseHistSt->SetBinEdges(0,ptBins);
+ NSparseHistSt->GetAxis(0)->SetName("pT_rec");
+ NSparseHistSt->GetAxis(1)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
+ NSparseHistSt->GetAxis(1)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
+ NSparseHistSt->GetAxis(2)->SetTitle("ID gen");
+ NSparseHistSt->GetAxis(2)->SetName("ID_gen");
+ NSparseHistSt->GetAxis(3)->SetTitle("isph");
+ NSparseHistSt->GetAxis(3)->SetName("isph");
+ NSparseHistSt->GetAxis(4)->SetTitle("y");
+ NSparseHistSt->GetAxis(4)->SetName("y");
+ fOutput->Add(NSparseHistSt);
+
+ //dimensions of THnSparse for the normalization
+ const Int_t nvarev=3;
+ // cent Q vec Nch
+ Int_t binsHistRealEv[nvarev] = { fnCentBins, fnQvecBins, fnNchBins};
+ Double_t xminHistRealEv[nvarev] = { 0., 0., 0.};
+ Double_t xmaxHistRealEv[nvarev] = { 100., 8., 2000.};
+ THnSparseF* NSparseHistEv = new THnSparseF("NSparseHistEv","NSparseHistEv",nvarev,binsHistRealEv,xminHistRealEv,xmaxHistRealEv);
+ NSparseHistEv->GetAxis(0)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
+ NSparseHistEv->GetAxis(0)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
+ NSparseHistEv->GetAxis(1)->SetTitle("Q vec");
+ NSparseHistEv->GetAxis(1)->SetName("Q_vec");
+ NSparseHistEv->GetAxis(2)->SetTitle("N charged");
+ NSparseHistEv->GetAxis(2)->SetName("N_ch");
+ fOutput->Add(NSparseHistEv);
+
+ PostData(1, fOutput );
+ PostData(2, fEventCuts);
+ PostData(3, fTrackCuts);
+ PostData(4, fHelperPID);
+}
+
+//________________________________________________________________________
+
+void AliAnalysisTaskSpectraAllChNanoAOD::UserExec(Option_t *)
+{
+ //Printf("An event");
+ // main event loop
+ fAOD = dynamic_cast<AliAODEvent*>(fInputEvent);
+ if (!fAOD) {
+ AliWarning("ERROR: AliAODEvent not available \n");
+ return;
+ }
+
+ if (strcmp(fAOD->ClassName(), "AliAODEvent"))
+ {
+ AliFatal("Not processing AODs");
+ }
+
+ AliNanoAODHeader * headNano = dynamic_cast<AliNanoAODHeader*>((TObject*)fAOD->GetHeader());
+
+ Bool_t isNano = (headNano != 0);
+
+ if(!isNano) {
+ if(!fEventCuts->IsSelected(fAOD,fTrackCuts))return;//event selection
+ }
+ //Default TPC priors
+ if(fHelperPID->GetPIDType()==kBayes)fHelperPID->GetPIDCombined()->SetDefaultTPCPriors();//FIXME maybe this can go in the UserCreateOutputObject?
+
+ Double_t Qvec=0.;//in case of MC we save space in the memory
+ if(!fIsMC){
+ if(!isNano) {
+ if(fVZEROside==0)Qvec=fEventCuts->GetqV0A();
+ else if (fVZEROside==1)Qvec=fEventCuts->GetqV0C();
+ } else {
+
+ if(fVZEROside==0) Qvec=headNano->GetVar(0); // FIXME: we need proper getters here
+ else if (fVZEROside==1)Qvec=headNano->GetVar(1);
+
+ }
+ }
+
+ Double_t Cent=isNano ? headNano->GetVar(2) : fEventCuts->GetCent();
+
+ // First do MC to fill up the MC particle array
+ TClonesArray *arrayMC = 0;
+ if (fIsMC)
+ {
+ arrayMC = (TClonesArray*) fAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName());
+ if (!arrayMC) {
+ AliFatal("Error: MC particles branch not found!\n");
+ }
+ Int_t nMC = arrayMC->GetEntries();
+ for (Int_t iMC = 0; iMC < nMC; iMC++)
+ {
+ AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(iMC);
+ if(!partMC->Charge()) continue;//Skip neutrals
+ if(fCharge != 0 && partMC->Charge()*fCharge < 0.) continue;//if fCharge != 0 only select fCharge
+
+ if(partMC->Eta() < fTrackCuts->GetEtaMin() || partMC->Eta() > fTrackCuts->GetEtaMax())continue;//ETA CUT ON GENERATED!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+ //pt cent IDgen isph y
+ Double_t varSt[5];
+ varSt[0]=partMC->Pt();
+ varSt[1]=Cent;
+ varSt[2]=(Double_t)fHelperPID->GetParticleSpecies(partMC);
+ varSt[3]=(Double_t)partMC->IsPhysicalPrimary();
+ varSt[4]=partMC->Y();
+ ((THnSparseF*)fOutput->FindObject("NSparseHistSt"))->Fill(varSt);//stack loop
+
+ //Printf("a particle");
+
+ }
+ }
+
+ //main loop on tracks
+
+ Int_t Nch = 0.;
+
+ for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
+ AliVTrack* track = (AliVTrack*) fAOD->GetTrack(iTracks);
+ if(fCharge != 0 && track->Charge() != fCharge) continue;//if fCharge != 0 only select fCharge
+ if(!isNano) {
+ if (!fTrackCuts->IsSelected((AliAODTrack*)track,kTRUE)) continue; //track selection (rapidity selection NOT in the standard cuts)
+ }
+
+ if(!fFillOnlyEvents){
+ Int_t IDrec=isNano ? GetNanoTrackID (track) : fHelperPID->GetParticleSpecies(track,kTRUE);//id from detector
+ Double_t y= 0;
+ if(isNano) y = ((AliNanoAODTrack*)track)->Y(fHelperPID->GetMass((AliHelperParticleSpecies_t)IDrec));
+ else y = ((AliAODTrack*)track)->Y(fHelperPID->GetMass((AliHelperParticleSpecies_t)IDrec));
+ Int_t IDgen=kSpUndefined;//set if MC
+ Int_t isph=-999;
+ Int_t iswd=-999;
+
+ if (arrayMC) {
+ AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(TMath::Abs(track->GetLabel()));
+ if (!partMC) {
+ AliError("Cannot get MC particle");
+ continue;
+ }
+ IDgen=fHelperPID->GetParticleSpecies(partMC);
+ isph=partMC->IsPhysicalPrimary();
+ iswd=partMC->IsSecondaryFromWeakDecay();//FIXME not working on old productions
+ }
+
+ //pt cent Q vec IDrec IDgen isph iswd y
+ Double_t varTrk[8];
+ varTrk[0]=track->Pt();
+ varTrk[1]=Cent;
+ varTrk[2]=Qvec;
+ varTrk[3]=(Double_t)IDrec;
+ varTrk[4]=(Double_t)IDgen;
+ varTrk[5]=(Double_t)isph;
+ varTrk[6]=(Double_t)iswd;
+ varTrk[7]=y;
+ ((THnSparseF*)fOutput->FindObject("NSparseHistTrk"))->Fill(varTrk);//track loop
+
+ //for nsigma PID fill double counting of ID
+ if(fHelperPID->GetPIDType()<kBayes && fDoDoubleCounting){//only nsigma
+ Bool_t *HasDC;
+ HasDC=fHelperPID->GetDoubleCounting(track,kTRUE);//get the array with double counting
+ for(Int_t ipart=0;ipart<kNSpecies;ipart++){
+ if(HasDC[ipart]==kTRUE){
+ varTrk[3]=(Double_t)ipart;
+ ((THnSparseF*)fOutput->FindObject("NSparseHistTrk"))->Fill(varTrk);//track loop
+ }
+ }
+ }
+
+ //fill all charged (3)
+ varTrk[3]=3.;
+ varTrk[4]=3.;
+ ((THnSparseF*)fOutput->FindObject("NSparseHistTrk"))->Fill(varTrk);//track loop
+ }//end if fFillOnlyEvents
+
+ //Printf("a track");
+
+ Nch++;
+ } // end loop on tracks
+
+ Double_t varEv[3];
+ varEv[0]=Cent;
+ varEv[1]=Qvec;
+ varEv[2]=Nch;
+ ((THnSparseF*)fOutput->FindObject("NSparseHistEv"))->Fill(varEv);//event loop
+
+ PostData(1, fOutput );
+ PostData(2, fEventCuts);
+ PostData(3, fTrackCuts);
+ PostData(4, fHelperPID);
+}
+
+//_________________________________________________________________
+void AliAnalysisTaskSpectraAllChNanoAOD::Terminate(Option_t *)
+{
+ // Terminate
+}
+
+
+Int_t AliAnalysisTaskSpectraAllChNanoAOD::GetNanoTrackID(AliVTrack * track) {
+ // Applies nsigma PID to nano tracks
+ AliNanoAODTrack * nanoTrack = dynamic_cast<AliNanoAODTrack*>(track);
+ if(!nanoTrack) AliFatal("Not a nano AOD track");
+
+ // Cache indexes
+ static const Int_t kcstNSigmaTPCPi = AliNanoAODTrackMapping::GetInstance()->GetVarIndex("cstNSigmaTPCPi");
+ static const Int_t kcstNSigmaTPCKa = AliNanoAODTrackMapping::GetInstance()->GetVarIndex("cstNSigmaTPCKa");
+ static const Int_t kcstNSigmaTPCPr = AliNanoAODTrackMapping::GetInstance()->GetVarIndex("cstNSigmaTPCPr");
+ static const Int_t kcstNSigmaTOFPi = AliNanoAODTrackMapping::GetInstance()->GetVarIndex("cstNSigmaTOFPi");
+ static const Int_t kcstNSigmaTOFKa = AliNanoAODTrackMapping::GetInstance()->GetVarIndex("cstNSigmaTOFKa");
+ static const Int_t kcstNSigmaTOFPr = AliNanoAODTrackMapping::GetInstance()->GetVarIndex("cstNSigmaTOFPr");
+
+ Double_t nSigmaPID = 3.0;
+
+
+
+ //get the identity of the particle with the minimum Nsigma
+ Double_t nsigmaPion=999., nsigmaKaon=999., nsigmaProton=999.;
+ if(nanoTrack->Pt() > fTrackCuts->GetPtTOFMatching()) {
+ nsigmaProton = TMath::Sqrt(nanoTrack->GetVar(kcstNSigmaTPCPr)*nanoTrack->GetVar(kcstNSigmaTPCPr)+nanoTrack->GetVar(kcstNSigmaTOFPr)*nanoTrack->GetVar(kcstNSigmaTOFPr));
+ nsigmaKaon = TMath::Sqrt(nanoTrack->GetVar(kcstNSigmaTPCKa)*nanoTrack->GetVar(kcstNSigmaTPCKa)+nanoTrack->GetVar(kcstNSigmaTOFKa)*nanoTrack->GetVar(kcstNSigmaTOFKa));
+ nsigmaPion = TMath::Sqrt(nanoTrack->GetVar(kcstNSigmaTPCPi)*nanoTrack->GetVar(kcstNSigmaTPCPi)+nanoTrack->GetVar(kcstNSigmaTOFPi)*nanoTrack->GetVar(kcstNSigmaTOFPi));
+ }
+ else {
+ nsigmaProton = TMath::Abs(nanoTrack->GetVar(kcstNSigmaTPCPr));
+ nsigmaKaon = TMath::Abs(nanoTrack->GetVar(kcstNSigmaTPCKa));
+ nsigmaPion = TMath::Abs(nanoTrack->GetVar(kcstNSigmaTPCPi));
+ }
+
+
+
+ // guess the particle based on the smaller nsigma (within nSigmaPID)
+ if( ( nsigmaKaon==nsigmaPion ) && ( nsigmaKaon==nsigmaProton )) return kSpUndefined;//if is the default value for the three
+
+ if( ( nsigmaKaon < nsigmaPion ) && ( nsigmaKaon < nsigmaProton ) && (nsigmaKaon < nSigmaPID)){
+ return kSpKaon;
+ }
+ if( ( nsigmaPion < nsigmaKaon ) && ( nsigmaPion < nsigmaProton ) && (nsigmaPion < nSigmaPID)){
+ return kSpPion;
+ }
+ if( ( nsigmaProton < nsigmaKaon ) && ( nsigmaProton < nsigmaPion ) && (nsigmaProton < nSigmaPID)){
+ return kSpProton;
+ }
+ // else, return undefined
+ return kSpUndefined;
+
+
+}
--- /dev/null
+#ifndef ALIANALYSISTASKSPECTRAALLCHNANOAOD_H
+#define ALIANALYSISTASKSPECTRAALLCHNANOAOD_H
+
+/* See cxx source for full Copyright notice */
+
+//-------------------------------------------------------------------------
+// AliAnalysisTaskSpectraAllChNanoAOD
+//
+//
+//
+//
+// Author: Leonardo Milano, CERN
+//-------------------------------------------------------------------------
+
+class AliAODEvent;
+class AliSpectraAODTrackCuts;
+class AliSpectraAODEventCuts;
+class AliHelperPID;
+
+#include "AliAnalysisTaskSE.h"
+
+class AliAnalysisTaskSpectraAllChNanoAOD : public AliAnalysisTaskSE
+{
+ public:
+ // constructors
+ AliAnalysisTaskSpectraAllChNanoAOD() : AliAnalysisTaskSE(),
+ fAOD(0x0),
+ fTrackCuts(0x0),
+ fEventCuts(0x0),
+ fHelperPID(0x0),
+ fIsMC(0),
+ fDoDoubleCounting(0),
+ fFillOnlyEvents(0),
+ fCharge(0),
+ fVZEROside(0),
+ fOutput(0x0),
+ fnCentBins(20),
+ fnQvecBins(40),
+ fnNchBins(200)
+ {}
+ AliAnalysisTaskSpectraAllChNanoAOD(const char *name);
+ virtual ~AliAnalysisTaskSpectraAllChNanoAOD() {
+ Printf("calling detructor of AliAnalysisTaskSpectraAllChNanoAOD - To be implemented");
+ }
+
+ void SetIsMC(Bool_t isMC = kFALSE) {fIsMC = isMC; };
+ Bool_t GetIsMC() const { return fIsMC;};
+
+ void SetDoDoubleCounting(Bool_t doDoubleCounting = kFALSE) {fDoDoubleCounting = doDoubleCounting; };
+ Bool_t GetDoDoubleCounting() const { return fDoDoubleCounting;};
+
+ void SetFillOnlyEvents(Bool_t fillOnlyEvents = kFALSE) {fFillOnlyEvents = fillOnlyEvents; };
+ Bool_t GetFillOnlyEvents() const { return fFillOnlyEvents;};
+
+ void SetCharge(Int_t charge = 0) {fCharge = charge; };
+ Int_t GetCharge() const { return fCharge;};
+
+ void SetVZEROside(Int_t side = 0) {fVZEROside = side; };
+ Int_t GetVZEROside() const { return fVZEROside;};
+
+ virtual void UserCreateOutputObjects();
+ virtual void UserExec(Option_t *option);
+ virtual void Terminate(Option_t *);
+
+ AliSpectraAODTrackCuts * GetTrackCuts() { return fTrackCuts; }
+ AliSpectraAODEventCuts * GetEventCuts() { return fEventCuts; }
+ AliHelperPID * GetHelperPID() { return fHelperPID; }
+ TList * GetOutputList() { return fOutput; }
+
+ void SetTrackCuts(AliSpectraAODTrackCuts * tc) { fTrackCuts = tc; }
+ void SetEventCuts(AliSpectraAODEventCuts * vc) { fEventCuts = vc; }
+ void SetHelperPID(AliHelperPID* pid) { fHelperPID = pid; }
+ void SetnCentBins(Int_t val) { fnCentBins = val; }
+ void SetnQvecBins(Int_t val) { fnQvecBins = val; }
+ void SetnNchBins(Int_t val) { fnNchBins = val; }
+
+
+ Int_t GetNanoTrackID(AliVTrack * track) ;
+
+ private:
+
+ AliAODEvent * fAOD; //! AOD object
+ AliSpectraAODTrackCuts * fTrackCuts; // Track Cuts
+ AliSpectraAODEventCuts * fEventCuts; // Event Cuts
+ AliHelperPID * fHelperPID; // points to class for PID
+ Bool_t fIsMC; // true if processing MC
+ Bool_t fDoDoubleCounting; // true is double counting for Nsigma accepetd
+ Bool_t fFillOnlyEvents; // if true fill only NSparseHistEv
+ Int_t fCharge; // charge to be selected
+ Int_t fVZEROside; // 0: VZERO-A 1: VZERO-C
+ TList * fOutput; // output list
+ Int_t fnCentBins; // number of bins for the centrality axis
+ Int_t fnQvecBins; // number of bins for the q vector axis
+ Int_t fnNchBins; // number of bins for the Nch axis
+ AliAnalysisTaskSpectraAllChNanoAOD(const AliAnalysisTaskSpectraAllChNanoAOD&);
+ AliAnalysisTaskSpectraAllChNanoAOD& operator=(const AliAnalysisTaskSpectraAllChNanoAOD&);
+
+ ClassDef(AliAnalysisTaskSpectraAllChNanoAOD, 6);
+};
+
+#endif
--- /dev/null
+// run.C
+//
+// Template run macro for AliBasicTask.cxx/.h with example layout of
+// physics selections and options, in macro and task.
+//
+// Author: Arvinder Palaha
+//
+class AliAnalysisGrid;
+
+//______________________________________________________________________________
+void runAllCh(
+ Bool_t isNano = 1,
+ const char* runtype = "local", // local, proof or grid
+ const char *gridmode = "test", // Set the run mode (can be "full", "test", "offline", "submit" or "terminate"). Full & Test work for proof
+ const bool bMCtruth = 1, // 1 = MCEvent handler is on (MC truth), 0 = MCEvent handler is off (MC reconstructed/real data)
+ const bool bMCphyssel = 0, // 1 = looking at MC truth or reconstructed, 0 = looking at real data
+ const Long64_t nentries = 2000, // for local and proof mode, ignored in grid mode. Set to 1234567890 for all events.
+ const Long64_t firstentry = 0, // for local and proof mode, ignored in grid mode
+ const char *proofdataset = "/alice/data/LHC10c_000120821_p1", // path to dataset on proof cluster, for proof analysis
+ const char *proofcluster = "alice-caf.cern.ch", // which proof cluster to use in proof mode
+ const char *taskname = "example_task" // sets name of grid generated macros
+ )
+{
+ // check run type
+ if(runtype != "local" && runtype != "proof" && runtype != "grid"){
+ Printf("\n\tIncorrect run option, check first argument of run macro");
+ Printf("\tint runtype = local, proof or grid\n");
+ return;
+ }
+ Printf("%s analysis chosen",runtype);
+
+ // load libraries
+ gSystem->Load("libCore.so");
+ gSystem->Load("libGeom.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libVMC");
+ gSystem->Load("libTree");
+ gSystem->Load("libProof");
+ gSystem->Load("libMatrix");
+ gSystem->Load("libMinuit");
+ gSystem->Load("libSTEERBase");
+ gSystem->Load("libESD");
+ gSystem->Load("libAOD");
+ gSystem->Load("libANALYSIS");
+ // return;
+ gSystem->Load("libOADB");
+ gSystem->Load("libANALYSISalice");
+ gSystem->Load("libTENDER");
+ gSystem->Load("libCORRFW");
+ // gSystem->Load("libPWG0base");
+ gSystem->Load("libMinuit");
+ gSystem->Load("libPWGTools");
+ gSystem->Load("libPWGLFspectra");
+ gSystem->Load("libPWGLFthermalfits");
+ gSystem->Load("libPWGDevNanoAOD.so");
+
+
+ // add aliroot indlude path
+ gROOT->ProcessLine(Form(".include %s/include",gSystem->ExpandPathName("$ALICE_ROOT")));
+ gROOT->SetStyle("Plain");
+
+ // analysis manager
+ AliAnalysisManager* mgr = new AliAnalysisManager(taskname);
+
+ // create the alien handler and attach it to the manager
+ AliAnalysisGrid *plugin = CreateAlienHandler(taskname, gridmode, proofcluster, proofdataset);
+ mgr->SetGridHandler(plugin);
+
+ AliInputEventHandler* iH = new AliAODInputHandler();
+ if(isNano) {
+ iH->SetEventSelection(new AliAnalysisNanoAODTrackCuts); // FIXME: we need this, otherwise we crash in AliAODInputHandler::BeginEvent where fIsSelectedResult = fEvent->GetHeader()->GetOfflineTrigger(). This is a temporary hack. In the future, it will be solved by using AliVHeader in the AliAODInputHandler.
+ }
+
+// AliAODInputHandler* iH = new AliAODInputHandler();
+// iH->SetInactiveBranches("tracks. vertices. v0s. cascades. jets. caloClusters. fmdClusters. pmdClusters. dimuons. AliAODZDC");
+// iH->SetInactiveBranches("*");
+// iH->SetCheckStatistics(kTRUE);
+ mgr->SetInputEventHandler(iH);
+
+ // mc event handlerrunEx01.C
+ if(bMCtruth) {
+ AliMCEventHandler* mchandler = new AliMCEventHandler();
+ // Not reading track references
+ mchandler->SetReadTR(kFALSE);
+ mgr->SetMCtruthEventHandler(mchandler);
+ }
+
+ // === Physics Selection Task ===
+ //
+ // In SelectCollisionCandidate(), default is kMB, so the task UserExec()
+ // function is only called for these events.
+ // Options are:
+ // kMB Minimum Bias trigger
+ // kMBNoTRD Minimum bias trigger where the TRD is not read out
+ // kMUON Muon trigger
+ // kHighMult High-Multiplicity Trigger
+ // kUserDefined For manually defined trigger selection
+ //
+ // Multiple options possible with the standard AND/OR operators && and ||
+ // These all have the usual offline SPD or V0 selections performed.
+ //
+ // With a pointer to the physics selection object using physSelTask->GetPhysicsSelection(),
+ // one can manually set the selected and background classes using:
+ // AddCollisionTriggerClass("+CINT1B-ABCE-NOPF-ALL")
+ // AddBGTriggerClass("+CINT1A-ABCE-NOPF-ALL");
+ //
+ // One can also specify multiple classes at once, or require a class to NOT
+ // trigger, for e.g.
+ // AddBGTriggerClass("+CSMBA-ABCE-NOPF-ALL -CSMBB-ABCE-NOPF-ALL");
+ //
+ // NOTE that manually setting the physics selection overrides the standard
+ // selection, so it must be done in completeness.
+ //
+ // ALTERNATIVELY, one can make the physics selection inside the task
+ // UserExec().
+ // For this case, comment out the task->SelectCol.... line,
+ // and see AliBasicTask.cxx UserExec() function for details on this.
+
+// gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
+// AliPhysicsSelectionTask *physSelTask = AddTaskPhysicsSelection(bMCphyssel);
+// if(!physSelTask) { Printf("no physSelTask"); return; }
+ //AliPhysicsSelection *physSel = physSelTask->GetPhysicsSelection();
+ //physSel->AddCollisionTriggerClass("+CINT1B-ABCE-NOPF-ALL");// #3119 #769");
+
+ // PID RESPONSE
+ if(!isNano){
+ gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
+ AliAnalysisTaskPIDResponse *taskPID=AddTaskPIDResponse(bMCtruth);
+ taskPID->SetUseTPCEtaCorrection(kTRUE);
+ taskPID->SetUserDataRecoPass(2);
+ }
+
+
+ // create task
+ gROOT->LoadMacro("AddTaskSpectraAllChNanoAOD.C");
+ //load calibration object for LHC10h
+ AliAnalysisTaskSpectraAllChNanoAOD* task = AddTaskSpectraAllChAOD(bMCtruth, // MC
+ 0, // cent min
+ 100, //cent max
+ -999999., // q min
+ 999999., // q max
+ -0.8, // eta min
+ 0.8, // eta max
+ 15.0, // pt max
+ 0.6, // pt tof
+ 1024, // trkbit
+ 100000, // dca
+ 70, // tpc cls
+ 3.0, // nsigma
+ 2, // pid type, 2 = nsigma TPC TOF
+ "_EvtOnly_F1024_VZEROA");
+ TFile * fCalib=TFile::Open("calibV0New.root");
+ task->GetEventCuts()->SetCalibFile(fCalib);
+ task->GetEventCuts()->SetIsLHC10h(kTRUE);
+ task->GetEventCuts()->SetCentralityMethod("V0M");
+
+ task->SetVZEROside(0); //0: VZERO-A 1: VZERO-C
+
+ task->SetFillOnlyEvents(0);
+ task->SetnCentBins(100);
+ task->SetnQvecBins(400);
+ task->SetnNchBins(1);
+
+ // task->SelectCollisionCandidates(AliVEvent::kMB); // if physics selection performed in UserExec(), this line should be commented
+ mgr->AddTask(task);
+
+ // set output root file name for different analysis
+ // TString outfilename = Form("list.%s.root",runtype);
+
+ // create containers for input/output
+ // AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ // AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("coutput1", TList::Class(), AliAnalysisManager::kOutputContainer, outfilename);
+
+ // connect input/output
+ // mgr->ConnectInput(task, 0, cinput);
+ // mgr->ConnectOutput(task, 1, coutput1);
+
+ // enable debug printouts
+ mgr->SetDebugLevel(10);
+ // mgr->SetNSysInfo(100);
+ if (!mgr->InitAnalysis()) return;
+ mgr->PrintStatus();
+
+ // start analysis
+ Printf("Starting Analysis....");
+ mgr->StartAnalysis(runtype,nentries,firstentry);
+}
+
+//______________________________________________________________________________
+AliAnalysisGrid* CreateAlienHandler(const char *taskname, const char *gridmode, const char *proofcluster, const char *proofdataset)
+{
+ AliAnalysisAlien *plugin = new AliAnalysisAlien();
+ // Set the run mode (can be "full", "test", "offline", "submit" or "terminate")
+ plugin->SetRunMode(gridmode);
+
+ // Set versions of used packages
+ plugin->SetAPIVersion("V1.1x");
+ plugin->SetROOTVersion("v5-34-08-5");
+ plugin->SetAliROOTVersion("vAN-20140317");
+
+ // Declare input data to be processed.
+// plugin->SetCheckCopy(kFALSE);
+
+ // Method 1: Create automatically XML collections using alien 'find' command.
+ // Define production directory LFN
+ plugin->SetGridDataDir("/alice/data/2010/LHC10b");
+ // On real reconstructed data:
+ // plugin->SetGridDataDir("/alice/data/2009/LHC09d");
+ // Set data search pattern
+ //plugin->SetDataPattern("*ESDs.root"); // THIS CHOOSES ALL PASSES
+ // Data pattern for reconstructed data
+ plugin->SetDataPattern("*ESDs/pass2/*ESDs.root"); // CHECK LATEST PASS OF DATA SET IN ALIENSH
+ // plugin->SetDataPattern("ESDs/pass2/AOD038/*AliAOD.root"); // CHECK LATEST PASS OF DATA SET IN ALIENSH
+ plugin->SetRunPrefix("000"); // real data
+ // ...then add run numbers to be considered
+ Int_t runlist[15]={117039, 146859, 146858, 146856, 146824, 146817, 146806, 146805, 146804, 146803, 146802, 146801, 146748, 146747, 146746};
+ for (Int_t ind=0; ind<1; ind++) {
+// plugin->AddRunNumber(138275);
+ plugin->AddRunNumber(runlist[ind]);
+ }
+ //plugin->SetRunRange(114917,115322);
+ plugin->SetNrunsPerMaster(10); // 1
+ plugin->SetOutputToRunNo();
+ // comment out the next line when using the "terminate" option, unless
+ // you want separate merged files for each run
+ plugin->SetMergeViaJDL();
+
+ // Method 2: Declare existing data files (raw collections, xml collections, root file)
+ // If no path mentioned data is supposed to be in the work directory (see SetGridWorkingDir())
+ // XML collections added via this method can be combined with the first method if
+ // the content is compatible (using or not tags)
+ // plugin->AddDataFile("tag.xml");
+ // plugin->AddDataFile("/alice/data/2008/LHC08c/000057657/raw/Run57657.Merged.RAW.tag.root");
+
+ // Define alien work directory where all files will be copied. Relative to alien $HOME.
+ plugin->SetGridWorkingDir(taskname);
+
+ // Declare alien output directory. Relative to working directory.
+ plugin->SetGridOutputDir("out"); // In this case will be $HOME/taskname/out
+
+ // Declare the analysis source files names separated by blancs. To be compiled runtime
+ // using ACLiC on the worker nodes.
+ plugin->SetAnalysisSource("AliAnalysisTaskEx01.cxx");
+
+ // Declare all libraries (other than the default ones for the framework. These will be
+ // loaded by the generated analysis macro. Add all extra files (task .cxx/.h) here.
+ plugin->SetAdditionalLibs("AliAnalysisTaskEx01.h AliAnalysisTaskEx01.cxx");
+
+ // Declare the output file names separated by blancs.
+ // (can be like: file.root or file.root@ALICE::Niham::File)
+ // To only save certain files, use SetDefaultOutputs(kFALSE), and then
+ // SetOutputFiles("list.root other.filename") to choose which files to save
+ plugin->SetDefaultOutputs();
+ //plugin->SetOutputFiles("list.root");
+
+ // Optionally set a name for the generated analysis macro (default MyAnalysis.C)
+ plugin->SetAnalysisMacro(Form("%s.C",taskname));
+
+ // Optionally set maximum number of input files/subjob (default 100, put 0 to ignore)
+ plugin->SetSplitMaxInputFileNumber(100);
+
+ // Optionally modify the executable name (default analysis.sh)
+ plugin->SetExecutable(Form("%s.sh",taskname));
+
+ // set number of test files to use in "test" mode
+ plugin->SetNtestFiles(10);
+
+ // Optionally resubmit threshold.
+ plugin->SetMasterResubmitThreshold(90);
+
+ // Optionally set time to live (default 30000 sec)
+ plugin->SetTTL(30000);
+
+ // Optionally set input format (default xml-single)
+ plugin->SetInputFormat("xml-single");
+
+ // Optionally modify the name of the generated JDL (default analysis.jdl)
+ plugin->SetJDLName(Form("%s.jdl",taskname));
+
+ // Optionally modify job price (default 1)
+ plugin->SetPrice(1);
+
+ // Optionally modify split mode (default 'se')
+ plugin->SetSplitMode("se");
+
+ //----------------------------------------------------------
+ //--- PROOF MODE SPECIFIC SETTINGS ------------
+ //----------------------------------------------------------
+ // Proof cluster
+ plugin->SetProofCluster(proofcluster);
+ // Dataset to be used
+ plugin->SetProofDataSet(proofdataset);
+ // May need to reset proof. Supported modes: 0-no reset, 1-soft, 2-hard
+ plugin->SetProofReset(0);
+ // May limit number of workers
+ plugin->SetNproofWorkers(0);
+ // May limit the number of workers per slave
+ plugin->SetNproofWorkersPerSlave(1);
+ // May use a specific version of root installed in proof
+ plugin->SetRootVersionForProof("current");
+ // May set the aliroot mode. Check http://aaf.cern.ch/node/83
+ plugin->SetAliRootMode("default"); // Loads AF libs by default
+ // May request ClearPackages (individual ClearPackage not supported)
+ plugin->SetClearPackages(kFALSE);
+ // Plugin test mode works only providing a file containing test file locations, used in "local" mode also
+ plugin->SetFileForTestMode("files.txt"); // file should contain path name to a local directory containg *ESDs.root etc
+ // Request connection to alien upon connection to grid
+ plugin->SetProofConnectGrid(kFALSE);
+ // Other PROOF specific parameters
+ plugin->SetProofParameter("PROOF_UseMergers","-1");
+ printf("Using: PROOF_UseMergers : %s\n", plugin->GetProofParameter("PROOF_UseMergers"));
+ return plugin;
+}
+