--- /dev/null
+AliAnalysisTaskSpectraAOD * AddTaskSpectraAOD(const char * outfilename)
+{
+ // TODO: add some parameters to set the centrality for this task, and maybe the name of the task
+ // TODO: shall I use the same file and different dirs for the different centralities?
+
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+ if (!mgr) {
+ ::Error("AddTaskSpectraAOD", "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("AddTaskSpectraAOD", "This task requires an input event handler");
+ return NULL;
+ }
+ TString inputDataType = mgr->GetInputEventHandler()->GetDataType(); // can be "ESD" or "AOD"
+
+ if (inputDataType != "AOD") {
+ Printf("ERROR! This task can only run on AODs!");
+ }
+
+ // Configure analysis
+ //===========================================================================
+
+
+ // Set I/O
+ AliAnalysisDataContainer *cinput0 = mgr->GetCommonInputContainer();
+ AliAnalysisTaskSpectraAOD *task = new AliAnalysisTaskSpectraAOD("TaskAODSpectra");
+ mgr->AddTask(task);
+
+ // Set the cuts
+ AliSpectraAODVertexCuts * vcuts = new AliSpectraAODVertexCuts("VertexCuts");
+ AliSpectraAODTrackCuts * tcuts = new AliSpectraAODTrackCuts ("TracksCuts");
+ tcuts->SetTrackType(6);
+ tcuts->SetEta(1.);
+ task->SetVertexCuts(vcuts);
+ task->SetTrackCuts (tcuts);
+
+
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ AliAnalysisDataContainer *coutputpt1 = mgr->CreateContainer("chistpt", AliSpectraAODHistoManager::Class(), AliAnalysisManager::kOutputContainer, outfilename);
+ AliAnalysisDataContainer *coutputpt2 = mgr->CreateContainer("cvcutpt", AliSpectraAODVertexCuts::Class(), AliAnalysisManager::kOutputContainer, outfilename);
+ AliAnalysisDataContainer *coutputpt3 = mgr->CreateContainer("ctcutpt", AliSpectraAODTrackCuts::Class(), AliAnalysisManager::kOutputContainer, outfilename);
+
+ mgr->ConnectInput(task, 0, cinput);
+ mgr->ConnectOutput(task, 1, coutputpt1);
+ mgr->ConnectOutput(task, 2, coutputpt2);
+ mgr->ConnectOutput(task, 3, coutputpt3);
+
+ return task;
+}
--- /dev/null
+/**************************************************************************\r
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: The ALICE Off-line Project. *\r
+ * Contributors are mentioned in the code where appropriate. *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+\r
+//-----------------------------------------------------------------\r
+// AliAnalysisTaskSpectraAOD class\r
+//-----------------------------------------------------------------\r
+\r
+#include "TChain.h"\r
+#include "TTree.h"\r
+#include "TLegend.h"\r
+#include "TH1F.h"\r
+#include "TH2F.h"\r
+#include "TCanvas.h"\r
+#include "AliAnalysisTask.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliAODTrack.h"\r
+#include "AliAODMCParticle.h"\r
+#include "AliAODEvent.h"\r
+#include "AliAODInputHandler.h"\r
+#include "AliAnalysisTaskSpectraAOD.h"\r
+#include "AliAnalysisTaskESDfilter.h"\r
+#include "AliAnalysisDataContainer.h"\r
+#include "AliSpectraAODHistoManager.h"\r
+#include "AliSpectraAODTrackCuts.h"\r
+#include "AliSpectraAODEventCuts.h"\r
+#include "AliCentrality.h"\r
+#include "TProof.h"\r
+#include "AliPID.h"\r
+#include "AliVEvent.h"\r
+#include "AliPIDResponse.h"\r
+#include <iostream>\r
+\r
+using namespace AliSpectraNameSpace;\r
+using namespace std;\r
+\r
+ClassImp(AliAnalysisTaskSpectraAOD) // EX1 This stuff tells root to implement the streamer, inspector methods etc (we discussed about it today)\r
+\r
+//________________________________________________________________________\r
+AliAnalysisTaskSpectraAOD::AliAnalysisTaskSpectraAOD(const char *name) : AliAnalysisTaskSE(name), fAOD(0), fHistMan(0), fTrackCuts(0), fEventCuts(0), fIsMC(0), fPIDResponse(0), fNSigmaPID(0), fYCut(0)\r
+{\r
+ // Default constructor\r
+ fNSigmaPID = 3.;\r
+ fYCut = .5;\r
+ DefineInput(0, TChain::Class());\r
+ DefineOutput(1, AliSpectraAODHistoManager::Class());\r
+ DefineOutput(2, AliSpectraAODEventCuts::Class());\r
+ DefineOutput(3, AliSpectraAODTrackCuts::Class());\r
+\r
+}\r
+//________________________________________________________________________\r
+Bool_t AliAnalysisTaskSpectraAOD::CheckYCut(AODParticleSpecies_t species, AliAODTrack* track) const\r
+{\r
+ // check if the rapidity is within the set range\r
+ // note: masses are hardcoded for now. we could look them up in the pdg database, but that would mean accecing it 100k+ times per run ...\r
+ Double_t y;\r
+ if (species == kProton) { y = track->Y(9.38271999999999995e-01); }\r
+ if ( species == kKaon ) { y = track->Y(4.93676999999999977e-01); }\r
+ if ( species == kPion) { y = track->Y(1.39570000000000000e-01); }\r
+ if (TMath::Abs(y) > fYCut || y < -998.) return kFALSE;\r
+ return kTRUE;\r
+}\r
+//____________________________________________________________________________\r
+Bool_t AliAnalysisTaskSpectraAOD::CheckYCut(AliAODMCParticle* particle) const\r
+{\r
+ // check if the rapidity is within the set range\r
+ Double_t y = particle->Y();\r
+ if (TMath::Abs(y) > fYCut || y < -998.) return kFALSE;\r
+ return kTRUE;\r
+}\r
+//________________________________________________________________________\r
+void AliAnalysisTaskSpectraAOD::UserCreateOutputObjects()\r
+{\r
+ // create output objects\r
+ fHistMan = new AliSpectraAODHistoManager("SpectraHistos");\r
+\r
+ if (!fTrackCuts) AliFatal("Track Cuts should be set in the steering macro");\r
+ if (!fEventCuts) AliFatal("Event Cuts should be set in the steering macro");\r
+\r
+ AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();\r
+ AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());\r
+ fPIDResponse = inputHandler->GetPIDResponse();\r
+\r
+ PostData(1, fHistMan);\r
+ PostData(2, fEventCuts);\r
+ PostData(3, fTrackCuts);\r
+\r
+}\r
+//________________________________________________________________________\r
+void AliAnalysisTaskSpectraAOD::UserExec(Option_t *)\r
+{\r
+ // main event loop\r
+ fAOD = dynamic_cast<AliAODEvent*>(fInputEvent);\r
+ if (strcmp(fAOD->ClassName(), "AliAODEvent"))\r
+ {\r
+ AliFatal("Not processing AODs");\r
+ }\r
+\r
+ if(!fPIDResponse) {\r
+ AliError("Cannot get pid response");\r
+ return;\r
+ }\r
+\r
+ if (!fEventCuts->IsSelected(fAOD)) return;\r
+\r
+ //AliCentrality fAliCentral*;\r
+// if ((fAOD->GetCentrality())->GetCentralityPercentile("V0M") > 5.) return;\r
+ \r
+ // First do MC to fill up the MC particle array, such that we can use it later\r
+ TClonesArray *arrayMC = 0;\r
+ if (fIsMC)\r
+ {\r
+ arrayMC = (TClonesArray*) fAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName());\r
+ if (!arrayMC) {\r
+ AliFatal("Error: MC particles branch not found!\n");\r
+ }\r
+ Int_t nMC = arrayMC->GetEntries();\r
+ for (Int_t iMC = 0; iMC < nMC; iMC++)\r
+ {\r
+ AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(iMC);\r
+ fHistMan->GetPtHistogram(kHistPtGen)->Fill(partMC->Pt());\r
+ \r
+ if(TMath::Abs(partMC->Eta()) > fTrackCuts->GetEta()) continue;\r
+ // check for true PID + and fill P_t histos \r
+ if (partMC->IsPhysicalPrimary() && CheckYCut(partMC) ) {// only primary vertices and y cut satisfied\r
+ if ( partMC->PdgCode() == 2212) { fHistMan->GetHistogram(kHistPtGenTruePrimaryProtonPlus)->Fill(partMC->Pt()); } \r
+ if ( partMC->PdgCode() == -2212) { fHistMan->GetHistogram(kHistPtGenTruePrimaryProtonMinus)->Fill(partMC->Pt()); } \r
+ if ( partMC->PdgCode() == 321) { fHistMan->GetHistogram(kHistPtGenTruePrimaryKaonPlus)->Fill(partMC->Pt()); } \r
+ if ( partMC->PdgCode() == -321) { fHistMan->GetHistogram(kHistPtGenTruePrimaryKaonMinus)->Fill(partMC->Pt()); } \r
+ if ( partMC->PdgCode() == 211) { fHistMan->GetHistogram(kHistPtGenTruePrimaryPionPlus)->Fill(partMC->Pt()); } \r
+ if ( partMC->PdgCode() == -211) { fHistMan->GetHistogram(kHistPtGenTruePrimaryPionMinus)->Fill(partMC->Pt()); } \r
+ }\r
+ }\r
+ }\r
+ \r
+ for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++)\r
+ {\r
+\r
+ AliAODTrack* track = fAOD->GetTrack(iTracks);\r
+ if (!fTrackCuts->IsSelected(track)) continue;\r
+\r
+ fHistMan->GetHistogram(kHistPtRec)->Fill(track->Pt()); // PT histo\r
+ fHistMan->GetPIDHistogram(kHistPIDTPC)->Fill(track->P(), track->GetTPCsignal()); // PID histo\r
+ fHistMan->GetPIDHistogram(kHistPIDTPCPt)->Fill(track->Pt(), track->GetTPCsignal());\r
+ \r
+ // Sigma identify particles (for both MC and regular data):\r
+ // note: as of now, very crude identification. how many sigmas are allowed? what to do with high pt?\r
+ AliVParticle *inEvHMain = dynamic_cast<AliVParticle *>(track);\r
+ Double_t nsigmaTPCkKaon = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kKaon)); \r
+ Double_t nsigmaTPCkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kProton));\r
+ Double_t nsigmaTPCkPion = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kPion)); \r
+ if( ( nsigmaTPCkKaon < nsigmaTPCkPion ) && ( nsigmaTPCkKaon < nsigmaTPCkProton )) { \r
+ \r
+ if ((nsigmaTPCkKaon > fNSigmaPID) || (!CheckYCut(kKaon, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaKaonPlus)->Fill(track->Pt()); } \r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaKaonMinus)->Fill(track->Pt()); } \r
+ }\r
+ if( ( nsigmaTPCkProton < nsigmaTPCkKaon ) && ( nsigmaTPCkProton < nsigmaTPCkPion ) ) {\r
+ if ( nsigmaTPCkProton > fNSigmaPID || (!CheckYCut(kProton, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaProtonPlus)->Fill(track->Pt()); }\r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaProtonMinus)->Fill(track->Pt()); }\r
+ }\r
+ if( (nsigmaTPCkPion < nsigmaTPCkProton ) && ( nsigmaTPCkPion < nsigmaTPCkKaon ) ) {\r
+ if (nsigmaTPCkPion > fNSigmaPID || (!CheckYCut(kPion, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaPionPlus)->Fill(track->Pt()); }\r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaPionMinus)->Fill(track->Pt()); }\r
+ }\r
+ /* MC Part */\r
+ if (arrayMC) {\r
+ AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(TMath::Abs(track->GetLabel()));\r
+ if (!partMC) { \r
+ AliError("Cannot get MC particle");\r
+ continue; }\r
+ if (CheckYCut(partMC)) {\r
+ // primaries, true pid\r
+ if ( partMC->PdgCode() == 2212) { fHistMan->GetHistogram(kHistPtRecTrueProtonPlus)->Fill(partMC->Pt()); \r
+ if (partMC->IsPhysicalPrimary()) {fHistMan->GetHistogram(kHistPtRecTruePrimaryProtonPlus)->Fill(partMC->Pt()); }}\r
+ if ( partMC->PdgCode() == -2212) { fHistMan->GetHistogram(kHistPtRecTrueProtonMinus)->Fill(partMC->Pt()); \r
+ if (partMC->IsPhysicalPrimary()) {fHistMan->GetHistogram(kHistPtRecTruePrimaryProtonMinus)->Fill(partMC->Pt()); }}\r
+ if ( partMC->PdgCode() == 321) { fHistMan->GetHistogram(kHistPtRecTrueKaonPlus)->Fill(partMC->Pt()); \r
+ if (partMC->IsPhysicalPrimary()) {fHistMan->GetHistogram(kHistPtRecTruePrimaryKaonPlus)->Fill(partMC->Pt()); }}\r
+ if ( partMC->PdgCode() == -321) { fHistMan->GetHistogram(kHistPtRecTrueKaonMinus)->Fill(partMC->Pt()); \r
+ if (partMC->IsPhysicalPrimary()) {fHistMan->GetHistogram(kHistPtRecTruePrimaryKaonMinus)->Fill(partMC->Pt()); }}\r
+ if ( partMC->PdgCode() == 211) { fHistMan->GetHistogram(kHistPtRecTruePionPlus)->Fill(partMC->Pt()); \r
+ if (partMC->IsPhysicalPrimary()) {fHistMan->GetHistogram(kHistPtRecTruePrimaryPionPlus)->Fill(partMC->Pt()); }}\r
+ if ( partMC->PdgCode() == -211) { fHistMan->GetHistogram(kHistPtRecTruePionMinus)->Fill(partMC->Pt()); \r
+ if (partMC->IsPhysicalPrimary()) {fHistMan->GetHistogram(kHistPtRecTruePrimaryPionMinus)->Fill(partMC->Pt()); }}\r
+ }\r
+\r
+ // primaries, sigma pid\r
+ if (partMC->IsPhysicalPrimary()) { \r
+ if( ( nsigmaTPCkKaon < nsigmaTPCkPion ) && ( nsigmaTPCkKaon < nsigmaTPCkProton ) ) {\r
+ if ( (nsigmaTPCkKaon > fNSigmaPID ) || (!CheckYCut(kKaon, track) ) ) continue; \r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaPrimaryKaonPlus)->Fill(track->Pt()); } \r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaPrimaryKaonMinus)->Fill(track->Pt()); } \r
+ }\r
+ if( ( nsigmaTPCkProton < nsigmaTPCkKaon ) && ( nsigmaTPCkProton < nsigmaTPCkPion ) ) {\r
+ if ( (nsigmaTPCkProton > fNSigmaPID ) || (!CheckYCut(kProton, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaPrimaryProtonPlus)->Fill(track->Pt()); }\r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaPrimaryProtonMinus)->Fill(track->Pt()); }\r
+ }\r
+ if( (nsigmaTPCkPion < nsigmaTPCkProton ) && ( nsigmaTPCkPion < nsigmaTPCkKaon ) ) {\r
+ if ( ( nsigmaTPCkPion > fNSigmaPID ) || (!CheckYCut(kPion, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaPrimaryPionPlus)->Fill(track->Pt()); }\r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaPrimaryPionMinus)->Fill(track->Pt()); }\r
+ }\r
+ }\r
+ //secondaries\r
+ else {\r
+ if( ( nsigmaTPCkKaon < nsigmaTPCkPion ) && ( nsigmaTPCkKaon < nsigmaTPCkProton ) ) { \r
+ if ( (nsigmaTPCkKaon > fNSigmaPID ) || (!CheckYCut(kKaon, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaSecondaryKaonPlus)->Fill(track->Pt()); } \r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaSecondaryKaonMinus)->Fill(track->Pt()); } \r
+ }\r
+ if( ( nsigmaTPCkProton < nsigmaTPCkKaon ) && ( nsigmaTPCkProton < nsigmaTPCkPion ) ) {\r
+ if ( (nsigmaTPCkProton > fNSigmaPID ) || (!CheckYCut(kProton, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaSecondaryProtonPlus)->Fill(track->Pt()); }\r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaSecondaryProtonMinus)->Fill(track->Pt()); }\r
+ }\r
+ if( (nsigmaTPCkPion < nsigmaTPCkProton ) && ( nsigmaTPCkPion < nsigmaTPCkKaon ) ) {\r
+ if ( ( nsigmaTPCkPion > fNSigmaPID ) || (!CheckYCut(kPion, track) ) ) continue;\r
+ if ( track->Charge() > 0 ) { fHistMan->GetHistogram(kHistPtRecSigmaSecondaryPionPlus)->Fill(track->Pt()); }\r
+ else { fHistMan->GetHistogram(kHistPtRecSigmaSecondaryPionMinus)->Fill(track->Pt()); }\r
+ }\r
+ }\r
+ }\r
+ }\r
+\r
+ PostData(1, fHistMan);\r
+ PostData(2, fEventCuts);\r
+ PostData(3, fTrackCuts);\r
+}\r
+ \r
+//_________________________________________________________________\r
+void AliAnalysisTaskSpectraAOD::Terminate(Option_t *)\r
+{\r
+ // Terminate\r
+}\r
--- /dev/null
+#ifndef ALIANALYSISTASKSPECTRAAOD_H\r
+#define ALIANALYSISTASKSPECTRAAOD_H\r
+\r
+/* See cxx source for full Copyright notice */\r
+\r
+//-------------------------------------------------------------------------\r
+// AliAnalysisTaskSpectraAOD\r
+//\r
+//\r
+//\r
+//\r
+// Author: Michele Floris, CERN\r
+//-------------------------------------------------------------------------\r
+\r
+class TH1F;\r
+class TH2F;\r
+class AliAODEvent;\r
+class AliSpectraAODHistoManager;\r
+class AliSpectraAODTrackCuts;\r
+class AliSpectraAODEventCuts;\r
+#include "AliSpectraAODHistoManager.h"\r
+#include "AliAnalysisTaskSE.h"\r
+\r
+\r
+class AliAnalysisTaskSpectraAOD : public AliAnalysisTaskSE\r
+{\r
+public:\r
+\r
+ // constructors\r
+ AliAnalysisTaskSpectraAOD() : AliAnalysisTaskSE(), fAOD(0), fHistMan(0), fTrackCuts(0), fEventCuts(0), fIsMC(0), fPIDResponse(0), fNSigmaPID(0), fYCut(0) {}\r
+ AliAnalysisTaskSpectraAOD(const char *name);\r
+ virtual ~AliAnalysisTaskSpectraAOD() {}\r
+\r
+ void SetIsMC(Bool_t isMC = kFALSE) {fIsMC = isMC; };\r
+ Bool_t GetIsMC() const { return fIsMC;};\r
+ void SetNSigmaForIdentification (Double_t sigma ) { fNSigmaPID = sigma; }\r
+ Double_t GetNSigmaForIdentification () const {return fNSigmaPID; }\r
+ void SetYCut (Double_t y ) { fYCut = y; }\r
+ Double_t GetYCut () const {return fYCut; }\r
+\r
+ virtual void UserCreateOutputObjects();\r
+ Bool_t CheckYCut(AliSpectraNameSpace::AODParticleSpecies_t species, AliAODTrack* track) const;\r
+ Bool_t CheckYCut(AliAODMCParticle* particle) const;\r
+ virtual void UserExec(Option_t *option);\r
+ virtual void Terminate(Option_t *);\r
+ void SetTrackCuts(AliSpectraAODTrackCuts * tc) { fTrackCuts = tc; }\r
+ void SetEventCuts(AliSpectraAODEventCuts * vc) { fEventCuts = vc; }\r
+\r
+private:\r
+\r
+ AliAODEvent * fAOD; //! AOD object\r
+ AliSpectraAODHistoManager * fHistMan; // Histogram Manager\r
+ AliSpectraAODTrackCuts * fTrackCuts; // Track Cuts\r
+ AliSpectraAODEventCuts * fEventCuts; // Event Cuts\r
+ Bool_t fIsMC;// true if processing MC\r
+ AliPIDResponse *fPIDResponse; // ! PID response object\r
+ Double_t fNSigmaPID; // Maximum number of sigmas allowed in particle identification\r
+ Double_t fYCut; // Maximum rapidity - calculated from identified particle's mass\r
+ AliAnalysisTaskSpectraAOD(const AliAnalysisTaskSpectraAOD&);\r
+ AliAnalysisTaskSpectraAOD& operator=(const AliAnalysisTaskSpectraAOD&);\r
+\r
+ ClassDef(AliAnalysisTaskSpectraAOD, 1);\r
+};\r
+\r
+#endif\r
--- /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. *
+ **************************************************************************/
+
+//-----------------------------------------------------------------
+// AliSpectraAODEventCuts class
+//-----------------------------------------------------------------
+
+#include "TChain.h"
+#include "TTree.h"
+#include "TLegend.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TCanvas.h"
+#include "AliAnalysisTask.h"
+#include "AliAnalysisManager.h"
+#include "AliAODTrack.h"
+#include "AliAODMCParticle.h"
+#include "AliAODEvent.h"
+#include "AliAODInputHandler.h"
+#include "AliAnalysisTaskESDfilter.h"
+#include "AliAnalysisDataContainer.h"
+#include "AliSpectraAODEventCuts.h"
+#include "AliSpectraAODHistoManager.h"
+#include <iostream>
+
+using namespace std;
+
+ClassImp(AliSpectraAODEventCuts)
+
+AliSpectraAODEventCuts::AliSpectraAODEventCuts(const char *name) : TNamed(name, "AOD Event Cuts"), fAOD(0), fIsSelected(0), fCentralityCutMin(0), fCentralityCutMax(0), fHistoCuts(0)
+{
+ // Constructor
+ fHistoCuts = new TH1I("fEventCuts", "Event Cuts", kNVtxCuts, -0.5, kNVtxCuts - 0.5);
+ fCentralityCutMin = 0.0; // default value of centrality cut minimum, 0 ~ no cut
+ fCentralityCutMax = 10000.0; // default value of centrality cut maximum, ~ no cut
+
+}
+
+//______________________________________________________
+Bool_t AliSpectraAODEventCuts::IsSelected(AliAODEvent * aod)
+{
+// Returns true if Event Cuts are selected and applied
+ fAOD = aod;
+ fIsSelected = (CheckVtxRange() && CheckCentralityCut());
+ if(fIsSelected) fHistoCuts->Fill(kAcceptedEvents);
+ return fIsSelected;
+}
+
+//______________________________________________________
+Bool_t AliSpectraAODEventCuts::CheckVtxRange()
+{
+ // reject events outside of range
+ AliAODVertex * vertex = fAOD->GetPrimaryVertex();
+ if (!vertex)
+ {
+ fHistoCuts->Fill(kVtxNoEvent);
+ return kFALSE;
+ }
+ if (TMath::Abs(vertex->GetZ()) < 10)
+ {
+ return kTRUE;
+ }
+ fHistoCuts->Fill(kVtxRange);
+ return kFALSE;
+}
+
+//______________________________________________________
+Bool_t AliSpectraAODEventCuts::CheckCentralityCut()
+{
+ // Check centrality cut
+ if ( (fAOD->GetCentrality()->GetCentralityPercentile("V0M") <= fCentralityCutMax) && (fAOD->GetCentrality()->GetCentralityPercentile("V0M") >= fCentralityCutMin) ) return kTRUE;
+ fHistoCuts->Fill(kVtxCentral);
+ return kFALSE;
+}
+
+//______________________________________________________
+void AliSpectraAODEventCuts::PrintCuts()
+{
+ // print info about event cuts
+ cout << "Event Stats" << endl;
+ cout << " > Number of accepted events: " << fHistoCuts->GetBinContent(kAcceptedEvents + 1) << endl;
+ cout << " > Vertex out of range: " << fHistoCuts->GetBinContent(kVtxRange + 1) << endl;
+ cout << " > Events cut by centrality: " << fHistoCuts->GetBinContent(kVtxCentral + 1) << endl;
+ cout << " > Events without vertex: " << fHistoCuts->GetBinContent(kVtxNoEvent + 1) << endl;
+ }
+//______________________________________________________
+
+Long64_t AliSpectraAODEventCuts::Merge(TCollection* list)
+{
+ // Merge a list of AliSpectraAODEventCuts objects with this.
+ // Returns the number of merged objects (including this).
+
+ // AliInfo("Merging");
+
+ if (!list)
+ return 0;
+
+ if (list->IsEmpty())
+ return 1;
+
+ TIterator* iter = list->MakeIterator();
+ TObject* obj;
+
+ // collections of all histograms
+ TList collections;
+
+ Int_t count = 0;
+
+ while ((obj = iter->Next())) {
+ AliSpectraAODEventCuts* entry = dynamic_cast<AliSpectraAODEventCuts*> (obj);
+ if (entry == 0)
+ continue;
+
+ TH1I * histo = entry->GetHistoCuts();
+ collections.Add(histo);
+ count++;
+ }
+
+ fHistoCuts->Merge(&collections);
+
+ delete iter;
+
+ return count+1;
+}
+
--- /dev/null
+#ifndef ALISPECTRAAODEVENTCUTS_H
+#define ALISPECTRAAODEVENTCUTS_H
+
+/* See cxx source for full Copyright notice */
+
+//-------------------------------------------------------------------------
+// AliSpectraAODEventCuts
+//
+//
+//
+//
+// Authors: Michele Floris, CERN, Philip Versteeg, UU, Redmer Bertens, UU
+//-------------------------------------------------------------------------
+
+class AliAODEvent;
+class AliSpectraAODHistoManager;
+
+#include "TNamed.h"
+
+class AliSpectraAODEventCuts : public TNamed
+{
+public:
+ enum { kAcceptedEvents = 0, kVtxRange, kVtxCentral, kVtxNoEvent, kNVtxCuts};
+
+ // Constructors
+ AliSpectraAODEventCuts() : TNamed(), fAOD(0), fIsSelected(0), fCentralityCutMin(0), fCentralityCutMax(0), fHistoCuts(0) {}
+ AliSpectraAODEventCuts(const char *name);
+ virtual ~AliSpectraAODEventCuts() {}
+
+ // Methods
+ Bool_t IsSelected(AliAODEvent * aod);
+ Bool_t CheckVtxRange();
+ Bool_t CheckCentralityCut();
+ void SetCentralityCutMin(Float_t cut) { fCentralityCutMin = cut; }
+ void SetCentralityCutMax(Float_t cut) { fCentralityCutMax = cut; }
+
+
+ TH1I * GetHistoCuts() { return fHistoCuts; }
+ Float_t GetCentralityMin() const { return fCentralityCutMin; }
+ Float_t GetCentralityMax() const { return fCentralityCutMax; }
+ void PrintCuts();
+ Float_t NumberOfEvents() { return fHistoCuts->GetBinContent(kAcceptedEvents+1); }
+
+ Long64_t Merge(TCollection* list);
+
+
+private:
+
+ AliAODEvent *fAOD; //! AOD event
+ Bool_t fIsSelected; // True if cuts are selected
+ Float_t fCentralityCutMin; // minimum centrality percentile
+ Float_t fCentralityCutMax; // maximum centrality percentile
+ TH1I *fHistoCuts; // Cuts statistics
+ AliSpectraAODEventCuts(const AliSpectraAODEventCuts&);
+ AliSpectraAODEventCuts& operator=(const AliSpectraAODEventCuts&);
+
+ ClassDef(AliSpectraAODEventCuts, 2);
+
+};
+#endif
+
--- /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. *
+ **************************************************************************/
+
+//-----------------------------------------------------------------
+// AliSpectraAODHistoManager class
+//-----------------------------------------------------------------
+
+#include "TChain.h"
+#include "TTree.h"
+#include "TLegend.h"
+#include "TH1.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TCanvas.h"
+#include "AliAnalysisTask.h"
+#include "AliAnalysisManager.h"
+#include "AliAODTrack.h"
+#include "AliAODMCParticle.h"
+#include "AliAODEvent.h"
+#include "AliAODInputHandler.h"
+#include "AliAnalysisTaskESDfilter.h"
+#include "AliAnalysisDataContainer.h"
+#include "AliSpectraAODHistoManager.h"
+#include <iostream>
+
+using namespace std;
+
+ClassImp(AliSpectraAODHistoManager)
+
+
+using namespace AliSpectraNameSpace;
+
+AliSpectraAODHistoManager::AliSpectraAODHistoManager(const char *name): TNamed(name, "AOD Spectra Histo Manager"), fOutputList(0)
+{
+ // ctor
+ fOutputList = new TList;
+ fOutputList->SetOwner();
+ Bool_t oldStatus = TH1::AddDirectoryStatus();
+ TH1::AddDirectory(kFALSE);
+ for (Int_t ihist = 0; ihist < kNHist ; ihist++)
+ {
+ if (ihist <= kNPtHist) BookPtHistogram(kHistName[ihist]); // PT histos
+ if (ihist > kNPtHist) BookPIDHistogram(kHistName[ihist]); // PID histos
+ }
+ TH1::AddDirectory(oldStatus);
+
+}
+//_______________________________________________________
+TH1F* AliSpectraAODHistoManager::BookPtHistogram(const char * name)
+{
+ // Return a pt histogram with predefined binning, set the ID and add it to the output list
+ AliInfo(Form("Booking pt histogram %s", name));
+
+ TH1F * hist = new TH1F(name, Form("P_{T} distribution (%s)", name), 40, 0, 1.2);
+ hist->GetXaxis()->SetTitle("P_{T} (GeV / c)");
+ hist->GetYaxis()->SetTitle("No. of events");
+ hist->SetMarkerStyle(kFullCircle);
+ hist->Sumw2();
+ fOutputList->Add(hist);
+
+ return hist;
+}
+//_____________________________________________________________________________
+
+TH2F* AliSpectraAODHistoManager::BookPIDHistogram(const char * name)
+{
+ // Return a pt histogram with predefined binning, set the ID and add it to the output list
+ AliInfo(Form("Booking pt histogram %s", name));
+
+ TH2F * hist = new TH2F(name, Form("Particle Identification (%s)", name), 100, 0, 1.2, 1000, 0, 1000);
+ hist->GetXaxis()->SetTitle("P (GeV / c)");
+ hist->GetYaxis()->SetTitle("TPC");
+// hist->Sumw2();
+ fOutputList->Add(hist);
+
+ return hist;
+}
+
+
+Long64_t AliSpectraAODHistoManager::Merge(TCollection* list)
+{
+ // Merge a list of AliSpectraAODHistoManager objects with this.
+ // Returns the number of merged objects (including this).
+
+ // AliInfo("Merging");
+
+ if (!list)
+ return 0;
+
+ if (list->IsEmpty())
+ return 1;
+
+ TIterator* iter = list->MakeIterator();
+ TObject* obj;
+
+ // collections of all histograms
+ TList collections;
+
+ Int_t count = 0;
+
+ while ((obj = iter->Next())) {
+ AliSpectraAODHistoManager* entry = dynamic_cast<AliSpectraAODHistoManager*> (obj);
+ if (entry == 0)
+ continue;
+
+ TList * hlist = entry->GetOutputList();
+ collections.Add(hlist);
+ count++;
+ }
+
+ fOutputList->Merge(&collections);
+
+ delete iter;
+
+ return count+1;
+}
+
--- /dev/null
+
+#ifndef ALISPECTRAAODHISTOMANAGER_H
+#define ALISPECTRAAODHISTOMANAGER_H
+
+/* See cxx source for full Copyright notice */
+
+//-------------------------------------------------------------------------
+// AliSpectraAODHistoManager
+//
+//
+//
+//
+// Authors: Michele Floris, CERN, Philip Versteeg, UU, Redmer Bertens, UU
+//-------------------------------------------------------------------------
+
+class AliAODEvent;
+class TH1F;
+class TH2F;
+class TList;
+#include "TNamed.h"
+
+namespace AliSpectraNameSpace
+{
+ enum AODPtHist_t
+ {
+ // 6 Pt Reconstructed Sigma
+ kHistPtRecSigmaProtonPlus = 0, // Pt histo for protons +, reconstructed tracks, sigma ID
+ kHistPtRecSigmaKaonPlus, // Pt histo for kaons +, reconsructed tracks, sigma ID
+ kHistPtRecSigmaPionPlus, // Pt histo for pions +, reconstructed tracks, sigma ID
+ kHistPtRecSigmaProtonMinus, // Pt histo for protons -, reconstructed tracks, sigma ID
+ kHistPtRecSigmaKaonMinus, // Pt histo for kaons -, reconstructed tracks, sigma ID
+ kHistPtRecSigmaPionMinus, // Pt histo for pions -, reconstructed tracks, sigma ID
+
+ // 6 Pt Reconstructed True
+ kHistPtRecTrueProtonPlus, // Pt histo for protons +, reconstructed tracks, true ID
+ kHistPtRecTrueKaonPlus, // Pt histo for kaons +, reconsructed tracks, true ID
+ kHistPtRecTruePionPlus, // Pt histo for pions +, reconstructed tracks, true ID
+ kHistPtRecTrueProtonMinus, // Pt histo for protons -, reconstructed tracks, true ID
+ kHistPtRecTrueKaonMinus, // Pt histo for kaons -, reconstructed tracks, true ID
+ kHistPtRecTruePionMinus, // Pt histo for pions -, reconstructed tracks, true ID
+
+ // 6 Pt Reconstructed Sigma Primary
+ kHistPtRecSigmaPrimaryProtonPlus, // Pt histo for protons +, reconstructed tracks, sigma ID, primary Event
+ kHistPtRecSigmaPrimaryKaonPlus, // Pt histo for kaons +, reconsructed tracks, sigma ID, primary Event
+ kHistPtRecSigmaPrimaryPionPlus, // Pt histo for pions +, reconstructed tracks, sigma ID, primary Event
+ kHistPtRecSigmaPrimaryProtonMinus, // Pt histo for protons -, reconstructed tracks, sigma ID, primary Event
+ kHistPtRecSigmaPrimaryKaonMinus, // Pt histo for kaons -, reconstructed tracks, sigma ID, primary Event
+ kHistPtRecSigmaPrimaryPionMinus, // Pt histo for pions -, reconstructed tracks, sigma ID, primary Event
+
+ // 6 Pt Reconstructed Sigma Secondary
+ kHistPtRecSigmaSecondaryProtonPlus, // Pt histo for protons +, reconstructed tracks, sigma ID, secondary Event
+ kHistPtRecSigmaSecondaryKaonPlus, // Pt histo for kaons +, reconsructed tracks, sigma ID, secondary Event
+ kHistPtRecSigmaSecondaryPionPlus, // Pt histo for pions +, reconstructed tracks, sigma ID, secondary Event
+ kHistPtRecSigmaSecondaryProtonMinus, // Pt histo for protons -, reconstructed tracks, sigma ID, secondary Event
+ kHistPtRecSigmaSecondaryKaonMinus, // Pt histo for kaons -, reconstructed tracks, sigma ID, secondary Event
+ kHistPtRecSigmaSecondaryPionMinus, // Pt histo for pions -, reconstructed tracks, sigma ID, secondary Event
+
+ // 6 Pt Generated True Primary
+ kHistPtGenTruePrimaryProtonPlus, // Pt histo for protons +, generated tracks, true ID, primary Event
+ kHistPtGenTruePrimaryKaonPlus, // Pt histo for kaons +, generated tracks, true ID, primary Event
+ kHistPtGenTruePrimaryPionPlus, // Pt histo for pions +, generated tracks, true ID, primary Event
+ kHistPtGenTruePrimaryProtonMinus, // Pt histo for protons -, generated tracks, true ID, primary Event
+ kHistPtGenTruePrimaryKaonMinus, // Pt histo for kaons -, generated tracks, true ID, primary Event
+ kHistPtGenTruePrimaryPionMinus, // Pt histo for pions -, generated tracks, true ID, primary Event
+ kNPtSpecies = kHistPtGenTruePrimaryPionMinus,
+
+ // 6 Pt Reconstructed True Primary
+ kHistPtRecTruePrimaryProtonPlus, // Pt histo for protons +, reconstructed tracks, true ID, primary event
+ kHistPtRecTruePrimaryKaonPlus, // Pt histo for kaons +, reconsructed tracks, true ID, primary event
+ kHistPtRecTruePrimaryPionPlus, // Pt histo for pions +, reconstructed tracks, true ID, primary event
+ kHistPtRecTruePrimaryProtonMinus, // Pt histo for protons -, reconstructed tracks, true ID, primary event
+ kHistPtRecTruePrimaryKaonMinus, // Pt histo for kaons -, reconstructed tracks, true ID, primary event
+ kHistPtRecTruePrimaryPionMinus, // Pt histo for pions -, reconstructed tracks, true ID, primary event
+
+ // Rest
+ kHistPtRec, // Pt histo for all particles, reconstructed tracks
+ kHistPtGen, // Pt histo for all particles, generated tracks
+ kNPtHist = kHistPtGen, // Number of pt-likehistos histos
+ kHistPIDTPC, // Particle Identification histo
+ kHistPIDTPCPt,
+ kNHist, // Total number of histos
+ }; // Type of events plotted in Pt Histogram
+
+ const char * kHistName[] =
+ {
+ // 6 Pt Reconstructed Sigma
+ "histPtRecSigmaProtonPlus", // Pt histo for protons +, reconstructed tracks, sigma ID
+ "histPtRecSigmaKaonPlus", // Pt histo for kaons +, reconsructed tracks, sigma ID
+ "histPtRecSigmaPionPlus", // Pt histo for pions +, reconstructed tracks, sigma ID
+ "histPtRecSigmaProtonMinus", // Pt histo for protons -, reconstructed tracks, sigma ID
+ "histPtRecSigmaKaonMinus", // Pt histo for kaons -, reconstructed tracks, sigma ID
+ "histPtRecSigmaPionMinus", // Pt histo for pions -, reconstructed tracks, sigma ID
+
+ // 6 Pt Reconstructed True
+ "histPtRecTrueProtonPlus", // Pt histo for protons +, reconstructed tracks, true ID
+ "histPtRecTrueKaonPlus", // Pt histo for kaons +, reconsructed tracks, true ID
+ "histPtRecTruePionPlus", // Pt histo for pions +, reconstructed tracks, true ID
+ "histPtRecTrueProtonMinus", // Pt histo for protons -, reconstructed tracks, true ID
+ "histPtRecTrueKaonMinus", // Pt histo for kaons -, reconstructed tracks, true ID
+ "histPtRecTruePionMinus", // Pt histo for pions -, reconstructed tracks, true ID
+
+ // 6 Pt Reconstructed Sigma Primary
+ "histPtRecSigmaPrimaryProtonPlus", // Pt histo for protons +, reconstructed tracks, sigma ID, primary Event
+ "histPtRecSigmaPrimaryKaonPlus", // Pt histo for kaons +, reconsructed tracks, sigma ID, primary Event
+ "histPtRecSigmaPrimaryPionPlus", // Pt histo for pions +, reconstructed tracks, sigma ID, primary Event
+ "histPtRecSigmaPrimaryProtonMinus", // Pt histo for protons -, reconstructed tracks, sigma ID, primary Event
+ "histPtRecSigmaPrimaryKaonMinus", // Pt histo for kaons -, reconstructed tracks, sigma ID, primary Event
+ "histPtRecSigmaPrimaryPionMinus", // Pt histo for pions -, reconstructed tracks, sigma ID, primary Event
+
+ // 6 Pt Reconstructed Sigma Seconday
+ "histPtRecSigmaSecondaryProtonPlus", // Pt histo for protons +, reconstructed tracks, sigma ID, secondary Event
+ "histPtRecSigmaSecondaryKaonPlus", // Pt histo for kaons +, reconsructed tracks, sigma ID, secondary Event
+ "histPtRecSigmaSecondaryPionPlus", // Pt histo for pions +, reconstructed tracks, sigma ID, secondary Event
+ "histPtRecSigmaSecondaryProtonMinus", // Pt histo for protons -, reconstructed tracks, sigma ID, secondary Event
+ "histPtRecSigmaSecondaryKaonMinus", // Pt histo for kaons -, reconstructed tracks, sigma ID, secondary Event
+ "histPtRecSigmaSecondaryPionMinus", // Pt histo for pions -, reconstructed tracks, sigma ID, secondary Event
+
+ // 6 Pt Reconstructed Sigma Primary
+ "histPtGenTruePrimaryProtonPlus", // Pt histo for protons +, generated tracks, sigma ID, primary Event
+ "histPtGenTruePrimaryKaonPlus", // Pt histo for kaons +, generated tracks, sigma ID, primary Event
+ "histPtGenTruePrimaryPionPlus", // Pt histo for pions +, generated tracks, sigma ID, primary Event
+ "histPtGenTruePrimaryProtonMinus", // Pt histo for protons -, generated tracks, sigma ID, primary Event
+ "histPtGenTruePrimaryKaonMinus", // Pt histo for kaons -, generated tracks, sigma ID, primary Event
+ "histPtGenTruePrimaryPionMinus", // Pt histo for pions -, generated tracks, sigma ID, primary Event
+
+ // 6 Pt Reconstructed True
+ "histPtRecTruePrimaryProtonPlus", // Pt histo for protons +, reconstructed tracks, true ID, primary event
+ "histPtRecTruePrimaryKaonPlus", // Pt histo for kaons +, reconsructed tracks, true ID, primary event
+ "histPtRecTruePrimaryPionPlus", // Pt histo for pions +, reconstructed tracks, true ID, primary event
+ "histPtRecTruePrimaryProtonMinus", // Pt histo for protons -, reconstructed tracks, true ID, primary event
+ "histPtRecTruePrimaryKaonMinus", // Pt histo for kaons -, reconstructed tracks, true ID, primary event
+ "histPtRecTruePrimaryPionMinus", // Pt histo for pions -, reconstructed tracks, true ID, primary event
+
+ // Rest
+ "histPtRec", // Pt histo for all particles, reconstructed tracks
+ "histPtGen", // Pt histo for all particles, generated tracks
+ "histPIDTPC", // Particle Identification histo
+ "histPIDTPCPt",
+ };
+
+ enum AODParticleSpecies_t
+ {
+ kProton = 0,
+ kKaon,
+ kPion,
+ kNSpecies = kPion,
+ }; // Particle species used in plotting
+
+ const char * kParticleSpecies[] =
+ {
+ "ProtonPlus",
+ "ProtonMinus",
+ "KaonPlus",
+ "KaonMinus",
+ "PionPlus",
+ "PionMinus",
+ };
+ enum AODHistoType_t
+ {
+ kHistSpectraRec = 0,
+ kHistSpectraGen,
+ kHNHistoTypes = kHistSpectraGen,
+ }; // Types of histos
+
+}
+
+using namespace AliSpectraNameSpace;
+
+class AliSpectraAODHistoManager : public TNamed
+{
+public:
+ AliSpectraAODHistoManager() : TNamed(), fOutputList(0) {}
+ AliSpectraAODHistoManager(const char *name);
+ virtual ~AliSpectraAODHistoManager() {}
+
+
+ TH1F* BookPtHistogram(const char * name);
+ TH2F* BookPIDHistogram(const char * name);
+ TH1* GetHistogram(UInt_t id) { return (TH1*) fOutputList->At(id); }
+ TH1* GetPtHistogram(UInt_t id) { return (TH1*) fOutputList->At(id); }
+ TH1* GetPtHistogram(const char * name) { return (TH1*) fOutputList->FindObject(name); }
+ TH1* GetPtHistogramByName(UInt_t id) { return (TH1*) fOutputList->FindObject(kHistName[id]); } // Use this if you want to read a file saved with a different histo list
+ TH2* GetPIDHistogram(UInt_t id) { return (TH2*) fOutputList->At(id); }
+ TH2* GetPIDHistogram(const char * name) { return (TH2*) fOutputList->FindObject(name); }
+ TH2* GetPIDHistogramByName(UInt_t id) { return (TH2*) fOutputList->FindObject(kHistName[id]); }// Use this if you want to read a file saved with a different histo list
+ TH1F* GetTH1F(UInt_t id) { return (TH1F*) GetPtHistogram(id); }
+ TH2F* GetTH2F(UInt_t id) { return (TH2F*) GetPIDHistogram(id); }
+
+ TList * GetOutputList() {return fOutputList;}
+
+ Long64_t Merge(TCollection* list);
+
+
+private:
+ TList *fOutputList; // List of Pt Histo's
+
+ AliSpectraAODHistoManager(const AliSpectraAODHistoManager&);
+ AliSpectraAODHistoManager& operator=(const AliSpectraAODHistoManager&);
+
+ ClassDef(AliSpectraAODHistoManager, 1);
+
+};
+#endif
+
--- /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. *
+ **************************************************************************/
+
+//-----------------------------------------------------------------
+// AliSpectraAODTrackCuts class
+//-----------------------------------------------------------------
+
+#include "TChain.h"
+#include "TTree.h"
+#include "TLegend.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TCanvas.h"
+#include "AliAnalysisTask.h"
+#include "AliAnalysisManager.h"
+#include "AliAODTrack.h"
+#include "AliAODMCParticle.h"
+#include "AliAODEvent.h"
+#include "AliAODInputHandler.h"
+#include "AliAnalysisTaskESDfilter.h"
+#include "AliAnalysisDataContainer.h"
+#include "AliSpectraAODTrackCuts.h"
+#include "AliSpectraAODHistoManager.h"
+#include <iostream>
+
+using namespace std;
+
+ClassImp(AliSpectraAODTrackCuts)
+
+
+AliSpectraAODTrackCuts::AliSpectraAODTrackCuts(const char *name) : TNamed(name, "AOD Track Cuts"), fIsSelected(0), fTrackBits(0), fEtaCut(0), fDCACut(0), fPCut(0), fPtCut(0), fHistoCuts(0), fTrack(0)
+
+{
+ // Constructor
+ fHistoCuts = new TH1I("fTrkCuts", "Track Cuts", kNTrkCuts, -0.5, kNTrkCuts - 0.5);
+ fEtaCut = 100000.0; // default value of eta cut ~ no cut
+ fDCACut = 100000.0; // default value of dca cut ~ no cut
+ fPCut = 100000.0; // default value of p cut ~ no cut
+ fPtCut = 100000.0; // default value of pt cut ~ no cut
+
+}
+
+//_______________________________________________________
+Bool_t AliSpectraAODTrackCuts::IsSelected(AliAODTrack * track)
+{
+// Returns true if Track Cuts are selected and applied
+ if (!track)
+ {
+ printf("ERROR: Could not receive track");
+ return kFALSE;
+ }
+ fTrack = track;
+ fIsSelected = (CheckTrackType() && CheckEtaCut() && CheckDCACut() && CheckPCut() && CheckPtCut());
+ return fIsSelected ;
+}
+//_________________________________________________________
+
+Bool_t AliSpectraAODTrackCuts::CheckTrackType()
+{
+ // Check track cuts
+ if (fTrack->TestFilterBit(fTrackBits)) return kTRUE;
+ fHistoCuts->Fill(kTrkBit);
+ return kFALSE;
+}
+//________________________________________________________
+Bool_t AliSpectraAODTrackCuts::CheckEtaCut()
+{
+ // Check eta cut
+ if (fTrack->Eta() < fEtaCut && fTrack->Eta() > - fEtaCut) return kTRUE;
+ fHistoCuts->Fill(kTrkEta);
+ return kFALSE;
+}
+//_______________________________________________________
+Bool_t AliSpectraAODTrackCuts::CheckDCACut()
+{
+ // Check DCA cut
+ if (fTrack->DCA() < fDCACut) return kTRUE;
+ fHistoCuts->Fill(kTrkDCA);
+ return kFALSE;
+}
+//________________________________________________________
+Bool_t AliSpectraAODTrackCuts::CheckPCut()
+{
+ // Check P cut
+ if (fTrack->P() < fPCut) return kTRUE;
+ fHistoCuts->Fill(kTrkP);
+ return kFALSE;
+}
+//_______________________________________________________
+Bool_t AliSpectraAODTrackCuts::CheckPtCut()
+{
+ // check Pt cut
+// if ((fTrack->Pt() < fPtCut) && (fTrack->Pt() > 0.3 )) return kTRUE;
+ if (fTrack->Pt() < fPtCut) return kTRUE;
+ fHistoCuts->Fill(kTrkPt);
+ return kFALSE;
+}
+//_______________________________________________________
+void AliSpectraAODTrackCuts::PrintCuts() const
+{
+ // Print cuts
+ cout << "Track Cuts" << endl;
+ cout << " > TrackBit\t" << fTrackBits << endl;
+ cout << " > Eta cut\t" << fEtaCut << endl;
+ cout << " > DCA cut\t" << fDCACut << endl;
+ cout << " > P cut\t" << fPCut << endl;
+ cout << " > Pt cut \t" << fPtCut << endl;
+}
+//_______________________________________________________
+void AliSpectraAODTrackCuts::SetTrackType(UInt_t bit)
+{
+ // Set the type of track to be used. The argument should be the bit number. The mask is produced automatically.
+ fTrackBits = (0x1 << (bit - 1));
+}
+//_______________________________________________________
+
+Long64_t AliSpectraAODTrackCuts::Merge(TCollection* list)
+{
+ // Merge a list of AliSpectraAODTrackCuts objects with this.
+ // Returns the number of merged objects (including this).
+
+ // AliInfo("Merging");
+
+ if (!list)
+ return 0;
+
+ if (list->IsEmpty())
+ return 1;
+
+ TIterator* iter = list->MakeIterator();
+ TObject* obj;
+
+ // collections of all histograms
+ TList collections;
+
+ Int_t count = 0;
+
+ while ((obj = iter->Next())) {
+ AliSpectraAODTrackCuts* entry = dynamic_cast<AliSpectraAODTrackCuts*> (obj);
+ if (entry == 0)
+ continue;
+
+ TH1I * histo = entry->GetHistoCuts();
+ collections.Add(histo);
+ count++;
+ }
+
+ fHistoCuts->Merge(&collections);
+
+ delete iter;
+
+ return count+1;
+}
+
--- /dev/null
+#ifndef ALISPECTRAAODTRACKCUTS_H
+#define ALISPECTRAAODTRACKCUTS_H
+
+/* See cxx source for full Copyright notice */
+
+//-------------------------------------------------------------------------
+// AliSpectraAODTrackCuts
+//
+//
+//
+//
+// Authors: Michele Floris, CERN, Philip Versteeg, UU, Redmer Bertens, UU
+//-------------------------------------------------------------------------
+
+class AliAODEvent;
+class AliSpectraAODHistoManager;
+
+#include "TNamed.h"
+
+class AliSpectraAODTrackCuts : public TNamed
+{
+public:
+
+ enum { kTrkBit = 0, kTrkEta, kTrkDCA, kTrkP, kTrkPt, kNTrkCuts};
+
+
+ AliSpectraAODTrackCuts() : TNamed(), fIsSelected(0), fTrackBits(0), fEtaCut(0), fPCut(0), fPtCut(0), fHistoCuts(0), fTrack(0) {}
+
+ AliSpectraAODTrackCuts(const char *name);
+ virtual ~AliSpectraAODTrackCuts() {} // To be implemented
+
+ Bool_t IsSelected(AliAODTrack * track);
+
+ void SetTrackType(UInt_t bit);
+ Bool_t CheckTrackType();
+ Bool_t CheckEtaCut();
+ Bool_t CheckDCACut();
+ Bool_t CheckPCut();
+ Bool_t CheckPtCut();
+ void PrintCuts() const;
+
+ UInt_t GetTrackType() const { return fTrackBits; }
+ TH1I * GetHistoCuts() { return fHistoCuts; }
+ void SetEta(Float_t eta) { fEtaCut = eta; }
+ void SetDCA(Float_t dca) { fDCACut = dca; }
+ void SetP(Float_t p) { fPCut = p; }
+ void SetPt(Float_t pt) { fPtCut = pt; }
+ Float_t GetEta() const { return fEtaCut; }
+ Float_t GetDCA() const { return fDCACut; }
+ Float_t GetP() const { return fPCut; }
+ Float_t GetPt() const { return fPtCut; }
+
+ Long64_t Merge(TCollection* list);
+
+
+private:
+
+ Bool_t fIsSelected; // True if cuts are selected
+ UInt_t fTrackBits; // Type of track to be used
+ Float_t fEtaCut; // Allowed absolute maximum value of Eta
+ Float_t fDCACut; // Maximum value of DCA
+ Float_t fPCut; // Maximum value of P
+ Float_t fPtCut; // Maximum value of Pt
+
+ TH1I * fHistoCuts; // Cuts statistics
+ AliAODTrack * fTrack; //! Track pointer
+
+ AliSpectraAODTrackCuts(const AliSpectraAODTrackCuts&);
+ AliSpectraAODTrackCuts& operator=(const AliSpectraAODTrackCuts&);
+
+ ClassDef(AliSpectraAODTrackCuts, 1);
+};
+#endif
+
--- /dev/null
+This code was developed in Summer 2011 by 2 summer students (Philip Versteeg, UU, Redmer Bertens, UU) under the supervision of Michele Floris. It can serve as a starting point for AOD analysis.
--- /dev/null
+#include "AliAODTrack.h"
+#include <iostream>
+#include "TFile.h"
+#include "TH1I.h"
+#include "TH2.h"
+#include "AliSpectraAODHistoManager.h"
+#include "AliSpectraAODEventCuts.h"
+#include "AliSpectraAODTrackCuts.h"
+#include "TCanvas.h"
+#include "TLegend.h"
+
+using namespace std;
+
+void analysis_macro(const char * dataFile = "Pt.AOD.1._data_ptcut.root", const char * mcFile ="Pt.AOD.1._MC.root") {
+
+ // gSystem->Load("libTree.so");
+ // gSystem->Load("libGeom.so");
+ // gSystem->Load("libVMC.so");
+ // gSystem->Load("libPhysics.so");
+ // gSystem->Load("libSTEERBase.so");
+ // gSystem->Load("libESD.so");
+ // gSystem->Load("libAOD.so");
+ // gSystem->Load("libANALYSIS.so");
+ // gSystem->Load("libANALYSISalice.so");
+ // gSystem->Load("libANALYSIS");
+ // gSystem->Load("libANALYSISalice");
+ // gROOT->LoadMacro("$ALICE_ROOT/PWG0/CreateESDChain.C");
+ // gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+ // gStyle->SetPalette(1);
+ // gStyle->SetFillColor(kWhite);
+
+ // gROOT->LoadMacro("AliSpectraAODTrackCuts.cxx+g");
+ // gROOT->LoadMacro("AliSpectraAODEventCuts.cxx+g");
+ // gROOT->LoadMacro("AliSpectraAODHistoManager.cxx+g");
+ // gROOT->LoadMacro("AliAnalysisTaskSpectraAOD.cxx+g");
+
+ const char * histoname[] =
+ {
+ // names of histos we want to draw. please observe order of histos, elsewise name will not make sense
+ "#frac{MC_{#sigma, rec}}{Data}, P^{+}",
+ "#frac{MC_{#sigma, rec}}{Data}, K^{+}",
+ "#frac{MC_{#sigma, rec}}{Data}, #pi^{+}",
+ "#frac{MC_{#sigma, rec}}{Data}, P^{-}",
+ "#frac{MC_{#sigma, rec}}{Data}, K^{-}",
+ "#frac{MC_{#sigma, rec}}{Data}, #pi^{-}",
+
+ "#frac{MC_{#sigma, rec, prim}}{Data}, P^{+}",
+ "#frac{MC_{#sigma, rec, prim}}{Data}, K^{+}",
+ "#frac{MC_{#sigma, rec, prim}}{Data}, #pi^{+}",
+ "#frac{MC_{#sigma, rec, prim}}{Data}, P^{-}",
+ "#frac{MC_{#sigma, rec, prim}}{Data}, K^{-}",
+ "#frac{MC_{#sigma, rec, prim}}{Data}, #pi^{-}",
+
+ "#frac{MC_{#sigma, rec, sec}}{Data}, P^{+}",
+ "#frac{MC_{#sigma, rec, sec}}{Data}, K^{+}",
+ "#frac{MC_{#sigma, rec, sec}}{Data}, #pi^{+}",
+ "#frac{MC_{#sigma, rec, sec}}{Data}, P^{-}",
+ "#frac{MC_{#sigma, rec, sec}}{Data}, K^{-}",
+ "#frac{MC_{#sigma, rec, sec}}{Data}, #pi^{-}",
+
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}}, P^{+}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}}, K^{+}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}}, #pi^{+}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}}, P^{-}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}}, K^{-}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}}, #pi^{-}",
+
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}} x Data, P^{+}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}} x Data, K^{+}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}} x Data, #pi^{+}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}} x Data, P^{-}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}} x Data, K^{-}",
+ "#frac{MC_{true, gen}}{MC_{#sigma, rec}} x Data, #pi^{-}",
+
+ "Data, P^{+}",
+ "Data, K^{+}",
+ "Data, #pi^{+}",
+ "Data, P^{-}",
+ "Data, K^{-}",
+ "Data, #pi^{-}",
+
+ "#frac{MC_{rec #sigma}}{MC_{true, rec}}, P^{+}",
+ "#frac{MC_{rec #sigma}}{MC_{true, rec}}, K^{+}",
+ "#frac{MC_{rec #sigma}}{MC_{true, rec}}, #pi^{+}",
+ "#frac{MC_{rec #sigma}}{MC_{true, rec}}, P^{-}",
+ "#frac{MC_{rec #sigma}}{MC_{true, rec}}, K^{-}",
+ "#frac{MC_{rec #sigma}}{MC_{true, rec}}, #pi^{-}",
+
+ "#frac{MC_{rec #sigma sec}}{MC_{#sigma, rec}}, P^{+}",
+ "#frac{MC_{rec #sigma sec}}{MC_{#sigma, rec}}, K^{+}",
+ "#frac{MC_{rec #sigma sec}}{MC_{#sigma, rec}}, #pi^{+}",
+ "#frac{MC_{rec #sigma sec}}{MC_{#sigma, rec}}, P^{-}",
+ "#frac{MC_{rec #sigma sec}}{MC_{#sigma, rec}}, K^{-}",
+ "#frac{MC_{rec #sigma sec}}{MC_{#sigma, rec}}, #pi^{-}",
+ };
+
+ const char * histotitle[] =
+ {
+ // titles of histos we want to draw. please observe order of histos, elsewise name will not make sense
+ "MC #sigma, rec / Data, P^{+}",
+ "MC #sigma, rec / Data, K^{+}",
+ "MC #sigma, rec / Data, #pi^{+}",
+ "MC #sigma, rec / Data, P^{-}",
+ "MC #sigma, rec / Data, K^{-}",
+ "MC #sigma, rec / Data, #pi^{-}",
+
+ "MC #sigma, rec, prim / Data, P^{+}",
+ "MC #sigma, rec, prim / Data, K^{+}",
+ "MC #sigma, rec, prim / Data, #pi^{+}",
+ "MC #sigma, rec, prim / Data, P^{-}",
+ "MC #sigma, rec, prim / Data, K^{-}",
+ "MC #sigma, rec, prim / Data, #pi^{-}",
+
+ "MC #sigma, rec, sec / Data, P^{+}",
+ "MC #sigma, rec, sec / Data, K^{+}",
+ "MC #sigma, rec, sec / Data, #pi^{+}",
+ "MC #sigma, rec, sec / Data, P^{-}",
+ "MC #sigma, rec, sec / Data, K^{-}",
+ "MC #sigma, rec, sec / Data, #pi^{-}",
+
+ "Correction factor, P^{+}",
+ "Correction factor, K^{+}",
+ "Correction factor, #pi^{+}",
+ "Correction factor, P^{-}",
+ "Correction factor, K^{-}",
+ "Correction factor, #pi^{-}",
+
+ "True Data, P^{+}",
+ "True Data, K^{+}",
+ "True Data, #pi^{+}",
+ "True Data, P^{-}",
+ "True Data, K^{-}",
+ "True Data, #pi^{-}",
+
+ "Raw Data, P^{+}",
+ "Raw Data, K^{+}",
+ "Raw Data, #pi^{+}",
+ "Raw Data, P^{-}",
+ "Raw Data, K^{-}",
+ "Raw Data, #pi^{-}",
+
+ "Identificaiton quality, P^{+}",
+ "Identificaiton quality, K^{+}",
+ "Identificaiton quality, #pi^{+}",
+ "Identificaiton quality, P^{-}",
+ "Identificaiton quality, K^{-}",
+ "Identificaiton quality, #pi^{-}",
+
+ "Contribution of secondaries, P^{+}",
+ "Contribution of secondaries, K^{+}",
+ "Contribution of secondaries, #pi^{+}",
+ "Contribution of secondaries, P^{-}",
+ "Contribution of secondaries, K^{-}",
+ "Contribution of secondaries, #pi^{-}",
+ };
+
+ // Open root MC file and get classes
+ cout << "Analysis Macro" << endl;
+ cout << " > Reading MC data" << endl;
+ TFile *_mc = TFile::Open(mcFile);
+ AliSpectraAODHistoManager* hman_mc = (AliSpectraAODHistoManager*) _mc->Get("SpectraHistos");
+ AliSpectraAODEventCuts* ecuts_mc = (AliSpectraAODEventCuts*) _mc->Get("Event Cuts");
+ AliSpectraAODTrackCuts* tcuts_mc = (AliSpectraAODTrackCuts*) _mc->Get("Track Cuts");
+ // print info about mc track and Event cuts
+ cout << " -- Info about MC -- "<< endl;
+ ecuts_mc->PrintCuts();
+ tcuts_mc->PrintCuts();
+ // get the mc histo's necessary to write the contamination histos from the rootfile
+ TH1F* spectrahistos_mc[30];
+ cout << " -- Reading and normalizing histograms -- " << endl;
+ for( Int_t i = 0 ; i <= AliSpectraNameSpace::kNPtSpecies ; i++ ) {
+ spectrahistos_mc[i] = dynamic_cast<TH1F*>(hman_mc->GetHistogram((AliSpectraNameSpace::AODPtHist_t)i));
+ Double_t events_mc = 1. / ecuts_mc->NumberOfEvents();
+ if(!spectrahistos_mc[i]) {
+ continue;
+ }
+ spectrahistos_mc[i]->Scale(events_mc, "width");
+ }
+ cout << " > Reading real data " << endl;
+ // proceed likewise for data
+ TFile *_data = TFile::Open(dataFile);
+ AliSpectraAODHistoManager* hman_data = (AliSpectraAODHistoManager*) _data->Get("SpectraHistos");
+ AliSpectraAODEventCuts* ecuts_data = (AliSpectraAODEventCuts*) _data->Get("Event Cuts");
+ AliSpectraAODTrackCuts* tcuts_data = (AliSpectraAODTrackCuts*) _data->Get("Track Cuts");
+ // print info about track and Event cuts
+ cout << " -- Info about data -- " << endl;
+ ecuts_data->PrintCuts();
+ tcuts_data->PrintCuts();
+ // get the histo's necessary to write the contamination histos from the rootfile
+ TH1F* spectrahistos_data[30];
+ cout << " -- Reading and normalizing histograms -- " << endl;
+ for( Int_t i = 0 ; i <= AliSpectraNameSpace::kNPtSpecies ; i++ ) {
+ spectrahistos_data[i] = dynamic_cast<TH1F*>(hman_data->GetHistogram((AliSpectraNameSpace::AODPtHist_t)i));
+ if(!spectrahistos_data[i]) {
+ continue;
+ }
+ Double_t events_data = 1. / ecuts_data->NumberOfEvents();
+ //normalize the histos
+ spectrahistos_data[i]->Scale(events_data, "width");
+ }
+ //create output file and select it as the active file
+ TFile * output = new TFile("analysis_output.root", "RECREATE");
+ output->cd();
+ // Write the scaled data and MC histos in the output file for convenience
+ output->mkdir("Data");
+ output->cd("Data");
+ for(Int_t ihist = 0; ihist <= AliSpectraNameSpace::kNPtSpecies; ihist++){
+ spectrahistos_data[ihist]->Write();
+ }
+ output->mkdir("MC");
+ output->cd("MC");
+ for(Int_t ihist = 0; ihist <= AliSpectraNameSpace::kNPtSpecies; ihist++){
+ spectrahistos_mc[ihist]->Write();
+ }
+ output->mkdir("ANALYSIS");
+ output->cd("ANALYSIS");
+ // perform some analysis!
+ // next loop: MC_Sigma_Rec over Rec_Data
+ TCanvas* c1 = new TCanvas("MC reconstructed over Data","MC reconstructed over Data");
+ c1->Divide(3,2);
+ TH1F* kRatioPtMCSigmaRecSigmaRec[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kRatioPtMCSigmaRecSigmaRec[i] = (TH1F*)(spectrahistos_mc[i])->Clone(histotitle[i]);
+ kRatioPtMCSigmaRecSigmaRec[i]->Divide(spectrahistos_data[i]);
+ c1->cd(1+i);
+ kRatioPtMCSigmaRecSigmaRec[i]->SetTitle(histoname[i]);
+ kRatioPtMCSigmaRecSigmaRec[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kRatioPtMCSigmaRecSigmaRec[i]->Draw();
+ kRatioPtMCSigmaRecSigmaRec[i]->Write();
+ }
+ // next loop: MC_sigma_primaries over Rec_data
+ TCanvas* c2 = new TCanvas("MC reconstructed primaries over Data", "MC reconstructed primaries over Data");
+ c2->Divide(3,2);
+ TH1F* kRatioPtMCSigmaRecPrimarySigmaRec[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kRatioPtMCSigmaRecPrimarySigmaRec[i] = (TH1F*)(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaPrimaryProtonPlus])->Clone(histotitle[i+6]);
+ kRatioPtMCSigmaRecPrimarySigmaRec[i]->Divide(spectrahistos_data[i]);
+ c2->cd(1+i);
+ kRatioPtMCSigmaRecPrimarySigmaRec[i]->SetTitle(histoname[i+6]);
+ kRatioPtMCSigmaRecPrimarySigmaRec[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kRatioPtMCSigmaRecPrimarySigmaRec[i]->Draw();
+ kRatioPtMCSigmaRecPrimarySigmaRec[i]->Write();
+ }
+ // next loop: MC_Sigma_secondaries over Rec_data
+ TCanvas* c3 = new TCanvas("MC reconstructed secondaries over Data","MC reconstructed secondaries over Data");
+ c3->Divide(3,2);
+ TH1F* kRatioPtMCSigmaRecSigmaRecSecondary[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kRatioPtMCSigmaRecSigmaRecSecondary[i] = (TH1F*)(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaSecondaryProtonPlus])->Clone(histotitle[i+12]);
+ kRatioPtMCSigmaRecSigmaRecSecondary[i]->Divide(spectrahistos_data[i]);
+ c3->cd(1+i);
+ kRatioPtMCSigmaRecSigmaRecSecondary[i]->SetTitle(histoname[12+i]);
+ kRatioPtMCSigmaRecSigmaRecSecondary[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kRatioPtMCSigmaRecSigmaRecSecondary[i]->Draw();
+ kRatioPtMCSigmaRecSigmaRecSecondary[i]->Write();
+ }
+ //
+ TCanvas* c4 = new TCanvas("MC true generated over MC true reconstructed","MC true generated over MC true reconstructed");
+ c4->Divide(3,2);
+ TH1F* kRatioPtMCTrueMCRec[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kRatioPtMCTrueMCRec[i] = (TH1F*)(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtGenTruePrimaryProtonPlus])->Clone(histotitle[i+18]);
+ kRatioPtMCTrueMCRec[i]->Divide(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaPrimaryProtonPlus]);
+ c4->cd(1+i);
+ kRatioPtMCTrueMCRec[i]->SetTitle(histoname[18+i]);
+ kRatioPtMCTrueMCRec[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kRatioPtMCTrueMCRec[i]->Draw();
+ kRatioPtMCTrueMCRec[i]->Write();
+ }
+
+ TCanvas* c5 = new TCanvas("True Data","True Data");
+ c5->Divide(3,2);
+ TH1F* kRatioPtTrueData[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kRatioPtTrueData[i] = (TH1F*)(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtGenTruePrimaryProtonPlus])->Clone(histotitle[i+24]);
+ kRatioPtTrueData[i]->Divide(spectrahistos_mc[i]); // Michele: Since we are correcting all contributions at once, this should be exaclty the same as in the data!
+ kRatioPtTrueData[i]->Multiply(spectrahistos_data[i]);
+ c5->cd(1+i);
+ kRatioPtTrueData[i]->SetTitle(histoname[24+i]);
+ kRatioPtTrueData[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kRatioPtTrueData[i]->Draw();
+ kRatioPtTrueData[i]->Write();
+ }
+
+ TCanvas* c6 = new TCanvas("Raw data", "Raw data");
+ c6->Divide(3,2);
+ TH1F* kRawData[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kRawData[i] = (TH1F*)(spectrahistos_data[i+AliSpectraNameSpace::kHistPtRecSigmaProtonPlus])->Clone(histotitle[i+30]);
+ c6->cd(1+i);
+ kRawData[i]->SetTitle(histoname[30+i]);
+ kRawData[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kRawData[i]->Draw();
+ kRawData[i]->Write();
+ }
+
+ TCanvas* c7 = new TCanvas("contamination", "contamination");
+ c7->Divide(3,2);
+ TH1F* kContamination[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kContamination[i] = (TH1F*)(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaProtonPlus])->Clone(histotitle[i+36]);
+ kContamination[i]->Divide(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecTrueProtonPlus]);
+ c7->cd(1+i);
+ kContamination[i]->SetTitle(histoname[36+i]);
+ kContamination[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kContamination[i]->Draw();
+ kContamination[i]->Write();
+ }
+
+ TCanvas* c8 = new TCanvas("contamination secondary", "contamination secondary");
+ c8->Divide(3,2);
+ TH1F* kSecContamination[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ kSecContamination[i] = (TH1F*)(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaSecondaryProtonPlus])->Clone(histotitle[i+42]);
+ kSecContamination[i]->Divide(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaProtonPlus]);
+ c8->cd(1+i);
+ kSecContamination[i]->SetTitle(histoname[42+i]);
+ kSecContamination[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kSecContamination[i]->Draw();
+ kSecContamination[i]->Write();
+ }
+
+ TCanvas* c9 = new TCanvas("true data and generated", "true data and generated");
+ c9->Divide(3,2);
+ TH1F* ktrueandgenerated[6];
+ for (Int_t i = 0 ; i < 6 ; i ++ )
+ {
+ // kSecContamination[i] = (TH1F*)(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaSecondaryProtonPlus])->Clone(histotitle[i+42]);
+ // kSecContamination[i]->Divide(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtRecSigmaProtonPlus]);
+ c9->cd(1+i);
+ // kSecContamination[i]->SetTitle(histoname[42+i]);
+ // kSecContamination[i]->GetYaxis()->SetTitle("#frac{d^{2} N}{dy dp_{T}} (c / GeV)");
+ kRatioPtTrueData[i]->Draw();
+ spectrahistos_mc[i+AliSpectraNameSpace::kHistPtGenTruePrimaryProtonPlus]->SetMarkerColor(kRed);
+ spectrahistos_mc[i+AliSpectraNameSpace::kHistPtGenTruePrimaryProtonPlus]->SetMarkerStyle(kOpenCircle);
+ spectrahistos_mc[i+AliSpectraNameSpace::kHistPtGenTruePrimaryProtonPlus]->SetLineColor(kRed);
+ spectrahistos_mc[i+AliSpectraNameSpace::kHistPtGenTruePrimaryProtonPlus]->Draw("same");
+
+ // legend
+ TLegend * leg = new TLegend(0.6,0.64,0.89,0.79); //coordinates are fractions
+ leg->AddEntry(spectrahistos_mc[i+AliSpectraNameSpace::kHistPtGenTruePrimaryProtonPlus],"Generated","p");
+ leg->AddEntry(kRatioPtTrueData[i],"Corrected","p");
+ leg->Draw();
+ }
+
+ cout << "Analysis complete" << endl;
+ cout << " - Please exit root to close and save output file " << endl;
+}
--- /dev/null
+void run()
+{
+
+ gSystem->Load("libTree.so");
+ gSystem->Load("libGeom.so");
+ gSystem->Load("libVMC.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libSTEERBase.so");
+ gSystem->Load("libESD.so");
+ gSystem->Load("libAOD.so");
+ gSystem->Load("libANALYSIS.so");
+ gSystem->Load("libANALYSISalice.so");
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libANALYSISalice");
+ gROOT->LoadMacro("$ALICE_ROOT/PWG0/CreateESDChain.C");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+ gStyle->SetPalette(1);
+ gStyle->SetFillColor(kWhite);
+
+ gROOT->LoadMacro("AliSpectraAODTrackCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODEventCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODHistoManager.cxx+g");
+ gROOT->LoadMacro("AliAnalysisTaskSpectraAOD.cxx+g");
+
+
+
+
+
+
+
+
+}
--- /dev/null
+void runAOD(Bool_t isMC = 0)
+{
+
+//gEnv->SetValue("XSec.GSI.DelegProxy", "2");
+//TProof::Open("pverstee@alice-caf.cern.ch");
+
+ gSystem->Load("libTree.so");
+ gSystem->Load("libGeom.so");
+ gSystem->Load("libVMC.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libSTEERBase.so");
+ gSystem->Load("libESD.so");
+ gSystem->Load("libAOD.so");
+ gSystem->Load("libANALYSIS.so");
+ gSystem->Load("libANALYSISalice.so");
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libANALYSISalice");
+ gROOT->LoadMacro("$ALICE_ROOT/PWG0/CreateESDChain.C");
+
+ ifstream flist;
+ flist.open("filelist.txt");
+ TChain * chain = new TChain("aodTree");
+ TString line;
+ while (line.ReadLine(flist))
+ {
+ gSystem->ExpandPathName(line);
+ chain->AddFile(line.Data());
+ cout << "Adding file " << line.Data() << endl;
+ }
+// gProof->UploadPackage("AF-v4-19-04-AN");
+// gProof->EnablePackage("AF-v4-19-04-AN");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+ gStyle->SetPalette(1);
+
+
+
+ AliAnalysisManager *mgr = new AliAnalysisManager("testAnalysis");
+ AliAODInputHandler* aodH = new AliAODInputHandler();
+ mgr->SetInputEventHandler(aodH);
+
+ gROOT->LoadMacro("AliSpectraAODTrackCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODVertexCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODHistoManager.cxx+g");
+ gROOT->LoadMacro("AliAnalysisTaskSpectraAOD.cxx+g");
+
+ gROOT->LoadMacro("$(ALICE_ROOT)/ANALYSIS/macros/AddTaskPIDResponse.C");
+ AliAnalysisTask * taskPID = AddTaskPIDResponse();
+ mgr->AddTask(taskPID);
+
+
+ AliAnalysisTaskSpectraAOD *task = new AliAnalysisTaskSpectraAOD("TaskAODExercise");
+ task->SetIsMC(isMC);
+ mgr->AddTask(task);
+
+ // Set the cuts
+ AliSpectraAODVertexCuts * vcuts = new AliSpectraAODVertexCuts("Vertex Cuts");
+ AliSpectraAODTrackCuts * tcuts = new AliSpectraAODTrackCuts("Tracks Cuts");
+ tcuts->SetTrackType(6);
+ tcuts->SetEta(1.);
+ task->SetVertexCuts(vcuts);
+ task->SetTrackCuts(tcuts);
+
+
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ AliAnalysisDataContainer *coutputpt1 = mgr->CreateContainer("chistpt", AliSpectraAODHistoManager::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1.root");
+ AliAnalysisDataContainer *coutputpt2 = mgr->CreateContainer("cvcutpt", AliSpectraAODVertexCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1.root");
+ AliAnalysisDataContainer *coutputpt3 = mgr->CreateContainer("ctcutpt", AliSpectraAODTrackCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1.root");
+
+ mgr->ConnectInput(task, 0, cinput);
+ mgr->ConnectOutput(task, 1, coutputpt1);
+ mgr->ConnectOutput(task, 2, coutputpt2);
+ mgr->ConnectOutput(task, 3, coutputpt3);
+ mgr->SetDebugLevel(2);
+
+ if (!mgr->InitAnalysis()) return;
+ mgr->PrintStatus();
+ mgr->StartAnalysis("local", chain, 100);
+}
--- /dev/null
+void runAODProof(const char * proofMode = "full")
+{
+
+ gEnv->SetValue("XSec.GSI.DelegProxy", "2");
+ // TProof::Open("rbertens@alice-caf.cern.ch");
+
+ gSystem->Load("libTree.so");
+ gSystem->Load("libGeom.so");
+ gSystem->Load("libVMC.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libSTEERBase.so");
+ gSystem->Load("libESD.so");
+ gSystem->Load("libAOD.so");
+ gSystem->Load("libANALYSIS.so");
+ gSystem->Load("libOADB.so");
+ gSystem->Load("libANALYSISalice.so");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+
+
+ AliAnalysisAlien * handler = new AliAnalysisAlien("test");
+ handler->SetOverwriteMode();
+ handler->SetRunMode(proofMode);
+ handler->SetProofReset(0);
+ handler->SetAliROOTVersion("v4-21-29-AN");
+ handler->SetProofCluster(Form("%s@alice-caf.cern.ch", gSystem->Getenv("CAFUSER")));
+ //handler->SetProofCluster("rbertens@skaf.saske.sk");
+ // handler->SetProofDataSet("/alice/data/LHC10h_000138653_p2_AOD049#aodTree");
+ handler->SetProofDataSet("/alice/sim/LHC11a10a_000138653_AOD048#aodTree");
+ handler->SetNproofWorkersPerSlave(1);
+ handler->SetAliRootMode("default");
+ handler->SetAdditionalLibs("AliSpectraAODHistoManager.cxx AliSpectraAODHistoManager.h AliSpectraAODEventCuts.cxx AliSpectraAODEventCuts.h AliSpectraAODTrackCuts.cxx AliSpectraAODTrackCuts.h AliAnalysisTaskSpectraAOD.cxx AliAnalysisTaskSpectraAOD.h");
+ handler->SetAnalysisSource("AliSpectraAODHistoManager.cxx+ AliSpectraAODEventCuts.cxx+ AliSpectraAODTrackCuts.cxx+ AliAnalysisTaskSpectraAOD.cxx+");
+ handler->SetFileForTestMode("filelist.txt"); // list of local files for testing
+
+ // handler->SetAliRootMode("");
+ handler->SetClearPackages();
+
+
+ AliAnalysisManager *mgr = new AliAnalysisManager("testAnalysis");
+ mgr->SetGridHandler(handler);
+ AliAODInputHandler* aodH = new AliAODInputHandler();
+ mgr->SetInputEventHandler(aodH);
+
+ gROOT->LoadMacro("AliSpectraAODTrackCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODEventCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODHistoManager.cxx+g");
+ gROOT->LoadMacro("AliAnalysisTaskSpectraAOD.cxx+g");
+
+ // Add PID task
+ gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
+ AliAnalysisTask * taskPID = AddTaskPIDResponse();
+ mgr->AddTask(taskPID);
+
+ AliAnalysisTaskSpectraAOD *task = new AliAnalysisTaskSpectraAOD("TaskAODExercise");
+ task->SetIsMC(1);
+ mgr->AddTask(task);
+
+ // Set the cuts
+ AliSpectraAODEventCuts * vcuts = new AliSpectraAODEventCuts("Event Cuts");
+ AliSpectraAODTrackCuts * tcuts = new AliSpectraAODTrackCuts("Track Cuts");
+ tcuts->SetTrackType(6);
+ tcuts->SetEta(1.);
+ tcuts->SetP(.7);
+ // vcuts->SetCentralityCutMin(0.0) // default
+ // vcuts->SetCentralityCutMax(0.0) // default
+ task->SetEventCuts(vcuts);
+ task->SetTrackCuts(tcuts);
+ task->SetNSigmaForIdentification(3.);
+
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ AliAnalysisDataContainer *coutputpt1 = mgr->CreateContainer("chistpt", AliSpectraAODHistoManager::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1.root");
+ AliAnalysisDataContainer *coutputpt2 = mgr->CreateContainer("cvcutpt", AliSpectraAODEventCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1.root");
+ AliAnalysisDataContainer *coutputpt3 = mgr->CreateContainer("ctcutpt", AliSpectraAODTrackCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1.root");
+
+ mgr->ConnectInput(task, 0, cinput);
+ mgr->ConnectOutput(task, 1, coutputpt1);
+ mgr->ConnectOutput(task, 2, coutputpt2);
+ mgr->ConnectOutput(task, 3, coutputpt3);
+ mgr->SetDebugLevel(2);
+
+ if (!mgr->InitAnalysis()) return;
+ mgr->PrintStatus();
+ mgr->StartAnalysis("proof");
+}
--- /dev/null
+void runAODProof(Int_t c, const char * proofMode = "full")
+{
+
+ gEnv->SetValue("XSec.GSI.DelegProxy", "2");
+
+ gSystem->Load("libTree.so");
+ gSystem->Load("libGeom.so");
+ gSystem->Load("libVMC.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libSTEERBase.so");
+ gSystem->Load("libESD.so");
+ gSystem->Load("libAOD.so");
+ gSystem->Load("libANALYSIS.so");
+ gSystem->Load("libOADB.so");
+ gSystem->Load("libANALYSISalice.so");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+
+
+ AliAnalysisAlien * handler = new AliAnalysisAlien("test");
+ handler->SetOverwriteMode();
+ handler->SetRunMode(proofMode);
+ handler->SetProofReset(0);
+ handler->SetAliROOTVersion("v4-21-29-AN");
+
+ handler->SetProofCluster(Form("%s@alice-caf.cern.ch", gSystem->Getenv("CAFUSER")));
+// handler->SetProofCluster(Form("%s@skaf.saske.sk",gSystem->Getenv("CAFUSER")));
+
+ // Set handler for Real DATA:
+ if (c == 1)
+ {
+ handler->SetProofDataSet("/alice/data/LHC10h_000138653_p2_AOD049#aodTree");
+ // handler->SetProofDataSet("/alice/sim/LHC11a10a_000139107_AOD048#aodTree|alice/sim/LHC11a10a_000138653_AOD048#aodTree");
+ }
+ if (c == 2)
+ {
+ handler->SetProofDataSet("/alice/sim/LHC11a10a_000139107_AOD048#aodTree|/alice/sim/LHC11a10a_000138653_AOD048#aodTree|/alice/sim/LHC11a10a_000139110_AOD048#aodTree|/alice/sim/LHC11a10a_000138662_AOD048#aodTree|/alice/sim/LHC11a10a_000138666_AOD048#aodTree|/alice/sim/LHC11a10a_000138795_AOD048#aodTree");
+ }
+ handler->SetNproofWorkersPerSlave(1);
+ handler->SetAliRootMode("default");
+ handler->SetAdditionalLibs("AliSpectraAODHistoManager.cxx AliSpectraAODHistoManager.h AliSpectraAODEventCuts.cxx AliSpectraAODEventCuts.h AliSpectraAODTrackCuts.cxx AliSpectraAODTrackCuts.h AliAnalysisTaskSpectraAOD.cxx AliAnalysisTaskSpectraAOD.h");
+ handler->SetAnalysisSource("AliSpectraAODHistoManager.cxx+ AliSpectraAODEventCuts.cxx+ AliSpectraAODTrackCuts.cxx+ AliAnalysisTaskSpectraAOD.cxx+");
+ handler->SetFileForTestMode("filelist.txt"); // list of local files for testing
+ // handler->SetAliRootMode("");
+ handler->SetClearPackages();
+
+
+ AliAnalysisManager *mgr = new AliAnalysisManager("testAnalysis");
+ mgr->SetGridHandler(handler);
+ AliAODInputHandler* aodH = new AliAODInputHandler();
+ mgr->SetInputEventHandler(aodH);
+
+ gROOT->LoadMacro("AliSpectraAODTrackCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODEventCuts.cxx+g");
+ gROOT->LoadMacro("AliSpectraAODHistoManager.cxx+g");
+ gROOT->LoadMacro("AliAnalysisTaskSpectraAOD.cxx+g");
+
+
+ // Add PID task
+ gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
+ Bool_t isMC = kFALSE;
+ if (c == 2 || c == 3) isMC = kTRUE;
+ AliAnalysisTask * taskPID = AddTaskPIDResponse(isMC);
+ mgr->AddTask(taskPID);
+
+
+ AliAnalysisTaskSpectraAOD *task = new AliAnalysisTaskSpectraAOD("TaskAODExercise");
+ mgr->AddTask(task);
+
+ // Set the cuts
+ AliSpectraAODEventCuts * vcuts = new AliSpectraAODEventCuts("Event Cuts");
+ AliSpectraAODTrackCuts * tcuts = new AliSpectraAODTrackCuts("Track Cuts");
+ tcuts->SetTrackType(6);
+ tcuts->SetEta(.8);
+ // tcuts->SetDCA(.1);
+ tcuts->SetPt(1.2);
+ vcuts->SetCentralityCutMax(5.); // example min max cuts
+ vcuts->SetCentralityCutMin(0.);
+ task->SetEventCuts(vcuts);
+ task->SetTrackCuts(tcuts);
+ task->SetNSigmaForIdentification(3.); // FIXME
+ task->SetYCut(.5);
+ // check for MC or real data
+ if (c == 2 || c == 3)
+ {
+ task->SetIsMC(kTRUE);
+
+
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ AliAnalysisDataContainer *coutputpt1 = mgr->CreateContainer("chistpt", AliSpectraAODHistoManager::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1._MC.root");
+ AliAnalysisDataContainer *coutputpt2 = mgr->CreateContainer("cvcutpt", AliSpectraAODEventCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1._MC.root");
+ AliAnalysisDataContainer *coutputpt3 = mgr->CreateContainer("ctcutpt", AliSpectraAODTrackCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1._MC.root");
+ }
+ if (c == 1)
+ {
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ AliAnalysisDataContainer *coutputpt1 = mgr->CreateContainer("chistpt", AliSpectraAODHistoManager::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1._data_ptcut.root");
+ AliAnalysisDataContainer *coutputpt2 = mgr->CreateContainer("cvcutpt", AliSpectraAODEventCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1._data_ptcut.root");
+ AliAnalysisDataContainer *coutputpt3 = mgr->CreateContainer("ctcutpt", AliSpectraAODTrackCuts::Class(), AliAnalysisManager::kOutputContainer, "Pt.AOD.1._data_ptcut.root");
+
+ }
+ mgr->ConnectInput(task, 0, cinput);
+ mgr->ConnectOutput(task, 1, coutputpt1);
+ mgr->ConnectOutput(task, 2, coutputpt2);
+ mgr->ConnectOutput(task, 3, coutputpt3);
+ mgr->SetDebugLevel(2);
+
+ if (!mgr->InitAnalysis()) return;
+ mgr->PrintStatus();
+ if (c == 3)
+ {
+ mgr->StartAnalysis("local");
+ }
+ if (c != 3)
+ {
+ mgr->StartAnalysis("proof");
+ }
+}
--- /dev/null
+// TODO:
+// 1. run with Many centrality bins at once
+#include <string.h>
+
+enum { kMyRunModeLocal = 0, kMyRunModeCAF, kMyRunModeGRID};
+
+TList * listToLoad = new TList();
+
+TChain * GetAnalysisChain(const char * incollection);
+
+void runProof(Char_t* data = "/alice/data/LHC10h_000138653_p2_AOD049", Long64_t nev = -1, Long64_t offset = 0, Bool_t debug = kFALSE, Int_t runMode = 1, Bool_t isMC = 0,
+ const char* option = "SAVE", Int_t workers = -1)
+{
+ // runMode:
+ //
+ // 0 local
+ // 1 proof
+ // 2 grid
+
+ if (nev < 0)
+ nev = 1234567890;
+
+ InitAndLoadLibs(runMode, workers, debug);
+
+ // Create the analysis manager
+ mgr = new AliAnalysisManager;
+
+ // Add AOD handler
+ AliAODInputHandler* AODH = new AliAODInputHandler;
+ mgr->SetInputEventHandler(AODH);
+ mgr->SetDebugLevel(10);
+
+ // Is this needed for AOD?
+ // if(isMC) {
+ // AliMCEventHandler* handler = new AliMCEventHandler;
+ // handler->SetPreReadMode(AliMCEventHandler::kLmPreRead);
+ // mgr->SetMCtruthEventHandler(handler);
+ // }
+
+
+ // If we are running on grid, we need the alien handler
+ if (runMode == kMyRunModeGRID)
+ {
+ // Create and configure the alien handler plugin
+ gROOT->LoadMacro("CreateAlienHandler.C");
+ AliAnalysisGrid *alienHandler = CreateAlienHandler(data, listToLoad, "full", isMC);
+ if (!alienHandler)
+ {
+ cout << "Cannot create alien handler" << endl;
+ exit(1);
+ }
+ mgr->SetGridHandler(alienHandler);
+ }
+
+ // Parse option strings
+ TString optionStr(option);
+
+ // remove SAVE option if set
+ // This is copied from a macro by Jan. The reason I kept it is that I may want to pass textual options to the new task at some point
+ Bool_t doSave = kFALSE;
+ TString optionStr(option);
+ if (optionStr.Contains("SAVE"))
+ {
+ optionStr = optionStr(0, optionStr.Index("SAVE")) + optionStr(optionStr.Index("SAVE") + 4, optionStr.Length());
+ doSave = kTRUE;
+ }
+
+ // AliLog::SetClassDebugLevel("AliESDtrackCuts", AliLog::kDebug);// FIXME
+
+ // load my task
+ gROOT->ProcessLine(".L AddTaskSpectraAOD.C");
+ AliAnalysisTaskSpectraAOD * task = AddTaskSpectraAOD("SpectraAOD.root");
+ task->SetIsMC(isMC);
+ // Init and run the analy
+ if (!mgr->InitAnalysis()) return;
+
+ mgr->PrintStatus();
+
+ if (runMode == kMyRunModeLocal)
+ {
+ // If running in local mode, create chain of ESD files
+ cout << "RUNNING LOCAL, CHAIN" << endl;
+ TChain * chain = GetAnalysisChain(data);
+ // chain->Print();
+ mgr->StartAnalysis("local", chain, nev);
+ }
+ else if (runMode == kMyRunModeCAF)
+ {
+ mgr->StartAnalysis("proof", TString(data) + "#aodTree", nev);
+ }
+ else if (runMode == kMyRunModeGRID)
+ {
+ mgr->StartAnalysis("grid");
+ }
+ else
+ {
+ cout << "ERROR: unknown run mode" << endl;
+ }
+
+ TString pathsuffix = "";
+ if (doSave) MoveOutput(data, pathsuffix.Data());
+
+}
+
+
+void MoveOutput(const char * data, const char * suffix = "")
+{
+
+ TString path("output/");
+ path = path + TString(data).Tokenize("/")->Last()->GetName() + suffix;
+
+ TString fileName = "SpectraAOD.root";
+ gSystem->mkdir(path, kTRUE);
+ gSystem->Rename(fileName, path + "/" + fileName);
+ // for(Int_t ibin = 0; ibin < 20; ibin++){
+ // TString fileBin = fileName;
+ // fileBin.ReplaceAll(".root",Form("_%2.2d.root",ibin));
+ // gSystem->Rename(fileBin, path + "/" + fileBin);
+ // }
+ Printf(">>>>> Moved files to %s", path.Data());
+}
+
+
+
+TChain * GetAnalysisChain(const char * incollection)
+{
+ // Builds a chain of esd files
+ // incollection can be
+ // - a single root file
+ // - an xml collection of files on alien
+ // - a ASCII containing a list of local root files
+ TChain* analysisChain = 0;
+ // chain
+ analysisChain = new TChain("aodTree");
+ if (TString(incollection).Contains(".root"))
+ {
+ analysisChain->Add(incollection);
+ }
+ else if (TString(incollection).Contains("xml"))
+ {
+ TGrid::Connect("alien://");
+ TAlienCollection * coll = TAlienCollection::Open(incollection);
+ while (coll->Next())
+ {
+ analysisChain->Add(TString("alien://") + coll->GetLFN());
+ }
+ }
+ else
+ {
+ ifstream file_collect(incollection);
+ TString line;
+ while (line.ReadLine(file_collect))
+ {
+ analysisChain->Add(line.Data());
+ }
+ }
+ analysisChain->GetListOfFiles()->Print();
+
+ return analysisChain;
+}
+
+
+void InitAndLoadLibs(Int_t runMode = kMyRunModeLocal, Int_t workers = 0, Bool_t debug = 0)
+{
+ // Loads libs and par files + custom task and classes
+
+ // Custom stuff to be loaded
+
+ listToLoad->Add(new TObjString("AliSpectraAODHistoManager.cxx+"));
+ listToLoad->Add(new TObjString("AliSpectraAODEventCuts.cxx+"));
+ listToLoad->Add(new TObjString("AliSpectraAODTrackCuts.cxx+"));
+ listToLoad->Add(new TObjString("AliAnalysisTaskSpectraAOD.cxx+"));
+
+ if (runMode == kMyRunModeCAF)
+ {
+ cout << "Init in CAF mode" << endl;
+
+ gEnv->SetValue("XSec.GSI.DelegProxy", "2");
+ TProof * p = TProof::Open("alice-caf.cern.ch", workers > 0 ? Form("workers=%d", workers) : "1x");
+ //TProof * p = TProof::Open("skaf.saske.sk", workers>0 ? Form("workers=%d",workers) : "");
+ p->Exec("TObject *o = gEnv->GetTable()->FindObject(\"Proof.UseMergers\"); gEnv->GetTable()->Remove(o);", kTRUE);
+
+ gProof->EnablePackage("VO_ALICE@AliRoot::v4-21-29-AN");
+ gSystem->Load("libCore.so");
+ gSystem->Load("libTree.so");
+ gSystem->Load("libGeom.so");
+ gSystem->Load("libVMC.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libSTEERBase");
+ gSystem->Load("libESD");
+ gSystem->Load("libAOD");
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libOADB");
+ gSystem->Load("libANALYSISalice");
+
+ // Enable the needed package
+ // gProof->UploadPackage("$ALICE_ROOT/obj/STEERBase");
+ // gProof->EnablePackage("$ALICE_ROOT/obj/STEERBase");
+ // gProof->UploadPackage("$ALICE_ROOT/obj/ESD");
+ // gProof->EnablePackage("$ALICE_ROOT/obj/ESD");
+ // gProof->UploadPackage("$ALICE_ROOT/obj/AOD");
+ // gProof->EnablePackage("$ALICE_ROOT/obj/AOD");
+ // gProof->UploadPackage("$ALICE_ROOT/obj/ANALYSIS");
+ // gProof->EnablePackage("$ALICE_ROOT/obj/ANALYSIS");
+ // gProof->UploadPackage("$ALICE_ROOT/obj/OADB");
+ // gProof->EnablePackage("$ALICE_ROOT/obj/OADB");
+ // gProof->UploadPackage("$ALICE_ROOT/obj/ANALYSISalice");
+ // gProof->EnablePackage("$ALICE_ROOT/obj/ANALYSISalice");
+ // gProof->UploadPackage("$ALICE_ROOT/obj/PWG0base");
+ // gProof->EnablePackage("$ALICE_ROOT/obj/PWG0base");
+ // gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWG0/multPb"));
+ // gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWG1/background"));
+ }
+ else
+ {
+ cout << "Init in Local or Grid mode" << endl;
+ gSystem->Load("libCore.so");
+ gSystem->Load("libTree.so");
+ gSystem->Load("libGeom.so");
+ gSystem->Load("libVMC.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libSTEERBase");
+ gSystem->Load("libESD");
+ gSystem->Load("libAOD");
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libOADB");
+ gSystem->Load("libANALYSISalice");
+ // Use AliRoot includes to compile our task
+ gROOT->ProcessLine(".include $ALICE_ROOT/include");
+
+ // gSystem->Load("libVMC");
+ // gSystem->Load("libTree");
+ // gSystem->Load("libSTEERBase");
+ // gSystem->Load("libESD");
+ // gSystem->Load("libAOD");
+ // gSystem->Load("libANALYSIS");
+ // gSystem->Load("libANALYSISalice");
+ // gSystem->Load("libPWG0base");
+
+ // gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWG0/multPb"));
+ // gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWG1/background"));
+ // gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWG1/background/"));
+ }
+ // Load helper classes
+ TIterator * iter = listToLoad->MakeIterator();
+ TObjString * name = 0;
+ while (name = (TObjString *)iter->Next())
+ {
+ gSystem->ExpandPathName(name->String());
+ cout << name->String().Data();
+ if (runMode == kMyRunModeCAF)
+ {
+ gProof->Load(name->String() + (debug ? "+g" : ""));
+ }
+ else
+ {
+ gROOT->LoadMacro(name->String() + (debug ? "+g" : ""));
+ }
+ }
+
+}