1 /**************************************************************************
\r
2 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
\r
4 * Author: The ALICE Off-line Project. *
\r
5 * Contributors are mentioned in the code where appropriate. *
\r
7 * Permission to use, copy, modify and distribute this software and its *
\r
8 * documentation strictly for non-commercial purposes is hereby granted *
\r
9 * without fee, provided that the above copyright notice appears in all *
\r
10 * copies and that both the copyright notice and this permission notice *
\r
11 * appear in the supporting documentation. The authors make no claims *
\r
12 * about the suitability of this software for any purpose. It is *
\r
13 * provided "as is" without express or implied warranty. *
\r
14 **************************************************************************/
\r
16 //-----------------------------------------------------------------
\r
17 // AliAnalysisTaskSpectraAllChAOD class
\r
18 //-----------------------------------------------------------------
\r
22 #include "TLegend.h"
\r
25 #include "THnSparse.h"
\r
26 #include "TCanvas.h"
\r
27 #include "AliAnalysisTask.h"
\r
28 #include "AliAODTrack.h"
\r
29 #include "AliAODMCParticle.h"
\r
30 #include "AliVParticle.h"
\r
31 #include "AliAODEvent.h"
\r
32 #include "AliAODInputHandler.h"
\r
33 #include "AliAnalysisTaskSpectraAllChAOD.h"
\r
34 #include "AliAnalysisTaskESDfilter.h"
\r
35 #include "AliAnalysisDataContainer.h"
\r
36 #include "AliSpectraAODTrackCuts.h"
\r
37 #include "AliSpectraAODEventCuts.h"
\r
38 #include "AliHelperPID.h"
\r
39 #include "AliPIDCombined.h"
\r
40 #include "AliCentrality.h"
\r
42 #include "AliVEvent.h"
\r
43 #include "AliStack.h"
\r
44 #include <TMCProcess.h>
\r
48 using namespace AliHelperPIDNameSpace;
\r
49 using namespace std;
\r
51 ClassImp(AliAnalysisTaskSpectraAllChAOD)
\r
53 //________________________________________________________________________
\r
54 AliAnalysisTaskSpectraAllChAOD::AliAnalysisTaskSpectraAllChAOD(const char *name) : AliAnalysisTaskSE(name),
\r
64 // Default constructor
\r
65 DefineInput(0, TChain::Class());
\r
66 DefineOutput(1, TList::Class());
\r
67 DefineOutput(2, AliSpectraAODEventCuts::Class());
\r
68 DefineOutput(3, AliSpectraAODTrackCuts::Class());
\r
69 DefineOutput(4, AliHelperPID::Class());
\r
72 //________________________________________________________________________
\r
73 void AliAnalysisTaskSpectraAllChAOD::UserCreateOutputObjects()
\r
75 // create output objects
\r
76 fOutput=new TList();
\r
77 fOutput->SetOwner();
\r
78 fOutput->SetName("fOutput");
\r
80 if (!fTrackCuts) AliFatal("Track Cuts should be set in the steering macro");
\r
81 if (!fEventCuts) AliFatal("Event Cuts should be set in the steering macro");
\r
82 if (!fHelperPID) AliFatal("HelperPID should be set in the steering macro");
\r
84 // binning common to all the THn
\r
85 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.};
\r
86 const Int_t nptBins=26;
\r
88 //dimensions of THnSparse for tracks
\r
89 const Int_t nvartrk=8;
\r
90 // pt cent Q vec IDrec IDgen isph iswd y
\r
91 Int_t binsHistRealTrk[nvartrk] = { nptBins, fnCentBins, fnQvecBins, 3, 3, 2, 2, 2};
\r
92 Double_t xminHistRealTrk[nvartrk] = { 0., 0., 0., -0.5, -0.5, -0.5, -0.5, -0.5};
\r
93 Double_t xmaxHistRealTrk[nvartrk] = { 10., 100., 10., 2.5, 2.5, 1.5, 1.5, 0.5};
\r
94 THnSparseF* NSparseHistTrk = new THnSparseF("NSparseHistTrk","NSparseHistTrk",nvartrk,binsHistRealTrk,xminHistRealTrk,xmaxHistRealTrk);
\r
95 NSparseHistTrk->GetAxis(0)->SetTitle("#it{p}_{T,rec}");
\r
96 NSparseHistTrk->GetAxis(0)->SetName("pT_rec");
\r
97 NSparseHistTrk->SetBinEdges(0,ptBins);
\r
98 NSparseHistTrk->GetAxis(1)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
\r
99 NSparseHistTrk->GetAxis(1)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
\r
100 NSparseHistTrk->GetAxis(2)->SetTitle("Q vec");
\r
101 NSparseHistTrk->GetAxis(2)->SetName("Q_vec");
\r
102 NSparseHistTrk->GetAxis(3)->SetTitle("ID rec");
\r
103 NSparseHistTrk->GetAxis(3)->SetName("ID_rec");
\r
104 NSparseHistTrk->GetAxis(4)->SetTitle("ID gen");
\r
105 NSparseHistTrk->GetAxis(4)->SetName("ID_gen");
\r
106 NSparseHistTrk->GetAxis(5)->SetTitle("isph");
\r
107 NSparseHistTrk->GetAxis(5)->SetName("isph");
\r
108 NSparseHistTrk->GetAxis(6)->SetTitle("iswd");
\r
109 NSparseHistTrk->GetAxis(6)->SetName("iswd");
\r
110 NSparseHistTrk->GetAxis(7)->SetTitle("y");
\r
111 NSparseHistTrk->GetAxis(7)->SetName("y");
\r
112 fOutput->Add(NSparseHistTrk);
\r
114 //dimensions of THnSparse for stack
\r
115 const Int_t nvarst=5;
\r
116 // pt cent IDgen isph y
\r
117 Int_t binsHistRealSt[nvarst] = { nptBins, fnCentBins, 3, 2, 2};
\r
118 Double_t xminHistRealSt[nvarst] = { 0., 0., -0.5, -0.5, -0.5};
\r
119 Double_t xmaxHistRealSt[nvarst] = { 10., 100., 2.5, 1.5, 0.5};
\r
120 THnSparseF* NSparseHistSt = new THnSparseF("NSparseHistSt","NSparseHistSt",nvarst,binsHistRealSt,xminHistRealSt,xmaxHistRealSt);
\r
121 NSparseHistSt->GetAxis(0)->SetTitle("#it{p}_{T,gen}");
\r
122 NSparseHistSt->SetBinEdges(0,ptBins);
\r
123 NSparseHistSt->GetAxis(0)->SetName("pT_rec");
\r
124 NSparseHistSt->GetAxis(1)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
\r
125 NSparseHistSt->GetAxis(1)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
\r
126 NSparseHistSt->GetAxis(2)->SetTitle("ID gen");
\r
127 NSparseHistSt->GetAxis(2)->SetName("ID_gen");
\r
128 NSparseHistSt->GetAxis(3)->SetTitle("isph");
\r
129 NSparseHistSt->GetAxis(3)->SetName("isph");
\r
130 NSparseHistSt->GetAxis(4)->SetTitle("y");
\r
131 NSparseHistSt->GetAxis(4)->SetName("y");
\r
132 fOutput->Add(NSparseHistSt);
\r
134 //dimensions of THnSparse for the normalization
\r
135 const Int_t nvarev=2;
\r
137 Int_t binsHistRealEv[nvarev] = { fnCentBins, fnQvecBins};
\r
138 Double_t xminHistRealEv[nvarev] = { 0., 0.};
\r
139 Double_t xmaxHistRealEv[nvarev] = { 100., 10.};
\r
140 THnSparseF* NSparseHistEv = new THnSparseF("NSparseHistEv","NSparseHistEv",nvarev,binsHistRealEv,xminHistRealEv,xmaxHistRealEv);
\r
141 NSparseHistEv->GetAxis(0)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
\r
142 NSparseHistEv->GetAxis(0)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
\r
143 NSparseHistEv->GetAxis(1)->SetTitle("Q vec");
\r
144 NSparseHistEv->GetAxis(1)->SetName("Q_vec");
\r
145 fOutput->Add(NSparseHistEv);
\r
147 PostData(1, fOutput );
\r
148 PostData(2, fEventCuts);
\r
149 PostData(3, fTrackCuts);
\r
150 PostData(4, fHelperPID);
\r
153 //________________________________________________________________________
\r
155 void AliAnalysisTaskSpectraAllChAOD::UserExec(Option_t *)
\r
157 //Printf("An event");
\r
159 fAOD = dynamic_cast<AliAODEvent*>(fInputEvent);
\r
161 AliWarning("ERROR: AliAODEvent not available \n");
\r
165 if (strcmp(fAOD->ClassName(), "AliAODEvent"))
\r
167 AliFatal("Not processing AODs");
\r
170 if(!fEventCuts->IsSelected(fAOD,fTrackCuts))return;//event selection
\r
172 //Default TPC priors
\r
173 if(fHelperPID->GetPIDType()==kBayes)fHelperPID->GetPIDCombined()->SetDefaultTPCPriors();//FIXME maybe this can go in the UserCreateOutputObject?
\r
175 Double_t Qvec=0.;//in case of MC we save space in the memory
\r
176 if(!fIsMC)Qvec=fEventCuts->GetqV0A();//FIXME we'll have to combine A and C
\r
177 Double_t Cent=fEventCuts->GetCent();
\r
179 // First do MC to fill up the MC particle array
\r
180 TClonesArray *arrayMC = 0;
\r
183 arrayMC = (TClonesArray*) fAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName());
\r
185 AliFatal("Error: MC particles branch not found!\n");
\r
187 Int_t nMC = arrayMC->GetEntries();
\r
188 for (Int_t iMC = 0; iMC < nMC; iMC++)
\r
190 AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(iMC);
\r
191 if(!partMC->Charge()) continue;//Skip neutrals
\r
192 if(partMC->Eta() < fTrackCuts->GetEtaMin() || partMC->Eta() > fTrackCuts->GetEtaMax())continue;//ETA CUT ON GENERATED!!!!!!!!!!!!!!!!!!!!!!!!!!
\r
194 //pt cent IDgen isph y
\r
196 varSt[0]=partMC->Pt();
\r
198 varSt[2]=(Double_t)fHelperPID->GetParticleSpecies(partMC);
\r
199 varSt[3]=(Double_t)partMC->IsPhysicalPrimary();
\r
200 varSt[4]=partMC->Y();
\r
201 ((THnSparseF*)fOutput->FindObject("NSparseHistSt"))->Fill(varSt);//stack loop
\r
203 //Printf("a particle");
\r
208 //main loop on tracks
\r
209 for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
\r
210 AliAODTrack* track = fAOD->GetTrack(iTracks);
\r
211 if (!fTrackCuts->IsSelected(track,kTRUE)) continue; //track selection (rapidity selection NOT in the standard cuts)
\r
212 Int_t IDrec=fHelperPID->GetParticleSpecies(track,kTRUE);//id from detector
\r
213 Double_t y= track->Y(fHelperPID->GetMass((AliHelperParticleSpecies_t)IDrec));
\r
214 Int_t IDgen=kSpUndefined;//set if MC
\r
219 AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(TMath::Abs(track->GetLabel()));
\r
221 AliError("Cannot get MC particle");
\r
224 IDgen=fHelperPID->GetParticleSpecies(partMC);
\r
225 isph=partMC->IsPhysicalPrimary();
\r
226 iswd=partMC->IsSecondaryFromWeakDecay();//FIXME not working on old productions
\r
229 //pt cent Q vec IDrec IDgen isph iswd y
\r
230 Double_t varTrk[8];
\r
231 varTrk[0]=track->Pt();
\r
234 varTrk[3]=(Double_t)IDrec;
\r
235 varTrk[4]=(Double_t)IDgen;
\r
236 varTrk[5]=(Double_t)isph;
\r
237 varTrk[6]=(Double_t)iswd;
\r
239 ((THnSparseF*)fOutput->FindObject("NSparseHistTrk"))->Fill(varTrk);//track loop
\r
241 //Printf("a track");
\r
243 } // end loop on tracks
\r
248 ((THnSparseF*)fOutput->FindObject("NSparseHistEv"))->Fill(varEv);//event loop
\r
250 PostData(1, fOutput );
\r
251 PostData(2, fEventCuts);
\r
252 PostData(3, fTrackCuts);
\r
253 PostData(4, fHelperPID);
\r
256 //_________________________________________________________________
\r
257 void AliAnalysisTaskSpectraAllChAOD::Terminate(Option_t *)
\r