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