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