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