]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EMCAL/AliEMCALJetFinderPlots.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
strdup replaced by malloc and strcpy
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALJetFinderPlots.cxx
CommitLineData
f7d5860b 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
ee6b678f 17/* $Id$ */
f7d5860b 18
185da5d3 19
f7d5860b 20//_________________________________________________________________________
21// Class for Filling JetFinder Plots
44f59d68 22// --
f7d5860b 23//*-- Author: Mark Horner (LBL/UCT)
44f59d68 24// --
25// --
f7d5860b 26
27
28#include "TMath.h"
29#include "AliEMCALJetFinderPlots.h"
30
31ClassImp(AliEMCALJetFinderPlots)
32
33AliEMCALJetFinderPlots::AliEMCALJetFinderPlots()
34{
44f59d68 35 // Constructor to initialise variables
f7d5860b 36 fInitialised = kFALSE;
37 fNominalEnergy = 0.0;
185da5d3 38 fConeRadius = 0.3;
39 fDebug = 0;
40 fOutput=0;
ab01dff2 41 fhFragmFcn=0;// = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
42 fhPartonFragmFcn=0;// = new TH1F("hPartonFragmFcn","Fragmentation Function",100,0,1);
43 fhPartonJT=0;// = new TH1F("hPartonJT","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
44 fhPartonPL=0;// = new TH1F("hPartonPL","Track Momentum Parallel to Parton Axis ",100,0.,100.);
45 fhJetJT=0;// = new TH1F("hJetJT","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
46 fhJetPL=0;// = new TH1F("hJetPL","Track Momentum Parallel to Jet Axis ",100,0.,100.);
47 fhJetEt=0;// = new TH1F("hJetEt","E_{T}^{reco}",250,0.,250.);
48 fhJetEta=0;// = new TH1F("hJetEta","#eta_{jet}^{reco}",180,-0.9,0.9);
49 fhJetPhi=0;// = new TH1F("hJetPhi","#phi_{jet}^{reco}",62,0.,3.1);
50 fhPartonEta=0;// = new TH1F("hPartonEta","#eta_{Parton}",180,-0.9,0.9);
51 fhPartonPhi=0;// = new TH1F("hPartonPhi","#phi_{Parton}",62,0.,3.1);
52 fhEtaDiff=0;// = new TH1F("hEtaDiff","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
53 fhPhiDiff=0;// = new TH1F("hPhiDiff","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
54 fhNJets=0;// = new TH1F("hNJets","N Reconstructed jets",11,-0.5,10.5);
55 fhEtaPhiSpread=0;
56
57fhFragmFcn2=0; // ("hFragmFcn2","Fragmentation Function",100,0,1);
58fhPartonFragmFcn2=0;// ("hFragmFcn2","Parton Fragmentation Function",100,0,1);
59fhPartonJT2=0; // ("hPartonJT2","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
60fhPartonPL2=0; // ("hPartonPL2","Track Momentum Parallel to Parton Axis ",100,0.,100.);
61fhJetJT2=0; // ("hJetJT2","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
62fhJetPL2=0; // ("hJetPL2","Track Momentum Parallel to Jet Axis ",100,0.,100.);
63fhJetEt2=0; // ("hJetEt2","E_{T}^{reco}",250,0.,250.);
64fhJetEta2=0; // ("hJetEta2","#eta_{jet}^{reco}",180,-0.9,0.9);
65fhJetPhi2=0; // ("hJetPhi2","#phi_{jet}^{reco}",62,0.,3.1);
66fhPartonEta2=0; // ("hPartonEta2","#eta_{Parton}",180,-0.9,0.9);
67fhPartonPhi2=0; // ("hPartonPhi2","#phi_{Parton}",62,0.,3.1);
68fhEtaDiff2=0; // ("hEtaDiff2","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
69fhPhiDiff2=0; // ("hPhiDiff2","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
70fhEtaPhiSpread2=0; // ("hEtaPhiSpread2","#eta - #phi Distribution
185da5d3 71 //of Reconstructed Jets",192,-0.7,0.7,288,pi/3,pi);
ab01dff2 72fhNJets2=0; // ("hNJets2","N Reconstructed jets",11,-0.5,10.5);
73fhJetEtSecond2=0; //("hJetEtSecond2","E_{T}^{reco}",250,0.,250.);
74fhJetEtRatio2=0; //("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
75fhEtaPhiDist2=0; //("hEtaPhiDist2","Angular Distance Between First and Second",100,0,3);
63131144 76fhInputOutput=0;
77// TH2F *fhInputOutput; //("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
78
79fhRecoBinPt=0; // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
80fhRecoBinPartonPt=0; // ("fhRecoBinPartonPt","Input Pt Distribution",100,0,1);
81fhRecoBinJetEt=0; // ("fhRecoJetEt","E_{T}^{reco}",250,0.,250.);
82fhRecoBinInputJetEt=0; // ("fhRecoInputJetEt","E_{T}^{reco}",250,0.,250.);
492f773f 83fhJetPT =0;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
84fhPartonPT =0;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
85fhJetPT2 =0;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
86fhPartonPT2 =0;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
87fhRecoBinFragmFcn =0;//new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",100,0,1);
88fhRecoBinPartonFragmFcn =0;// new TH1F("fhRecoBinPartonFragmFcn","Input Bin Fragm Fcn Distribution",100,0,1);
da2f0870 89fhJetInvE=0;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
90fhJetInvE2=0;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
1515a8af 91
92
93
f7d5860b 94}
95
96void AliEMCALJetFinderPlots::InitPlots()
97{
185da5d3 98//========================= CASE 1 =======================================
ab01dff2 99 fhFragmFcn = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
100 fhFragmFcn->Sumw2();
492f773f 101 fhJetPT = new TH1F("hJetPT","P_{T} Distribution",200,0,200);
102 fhJetPT->Sumw2();
103 fhPartonPT = new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
104 fhPartonPT->Sumw2();
ab01dff2 105 fhPartonFragmFcn = new TH1F("hPartonFragmFcn","Parton Fragmentation Function",100,0,1);
106 fhPartonFragmFcn->Sumw2();
107 fhPartonJT = new TH1F("hPartonJT","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
108 fhPartonJT->Sumw2();
109 fhPartonPL = new TH1F("hPartonPL","Track Momentum Parallel to Parton Axis ",100,0.,100.);
110 fhPartonPL->Sumw2();
111 fhJetJT = new TH1F("hJetJT","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
112 fhJetJT->Sumw2();
113 fhJetPL = new TH1F("hJetPL","Track Momentum Parallel to Jet Axis ",100,0.,100.);
114 fhJetPL->Sumw2();
115 fhJetEt = new TH1F("hJetEt","E_{T}^{reco}",250,0.,250.);
116 fhJetEt->Sumw2();
117 fhJetEta = new TH1F("hJetEta","#eta_{jet}^{reco}",180,-0.9,0.9);
118 fhJetEta->Sumw2();
119 fhJetPhi = new TH1F("hJetPhi","#phi_{jet}^{reco}",62,0.,3.1);
120 fhJetPhi->Sumw2();
121 fhPartonEta = new TH1F("hPartonEta","#eta_{Parton}",180,-0.9,0.9);
122 fhPartonEta->Sumw2();
123 fhPartonPhi = new TH1F("hPartonPhi","#phi_{Parton}",62,0.,3.1);
124 fhPartonPhi->Sumw2();
125 fhEtaDiff = new TH1F("hEtaDiff","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
126 fhEtaDiff->Sumw2();
127 fhPhiDiff = new TH1F("hPhiDiff","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
128 fhPhiDiff->Sumw2();
129 fhNJets = new TH1F("hNJets","N Reconstructed jets",11,-0.5,10.5);
130 fhNJets->Sumw2();
131 fhEtaPhiSpread = new TH2F("hEtaPhiSpread","#eta - #phi Distribution of Reconstructed Jets",100,-0.5,0.5,100,-0.5,0.5);
132 fhEtaPhiSpread->Sumw2();
133 fhNJets->SetXTitle("N_{jets}^{reco}/event");
134 fhNJets->SetYTitle("N_{events}");
f7d5860b 135
136 //Jet properties
ab01dff2 137 fhJetEt->SetFillColor(16);
138 fhJetEt->SetXTitle("E_{T}^{reco}");
f7d5860b 139
ab01dff2 140 fhJetEta->SetFillColor(16);
141 fhJetEta->SetXTitle("#eta_{jet}^{reco}");
f7d5860b 142
ab01dff2 143 fhJetPhi->SetFillColor(16);
144 fhJetPhi->SetXTitle("#phi_{jet}^{reco}");
f7d5860b 145
ab01dff2 146 fhPartonEta->SetFillColor(16);
147 fhPartonEta->SetXTitle("#eta_{parton}");
f7d5860b 148
ab01dff2 149 fhPartonPhi->SetFillColor(16);
150 fhPartonPhi->SetXTitle("#phi_{parton}");
f7d5860b 151
ab01dff2 152 fhPartonPL->SetXTitle("p (GeV/c)");
153 fhPartonJT->SetXTitle("p (GeV/c)");
f7d5860b 154
ab01dff2 155 fhPartonFragmFcn->SetXTitle("Z = p_{T}^{Chg}/E_{T}^{parton}");
f7d5860b 156
157 //Jet component properties
158
ab01dff2 159 fhJetPL->SetXTitle("p (GeV/c)");
160 fhJetJT->SetXTitle("p (GeV/c)");
161 fhFragmFcn->SetXTitle("Z = p_{T}^{Chg}/E_{T}^{reco}");
162 fhPartonFragmFcn->SetXTitle("Z = p_{T}^{Chg}/E_{T}^{reco}");
f7d5860b 163
ab01dff2 164 fhEtaDiff->SetXTitle("#eta_{jet}^{reco}-#eta_{jet}^{input}");
165 fhPhiDiff->SetXTitle("#phi_{jet}^{reco}-#phi_{jet}^{input}");
166 fhEtaPhiSpread->SetXTitle("#eta");
167 fhEtaPhiSpread->SetYTitle("#phi");
185da5d3 168
169//======================= CASE 2 ======================================
170
171
ab01dff2 172fhFragmFcn2 = new TH1F("hFragmFcn2","Fragmentation Function",100,0,1);
173fhFragmFcn2->Sumw2();
492f773f 174 fhJetPT2 = new TH1F("hJetPT2","P_{T} Distribution",200,0,200);
175 fhJetPT2->Sumw2();
176 fhPartonPT2 = new TH1F("hPartonPT2","Parton P_{T} Distribution",200,0,1);
177 fhPartonPT2->Sumw2();
ab01dff2 178fhPartonFragmFcn2 = new TH1F("hPartonFragmFcn2","Parton Fragmentation Function",100,0,1);
179fhPartonFragmFcn2->Sumw2();
180fhPartonJT2 = new TH1F("hPartonJT2","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
181fhPartonJT2->Sumw2();
182fhPartonPL2 = new TH1F("hPartonPL2","Track Momentum Parallel to Parton Axis ",100,0.,100.);
183fhPartonPL2->Sumw2();
184fhJetJT2 = new TH1F("hJetJT2","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
185fhJetJT2->Sumw2();
186fhJetPL2 = new TH1F("hJetPL2","Track Momentum Parallel to Jet Axis ",100,0.,100.);
187fhJetPL2->Sumw2();
188fhJetEt2 = new TH1F("hJetEt2","E_{T}^{reco}",250,0.,250.);
189fhJetEt2->Sumw2();
190fhJetEta2 = new TH1F("hJetEta2","#eta_{jet}^{reco}",180,-0.9,0.9);
191fhJetEta2->Sumw2();
192fhJetPhi2 = new TH1F("hJetPhi2","#phi_{jet}^{reco}",62,0.,3.1);
193fhJetPhi2->Sumw2();
194fhPartonEta2 = new TH1F("hPartonEta2","#eta_{Parton}",180,-0.9,0.9);
195fhPartonEta2->Sumw2();
196fhPartonPhi2 = new TH1F("hPartonPhi2","#phi_{Parton}",62,0.,3.1);
197fhPartonPhi2->Sumw2();
198fhEtaDiff2 = new TH1F("hEtaDiff2","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
199fhEtaDiff2->Sumw2();
200fhPhiDiff2 = new TH1F("hPhiDiff2","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
201fhPhiDiff2->Sumw2();
202fhEtaPhiSpread2 = new TH2F("hEtaPhiSpread2","#eta - #phi Distribution of Reconstructed Jets",100,-0.5,0.5,100,-0.5,0.5);
203fhEtaPhiSpread2->Sumw2();
204fhNJets2 = new TH1F("hNJets2","N Reconstructed jets",11,-0.5,10.5);
205fhNJets2->Sumw2();
206fhJetEtSecond2 = new TH1F("hJetEtSecond2","E_{T}^{reco}",250,0.,250.);
207fhJetEtSecond2->Sumw2();
208fhJetEtRatio2 = new TH1F("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
209fhJetEtRatio2->Sumw2();
210fhEtaPhiDist2 = new TH1F("hEtaPhiDist2","Angular Distance Between First and Second",100,0,3);
211fhEtaPhiDist2->Sumw2();
63131144 212
213fhInputOutput= new TH2F("hInputOutput","Input and Reconstruction Correlations;Input;Output",200,0,200,200,0,200); //("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
214
215//============================== Reconstruction Bin Comparison ============================================
216
492f773f 217fhRecoBinPt =new TH1F("fhRecoBinPt","Reconstructed Pt Distribution",200,0,200);
63131144 218fhRecoBinPt->Sumw2();
492f773f 219fhRecoBinPartonPt = new TH1F("fhRecoBinPartonPt","Input Pt Distribution",200,0,200);
63131144 220fhRecoBinPartonPt->Sumw2();
492f773f 221fhRecoBinFragmFcn =new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",100,0,1);
222fhRecoBinFragmFcn->Sumw2();
223fhRecoBinPartonFragmFcn = new TH1F("fhRecoBinPartonFragmFcn","Input Bin Fragm Fcn Distribution",100,0,1);
224fhRecoBinPartonFragmFcn->Sumw2();
63131144 225fhRecoBinJetEt = new TH1F("fhRecoJetEt","E_{T}^{reco}",250,0.,250.);
226fhRecoBinJetEt->Sumw2();
227fhRecoBinInputJetEt = new TH1F("fhRecoInputJetEt","E_{T}^{reco}",250,0.,250.);
228fhRecoBinInputJetEt->Sumw2();
da2f0870 229
230
231fhJetInvE = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
232fhJetInvE->Sumw2();
233fhJetInvE2 = new TH1F("fhJetInvE2","#frac{1}{E_{R}}",100,0,1);
234fhJetInvE2->Sumw2();
da2f0870 235
236
237
f7d5860b 238 fInitialised = kTRUE;
239
240}
241
44f59d68 242AliEMCALJetFinderPlots::~AliEMCALJetFinderPlots()
243{
244 // To ensure that all requested memory is returned
ab01dff2 245delete fhFragmFcn;// = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
246delete fhPartonFragmFcn;// = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
247delete fhPartonJT;// = new TH1F("hPartonJT","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
248delete fhPartonPL;// = new TH1F("hPartonPL","Track Momentum Parallel to Parton Axis ",100,0.,100.);
249delete fhJetJT;// = new TH1F("hJetJT","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
250delete fhJetPL;// = new TH1F("hJetPL","Track Momentum Parallel to Jet Axis ",100,0.,100.);
251delete fhJetEt;// = new TH1F("hJetEt","E_{T}^{reco}",250,0.,250.);
252delete fhJetEta;// = new TH1F("hJetEta","#eta_{jet}^{reco}",180,-0.9,0.9);
253delete fhJetPhi;// = new TH1F("hJetPhi","#phi_{jet}^{reco}",62,0.,3.1);
254delete fhPartonEta;// = new TH1F("hPartonEta","#eta_{Parton}",180,-0.9,0.9);
255delete fhPartonPhi;// = new TH1F("hPartonPhi","#phi_{Parton}",62,0.,3.1);
256delete fhEtaDiff;// = new TH1F("hEtaDiff","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
257delete fhPhiDiff;// = new TH1F("hPhiDiff","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
258delete fhNJets;// = new TH1F("hNJets","N Reconstructed jets",11,-0.5,10.5);
259delete fhEtaPhiSpread;
260
261 delete fhFragmFcn2; // ("hFragmFcn2","Fragmentation Function",100,0,1);
262 delete fhPartonFragmFcn2;// ("hFragmFcn2","Parton Fragmentation Function",100,0,1);
263 delete fhPartonJT2; // ("hPartonJT2","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
264 delete fhPartonPL2; // ("hPartonPL2","Track Momentum Parallel to Parton Axis ",100,0.,100.);
265 delete fhJetJT2; // ("hJetJT2","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
266 delete fhJetPL2; // ("hJetPL2","Track Momentum Parallel to Jet Axis ",100,0.,100.);
267 delete fhJetEt2; // ("hJetEt2","E_{T}^{reco}",250,0.,250.);
268 delete fhJetEta2; // ("hJetEta2","#eta_{jet}^{reco}",180,-0.9,0.9);
269 delete fhJetPhi2; // ("hJetPhi2","#phi_{jet}^{reco}",62,0.,3.1);
270 delete fhPartonEta2; // ("hPartonEta2","#eta_{Parton}",180,-0.9,0.9);
271 delete fhPartonPhi2; // ("hPartonPhi2","#phi_{Parton}",62,0.,3.1);
272 delete fhEtaDiff2; // ("hEtaDiff2","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
273 delete fhPhiDiff2; // ("hPhiDiff2","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
274 delete fhEtaPhiSpread2; // ("hEtaPhiSpread2","#eta - #phi Distribution
185da5d3 275 //of Reconstructed Jets",192,-0.7,0.7,288,pi/3,pi);
ab01dff2 276 delete fhNJets2; // ("hNJets2","N Reconstructed jets",11,-0.5,10.5);
277 delete fhJetEtSecond2; //("hJetEtSecond2","E_{T}^{reco}",250,0.,250.);
278 delete fhJetEtRatio2; //("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
279 delete fhEtaPhiDist2; //("hEtaPhiDist2","Angular Distance Between First and Second",100,0,3);
185da5d3 280
63131144 281 delete fhRecoBinPt; // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
282 delete fhRecoBinPartonPt; // ("fhRecoBinPartonPt","Input Pt Distribution",100,0,1);
283 delete fhRecoBinJetEt; // ("fhRecoJetEt","E_{T}^{reco}",250,0.,250.);
284 delete fhRecoBinInputJetEt; // ("fhRecoInputJetEt","E_{T}^{reco}",250,0.,250.);
285
492f773f 286 delete fhJetPT ;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
287 delete fhPartonPT ;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
288 delete fhJetPT2 ;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
289 delete fhPartonPT2 ;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
290 delete fhRecoBinFragmFcn;//new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",100,0,1);
291 delete fhRecoBinPartonFragmFcn;// new TH1F("fhRecoBinPartonFragmFcn","Input Bin Fragm Fcn Distribution",100,0,1);
1515a8af 292 delete fhJetInvE;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
293 delete fhJetInvE2;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
185da5d3 294
f7d5860b 295}
296
297void AliEMCALJetFinderPlots::FillFromOutput(AliEMCALJetFinderOutput* output)
298{
44f59d68 299 // Fill histograms from an output object
f7d5860b 300if (!fInitialised) InitPlots();
301 fOutput = output;
185da5d3 302 if (!fOutput) return;
ab01dff2 303 fhNJets->Fill(fOutput->GetNJets());
63131144 304 Bool_t doesJetMeetBinCriteria = 0;
305 AliEMCALJet* jethighest=0;
306 AliEMCALJet* jetsecond=0;
4c162f44 307 // Find Highest and Second Highest Jet
308 // (NB!!!!!!!) Pointing into the EMCAL!!!!!!
63131144 309
310// =========================== All cases ===================================
1050d19e 311
312
4c162f44 313 // I will make a little array of jet indices for which jets are in
314 // the EMCAL then my counter can loop from 0 below - but it will
315 // be the index of the array of applicable jets
1050d19e 316
317 Int_t appjet[4];
318 Int_t numappjet=0;
319
320 for (Int_t appc=0;appc<fOutput->GetNJets();appc++)
321 { // Check all jets for applicability
322 Float_t eta = fOutput->GetJet(appc)->Eta();
323 Float_t phi = fOutput->GetJet(appc)->Phi();
324 if (eta > -0.7 && eta < 0.7 && phi > 1./3.*TMath::Pi() && phi < TMath::Pi())
325 { // Then jet is applicable
326 appjet[numappjet]=appc;
327 numappjet++;
328 }
329 }
4c162f44 330
1050d19e 331
332if (numappjet >=1)
333{
334 Float_t theta = 2.0*atan(exp(-fOutput->GetParton(0)->Eta()));
335 Float_t et = fOutput->GetParton(0)->Energy() * TMath::Sin(theta);
63131144 336 if (fOutput->GetNJets()>1)
337 {
4c162f44 338 for (Int_t counter = 0; counter<numappjet;counter++)
63131144 339 {
340 if (counter==0)
341 {
4c162f44 342 jethighest = fOutput->GetJet(appjet[0]);
343 jetsecond = fOutput->GetJet(appjet[1]);
63131144 344 }
345 if (counter>0)
346 {
347 Float_t energyhighest = jethighest->Energy();
348 Float_t energysecond = jetsecond->Energy();
349
4c162f44 350 if ((fOutput->GetJet(appjet[counter]))->Energy()>energyhighest)
63131144 351 {
352 jetsecond=jethighest;
4c162f44 353 jethighest=fOutput->GetJet(appjet[counter]);
354 }else if ((fOutput->GetJet(appjet[counter]))->Energy()>energysecond)
63131144 355 {
4c162f44 356 jetsecond=fOutput->GetJet(appjet[counter]);
63131144 357 }
358
359 }
360 }
361 }else
362 {
4c162f44 363 Float_t eta = fOutput->GetJet(0)->Eta();
364 Float_t phi = fOutput->GetJet(0)->Phi();
365 if (eta > -0.7 && eta < 0.7 && phi > 1./3.*TMath::Pi() && phi < TMath::Pi())
366 { // Then jet is applicable
367 jethighest=fOutput->GetJet(0);
368 jetsecond=0;
369 }else
370 {
371 Error("FillFromOutput","There is only one jet and it isn't in the area of applicability");
372 }
373
63131144 374 }
375 if ( 95.0 < jethighest->Energy() && jethighest->Energy() < 105.0 )
376 {
377 doesJetMeetBinCriteria = 1;
378 fhRecoBinJetEt->Fill(jethighest->Energy());
379 fhRecoBinInputJetEt->Fill(et);
380 }
4c162f44 381 fhInputOutput->Fill(et,jethighest->Energy());
63131144 382
383}
384
1050d19e 385if (numappjet > 1)
185da5d3 386{
387//========================= CASE 2 ===========================
388 Int_t nPartons = fOutput->GetNPartons();
ab01dff2 389 fhNJets2->Fill(fOutput->GetNJets());
185da5d3 390 AliEMCALParton* parton;
185da5d3 391
392 // End finding highest and second highest and continue
ab01dff2 393 fhJetEt2->Fill(jethighest->Energy());
da2f0870 394 fhJetInvE2->Fill(1.0/jethighest->Energy());
ab01dff2 395 fhJetEta2->Fill(jethighest->Eta() );
396 fhJetPhi2->Fill(jethighest->Phi() );
185da5d3 397 if (nPartons ==0) return;
398 parton = fOutput->GetParton(0);
399
ab01dff2 400 fhPartonEta2->Fill( parton->Eta() );
401 fhPartonPhi2->Fill( parton->Phi() );
185da5d3 402
403 //hJetEtDiff->Fill( jet->Energy() - parton->Energy() );
ab01dff2 404 fhEtaDiff2->Fill( jethighest->Eta() - parton->Eta() );
405 fhPhiDiff2->Fill( jethighest->Phi() - parton->Phi() );
406 fhEtaPhiSpread2->Fill(jethighest->Eta()-parton->Eta(),jethighest->Phi() - parton->Phi());
407 fhJetEtSecond2->Fill(jetsecond->Energy());
408 fhJetEtRatio2->Fill(jetsecond->Energy()/jethighest->Energy());
409 fhEtaPhiDist2->Fill( TMath::Sqrt((jethighest->Eta() - jetsecond->Eta())*(jethighest->Eta() - jetsecond->Eta())
185da5d3 410 + (jethighest->Phi() - jetsecond->Phi())*(jethighest->Phi() - jetsecond->Phi()) ));
411 /*
412 Float_t *pt,*phi,*eta;
413 Int_t *pdg;
414 pt = new Float_t[parton->GetNTracks()];
415 eta = new Float_t[parton->GetNTracks()];
416 phi = new Float_t[parton->GetNTracks()];
417 pdg = new Int_t[parton->GetNTracks()];*/
418
419
420
421 Float_t pt[2000];
422 Float_t eta[2000];
423 Float_t phi[2000];
424 Int_t pdg[2000];
425
426 parton->GetTrackList(pt,eta,phi,pdg);
427 for(Int_t iT=0; iT< parton->GetNTracks() ; iT++ )
428 {
429 if ( (eta[iT]-parton->Eta())*(eta[iT]-parton->Eta())+
430 (phi[iT]-parton->Phi())*(phi[iT]-parton->Phi()) >fConeRadius * fConeRadius ) continue;
431 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
432 Double_t rt = 2.0*atan(exp(-parton->Eta()));
433 Double_t ctt = cos(tt);
434 Double_t crt = cos(rt);
435 Double_t stt = sin(tt);
436 Double_t srt = sin(rt);
437 Double_t ctp = cos(phi[iT]);
438 Double_t crp = cos(parton->Phi());
439 Double_t stp = sin(phi[iT]);
440 Double_t srp = sin(parton->Phi());
63131144 441 //Double_t alpha = acos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
442 Double_t alpha;
443 if (TMath::Abs(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt) > 0.9990)
444 {
445 alpha = 0.0;
446 }else
447 {
448 alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
449 }
185da5d3 450 Double_t correctp = pt[iT]/stt;
ab01dff2 451 fhPartonPL2->Fill( correctp*cos(alpha));
185da5d3 452 if ( (parton->Eta()-eta[iT])*(parton->Eta()-eta[iT]) +
453 (parton->Phi()-phi[iT])*(parton->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 454 fhPartonJT2->Fill( correctp*sin(alpha));
1515a8af 455 fhPartonPT2->Fill(correctp*sin(tt));
185da5d3 456 if (fNominalEnergy == 0.0) {
ab01dff2 457 fhPartonFragmFcn2->Fill( correctp*sin(tt)/parton->Energy() );
185da5d3 458 }else
459 {
ab01dff2 460 fhPartonFragmFcn2->Fill(correctp*sin(tt)/fNominalEnergy);
63131144 461 }
462 if (doesJetMeetBinCriteria)
463 {
464 fhRecoBinPartonPt->Fill(correctp*sin(tt)); // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
185da5d3 465 }
466 }// loop over tracks
467
468/*
469 pt = new Float_t[jet->NTracks()];
470 eta = new Float_t[jet->NTracks()];
471 phi = new Float_t[jet->NTracks()];
472 pdg = new Int_t[jet->NTracks()];*/
473 jethighest->TrackList(pt,eta,phi,pdg);
474 for(Int_t iT=0; iT< jethighest->NTracks() ; iT++ )
475 {
476 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
477 Double_t rt = 2.0*atan(exp(-jethighest->Eta()));
478 Double_t ctt = cos(tt);
479 Double_t crt = cos(rt);
480 Double_t stt = sin(tt);
481 Double_t srt = sin(rt);
482 Double_t ctp = cos(phi[iT]);
483 Double_t crp = cos(jethighest->Phi());
484 Double_t stp = sin(phi[iT]);
485 Double_t srp = sin(jethighest->Phi());
63131144 486 // Double_t alpha = acos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
487 Double_t alpha;
488 if (TMath::Abs(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt) > 0.9990)
489 {
490 alpha = 0.0;
491 }else
492 {
493 alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
494 }
185da5d3 495 Double_t correctp = pt[iT]/stt;
ab01dff2 496 fhJetPL2->Fill( correctp*cos(alpha));
185da5d3 497 if ( (jethighest->Eta()-eta[iT])*(jethighest->Eta()-eta[iT]) +
498 (jethighest->Phi()-phi[iT])*(jethighest->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 499 fhJetJT2->Fill( correctp*sin(alpha));
1515a8af 500 fhJetPT2->Fill(correctp*sin(tt));
185da5d3 501 if (fNominalEnergy==0.0){
da2f0870 502 fhFragmFcn2->Fill( correctp*sin(tt)/jethighest->Energy() );
185da5d3 503 } else
504 {
ab01dff2 505 fhFragmFcn2->Fill( correctp*sin(tt)/fNominalEnergy );
185da5d3 506 }
63131144 507 if (doesJetMeetBinCriteria)
508 {
509 fhRecoBinPt->Fill(correctp*sin(tt)); // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
510 }
185da5d3 511 }// loop over tracks
512 }
513
1050d19e 514 if (numappjet == 1)
185da5d3 515 {
516
517//========================= CASE 1 ===========================
518 Int_t nPartons = fOutput->GetNPartons();
f7d5860b 519 if (fOutput->GetNJets()!=1) return;
520 AliEMCALParton* parton;
521 AliEMCALJet* jet;
1050d19e 522 jet = jethighest;//fOutput->GetJet(0);
ab01dff2 523 fhJetEt->Fill(jet->Energy());
da2f0870 524 fhJetInvE->Fill(1.0/jethighest->Energy());
ab01dff2 525 fhJetEta->Fill(jet->Eta() );
526 fhJetPhi->Fill(jet->Phi() );
f7d5860b 527 if (nPartons ==0) return;
528 parton = fOutput->GetParton(0);
185da5d3 529
ab01dff2 530 fhPartonEta->Fill( parton->Eta() );
531 fhPartonPhi->Fill( parton->Phi() );
f7d5860b 532
533 //hJetEtDiff->Fill( jet->Energy() - parton->Energy() );
ab01dff2 534 fhEtaDiff->Fill( jet->Eta() - parton->Eta() );
535 fhPhiDiff->Fill( jet->Phi() - parton->Phi() );
536 fhEtaPhiSpread->Fill(jet->Eta()-parton->Eta(),jet->Phi() - parton->Phi());
f7d5860b 537 /*
538 Float_t *pt,*phi,*eta;
539 Int_t *pdg;
540 pt = new Float_t[parton->GetNTracks()];
541 eta = new Float_t[parton->GetNTracks()];
542 phi = new Float_t[parton->GetNTracks()];
543 pdg = new Int_t[parton->GetNTracks()];*/
544
545
546
547 Float_t pt[2000];
548 Float_t eta[2000];
549 Float_t phi[2000];
550 Int_t pdg[2000];
551
552 parton->GetTrackList(pt,eta,phi,pdg);
553 for(Int_t iT=0; iT< parton->GetNTracks() ; iT++ )
554 {
555 if ( (eta[iT]-parton->Eta())*(eta[iT]-parton->Eta())+
556 (phi[iT]-parton->Phi())*(phi[iT]-parton->Phi()) >fConeRadius * fConeRadius ) continue;
557 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
558 Double_t rt = 2.0*atan(exp(-parton->Eta()));
559 Double_t ctt = cos(tt);
560 Double_t crt = cos(rt);
561 Double_t stt = sin(tt);
562 Double_t srt = sin(rt);
563 Double_t ctp = cos(phi[iT]);
564 Double_t crp = cos(parton->Phi());
565 Double_t stp = sin(phi[iT]);
566 Double_t srp = sin(parton->Phi());
63131144 567 // Double_t alpha = acos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
568 Double_t alpha;
569 if (TMath::Abs(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt) > 0.9990)
570 {
571 alpha = 0.0;
572 }else
573 {
574 alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt); }
f7d5860b 575 Double_t correctp = pt[iT]/stt;
ab01dff2 576 fhPartonPL->Fill( correctp*cos(alpha));
f7d5860b 577 if ( (parton->Eta()-eta[iT])*(parton->Eta()-eta[iT]) +
578 (parton->Phi()-phi[iT])*(parton->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 579 fhPartonJT->Fill( correctp*sin(alpha));
f7d5860b 580 if (fNominalEnergy == 0.0) {
ab01dff2 581 fhPartonFragmFcn->Fill( correctp*sin(tt)/parton->Energy() );
1515a8af 582 fhPartonPT->Fill(correctp*sin(tt));
f7d5860b 583 }else
584 {
ab01dff2 585 fhPartonFragmFcn->Fill(correctp*sin(tt)/fNominalEnergy);
f7d5860b 586 }
63131144 587 if (doesJetMeetBinCriteria)
588 {
589 fhRecoBinPartonPt->Fill(correctp*sin(tt)); // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
590 }
f7d5860b 591 }// loop over tracks
592
593/*
594 pt = new Float_t[jet->NTracks()];
595 eta = new Float_t[jet->NTracks()];
596 phi = new Float_t[jet->NTracks()];
597 pdg = new Int_t[jet->NTracks()];*/
598 jet->TrackList(pt,eta,phi,pdg);
599 for(Int_t iT=0; iT< jet->NTracks() ; iT++ )
600 {
601 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
602 Double_t rt = 2.0*atan(exp(-jet->Eta()));
603 Double_t ctt = cos(tt);
604 Double_t crt = cos(rt);
605 Double_t stt = sin(tt);
606 Double_t srt = sin(rt);
607 Double_t ctp = cos(phi[iT]);
608 Double_t crp = cos(jet->Phi());
609 Double_t stp = sin(phi[iT]);
610 Double_t srp = sin(jet->Phi());
63131144 611 //Info("plots","acos(%1.16f)\nstt=%f\npt=%f",crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt,stt,pt[iT]);
612 //Info("plots","diff to 1 %f",1.0-crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
613 Double_t alpha;
614 if (TMath::Abs(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt) > 0.9990)
615 {
616 alpha = 0.0;
617 }else
618 {
619 alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
620 }
f7d5860b 621 Double_t correctp = pt[iT]/stt;
ab01dff2 622 fhJetPL->Fill( correctp*cos(alpha));
f7d5860b 623 if ( (jet->Eta()-eta[iT])*(jet->Eta()-eta[iT]) +
624 (jet->Phi()-phi[iT])*(jet->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 625 fhJetJT->Fill( correctp*sin(alpha));
1515a8af 626 fhJetPT->Fill(correctp*sin(tt));
f7d5860b 627 if (fNominalEnergy==0.0){
da2f0870 628 fhFragmFcn->Fill( correctp*sin(tt)/jet->Energy() ); // This is the jet fragmentation function
f7d5860b 629 } else
630 {
ab01dff2 631 fhFragmFcn->Fill( correctp*sin(tt)/fNominalEnergy );
f7d5860b 632 }
63131144 633 if (doesJetMeetBinCriteria)
634 {
635 fhRecoBinPt->Fill(correctp*sin(tt)); // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
636 }
f7d5860b 637 }// loop over tracks
185da5d3 638 }
f7d5860b 639}
640
641