]>
Commit | Line | Data |
---|---|---|
62567a80 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
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 | // JetAnalysis class | |
83a444b1 | 18 | // Analyse Jets (already found jets) |
19 | // Authors: andreas.morsch@cern.ch, jgcn@mail.cern.ch | |
20 | // mercedes.lopez.noriega@cern.ch | |
62567a80 | 21 | //--------------------------------------------------------------------- |
22 | ||
83a444b1 | 23 | #include <Riostream.h> |
62567a80 | 24 | #include "AliJetAnalysis.h" |
25 | ClassImp(AliJetAnalysis) | |
83a444b1 | 26 | |
27 | // root | |
28 | #include <Riostream.h> | |
62567a80 | 29 | #include <TH1F.h> |
83a444b1 | 30 | #include <TH2F.h> |
62567a80 | 31 | #include <TProfile.h> |
62567a80 | 32 | #include <TFile.h> |
33 | #include <TTree.h> | |
34 | #include <TStyle.h> | |
35 | #include <TSystem.h> | |
36 | #include <TLorentzVector.h> | |
83a444b1 | 37 | #include <TMath.h> |
38 | // aliroot | |
62567a80 | 39 | #include "AliJetProductionDataPDC2004.h" |
40 | #include "AliJet.h" | |
83a444b1 | 41 | #include "AliJetKineReaderHeader.h" |
62567a80 | 42 | #include "AliJetESDReaderHeader.h" |
43 | #include "AliUA1JetHeader.h" | |
44 | #include "AliLeading.h" | |
83a444b1 | 45 | |
46 | ||
47 | AliJetAnalysis::AliJetAnalysis() | |
48 | { | |
49 | // Constructor | |
50 | fDirectory = 0x0; | |
51 | fBkgdDirectory = 0x0; | |
52 | fFile = "anajets.root"; | |
53 | fEventMin = 0; | |
54 | fEventMax = -1; | |
55 | fRunMin = 0; | |
56 | fRunMax = 11; | |
57 | fminMult = 0; | |
58 | fPercentage = -1.0; | |
59 | fEfactor = 1.0; | |
60 | SetReaderHeader(); | |
61 | ||
62 | // for Analize() | |
63 | fPythia = kFALSE; | |
64 | fDoPart = kTRUE; | |
65 | fDoGenJ = kTRUE; | |
66 | fDoRecJ = kTRUE; | |
67 | fDoKine = kTRUE; | |
68 | fDoCorr = kTRUE; | |
69 | fDoCorr50 = kFALSE; | |
70 | fDoShap = kTRUE; | |
71 | fDoFrag = kTRUE; | |
72 | fDoTrig = kTRUE; | |
73 | fDoJt = kTRUE; | |
74 | fDodNdxi = kTRUE; | |
75 | fDoBkgd = kFALSE; | |
76 | ||
77 | fPart = 0; | |
78 | fGenJ = 0; | |
79 | fRecJ = 0; | |
80 | fRecB = 0; | |
81 | fWeight = 1.0; | |
82 | fWdEdr = 0.0; | |
83 | fPartPtCut = 0.0; | |
84 | fWJt = 0.0; | |
85 | fWdNdxi = 0.0; | |
86 | ||
87 | // for bkgd | |
88 | fp0 = 0.0; | |
89 | fPtJ = 0.0; | |
90 | fEJ = 0.0; | |
91 | fEtaJ = 0.0; | |
92 | fPhiJ = 0.0; | |
93 | fjv3X = 0.0; | |
94 | fjv3Y = 0.0; | |
95 | fjv3Z = 0.0; | |
96 | ||
97 | // kine histos | |
98 | fRKineEneH = 0; | |
99 | fRKinePtH = 0; | |
100 | fRKinePhiH = 0; | |
101 | fRKineEtaH = 0; | |
102 | fGKineEneH = 0; | |
103 | fGKinePtH = 0; | |
104 | fGKinePhiH = 0; | |
105 | fGKineEtaH = 0; | |
106 | fPKineEneH = 0; | |
107 | fPKinePtH = 0; | |
108 | fPKinePhiH = 0; | |
109 | fPKineEtaH = 0; | |
110 | ||
111 | // correlation histograms | |
112 | fPGCorrEneH = 0; | |
113 | fPGCorrPtH = 0; | |
114 | fPGCorrEtaH = 0; | |
115 | fPGCorrPhiH = 0; | |
116 | fPRCorrEneH = 0; | |
117 | fPRCorrPtH = 0; | |
118 | fPRCorrEtaH = 0; | |
119 | fPRCorrPhiH = 0; | |
120 | fRGCorrEneH = 0; | |
121 | fRGCorrPtH = 0; | |
122 | fRGCorrEtaH = 0; | |
123 | fRGCorrPhiH = 0; | |
124 | ||
125 | // correlation histograms when one particle | |
126 | // has more than 50% of the energy of the jet | |
127 | fPRCorr50EneH = 0; | |
128 | fPRCorr50PtH = 0; | |
129 | fPRCorr50EtaH = 0; | |
130 | fPRCorr50PhiH = 0; | |
131 | fRGCorr50EneH = 0; | |
132 | fRGCorr50PtH = 0; | |
133 | fRGCorr50EtaH = 0; | |
134 | fRGCorr50PhiH = 0; | |
135 | ||
136 | // shape histos | |
137 | fWShapR = 0.0; | |
138 | fRShapSelH = 0; | |
139 | fRShapRejH = 0; | |
140 | fRShapAllH = 0; | |
141 | ||
142 | // fragmentation function histos | |
143 | fWFragR = 0.0; | |
144 | fRFragSelH = 0; | |
145 | fRFragRejH = 0; | |
146 | fRFragAllH = 0; | |
147 | ||
148 | // trigger bias histos | |
149 | fGTriggerEneH = 0; | |
150 | fRTriggerEneH = 0; | |
151 | fGPTriggerEneH = 0; | |
152 | fPTriggerEneH = 0; | |
153 | ||
154 | // dE/dr histo and its dr | |
155 | fdEdrH = 0; | |
156 | fdEdrB = 0; | |
157 | fPtEneH2 = 0; | |
158 | fdEdrW = 0; | |
159 | fdrdEdr = 1.0; | |
160 | ||
161 | // jt histo and its dr | |
162 | fJtH = 0; | |
163 | fJtB = 0; | |
164 | fJtW = 0; | |
165 | fdrJt = 1.0; | |
166 | ||
167 | // dN/dxi histo and its dr | |
168 | fdNdxiH = 0; | |
169 | fdNdxiB = 0; | |
170 | fdNdxiW = 0; | |
171 | fPtEneH = 0; | |
172 | fdrdNdxi = 0.7; | |
173 | ||
174 | // initialize weight for dE/dr histo | |
175 | SetdEdrWeight(); | |
176 | } | |
177 | ||
178 | AliJetAnalysis::~AliJetAnalysis() | |
179 | { | |
180 | // Destructor | |
181 | } | |
182 | ||
183 | ||
184 | //////////////////////////////////////////////////////////////////////// | |
185 | // define histogrames | |
186 | ||
187 | void AliJetAnalysis::DefineHistograms() | |
188 | { | |
189 | // Define the histograms to be filled | |
190 | if (fDoKine) DefineKineH(); | |
191 | if (fDoCorr) DefineCorrH(); | |
192 | if (fDoCorr50) DefineCorr50H(); | |
193 | if (fDoShap) DefineShapH(); | |
194 | if (fDoFrag) DefineFragH(); | |
195 | if (fDoTrig) DefineTrigH(); | |
196 | if (fDoJt) DefineJtH(); | |
197 | if (fDodNdxi) DefinedNdxiH(); | |
198 | } | |
199 | ||
200 | void AliJetAnalysis::DefineKineH() | |
201 | { | |
202 | // leading particle | |
203 | if (fDoPart) { | |
204 | fPKineEneH = new TH1F("PKineEne","Energy of leading particle",50,0,200); | |
205 | SetProperties(fPKineEneH,"Energy (GeV)","Entries"); | |
206 | fPKinePtH = new TH1F("PKinePt","Pt of leading particle",50,0,200); | |
207 | SetProperties(fPKinePtH,"P_{T} (GeV)","Entries"); | |
208 | fPKinePhiH = new TH1F("PKinePhiH","Azimuthal angle of leading particle", | |
209 | 90,0.,2.0*TMath::Pi()); | |
210 | SetProperties(fPKinePhiH,"#phi","Entries"); | |
211 | fPKineEtaH = new TH1F("PKineEtaH","Pseudorapidity of leading particle", | |
212 | 40,-1.0,1.0); | |
213 | SetProperties(fPKineEtaH,"#eta","Entries"); | |
214 | } | |
215 | // leading generated jet | |
216 | if (fDoGenJ) { | |
217 | fGKineEneH = new TH1F("GKineEne","Energy of generated jet",50,0,200); | |
218 | SetProperties(fGKineEneH,"Energy (GeV)","Entries"); | |
219 | fGKinePtH = new TH1F("GKinePt","Pt of generated jet",50,0,200); | |
220 | SetProperties(fGKinePtH,"P_{T} (GeV)","Entries"); | |
221 | fGKinePhiH = new TH1F("GKinePhiH","Azimuthal angle of generated jet", | |
222 | 90,0.,2.0*TMath::Pi()); | |
223 | SetProperties(fGKinePhiH,"#phi","Entries"); | |
224 | fGKineEtaH = new TH1F("GKineEtaH","Pseudorapidity of generated jet", | |
225 | 40,-1.0,1.0); | |
226 | SetProperties(fGKineEtaH,"#eta","Entries"); | |
227 | } | |
228 | // leading reconstructed jet | |
229 | if (fDoRecJ) { | |
230 | fRKineEneH = new TH1F("RKineEne","Energy of reconstructed jet",50,0,200); | |
231 | SetProperties(fRKineEneH,"Energy (GeV)","Entries"); | |
232 | fRKinePtH = new TH1F("RKinePt","Pt of reconstructed jet",50,0,200); | |
233 | SetProperties(fRKinePtH,"P_{T} (GeV)","Entries"); | |
234 | fRKinePhiH = new TH1F("RKinePhiH","Azimuthal angle of reconstructed jet", | |
235 | 90,0.,2.0*TMath::Pi()); | |
236 | SetProperties(fRKinePhiH,"#phi","Entries"); | |
237 | fRKineEtaH = new TH1F("RKineEtaH","Pseudorapidity of reconstructed jet", | |
238 | 40,-1.0,1.0); | |
239 | SetProperties(fRKineEtaH,"#eta","Entries"); | |
240 | } | |
241 | } | |
62567a80 | 242 | |
83a444b1 | 243 | void AliJetAnalysis::DefineCorrH() |
62567a80 | 244 | { |
83a444b1 | 245 | // correlation |
246 | if (fDoPart && fDoGenJ) { | |
247 | fPGCorrEneH = new TH2F("PGCorrEne","Energy correlation part-gen jet", | |
248 | 40,0,200,40,0,200); | |
249 | SetProperties(fPGCorrEneH,"Part Energy (GeV)","Gen Jet Energy (GeV)"); | |
250 | fPGCorrPtH = new TH2F("PGCorrPt","Pt correlation part-gen jet", | |
251 | 40,0,200,40,0,200); | |
252 | SetProperties(fPGCorrPtH,"Part P_{T} (GeV)","Gen Jet P_{T} (GeV)"); | |
253 | fPGCorrEtaH = new TH2F("PGCorrEta","Pseudorapidity correlation part-gen jet", | |
254 | 40,-1.0,1.0,40,-1.0,1.0); | |
255 | SetProperties(fPGCorrEtaH,"Part #eta","Gen Jet #eta"); | |
256 | fPGCorrPhiH = new TH2F("PGCorrPhi","Azimuthal angle correlation part-gen jet", | |
257 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
258 | SetProperties(fPGCorrPhiH,"Part #phi","Gen Jet #phi"); | |
259 | } | |
260 | if (fDoPart && fDoRecJ) { | |
261 | fPRCorrEneH = new TH2F("PRCorrEne","Energy correlation part-rec jet", | |
262 | 40,0,200,40,0,200); | |
263 | SetProperties(fPRCorrEneH,"Part Jet Energy (GeV)","Rec Jet Energy (GeV)"); | |
264 | fPRCorrPtH = new TH2F("PRCorrPt","Pt correlation part-rec jet", | |
265 | 40,0,200,40,0,200); | |
266 | SetProperties(fPRCorrPtH,"Part Jet P_{T} (GeV)","Rec Jet P_{T} (GeV)"); | |
267 | fPRCorrEtaH = new TH2F("PRCorrEta","Pseudorapidity correlation part-rec jet", | |
268 | 40,-1.0,1.0,40,-1.0,1.0); | |
269 | SetProperties(fPRCorrEtaH,"part #eta","Rec Jet #eta"); | |
270 | fPRCorrPhiH = new TH2F("PRCorrPhi","Azimuthal angle correlation part-rec jet", | |
271 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
272 | SetProperties(fPRCorrPhiH,"Part #phi","Rec Jet #phi"); | |
273 | } | |
274 | if (fDoGenJ && fDoRecJ) { | |
275 | fRGCorrEneH = new TH2F("RGCorrEne","Energy correlation rec jet-gen jet", | |
276 | 40,0,200,40,0,200); | |
277 | SetProperties(fRGCorrEneH,"Rec Jet Energy (GeV)","Gen Jet Energy (GeV)"); | |
278 | fRGCorrPtH = new TH2F("RGCorrPt","Pt correlation rec jet-gen jet", | |
279 | 40,0,200,40,0,200); | |
280 | SetProperties(fRGCorrPtH,"Rec Jet P_{T} (GeV)","Gen Jet P_{T} (GeV)"); | |
281 | fRGCorrEtaH = new TH2F("RGCorrEta","Pseudorapidity correlation rec jet-gen jet", | |
282 | 40,-1.0,1.0,40,-1.0,1.0); | |
283 | SetProperties(fRGCorrEtaH,"Rec Jet #eta","Gen Jet #eta"); | |
284 | fRGCorrPhiH = new TH2F("RGCorrPhi","Azimuthal angle correlation rec jet-gen jet", | |
285 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
286 | SetProperties(fRGCorrPhiH,"Rec Jet #phi","Gen Jet #phi"); | |
287 | } | |
62567a80 | 288 | } |
289 | ||
83a444b1 | 290 | void AliJetAnalysis::DefineCorr50H() |
62567a80 | 291 | { |
83a444b1 | 292 | // correlation |
293 | if (fDoPart && fDoRecJ) { | |
294 | fPRCorr50EneH = new TH2F("PRCorr50Ene","Energy correlation part-rec jet", | |
295 | 40,0,200,40,0,200); | |
296 | SetProperties(fPRCorr50EneH,"Part Jet Energy (GeV)","Rec Jet Energy (GeV)"); | |
297 | fPRCorr50PtH = new TH2F("PRCorr50Pt","Pt correlation part-rec jet", | |
298 | 40,0,200,40,0,200); | |
299 | SetProperties(fPRCorr50PtH,"Part Jet P_{T} (GeV)","Rec Jet P_{T} (GeV)"); | |
300 | fPRCorr50EtaH = new TH2F("PRCorr50Eta","Pseudorapidity correlation part-rec jet", | |
301 | 40,-1.0,1.0,40,-1.0,1.0); | |
302 | SetProperties(fPRCorr50EtaH,"part #eta","Rec Jet #eta"); | |
303 | fPRCorr50PhiH = new TH2F("PRCorr50Phi","Azimuthal angle correlation part-rec jet", | |
304 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
305 | SetProperties(fPRCorr50PhiH,"Part #phi","Rec Jet #phi"); | |
306 | } | |
307 | ||
308 | if (fDoGenJ && fDoRecJ) { | |
309 | fRGCorr50EneH = new TH2F("RGCorr50Ene","Energy correlation rec jet-gen jet", | |
310 | 40,0,200,40,0,200); | |
311 | SetProperties(fRGCorr50EneH,"Rec Jet Energy (GeV)","Gen Jet Energy (GeV)"); | |
312 | fRGCorr50PtH = new TH2F("RGCorr50Pt","Pt correlation rec jet-gen jet", | |
313 | 40,0,200,40,0,200); | |
314 | SetProperties(fRGCorr50PtH,"Rec Jet P_{T} (GeV)","Gen Jet P_{T} (GeV)"); | |
315 | fRGCorr50EtaH = new TH2F("RGCorr50Eta","Pseudorapidity correlation rec jet-gen jet", | |
316 | 40,-1.0,1.0,40,-1.0,1.0); | |
317 | SetProperties(fRGCorr50EtaH,"Rec Jet #eta","Gen Jet #eta"); | |
318 | fRGCorr50PhiH = new TH2F("RGCorr50Phi","Azimuthal angle correlation rec jet-gen jet", | |
319 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
320 | SetProperties(fRGCorr50PhiH,"Rec Jet #phi","Gen Jet #phi"); | |
321 | } | |
322 | } | |
323 | ||
324 | void AliJetAnalysis::DefineShapH() | |
325 | { | |
326 | // leading reconstructed jet | |
327 | if (fDoRecJ) { | |
328 | fdEdrH = new TH2F("fdEdrH","dE/dr histo",20,0,1,40,0,200); | |
329 | SetProperties(fdEdrH,"r","Rec Jet P_{T}"); | |
330 | fdEdrB = new TH2F("fdEdrB","dE/dr bkgdhisto",20,0,1,40,0,200); | |
331 | SetProperties(fdEdrB,"r","Rec P_{T}"); | |
332 | fPtEneH2 = new TH2F("fPtEneH2","fPtEneH2",40,0,200,40,0,200); | |
333 | SetProperties(fPtEneH2,"Rec Jet E","Rec Jet P_{T}"); | |
334 | fdEdrW = new TH1F("fdEdrW","weights for dE/dr",40,0,200); | |
335 | SetProperties(fdEdrW,"Rec Jet P_{T}","weights"); | |
62567a80 | 336 | |
83a444b1 | 337 | fRShapSelH = new TH1F("RShapSel","Shape of generated jets (sel part)",20,0.,1.); |
338 | SetProperties(fRShapSelH,"r","1/N_{JET}#Sigma P_{T}(0,r)/P_{T}(0,R_{JET})"); | |
339 | fRShapRejH = new TH1F("RShapRej","Shape of generated jets (rej part)",20,0.,1.); | |
340 | SetProperties(fRShapRejH,"r","1/N_{JET}#Sigma P_{T}(0,r)/P_{T}(0,R_{JET})"); | |
341 | fRShapAllH = new TH1F("RShapAll","Shape of generated jets (all part)",20,0.,1.); | |
342 | SetProperties(fRShapAllH,"r","1/N_{JET}#Sigma P_{T}(0,r)/P_{T}(0,R_{JET})"); | |
343 | } | |
344 | } | |
62567a80 | 345 | |
83a444b1 | 346 | void AliJetAnalysis::DefineJtH() |
347 | { | |
348 | // Define the histogram for J_T | |
349 | if (fDoRecJ) { | |
350 | fJtH = new TH2F("fJtH","J_{T} histogram",80,0.,4.,40,0.,200.); | |
351 | SetProperties(fJtH,"J_{T}","Rec Jet P_{T}"); | |
352 | fJtB = new TH2F("fJtB","J_{T} bkgd histogram",80,0.,4.,40,0.,200.); | |
353 | SetProperties(fJtB,"J_{T}","Rec P_{T}"); | |
354 | fJtW = new TH1F("fJtW","J_{T} weight",40,0,200); | |
355 | SetProperties(fJtW,"J_{T}W","weight"); | |
356 | } | |
357 | } | |
62567a80 | 358 | |
83a444b1 | 359 | void AliJetAnalysis::DefinedNdxiH() |
360 | { | |
361 | // Define the histogram for dN/dxi | |
362 | if (fDoRecJ) { | |
363 | fdNdxiH = new TH2F("fdNdxiH","dN/d#xi histo",200,0,10,40,0,200); | |
364 | SetProperties(fdNdxiH,"xi","Rec Jet P_{T}"); | |
365 | fdNdxiB = new TH2F("fdNdxiB","dN/d#xi bkgd histo",200,0,10,40,0,200); | |
366 | SetProperties(fdNdxiB,"xi","Rec P_{T}"); | |
367 | fdNdxiW = new TH1F("fdNdxiW","dN/d#xi histo",40,0,200); | |
368 | SetProperties(fdNdxiW,"Rec Jet P_{T}","weights"); | |
369 | fPtEneH = new TH2F("fPtEneH","fPtEneH",40,0,200,40,0,200); | |
370 | SetProperties(fPtEneH,"Rec Jet E","Rec Jet P_{T}"); | |
371 | } | |
372 | } | |
62567a80 | 373 | |
83a444b1 | 374 | void AliJetAnalysis::DefineFragH() |
375 | { | |
376 | // leading reconstructed jet | |
377 | if (fDoRecJ) { | |
378 | fRFragSelH = new TH1F("RFragSel","Frag Fun of reconstructed jets (sel part)",20,0.,10.); | |
379 | SetProperties(fRFragSelH,"#xi=ln(E_{JET}/E_{i})","1/N_{JET}dN_{ch}/d#xi"); | |
380 | fRFragRejH = new TH1F("RFragRej","Frag Fun of reconstructed jets (rej part)",20,0.,10.); | |
381 | SetProperties(fRFragRejH,"#xi=ln(E_{JET}/E_{i})","1/N_{JET}dN_{ch}/d#xi"); | |
382 | fRFragAllH = new TH1F("RFragAll","Frag Fun of reconstructed jets (all part)",20,0.,10.); | |
383 | SetProperties(fRFragAllH,"#xi=ln(E_{JET}/E_{i})","1/N_{JET}dN_{ch}/d#xi"); | |
384 | } | |
385 | } | |
62567a80 | 386 | |
83a444b1 | 387 | void AliJetAnalysis::DefineTrigH() |
388 | { | |
389 | // generated energy | |
390 | fGTriggerEneH = new TProfile("GTriggerEne","Generated energy (trigger bias)", | |
391 | 20, 0., 200., 0., 1., "S"); | |
392 | fGTriggerEneH->SetXTitle("E_{gen}"); | |
393 | fGTriggerEneH->SetYTitle("E_{rec}/E_{gen}"); | |
394 | // reconstructed energy | |
395 | fRTriggerEneH = new TProfile("RTriggerEne","Reconstructed energy (trigger bias)", | |
396 | 20, 0., 200., 0., 1., "S"); | |
397 | fRTriggerEneH->SetXTitle("E_{rec}"); | |
398 | fRTriggerEneH->SetYTitle("E_{rec}/E_{gen}"); | |
399 | // generated energy vs generated/leading | |
400 | fGPTriggerEneH = new TProfile("GPTriggerEne","Generated energy (trigger bias)", | |
401 | 20, 0., 200., 0., 1., "S"); | |
402 | fGPTriggerEneH->SetXTitle("E_{gen}"); | |
403 | fGPTriggerEneH->SetYTitle("E_{L}/E_{gen}"); | |
404 | // leading particle energy | |
405 | fPTriggerEneH = new TProfile("PTriggerEne","Leading particle energy (trigger bias)", | |
406 | 20, 0., 200., 0., 1., "S"); | |
407 | fPTriggerEneH->SetXTitle("E_{L}/0.2"); | |
408 | fPTriggerEneH->SetYTitle("E_{L}/E_{gen}"); | |
62567a80 | 409 | |
83a444b1 | 410 | } |
62567a80 | 411 | |
83a444b1 | 412 | void AliJetAnalysis::SetProperties(TH1* h,const char* x, const char* y) const |
413 | { | |
414 | //Set properties of histos (x and y title and error propagation) | |
415 | h->SetXTitle(x); | |
416 | h->SetYTitle(y); | |
417 | h->Sumw2(); | |
418 | } | |
419 | ||
420 | //////////////////////////////////////////////////////////////////////// | |
421 | // compute weights for dE/dr histogram | |
62567a80 | 422 | |
83a444b1 | 423 | void AliJetAnalysis::SetdEdrWeight() |
424 | { | |
425 | // Due to the limited acceptance, each bin in the dE/dr has a different | |
426 | // weight given by the ratio of the area of a disk dR to the area | |
427 | // within the eta acceptance. The weight depends on the eta position | |
428 | // of the jet. Here a look up table for the weights is computed. It | |
429 | // assumes |etaJet|<0.5 and |eta_lego|<0.9. It makes bins of 0.05 | |
430 | // units in eta. Note that this is fixed by the bin width chosen for | |
431 | // the histogram --> this weights are tailored for the specific | |
432 | // histogram definition used here! | |
433 | ||
434 | // two dimensional table: first index, bin in eta of jet, second | |
435 | // index bin of dE/dr histo | |
436 | ||
437 | Int_t nBins = 20; | |
438 | Float_t xMin = 0.0; | |
439 | Float_t xMax = 1.0; | |
440 | Float_t binSize = (xMax-xMin)/nBins; | |
62567a80 | 441 | |
83a444b1 | 442 | Float_t da,ds,r1,r2,h1,h2,a1,a2,theta1,theta2,area1,area2; |
443 | Int_t ji; | |
444 | for (Int_t i=0;i<(nBins/2);i++) { | |
445 | // index of first histo bin needing a scale factor | |
446 | ji=(nBins-2)-i; | |
447 | for(Int_t j=0;j<nBins;j++) { | |
448 | // area of ring. | |
449 | da = TMath::Pi()*(binSize*binSize)*(1.0+2.0*j); | |
450 | ds = 1.0; | |
451 | if (j>=ji) { // compute missing area using segments of circle | |
452 | r1=j*binSize; | |
453 | r2=(j+1)*binSize; | |
454 | h1=(j-ji)*binSize; | |
455 | h2=(j+1-ji)*binSize; | |
456 | a1=2.0*TMath::Sqrt(2.0*h1*r1-h1*h1); | |
457 | a2=2.0*TMath::Sqrt(2.0*h2*r2-h2*h2); | |
458 | theta1=2*TMath::ACos((r1-h1)/r1); | |
459 | theta2=2*TMath::ACos((r2-h2)/r2); | |
460 | area1=binSize*(r1*r1*theta1-a1*(r1-h1)); | |
461 | area2=binSize*(r2*r2*theta2-a2*(r2-h2)); | |
462 | ds = (da-(area2-area1))/da; | |
463 | } | |
464 | fWeightdEdr[i][j]=ds/da; | |
465 | } | |
466 | } | |
467 | } | |
468 | ||
469 | // get weight for dE/dr histogram | |
470 | Float_t AliJetAnalysis::GetdEdrWeight(Float_t etaJet, Float_t r) | |
471 | { | |
472 | // Return the correponding weight for the dE/dr plot | |
473 | Int_t nBins = 20; | |
474 | Float_t xMin = 0.0; | |
475 | Float_t xMax = 1.0; | |
476 | Float_t binSize = (xMax-xMin)/nBins; | |
477 | ||
478 | Float_t eta = TMath::Abs(etaJet); | |
479 | if ((eta > 0.5) || (r > fdrdEdr)) return 0.0; | |
480 | Int_t binJet = (Int_t) TMath::Floor(eta/binSize); | |
481 | Int_t binR = (Int_t) TMath::Floor(r/binSize); | |
482 | Float_t w = fWeightdEdr[binJet][binR]; | |
483 | return w; | |
484 | } | |
485 | ||
486 | ||
487 | //////////////////////////////////////////////////////////////////////// | |
488 | // fill histograms | |
489 | ||
490 | void AliJetAnalysis::FillHistograms() | |
491 | { | |
492 | // fill histograms | |
62567a80 | 493 | |
494 | // Run data | |
83a444b1 | 495 | AliJetProductionDataPDC2004* runData = new AliJetProductionDataPDC2004(); |
496 | ||
497 | // Loop over runs | |
498 | TFile* jFile = 0x0; | |
499 | ||
500 | for (Int_t iRun = fRunMin; iRun <= fRunMax; iRun++) { | |
501 | // Open file | |
502 | char fn[256]; | |
503 | sprintf(fn, "%s/%s.root", fDirectory, (runData->GetRunTitle(iRun)).Data()); | |
504 | jFile = new TFile(fn); | |
505 | printf(" Analyzing run: %d %s\n", iRun,fn); | |
506 | ||
507 | // Get reader header and events to be looped over | |
508 | AliJetReaderHeader *jReaderH = (AliJetReaderHeader*)(jFile->Get(fReaderHeader)); | |
509 | if (fEventMin == -1) fEventMin = jReaderH->GetFirstEvent(); | |
510 | if (fEventMax == -1) { | |
511 | fEventMax = jReaderH->GetLastEvent(); | |
512 | } else { | |
513 | fEventMax = TMath::Min(fEventMax, jReaderH->GetLastEvent()); | |
514 | } | |
515 | ||
516 | AliUA1JetHeader *jH = (AliUA1JetHeader *) (jFile->Get("AliUA1JetHeader")); | |
62567a80 | 517 | |
83a444b1 | 518 | // Calculate weight |
519 | fWeight = runData->GetWeight(iRun)/ ( (Float_t) (fEventMax - fEventMin + 1)); | |
62567a80 | 520 | |
83a444b1 | 521 | // Loop over events |
522 | for (Int_t i = fEventMin; i < fEventMax; i++) { | |
523 | if (i%100 == 0) printf(" Analyzing event %d / %d \n",i,fEventMax); | |
524 | ||
525 | // Get next tree with AliJet | |
526 | char nameT[100]; | |
527 | sprintf(nameT, "TreeJ%d",i); | |
528 | TTree *jetT =(TTree *)(jFile->Get(nameT)); | |
529 | if (fDoRecJ) jetT->SetBranchAddress("FoundJet", &fRecJ); | |
530 | if (fDoGenJ) jetT->SetBranchAddress("GenJet", &fGenJ); | |
531 | if (fDoPart) jetT->SetBranchAddress("LeadingPart", &fPart); | |
532 | jetT->GetEntry(0); | |
533 | ||
534 | int nJets = fRecJ->GetNJets(); | |
535 | ||
536 | TArrayI inJet = fRecJ->GetInJet(); | |
537 | if(inJet.GetSize()>fminMult){ // removes events with low multiplicity | |
538 | if (fDoKine) FillKineH(); | |
539 | if (fDoCorr) FillCorrH(); | |
540 | if (fDoCorr50) FillCorr50H(); | |
541 | if (fDoShap) FillShapH(jH->GetRadius()); | |
542 | if (fDoFrag) FillFragH(); | |
543 | if (fDoTrig) FillTrigH(); | |
544 | if (fDoJt) FillJtH(); | |
545 | if (fDodNdxi) FilldNdxiH(); | |
546 | } | |
547 | delete jetT; // jet should be deleted before creating a new one | |
548 | if(inJet.GetSize()>fminMult){ // removes events with low multiplicity | |
549 | if (fDoBkgd && nJets>0) FillBkgd(i,iRun); | |
550 | } | |
551 | } // end loop over events in one file | |
552 | if (jFile) jFile->Close(); | |
553 | delete jFile; | |
554 | } // end loop over files | |
555 | } | |
556 | ||
557 | void AliJetAnalysis::FillKineH() | |
558 | { | |
559 | // leading particle | |
560 | if (fDoPart && fPart->LeadingFound()){ | |
561 | fPKineEneH->Fill(fPart->GetE(),fWeight); | |
562 | fPKinePtH->Fill(fPart->GetPt(),fWeight); | |
563 | fPKinePhiH->Fill(fPart->GetPhi(),fWeight); | |
564 | fPKineEtaH->Fill(fPart->GetEta(),fWeight); | |
565 | } | |
566 | // leading generated jet | |
567 | if (fDoGenJ && fGenJ->GetNJets()> 0){ | |
568 | fGKineEneH->Fill(fGenJ->GetE(0),fWeight); | |
569 | fGKinePtH->Fill(fGenJ->GetPt(0),fWeight); | |
570 | fGKinePhiH->Fill(fGenJ->GetPhi(0),fWeight); | |
571 | fGKineEtaH->Fill(fGenJ->GetEta(0),fWeight); | |
572 | } | |
573 | // leading reconstructed jet | |
574 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
575 | TArrayF p=fRecJ->GetPtFromSignal(); | |
576 | if (p[0]>fPercentage) { | |
577 | fRKineEneH->Fill(fRecJ->GetE(0)/fEfactor,fWeight); | |
578 | fRKinePtH->Fill(fRecJ->GetPt(0),fWeight); | |
579 | fRKinePhiH->Fill(fRecJ->GetPhi(0),fWeight); | |
580 | fRKineEtaH->Fill(fRecJ->GetEta(0),fWeight); | |
581 | } | |
582 | } | |
583 | } | |
584 | ||
585 | void AliJetAnalysis::FillCorrH() | |
586 | { | |
587 | // Fill correlation histograms | |
588 | if (fDoPart && fPart->LeadingFound() && fDoGenJ && fGenJ->GetNJets()> 0) | |
589 | Correlation(fPart->GetLeading(),fGenJ->GetJet(0), | |
590 | fPGCorrEneH,fPGCorrPtH,fPGCorrEtaH,fPGCorrPhiH); | |
591 | if (fDoPart && fPart->LeadingFound() && fDoRecJ && fRecJ->GetNJets()> 0) { | |
592 | TArrayF p=fRecJ->GetPtFromSignal(); | |
593 | if (p[0]>fPercentage) | |
594 | Correlation(fPart->GetLeading(),fRecJ->GetJet(0), | |
595 | fPRCorrEneH,fPRCorrPtH,fPRCorrEtaH,fPRCorrPhiH); | |
596 | } | |
597 | if (fDoGenJ && fGenJ->GetNJets()> 0 && fDoRecJ && fRecJ->GetNJets()> 0) { | |
598 | TArrayF p=fRecJ->GetPtFromSignal(); | |
599 | if (p[0]>fPercentage) | |
600 | Correlation(fRecJ->GetJet(0), fGenJ->GetJet(0), | |
601 | fRGCorrEneH,fRGCorrPtH,fRGCorrEtaH,fRGCorrPhiH); | |
602 | } | |
603 | } | |
604 | ||
605 | void AliJetAnalysis::Correlation(TLorentzVector *lv1,TLorentzVector *lv2, | |
606 | TH2F *h1, TH2F *h2, TH2F *h3, TH2F *h4) | |
607 | { | |
608 | // Correlation histograms | |
609 | h1->Fill(lv1->E(),lv2->E(),fWeight); | |
610 | h2->Fill(lv1->Pt(),lv2->Pt(),fWeight); | |
611 | h3->Fill(lv1->Eta(),lv2->Eta(),fWeight); | |
612 | Float_t p1, p2; | |
613 | p1 = ((lv1->Phi() < 0) ? (lv1->Phi()) + 2. * TMath::Pi() : lv1->Phi()); | |
614 | p2 = ((lv2->Phi() < 0) ? (lv2->Phi()) + 2. * TMath::Pi() : lv2->Phi()); | |
615 | h4->Fill(p1,p2,fWeight); | |
616 | } | |
617 | ||
618 | void AliJetAnalysis::FillCorr50H() | |
619 | { | |
620 | // Fill correlation histograms when one particle has > 50% of jet energy | |
621 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
622 | TArrayF p = fRecJ->GetPtFromSignal(); | |
623 | if (p[0]>fPercentage) { | |
624 | if (fDoPart && fPart->LeadingFound()) | |
625 | Correlation50(fRecJ, fPart->GetLeading(),fRecJ->GetJet(0), | |
626 | fPRCorr50EneH,fPRCorr50PtH,fPRCorr50EtaH,fPRCorr50PhiH); | |
627 | if (fDoGenJ && fGenJ->GetNJets()> 0) | |
628 | Correlation50(fRecJ, fRecJ->GetJet(0), fGenJ->GetJet(0), | |
629 | fRGCorr50EneH,fRGCorr50PtH,fRGCorr50EtaH,fRGCorr50PhiH); | |
630 | } | |
631 | } | |
632 | } | |
633 | ||
634 | void AliJetAnalysis::Correlation50(AliJet *j,TLorentzVector *lv1,TLorentzVector *lv2, | |
635 | TH2F *h1, TH2F *h2, TH2F *h3, TH2F *h4) | |
636 | { | |
637 | // Correlation histograms when one particle has > 50% of jet energy | |
638 | TArrayF ptin = j->GetPtIn(); | |
639 | TArrayF etain = j->GetEtaIn(); | |
640 | TArrayI inJet = j->GetInJet(); | |
641 | ||
642 | Int_t flag = 0; | |
643 | for(Int_t i=0;i<(inJet.GetSize());i++) { | |
644 | if (inJet[i] == 1) { | |
645 | Float_t t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
646 | Float_t x = (ptin[i]*TMath::Sqrt(1.+1./(t1*t1)))/(j->GetE(0)); | |
647 | if (x>0.5) flag++; | |
648 | } | |
649 | } | |
650 | if (flag>1) cout << " Flag = " << flag << endl; | |
651 | if (flag == 1) { | |
652 | h1->Fill(lv1->E(),lv2->E(),fWeight); | |
653 | h2->Fill(lv1->Pt(),lv2->Pt(),fWeight); | |
654 | h3->Fill(lv1->Eta(),lv2->Eta(),fWeight); | |
655 | Float_t p1, p2; | |
656 | p1 = ((lv1->Phi() < 0) ? | |
657 | (lv1->Phi()) + 2. * TMath::Pi() : lv1->Phi()); | |
658 | p2 = ((lv2->Phi() < 0) ? | |
659 | (lv2->Phi()) + 2. * TMath::Pi() : lv2->Phi()); | |
660 | h4->Fill(p1,p2,fWeight); | |
661 | } | |
662 | } | |
663 | ||
664 | void AliJetAnalysis::FillJtH() | |
665 | { | |
666 | // Fill J_T histogram | |
667 | if (fRecJ->GetNJets()> 0) { | |
668 | fjv3X = 0.0; fjv3Y = 0.0; fjv3Z = 0.0; | |
669 | TArrayF p=fRecJ->GetPtFromSignal(); | |
670 | if (p[0]>fPercentage) { | |
671 | // initialize | |
672 | const TVector3 kJv3 = fRecJ->GetJet(0)->Vect(); | |
673 | fjv3X = kJv3.X(); fjv3Y = kJv3.Y(); fjv3Z = kJv3.Z(); | |
674 | TVector3 trk; | |
675 | TArrayI inJet = fRecJ->GetInJet(); | |
676 | TArrayF etain = fRecJ->GetEtaIn(); | |
677 | TArrayF ptin = fRecJ->GetPtIn(); | |
678 | TArrayF phiin = fRecJ->GetPhiIn(); | |
679 | Float_t deta, dphi,jt, dr; | |
680 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
681 | deta = etain[i] - fRecJ->GetEta(0); | |
682 | dphi = phiin[i] - fRecJ->GetPhi(0); | |
683 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
684 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
685 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
686 | if ((dr<fdrJt) && (ptin[i] > fPartPtCut)) { | |
687 | trk.SetPtEtaPhi(ptin[i],etain[i],phiin[i]); | |
688 | jt = TMath::Sin(trk.Angle(kJv3))*trk.Mag(); | |
689 | fJtH->Fill(jt,fRecJ->GetPt(0),fWeight); | |
62567a80 | 690 | } |
83a444b1 | 691 | } |
692 | fJtW->Fill(fRecJ->GetPt(0),fWeight); | |
693 | fWJt+=fWeight; | |
694 | } | |
695 | } | |
696 | } | |
697 | ||
698 | void AliJetAnalysis::FilldNdxiH() | |
699 | { | |
700 | // Fill dN/d#xi histograms | |
701 | if (fRecJ->GetNJets()> 0) { | |
702 | TArrayF p=fRecJ->GetPtFromSignal(); | |
703 | if (p[0]>fPercentage) { | |
704 | fp0 = p[0]; | |
705 | TArrayI inJet = fRecJ->GetInJet(); | |
706 | TArrayF etain = fRecJ->GetEtaIn(); | |
707 | TArrayF ptin = fRecJ->GetPtIn(); | |
708 | TArrayF phiin = fRecJ->GetPhiIn(); | |
709 | Float_t xi,t1,ene,dr,deta,dphi; | |
710 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
711 | deta = etain[i] - fRecJ->GetEta(0); | |
712 | dphi = phiin[i] - fRecJ->GetPhi(0); | |
713 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
714 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
715 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
716 | if ((dr<fdrdNdxi) && (ptin[i] > fPartPtCut)) { | |
717 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
718 | ene = ptin[i]*TMath::Sqrt(1.+1./(t1*t1)); | |
719 | xi = (Float_t) TMath::Log((fRecJ->GetE(0)/fEfactor)/ene); | |
720 | fdNdxiH->Fill(xi,fRecJ->GetPt(0),fWeight); | |
721 | } | |
722 | } | |
723 | fdNdxiW->Fill(fRecJ->GetPt(0),fWeight); | |
724 | fPtEneH->Fill(fRecJ->GetE(0)/fEfactor,fRecJ->GetPt(0),fWeight); | |
725 | fWdNdxi+=fWeight; | |
726 | } | |
727 | } | |
728 | } | |
729 | ||
730 | void AliJetAnalysis::FillBkgd(Int_t eventN, Int_t runN) | |
731 | { | |
732 | // Background calculated with hijing events (no pythia) | |
733 | if (fp0>fPercentage) { | |
734 | fPtJ=0.,fEJ=0.,fEtaJ=0.,fPhiJ=0.; | |
735 | // Background calculated with hijing events (no pythia) | |
736 | AliJetProductionDataPDC2004* runDataB = new AliJetProductionDataPDC2004(); | |
737 | TFile* jFileB =0x0;; | |
738 | char fnB[256]; | |
62567a80 | 739 | |
83a444b1 | 740 | sprintf(fnB, "%s/%s.root", fBkgdDirectory, (runDataB->GetRunTitle(runN)).Data()); |
741 | jFileB = new TFile(fnB); | |
742 | ||
743 | char nameB[100]; | |
744 | sprintf(nameB, "TreeJ%d",eventN); | |
745 | TTree *jetB =(TTree *)(jFileB->Get(nameB)); | |
746 | jetB->SetBranchAddress("FoundJet", &fRecB); | |
747 | jetB->GetEntry(0); | |
62567a80 | 748 | |
83a444b1 | 749 | TArrayI inJetB = fRecB->GetInJet(); |
750 | TArrayF etainB = fRecB->GetEtaIn(); | |
751 | TArrayF ptinB = fRecB->GetPtIn(); | |
752 | TArrayF phiinB = fRecB->GetPhiIn(); | |
753 | fPtJ = fRecJ->GetPt(0); | |
754 | fEJ = fRecJ->GetE(0); | |
755 | fEtaJ = fRecJ->GetEta(0); | |
756 | fPhiJ = fRecJ->GetPhi(0); | |
757 | Float_t t1,ene,xi,detaB,dphiB,drB,jt,wB; | |
758 | TVector3 trkB; | |
759 | TVector3 jv3B; | |
760 | jv3B.SetX(fjv3X); jv3B.SetY(fjv3Y); jv3B.SetZ(fjv3Z); | |
761 | ||
762 | for(Int_t k=0;k<inJetB.GetSize();k++){ | |
763 | if(ptinB[k] > fPartPtCut){ | |
764 | detaB = etainB[k] - fEtaJ; | |
765 | dphiB = phiinB[k] - fPhiJ; | |
766 | if (dphiB < -TMath::Pi()) dphiB= -dphiB - 2.0 * TMath::Pi(); | |
767 | if (dphiB > TMath::Pi()) dphiB = 2.0 * TMath::Pi() - dphiB; | |
768 | drB = TMath::Sqrt(detaB * detaB + dphiB * dphiB); | |
769 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etainB[k]))); | |
770 | ene = ptinB[k]*TMath::Sqrt(1.+1./(t1*t1)); | |
771 | trkB.SetPtEtaPhi(ptinB[k],etainB[k],phiinB[k]); | |
772 | // --- dN/dxi | |
773 | if (drB<fdrdNdxi) { | |
774 | xi = (Float_t) TMath::Log(fEJ/ene); | |
775 | fdNdxiB->Fill(xi,fPtJ,fWeight); | |
776 | } | |
777 | // --- Jt | |
778 | if (drB<fdrJt) { | |
779 | jt = TMath::Sin(trkB.Angle(jv3B))*(trkB.Mag()); | |
780 | fJtB->Fill(jt,fPtJ,fWeight); | |
781 | } | |
782 | // --- dE/dr | |
783 | if (drB<fdrdEdr) { | |
784 | wB = GetdEdrWeight(fEtaJ,drB)*fWeight*ene; | |
785 | fdEdrB->Fill(drB,fPtJ,wB); | |
786 | } | |
787 | } | |
788 | } | |
789 | delete jetB; | |
790 | if (jFileB) jFileB->Close(); | |
791 | delete jFileB; | |
792 | } | |
793 | } | |
62567a80 | 794 | |
83a444b1 | 795 | void AliJetAnalysis::FillShapH(Float_t r) |
796 | { | |
797 | // Fill jet shape histograms | |
798 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
799 | TArrayF p=fRecJ->GetPtFromSignal(); | |
800 | if (p[0]>fPercentage) { | |
801 | Shape(fRecJ,fRShapSelH,fRShapRejH,fRShapAllH,fdEdrH,fPtEneH2,fdEdrW,r); | |
802 | fWdEdr+=fWeight; | |
803 | fWShapR+=fWeight; | |
804 | } | |
805 | } | |
806 | } | |
62567a80 | 807 | |
83a444b1 | 808 | void AliJetAnalysis::Shape(AliJet *j,TH1F* hs, TH1F* hr, TH1F* ha, |
809 | TH2F *hdedr, TH2F *hptene, TH1F* wdedr, Float_t radius) | |
810 | { | |
811 | // initialize | |
812 | TArrayI inJet = j->GetInJet(); | |
813 | TArrayF etain = j->GetEtaIn(); | |
814 | TArrayF ptin = j->GetPtIn(); | |
815 | TArrayF phiin = j->GetPhiIn(); | |
62567a80 | 816 | |
83a444b1 | 817 | // first compute dE/dr histo |
818 | Float_t etaj = j->GetEta(0); | |
819 | Float_t ene,w,deta,dphi,dr,t1; | |
820 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
821 | deta = etain[i] - j->GetEta(0); | |
822 | dphi = phiin[i] - j->GetPhi(0); | |
823 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
824 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
825 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
826 | if ((dr<fdrdEdr) && (ptin[i] > fPartPtCut)) { | |
827 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
828 | ene = ptin[i]*TMath::Sqrt(1.+1./(t1*t1)); | |
829 | w = GetdEdrWeight(etaj,dr)*fWeight*ene; | |
830 | hdedr->Fill(dr,j->GetPt(0),w); | |
831 | } | |
832 | } | |
833 | hptene->Fill(fRecJ->GetE(0),fRecJ->GetPt(0),fWeight); | |
834 | wdedr->Fill(j->GetPt(0),fWeight); | |
62567a80 | 835 | |
83a444b1 | 836 | // now compute shape histos |
837 | Int_t nBins = ha->GetNbinsX(); | |
838 | Float_t xMin = ha->GetXaxis()->GetXmin(); | |
839 | Float_t xMax = ha->GetXaxis()->GetXmax(); | |
840 | Float_t binSize,halfBin; | |
841 | binSize = (xMax-xMin)/nBins; | |
842 | halfBin = binSize/2; | |
843 | Float_t rptS[20], rptR[20], rptA[20]; | |
844 | for(Int_t i=0;i<nBins;i++) rptS[i]=rptR[i]=rptA[i]=0.0; | |
845 | // fill bins in r for leading jet | |
846 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
847 | if (inJet[i] == 1 || (inJet[i] == -1 && fPythia)) { | |
848 | deta = etain[i] - j->GetEta(0); | |
849 | dphi = phiin[i] - j->GetPhi(0); | |
850 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
851 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
852 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
853 | Float_t rR = dr/radius; | |
854 | Int_t bin = (Int_t) TMath::Floor(rR/binSize); | |
855 | rptA[bin]+=ptin[i]/(j->GetPt(0)); | |
856 | if (inJet[i] == 1) rptS[bin]+=ptin[i]/(j->GetPt(0)); | |
857 | if (fPythia && inJet[i] == -1) | |
858 | rptR[bin]+=ptin[i]/(j->GetPt(0)); | |
859 | } | |
860 | } | |
62567a80 | 861 | |
83a444b1 | 862 | // compute shape and fill histogram |
863 | Float_t ptS,ptR,ptA,r; | |
864 | ptS=ptR=ptA=0.0; | |
865 | for (Int_t i=0;i<nBins;i++) { | |
866 | ptS+=rptS[i]; | |
867 | ptR+=rptR[i]; | |
868 | ptA+=rptA[i]; | |
869 | r=(i+1)*binSize-halfBin; | |
870 | hs->Fill(r,ptS*fWeight); | |
871 | if(fPythia) { | |
872 | hr->Fill(r,ptR*fWeight); | |
873 | ha->Fill(r,ptA*fWeight); | |
874 | } | |
875 | } | |
876 | } | |
877 | ||
878 | void AliJetAnalysis::FillFragH() | |
879 | { | |
880 | // Fill fragmentation histogram | |
881 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
882 | TArrayF p=fRecJ->GetPtFromSignal(); | |
883 | if (p[0]>fPercentage) { | |
884 | FragFun(fRecJ,fRFragSelH,fRFragRejH,fRFragAllH); | |
885 | fWFragR+=fWeight; | |
886 | } | |
887 | } | |
888 | } | |
62567a80 | 889 | |
83a444b1 | 890 | void AliJetAnalysis::FragFun(AliJet *j,TH1F* hs, TH1F* hr, TH1F* ha) |
891 | { | |
892 | // Calculate de fragmentation function | |
893 | TArrayI inJet = j->GetInJet(); | |
894 | TArrayF etain = j->GetEtaIn(); | |
895 | TArrayF ptin = j->GetPtIn(); | |
62567a80 | 896 | |
83a444b1 | 897 | Float_t xi,t1,ene; |
62567a80 | 898 | |
83a444b1 | 899 | for(Int_t i=0;i<(inJet.GetSize());i++) { |
900 | if (inJet[i] == 1 || (inJet[i] == -1 && fPythia)) { | |
901 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
902 | ene = ptin[i]*TMath::Sqrt(1.+1./(t1*t1)); | |
903 | xi = (Float_t) TMath::Log((j->GetE(0))/ene); | |
904 | if (fPythia) ha->Fill(xi,fWeight); | |
905 | if (inJet[i] == 1) hs->Fill(xi,fWeight); | |
906 | if (fPythia && inJet[i] == -1) hr->Fill(xi,fWeight); | |
907 | } | |
908 | } | |
909 | } | |
62567a80 | 910 | |
83a444b1 | 911 | void AliJetAnalysis:: FillTrigH() |
912 | { | |
913 | // Fill trigger bias histograms | |
914 | if(fDoTrig && fGenJ->GetNJets()>0 && fRecJ->GetNJets()>0){ | |
915 | TArrayF p=fRecJ->GetPtFromSignal(); | |
916 | if (p[0]>fPercentage) { | |
917 | float genEne = fGenJ->GetE(0); | |
918 | float recEne = fRecJ->GetE(0); | |
919 | float eneRatio = (recEne/fEfactor)/genEne; | |
920 | ||
921 | fGTriggerEneH->Fill(genEne, eneRatio, fWeight); | |
922 | fRTriggerEneH->Fill(recEne/fEfactor, eneRatio, fWeight); | |
923 | ||
924 | if (fPart->LeadingFound()){ | |
925 | float leaEne = fPart->GetE(); | |
926 | float eneRatio2 = leaEne/genEne; | |
927 | fGPTriggerEneH->Fill(genEne, eneRatio2, fWeight); | |
928 | fPTriggerEneH->Fill(leaEne/0.2, eneRatio2, fWeight); | |
929 | } | |
930 | } | |
931 | } | |
932 | } | |
933 | ||
934 | ||
935 | //////////////////////////////////////////////////////////////////////// | |
936 | // Normalize histogrames | |
937 | ||
938 | void AliJetAnalysis::NormHistograms() | |
939 | { | |
940 | // normalize shape histograms | |
941 | if (fDoShap) { | |
942 | if (fDoRecJ && fWShapR>0.) { // leading reconstructed jet | |
943 | fRShapSelH->Scale(1.0/fWShapR); | |
944 | fRShapRejH->Scale(1.0/fWShapR); | |
945 | fRShapAllH->Scale(1.0/fWShapR); | |
946 | } | |
947 | } | |
62567a80 | 948 | |
83a444b1 | 949 | // normalize fragmentation function histograms |
950 | if (fDoFrag) { | |
951 | if (fDoRecJ && fWFragR>0.) { // leading reconstructed jet | |
952 | fRFragSelH->Scale(2.0/fWFragR); | |
953 | fRFragRejH->Scale(2.0/fWFragR); | |
954 | fRFragAllH->Scale(2.0/fWFragR); | |
955 | } | |
956 | } | |
957 | } | |
958 | ||
959 | //////////////////////////////////////////////////////////////////////// | |
960 | ||
961 | void AliJetAnalysis::PlotHistograms() | |
962 | { | |
963 | // Plot histogramas (to be done...) | |
964 | if (fDoKine) PlotKineH(); | |
965 | if (fDoCorr) PlotCorrH(); | |
966 | if (fDoCorr50) PlotCorr50H(); | |
967 | if (fDoShap) PlotShapH(); | |
968 | if (fDoFrag) PlotFragH(); | |
969 | if (fDoTrig) PlotTrigH(); | |
970 | } | |
971 | ||
972 | void AliJetAnalysis::PlotKineH() const | |
973 | { | |
974 | // missing | |
975 | if (fDoPart) ; | |
976 | if (fDoGenJ) ; | |
977 | if (fDoRecJ) ; | |
978 | } | |
979 | ||
980 | void AliJetAnalysis::PlotCorrH() const | |
981 | { | |
982 | // missing | |
983 | if (fDoPart && fDoGenJ) ; | |
984 | if (fDoPart && fDoRecJ) ; | |
985 | if (fDoGenJ && fDoRecJ) ; | |
986 | } | |
987 | void AliJetAnalysis::PlotCorr50H() const | |
988 | { | |
989 | // missing | |
990 | if (fDoPart && fDoGenJ) ; | |
991 | if (fDoPart && fDoRecJ) ; | |
992 | if (fDoGenJ && fDoRecJ) ; | |
993 | } | |
994 | ||
995 | void AliJetAnalysis::PlotShapH() const | |
996 | { | |
997 | // missing | |
998 | if (fDoGenJ) ; | |
999 | if (fDoRecJ) ; | |
1000 | } | |
1001 | ||
1002 | void AliJetAnalysis::PlotFragH() const | |
1003 | { | |
1004 | // missing | |
1005 | if (fDoGenJ) ; | |
1006 | if (fDoRecJ) ; | |
1007 | } | |
1008 | ||
1009 | void AliJetAnalysis::PlotTrigH() | |
1010 | { | |
1011 | // missing | |
1012 | ||
1013 | } | |
1014 | ||
1015 | //////////////////////////////////////////////////////////////////////// | |
1016 | ||
1017 | void AliJetAnalysis::SaveHistograms() | |
1018 | { | |
1019 | // Save histograms | |
1020 | TFile *fOut = new TFile(fFile,"recreate"); | |
1021 | fOut->cd(); | |
1022 | if (fDoKine) SaveKineH(); | |
1023 | if (fDoCorr) SaveCorrH(); | |
1024 | if (fDoCorr50) SaveCorr50H(); | |
1025 | if (fDoShap) SaveShapH(); | |
1026 | if (fDoFrag) SaveFragH(); | |
1027 | if (fDoTrig) SaveTrigH(); | |
1028 | if (fDoJt) SaveJtH(); | |
1029 | if (fDodNdxi) SavedNdxiH(); | |
1030 | fOut->Close(); | |
1031 | } | |
1032 | ||
1033 | void AliJetAnalysis::SaveKineH() | |
1034 | { | |
1035 | // Save kinematic histograms | |
1036 | if (fDoPart) { | |
1037 | fPKineEneH->Write(); | |
1038 | fPKinePtH->Write(); | |
1039 | fPKinePhiH->Write(); | |
1040 | fPKineEtaH->Write(); | |
1041 | } | |
62567a80 | 1042 | |
83a444b1 | 1043 | if (fDoGenJ) { |
1044 | fGKineEneH->Write(); | |
1045 | fGKinePtH->Write(); | |
1046 | fGKinePhiH->Write(); | |
1047 | fGKineEtaH->Write(); | |
1048 | } | |
62567a80 | 1049 | |
83a444b1 | 1050 | if (fDoRecJ) { |
1051 | fRKineEneH->Write(); | |
1052 | fRKinePtH->Write(); | |
1053 | fRKinePhiH->Write(); | |
1054 | fRKineEtaH->Write(); | |
1055 | } | |
1056 | } | |
1057 | ||
1058 | void AliJetAnalysis::SaveCorrH() | |
1059 | { | |
1060 | // Save correlation histograms | |
1061 | if (fDoPart && fDoGenJ) { | |
1062 | fPGCorrEneH->Write(); | |
1063 | fPGCorrPtH->Write(); | |
1064 | fPGCorrEtaH->Write(); | |
1065 | fPGCorrPhiH->Write(); | |
1066 | } | |
62567a80 | 1067 | |
83a444b1 | 1068 | if (fDoPart && fDoRecJ) { |
1069 | fPRCorrEneH->Write(); | |
1070 | fPRCorrPtH->Write(); | |
1071 | fPRCorrEtaH->Write(); | |
1072 | fPRCorrPhiH->Write(); | |
1073 | } | |
62567a80 | 1074 | |
83a444b1 | 1075 | if (fDoGenJ && fDoRecJ) { |
1076 | fRGCorrEneH->Write(); | |
1077 | fRGCorrPtH->Write(); | |
1078 | fRGCorrEtaH->Write(); | |
1079 | fRGCorrPhiH->Write(); | |
1080 | } | |
1081 | } | |
1082 | ||
1083 | void AliJetAnalysis::SaveCorr50H() | |
1084 | { | |
1085 | // Save correlation histograms (special case) | |
1086 | if (fDoPart && fDoRecJ) { | |
1087 | fPRCorr50EneH->Write(); | |
1088 | fPRCorr50PtH->Write(); | |
1089 | fPRCorr50EtaH->Write(); | |
1090 | fPRCorr50PhiH->Write(); | |
1091 | } | |
1092 | if (fDoGenJ && fDoRecJ) { | |
1093 | fRGCorr50EneH->Write(); | |
1094 | fRGCorr50PtH->Write(); | |
1095 | fRGCorr50EtaH->Write(); | |
1096 | fRGCorr50PhiH->Write(); | |
1097 | } | |
1098 | } | |
1099 | ||
1100 | void AliJetAnalysis::SaveShapH() | |
1101 | { | |
1102 | // Save jet shape histograms | |
1103 | if (fDoRecJ) { | |
1104 | fRShapSelH->Write(); | |
1105 | fdEdrH->Write(); | |
1106 | if(fDoBkgd) fdEdrB->Write(); | |
1107 | fPtEneH2->Write(); | |
1108 | fdEdrW->Write(); | |
1109 | if (fPythia){ | |
1110 | fRShapRejH->Write(); | |
1111 | fRShapAllH->Write(); | |
1112 | } | |
1113 | } | |
1114 | } | |
1115 | ||
1116 | void AliJetAnalysis::SaveJtH() | |
1117 | { | |
1118 | // Save J_T histograms | |
1119 | if (fDoRecJ) { | |
1120 | fJtH->Write(); | |
1121 | if(fDoBkgd) fJtB->Write(); | |
1122 | fJtW->Write(); | |
1123 | } | |
1124 | } | |
1125 | ||
1126 | void AliJetAnalysis::SavedNdxiH() | |
1127 | { | |
1128 | // Save dN/d#xi histograms | |
1129 | if (fDoRecJ) { | |
1130 | fdNdxiH->Write(); | |
1131 | if(fDoBkgd) fdNdxiB->Write(); | |
1132 | fPtEneH->Write(); | |
1133 | fdNdxiW->Write(); | |
1134 | } | |
1135 | } | |
1136 | ||
1137 | void AliJetAnalysis::SaveFragH() | |
1138 | { | |
1139 | // Save fragmentation histograms | |
1140 | if (fDoRecJ) { | |
1141 | fRFragSelH->Write(); | |
1142 | if(fPythia){ | |
1143 | fRFragRejH->Write(); | |
1144 | fRFragAllH->Write(); | |
1145 | } | |
1146 | } | |
1147 | } | |
1148 | ||
1149 | void AliJetAnalysis::SaveTrigH() | |
1150 | { | |
1151 | // Save trigger bias histograms | |
1152 | if(fDoTrig){ | |
1153 | fGTriggerEneH->Write(); | |
1154 | fRTriggerEneH->Write(); | |
1155 | fGPTriggerEneH->Write(); | |
1156 | fPTriggerEneH->Write(); | |
1157 | } | |
1158 | } | |
1159 | ||
1160 | //////////////////////////////////////////////////////////////////////// | |
1161 | // main Analysis function | |
1162 | ||
1163 | void AliJetAnalysis::Analyze() | |
1164 | ||
1165 | { | |
1166 | // Kinematics for | |
1167 | // leading particle | |
1168 | // leading generated jet | |
1169 | // leading reconstructed jet | |
1170 | ||
1171 | // Correlations amd resolutions | |
1172 | // a) correlations in energy, pt, phi, eta | |
1173 | // b) resolutions in energy, pt, phi, eta, r | |
1174 | // leading particle and leading generated jet | |
1175 | // leading particle and leading reconstructed jet | |
1176 | // leading generated jet and leading reconstructed jet | |
1177 | ||
1178 | // Fragmentation functions and Shapes | |
1179 | // a) integrated over all pt | |
1180 | // b) in 3 flavors: | |
1181 | // b.1) only for user selected particles in jet | |
1182 | // b.2) only for user rejected particles in jet | |
1183 | // b.3) for all particles in jet | |
1184 | ||
1185 | DefineHistograms(); | |
1186 | FillHistograms(); | |
1187 | NormHistograms(); | |
1188 | PlotHistograms(); | |
1189 | SaveHistograms(); | |
1190 | } | |
62567a80 | 1191 |