]>
Commit | Line | Data |
---|---|---|
1 | /************************************************************************** | |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: Satyajit Jena. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | ||
17 | //=========================================================================// | |
18 | // AliEbyE Fluctuaion Analysis for PID // | |
19 | // Deepika Jena | drathee@cern.ch // | |
20 | // Satyajit Jena | sjena@cern.ch // | |
21 | //=========================================================================// | |
22 | ||
23 | #include "TChain.h" | |
24 | #include "TList.h" | |
25 | #include "TFile.h" | |
26 | #include "TTree.h" | |
27 | #include "TH1D.h" | |
28 | #include "TH2F.h" | |
29 | #include "TH3F.h" | |
30 | #include "TProfile.h" | |
31 | ||
32 | #include "AliAnalysisTask.h" | |
33 | #include "AliAnalysisManager.h" | |
34 | ||
35 | #include "AliVEvent.h" | |
36 | #include "AliESDEvent.h" | |
37 | #include "AliMCEvent.h" | |
38 | #include "AliAODEvent.h" | |
39 | ||
40 | #include "AliStack.h" | |
41 | #include "AliGenEventHeader.h" | |
42 | #include "AliGenPythiaEventHeader.h" | |
43 | #include "AliGenHijingEventHeader.h" | |
44 | #include "AliGenDPMjetEventHeader.h" | |
45 | #include "TDatabasePDG.h" | |
46 | ||
47 | #include "AliEbyEPidRatioTaskOnFly.h" | |
48 | ||
49 | ClassImp(AliEbyEPidRatioTaskOnFly) | |
50 | ||
51 | //----------------------------------------------------------------------- | |
52 | AliEbyEPidRatioTaskOnFly::AliEbyEPidRatioTaskOnFly( const char *name ) | |
53 | : AliAnalysisTaskSE( name ), | |
54 | fThnList(0), | |
55 | fPtLowerLimit(0.2), | |
56 | fPtHigherLimit(5.), | |
57 | fEtaLowerLimit(-1.), | |
58 | fEtaHigherLimit(1.), | |
59 | fCentrality(-1), | |
60 | fOrder(8), | |
61 | fRedFactp(NULL) | |
62 | ||
63 | ||
64 | { | |
65 | ||
66 | DefineOutput(1, TList::Class()); | |
67 | } | |
68 | ||
69 | AliEbyEPidRatioTaskOnFly::~AliEbyEPidRatioTaskOnFly() | |
70 | { | |
71 | if (fThnList) delete fThnList; | |
72 | ||
73 | for (Int_t ii = 0; ii <= fOrder; ++ii) | |
74 | if (fRedFactp[ii]) delete[] fRedFactp[ii]; | |
75 | if (fRedFactp) delete[] fRedFactp; | |
76 | ||
77 | } | |
78 | ||
79 | ||
80 | ||
81 | ||
82 | //--------------------------------------------------------------------------------- | |
83 | void AliEbyEPidRatioTaskOnFly::UserCreateOutputObjects() { | |
84 | fThnList = new TList(); | |
85 | fThnList->SetOwner(kTRUE); | |
86 | ||
87 | const Char_t *name = "Mc"; | |
88 | const Char_t *title = Form(" #eta [%2.1f-%2.1f]",fEtaLowerLimit,fEtaHigherLimit); | |
89 | ||
90 | ||
91 | TString sName(name); | |
92 | TString sTitle(title); | |
93 | ||
94 | ||
95 | fRedFactp = new Double_t*[fOrder+1]; | |
96 | for (Int_t ii = 0 ; ii <= fOrder; ++ii) | |
97 | fRedFactp[ii] = new Double_t[2]; | |
98 | ||
99 | //TList *list[4]; | |
100 | fThnList->Add(new TList); | |
101 | TList *list = static_cast<TList*>(fThnList->Last()); | |
102 | list->SetName(Form("f%s",name)); | |
103 | list->SetOwner(kTRUE); | |
104 | ||
105 | for (Int_t iPid = 0; iPid < 4; ++iPid) { | |
106 | TString sNetTitle(Form("%s - %s",fgkPidLatex[iPid][1],fgkPidLatex[iPid][0])); | |
107 | ||
108 | list->Add(new TProfile(Form("fProfTot%sPlus%s",fgkPidName[iPid],name), | |
109 | Form("(%s) : %s;Centrality(100);(%s)",fgkPidName[iPid], sTitle.Data(), sNetTitle.Data()), | |
110 | 20,-0.5,19.5)); | |
111 | ||
112 | list->Add(new TProfile(Form("fProfTot%sMinus%s",fgkPidName[iPid],name), | |
113 | Form("(%s) : %s;Centrality(100);(%s)",fgkPidName[iPid], sTitle.Data(), sNetTitle.Data()), | |
114 | 20,-0.5,19.5)); | |
115 | ||
116 | for (Int_t idx = 1; idx <= fOrder; ++idx) { | |
117 | list->Add(new TProfile(Form("fProf%s%sNet%dM",fgkPidName[iPid],name, idx), | |
118 | Form("(%s)^{%d} : %s;Centrality(100);(%s)^{%d}",sNetTitle.Data(), idx, sTitle.Data(), sNetTitle.Data(), idx), | |
119 | 20,-0.5,19.5)); | |
120 | } | |
121 | ||
122 | for (Int_t ii = 0; ii <= fOrder; ++ii) { | |
123 | for (Int_t kk = 0; kk <= fOrder; ++kk) { | |
124 | list->Add(new TProfile(Form("fProf%s%sNetF%02d%02d",fgkPidName[iPid], name, ii, kk), | |
125 | Form("f_{%02d%02d} : %s;Centrality(100);f_{%02d%02d}", ii, kk, sTitle.Data(), ii, kk), | |
126 | 20,-0.5,19.5)); | |
127 | } | |
128 | } | |
129 | ||
130 | } | |
131 | ||
132 | ||
133 | for (Int_t iPhy = 0; iPhy < 46; ++iPhy) { | |
134 | list->Add(new TProfile(Form("fProf%sNu%02d",name,iPhy),Form("Physics Variable for index %d | %s ; Centrality;",iPhy,name),20,-0.5,19.5)); | |
135 | } | |
136 | ||
137 | ||
138 | PostData(1, fThnList); | |
139 | } | |
140 | ||
141 | ||
142 | //---------------------------------------------------------------------------------- | |
143 | void AliEbyEPidRatioTaskOnFly::UserExec( Option_t * ){ | |
144 | AliMCEvent* mcEvent = MCEvent(); | |
145 | if (!mcEvent) { | |
146 | Printf("ERROR: Could not retrieve MC event"); | |
147 | return; | |
148 | } | |
149 | ||
150 | for (Int_t ii = 0 ; ii < 4; ++ii) | |
151 | for (Int_t kk = 0 ; kk < 2; ++kk) | |
152 | fNp[ii][kk] = 0; | |
153 | ||
154 | ||
155 | const AliVVertex *vtxMC = mcEvent->GetPrimaryVertex(); | |
156 | if (vtxMC->GetZ() > 10.) return; | |
157 | ||
158 | AliStack *stack = mcEvent->Stack(); | |
159 | ||
160 | AliGenEventHeader* genHeader = mcEvent->GenEventHeader(); | |
161 | if(!genHeader){ | |
162 | printf(" Event generator header not available!!!\n"); | |
163 | return; | |
164 | } | |
165 | ||
166 | Double_t imp = 21; | |
167 | ||
168 | if(genHeader->InheritsFrom(AliGenHijingEventHeader::Class())){ | |
169 | imp = ((AliGenHijingEventHeader*) genHeader)->ImpactParameter(); | |
170 | } | |
171 | if (imp < 20) fCentrality = Int_t(imp); | |
172 | else return; | |
173 | ||
174 | Int_t nParticles = stack->GetNtrack(); | |
175 | for (Int_t iParticle = 0; iParticle < nParticles; iParticle++) { | |
176 | TParticle* part = stack->Particle(iParticle); | |
177 | if (!part) { | |
178 | Printf(" No Particle Available "); | |
179 | continue; | |
180 | } | |
181 | ||
182 | Float_t pt = part->Pt(); | |
183 | if(pt < fPtLowerLimit || pt > fPtHigherLimit ) continue; | |
184 | Float_t gEta = part->Eta(); | |
185 | if(gEta < fEtaLowerLimit || gEta > fEtaHigherLimit ) continue; | |
186 | ||
187 | // Float_t gCharge = (pdgPart ? pdgPart->Charge() : 0); | |
188 | ||
189 | Int_t pid = part->GetPdgCode(); | |
190 | if(pid == 211) { fNp[1][1]++; fNp[0][1]++; } | |
191 | else if(pid == -211) { fNp[1][0]++; fNp[0][0]++; } | |
192 | else if(pid == 321) { fNp[2][1]++; fNp[0][1]++; } | |
193 | else if(pid == -321) { fNp[2][0]++; fNp[0][0]++; } | |
194 | else if(pid == 2212) { fNp[3][1]++; fNp[0][1]++; } | |
195 | else if(pid == -2212) { fNp[3][0]++; fNp[0][0]++; } | |
196 | ||
197 | } | |
198 | ||
199 | /*Printf("%6d %6d %6d %6d %6d %6d %6d %6d %6d", fCentrality, | |
200 | fNp[0][0], fNp[0][1], | |
201 | fNp[1][0], fNp[1][1], | |
202 | fNp[2][0], fNp[2][1], | |
203 | fNp[3][0], fNp[3][1]);*/ | |
204 | ||
205 | ||
206 | FillHistSetCent(); | |
207 | PostData(1, fThnList); | |
208 | } | |
209 | ||
210 | void AliEbyEPidRatioTaskOnFly::Terminate( Option_t * ){ | |
211 | Info("AliEbyEPidRatioTaskOnFly"," Task Successfully finished"); | |
212 | } | |
213 | //________________________________________________________________________ | |
214 | const Char_t* AliEbyEPidRatioTaskOnFly::fgkPidName[4] = {"Nch","Npi","Nka","Npr"}; | |
215 | //________________________________________________________________________ | |
216 | const Char_t* AliEbyEPidRatioTaskOnFly::fgkPidLatex[4][2] = {{"N_{-}","N_{+}"}, {"N_{#pi^{-}}","N_{#pi^{+}}"},{"N_{K^{-}}","N_{K^{+}}"}, {"N_{#bar{p}}","N_{p}"}}; | |
217 | //________________________________________________________________________ | |
218 | const Char_t* AliEbyEPidRatioTaskOnFly::fgkPidTitles[4][2] = {{"Negative","Positive"},{"Anti-Pions","Pions"},{"Anti-Kaons","Kaons"}, {"Anti-Protons","Protons"}}; | |
219 | //________________________________________________________________________ | |
220 | void AliEbyEPidRatioTaskOnFly::FillHistSetCent() { | |
221 | ||
222 | const Char_t *name = "Mc"; | |
223 | ||
224 | TList *list = static_cast<TList*>(fThnList->FindObject(Form("f%s",name))); | |
225 | ||
226 | for (Int_t iPid = 0; iPid < 4; ++iPid) { | |
227 | Int_t deltaNp = fNp[iPid][1]-fNp[iPid][0]; | |
228 | Double_t delta = 1.; | |
229 | (static_cast<TProfile*>(list->FindObject(Form("fProfTot%sPlus%s", fgkPidName[iPid], name))))->Fill(fCentrality, fNp[iPid][1]); | |
230 | (static_cast<TProfile*>(list->FindObject(Form("fProfTot%sMinus%s", fgkPidName[iPid], name))))->Fill(fCentrality, fNp[iPid][0]); | |
231 | for (Int_t idxOrder = 1; idxOrder <= fOrder; ++idxOrder) { | |
232 | delta *= deltaNp; | |
233 | (static_cast<TProfile*>(list->FindObject(Form("fProf%s%sNet%dM", fgkPidName[iPid], name, idxOrder))))->Fill(fCentrality, delta); | |
234 | } | |
235 | ||
236 | for (Int_t idxOrder = 0; idxOrder <= fOrder; ++ idxOrder) { | |
237 | fRedFactp[idxOrder][0] = 1.; | |
238 | fRedFactp[idxOrder][1] = 1.; | |
239 | } | |
240 | ||
241 | for (Int_t idxOrder = 1; idxOrder <= fOrder; ++ idxOrder) { | |
242 | fRedFactp[idxOrder][0] = fRedFactp[idxOrder-1][0] * Double_t(fNp[iPid][0]-(idxOrder-1)); | |
243 | fRedFactp[idxOrder][1] = fRedFactp[idxOrder-1][1] * Double_t(fNp[iPid][1]-(idxOrder-1)); | |
244 | } | |
245 | ||
246 | for (Int_t ii = 0; ii <= fOrder; ++ii) { | |
247 | for (Int_t kk = 0; kk <= fOrder; ++kk) { | |
248 | Double_t fik = fRedFactp[ii][1] * fRedFactp[kk][0]; | |
249 | (static_cast<TProfile*>(list->FindObject(Form("fProf%s%sNetF%02d%02d", fgkPidName[iPid], name, ii, kk))))->Fill(fCentrality, fik); | |
250 | } | |
251 | } | |
252 | } | |
253 | ||
254 | //Printf("%6d %20s %6.2f %6d %6d %6d %6d %6d %6d %6d %6d", idx, name, centralityBin, | |
255 | // fNp[0][1], fNp[0][0], | |
256 | // fNp[1][1], fNp[1][0], | |
257 | /// fNp[2][1], fNp[2][0], | |
258 | // fNp[3][1], fNp[3][0]); | |
259 | // | |
260 | ||
261 | Int_t a[6][4]; Int_t b[22]; | |
262 | for (Int_t iPid = 0; iPid < 4; ++iPid) { | |
263 | a[0][iPid] = fNp[iPid][1]+fNp[iPid][0]; // 0 n+ + n- | |
264 | a[1][iPid] = fNp[iPid][1]; // 1 n+ | |
265 | a[2][iPid] = fNp[iPid][0]; // 2 n- | |
266 | a[3][iPid] = fNp[iPid][1]*fNp[iPid][0]; // 3 n+ . n- | |
267 | a[4][iPid] = fNp[iPid][1]*(fNp[iPid][1]-1); // 4 n+ (n+ - 1) | |
268 | a[5][iPid] = fNp[iPid][0]*(fNp[iPid][0]-1); // 5 n- (n- - 1) | |
269 | ||
270 | // Printf("%6d %20s %6.2f %6d %6d %6d ", idx, name, centralityBin, | |
271 | // a[0][iPid], a[1][iPid], a[2][iPid]); | |
272 | ||
273 | } | |
274 | ||
275 | b[0] = a[0][0]*a[0][2]; // 24 N K | |
276 | b[1] = a[0][1]*a[0][2]; // 25 Pi K | |
277 | b[2] = a[1][1]*a[1][2]; // 26 pi+ k+ | |
278 | b[3] = a[1][1]*a[2][2]; // 27 pi+ k- | |
279 | b[4] = a[2][1]*a[1][2]; // 28 pi- k+ | |
280 | b[5] = a[2][1]*a[2][2]; // 29 pi- k- | |
281 | ||
282 | b[6] = a[0][0]*a[0][3]; // 30 N P | |
283 | b[7] = a[0][2]*a[0][3]; // 31 K P | |
284 | b[8] = a[1][2]*a[1][3]; // 32 k+ p+ | |
285 | b[9] = a[1][2]*a[2][3]; // 33 k+ p- | |
286 | b[10] = a[2][2]*a[1][3]; // 34 k- p+ | |
287 | b[11] = a[2][2]*a[2][3]; // 35 k- p- | |
288 | ||
289 | b[12] = a[0][0]*a[0][1]; // 36 N Pi | |
290 | b[13] = a[0][3]*a[0][1]; // 37 P Pi | |
291 | b[14] = a[1][3]*a[1][1]; // 38 p+ pi+ | |
292 | b[15] = a[1][3]*a[2][1]; // 39 p+ pi- | |
293 | b[16] = a[2][3]*a[1][1]; // 40 p- pi+ | |
294 | b[17] = a[2][3]*a[2][1]; // 41 p- pi- | |
295 | ||
296 | b[18] = a[0][0]*(a[0][0] - 1); // 42 N ( N - 1 ) | |
297 | b[19] = a[0][1]*(a[0][1] - 1); // 43 Pi( Pi- 1 ) | |
298 | b[20] = a[0][2]*(a[0][1] - 1); // 44 K ( K - 1 ) | |
299 | b[21] = a[0][3]*(a[0][3] - 1); // 45 P ( P - 1 ) | |
300 | // TList *list_nu = static_cast<TList*>(fOutList->FindObject(Form("f%s_nu",name))); | |
301 | Int_t k = 0; | |
302 | for (Int_t j = 0; j < 4; j++) { | |
303 | for (Int_t i = 0; i < 6; i++) { | |
304 | (static_cast<TProfile*>(list->FindObject(Form("fProf%sNu%02d", name,k))))->Fill(fCentrality,a[i][j]); | |
305 | k++; | |
306 | } | |
307 | } | |
308 | ||
309 | for (Int_t j = 0; j < 22; j++) { | |
310 | (static_cast<TProfile*>(list->FindObject(Form("fProf%sNu%02d", name,j+24))))->Fill(fCentrality,b[j]); | |
311 | } | |
312 | ||
313 | /* | |
314 | Printf("%6d %6d %6d %6d %6d %6d %6d %6d %6d", fCentrality, | |
315 | fNp[0][0], fNp[0][1], | |
316 | fNp[1][0], fNp[1][1], | |
317 | fNp[2][0], fNp[2][1], | |
318 | fNp[3][0], fNp[3][1]); | |
319 | */ | |
320 | ||
321 | return; | |
322 | } | |
323 |