Made more robust
[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);
185da5d3 76
f7d5860b 77}
78
79void AliEMCALJetFinderPlots::InitPlots()
80{
185da5d3 81//========================= CASE 1 =======================================
ab01dff2 82 fhFragmFcn = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
83 fhFragmFcn->Sumw2();
84 fhPartonFragmFcn = new TH1F("hPartonFragmFcn","Parton Fragmentation Function",100,0,1);
85 fhPartonFragmFcn->Sumw2();
86 fhPartonJT = new TH1F("hPartonJT","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
87 fhPartonJT->Sumw2();
88 fhPartonPL = new TH1F("hPartonPL","Track Momentum Parallel to Parton Axis ",100,0.,100.);
89 fhPartonPL->Sumw2();
90 fhJetJT = new TH1F("hJetJT","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
91 fhJetJT->Sumw2();
92 fhJetPL = new TH1F("hJetPL","Track Momentum Parallel to Jet Axis ",100,0.,100.);
93 fhJetPL->Sumw2();
94 fhJetEt = new TH1F("hJetEt","E_{T}^{reco}",250,0.,250.);
95 fhJetEt->Sumw2();
96 fhJetEta = new TH1F("hJetEta","#eta_{jet}^{reco}",180,-0.9,0.9);
97 fhJetEta->Sumw2();
98 fhJetPhi = new TH1F("hJetPhi","#phi_{jet}^{reco}",62,0.,3.1);
99 fhJetPhi->Sumw2();
100 fhPartonEta = new TH1F("hPartonEta","#eta_{Parton}",180,-0.9,0.9);
101 fhPartonEta->Sumw2();
102 fhPartonPhi = new TH1F("hPartonPhi","#phi_{Parton}",62,0.,3.1);
103 fhPartonPhi->Sumw2();
104 fhEtaDiff = new TH1F("hEtaDiff","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
105 fhEtaDiff->Sumw2();
106 fhPhiDiff = new TH1F("hPhiDiff","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
107 fhPhiDiff->Sumw2();
108 fhNJets = new TH1F("hNJets","N Reconstructed jets",11,-0.5,10.5);
109 fhNJets->Sumw2();
110 fhEtaPhiSpread = new TH2F("hEtaPhiSpread","#eta - #phi Distribution of Reconstructed Jets",100,-0.5,0.5,100,-0.5,0.5);
111 fhEtaPhiSpread->Sumw2();
112 fhNJets->SetXTitle("N_{jets}^{reco}/event");
113 fhNJets->SetYTitle("N_{events}");
f7d5860b 114
115 //Jet properties
ab01dff2 116 fhJetEt->SetFillColor(16);
117 fhJetEt->SetXTitle("E_{T}^{reco}");
f7d5860b 118
ab01dff2 119 fhJetEta->SetFillColor(16);
120 fhJetEta->SetXTitle("#eta_{jet}^{reco}");
f7d5860b 121
ab01dff2 122 fhJetPhi->SetFillColor(16);
123 fhJetPhi->SetXTitle("#phi_{jet}^{reco}");
f7d5860b 124
ab01dff2 125 fhPartonEta->SetFillColor(16);
126 fhPartonEta->SetXTitle("#eta_{parton}");
f7d5860b 127
ab01dff2 128 fhPartonPhi->SetFillColor(16);
129 fhPartonPhi->SetXTitle("#phi_{parton}");
f7d5860b 130
ab01dff2 131 fhPartonPL->SetXTitle("p (GeV/c)");
132 fhPartonJT->SetXTitle("p (GeV/c)");
f7d5860b 133
ab01dff2 134 fhPartonFragmFcn->SetXTitle("Z = p_{T}^{Chg}/E_{T}^{parton}");
f7d5860b 135
136 //Jet component properties
137
ab01dff2 138 fhJetPL->SetXTitle("p (GeV/c)");
139 fhJetJT->SetXTitle("p (GeV/c)");
140 fhFragmFcn->SetXTitle("Z = p_{T}^{Chg}/E_{T}^{reco}");
141 fhPartonFragmFcn->SetXTitle("Z = p_{T}^{Chg}/E_{T}^{reco}");
f7d5860b 142
ab01dff2 143 fhEtaDiff->SetXTitle("#eta_{jet}^{reco}-#eta_{jet}^{input}");
144 fhPhiDiff->SetXTitle("#phi_{jet}^{reco}-#phi_{jet}^{input}");
145 fhEtaPhiSpread->SetXTitle("#eta");
146 fhEtaPhiSpread->SetYTitle("#phi");
185da5d3 147
148//======================= CASE 2 ======================================
149
150
ab01dff2 151fhFragmFcn2 = new TH1F("hFragmFcn2","Fragmentation Function",100,0,1);
152fhFragmFcn2->Sumw2();
153fhPartonFragmFcn2 = new TH1F("hPartonFragmFcn2","Parton Fragmentation Function",100,0,1);
154fhPartonFragmFcn2->Sumw2();
155fhPartonJT2 = new TH1F("hPartonJT2","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
156fhPartonJT2->Sumw2();
157fhPartonPL2 = new TH1F("hPartonPL2","Track Momentum Parallel to Parton Axis ",100,0.,100.);
158fhPartonPL2->Sumw2();
159fhJetJT2 = new TH1F("hJetJT2","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
160fhJetJT2->Sumw2();
161fhJetPL2 = new TH1F("hJetPL2","Track Momentum Parallel to Jet Axis ",100,0.,100.);
162fhJetPL2->Sumw2();
163fhJetEt2 = new TH1F("hJetEt2","E_{T}^{reco}",250,0.,250.);
164fhJetEt2->Sumw2();
165fhJetEta2 = new TH1F("hJetEta2","#eta_{jet}^{reco}",180,-0.9,0.9);
166fhJetEta2->Sumw2();
167fhJetPhi2 = new TH1F("hJetPhi2","#phi_{jet}^{reco}",62,0.,3.1);
168fhJetPhi2->Sumw2();
169fhPartonEta2 = new TH1F("hPartonEta2","#eta_{Parton}",180,-0.9,0.9);
170fhPartonEta2->Sumw2();
171fhPartonPhi2 = new TH1F("hPartonPhi2","#phi_{Parton}",62,0.,3.1);
172fhPartonPhi2->Sumw2();
173fhEtaDiff2 = new TH1F("hEtaDiff2","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
174fhEtaDiff2->Sumw2();
175fhPhiDiff2 = new TH1F("hPhiDiff2","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
176fhPhiDiff2->Sumw2();
177fhEtaPhiSpread2 = new TH2F("hEtaPhiSpread2","#eta - #phi Distribution of Reconstructed Jets",100,-0.5,0.5,100,-0.5,0.5);
178fhEtaPhiSpread2->Sumw2();
179fhNJets2 = new TH1F("hNJets2","N Reconstructed jets",11,-0.5,10.5);
180fhNJets2->Sumw2();
181fhJetEtSecond2 = new TH1F("hJetEtSecond2","E_{T}^{reco}",250,0.,250.);
182fhJetEtSecond2->Sumw2();
183fhJetEtRatio2 = new TH1F("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
184fhJetEtRatio2->Sumw2();
185fhEtaPhiDist2 = new TH1F("hEtaPhiDist2","Angular Distance Between First and Second",100,0,3);
186fhEtaPhiDist2->Sumw2();
185da5d3 187
f7d5860b 188 fInitialised = kTRUE;
189
190}
191
44f59d68 192AliEMCALJetFinderPlots::~AliEMCALJetFinderPlots()
193{
194 // To ensure that all requested memory is returned
ab01dff2 195delete fhFragmFcn;// = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
196delete fhPartonFragmFcn;// = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
197delete fhPartonJT;// = new TH1F("hPartonJT","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
198delete fhPartonPL;// = new TH1F("hPartonPL","Track Momentum Parallel to Parton Axis ",100,0.,100.);
199delete fhJetJT;// = new TH1F("hJetJT","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
200delete fhJetPL;// = new TH1F("hJetPL","Track Momentum Parallel to Jet Axis ",100,0.,100.);
201delete fhJetEt;// = new TH1F("hJetEt","E_{T}^{reco}",250,0.,250.);
202delete fhJetEta;// = new TH1F("hJetEta","#eta_{jet}^{reco}",180,-0.9,0.9);
203delete fhJetPhi;// = new TH1F("hJetPhi","#phi_{jet}^{reco}",62,0.,3.1);
204delete fhPartonEta;// = new TH1F("hPartonEta","#eta_{Parton}",180,-0.9,0.9);
205delete fhPartonPhi;// = new TH1F("hPartonPhi","#phi_{Parton}",62,0.,3.1);
206delete fhEtaDiff;// = new TH1F("hEtaDiff","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
207delete fhPhiDiff;// = new TH1F("hPhiDiff","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
208delete fhNJets;// = new TH1F("hNJets","N Reconstructed jets",11,-0.5,10.5);
209delete fhEtaPhiSpread;
210
211 delete fhFragmFcn2; // ("hFragmFcn2","Fragmentation Function",100,0,1);
212 delete fhPartonFragmFcn2;// ("hFragmFcn2","Parton Fragmentation Function",100,0,1);
213 delete fhPartonJT2; // ("hPartonJT2","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
214 delete fhPartonPL2; // ("hPartonPL2","Track Momentum Parallel to Parton Axis ",100,0.,100.);
215 delete fhJetJT2; // ("hJetJT2","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
216 delete fhJetPL2; // ("hJetPL2","Track Momentum Parallel to Jet Axis ",100,0.,100.);
217 delete fhJetEt2; // ("hJetEt2","E_{T}^{reco}",250,0.,250.);
218 delete fhJetEta2; // ("hJetEta2","#eta_{jet}^{reco}",180,-0.9,0.9);
219 delete fhJetPhi2; // ("hJetPhi2","#phi_{jet}^{reco}",62,0.,3.1);
220 delete fhPartonEta2; // ("hPartonEta2","#eta_{Parton}",180,-0.9,0.9);
221 delete fhPartonPhi2; // ("hPartonPhi2","#phi_{Parton}",62,0.,3.1);
222 delete fhEtaDiff2; // ("hEtaDiff2","#eta_{jet}^{reco}-#eta_{jet}^{input}",100,-0.5,0.5);
223 delete fhPhiDiff2; // ("hPhiDiff2","#phi_{jet}^{reco}-#phi_{jet}^{input}",100,-0.5,0.5);
224 delete fhEtaPhiSpread2; // ("hEtaPhiSpread2","#eta - #phi Distribution
185da5d3 225 //of Reconstructed Jets",192,-0.7,0.7,288,pi/3,pi);
ab01dff2 226 delete fhNJets2; // ("hNJets2","N Reconstructed jets",11,-0.5,10.5);
227 delete fhJetEtSecond2; //("hJetEtSecond2","E_{T}^{reco}",250,0.,250.);
228 delete fhJetEtRatio2; //("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
229 delete fhEtaPhiDist2; //("hEtaPhiDist2","Angular Distance Between First and Second",100,0,3);
185da5d3 230
231
232
f7d5860b 233}
234
235void AliEMCALJetFinderPlots::FillFromOutput(AliEMCALJetFinderOutput* output)
236{
44f59d68 237 // Fill histograms from an output object
f7d5860b 238if (!fInitialised) InitPlots();
239 fOutput = output;
185da5d3 240 if (!fOutput) return;
ab01dff2 241 fhNJets->Fill(fOutput->GetNJets());
185da5d3 242if (fOutput->GetNJets()>1)
243{
244//========================= CASE 2 ===========================
245 Int_t nPartons = fOutput->GetNPartons();
ab01dff2 246 fhNJets2->Fill(fOutput->GetNJets());
185da5d3 247 AliEMCALParton* parton;
248 AliEMCALJet* jethighest=0;
249 AliEMCALJet* jetsecond=0;
250 // Find Highest and Second Highest Jet
251 for (Int_t counter = 0; counter<fOutput->GetNJets();counter++)
252 {
253 if (counter==0){
254 jethighest = fOutput->GetJet(0);
255 jetsecond = fOutput->GetJet(1);
256 }
257 if (counter>0)
258 {
259 Float_t energyhighest = jethighest->Energy();
260 Float_t energysecond = jetsecond->Energy();
261
262 if ((fOutput->GetJet(counter))->Energy()>energyhighest)
263 {
264 jetsecond=jethighest;
265 jethighest=fOutput->GetJet(counter);
266 }else if ((fOutput->GetJet(counter))->Energy()>energysecond)
267 {
268 jetsecond=fOutput->GetJet(counter);
269 }
270 }
271 }
272
273 // End finding highest and second highest and continue
ab01dff2 274 fhJetEt2->Fill(jethighest->Energy());
275 fhJetEta2->Fill(jethighest->Eta() );
276 fhJetPhi2->Fill(jethighest->Phi() );
185da5d3 277 if (nPartons ==0) return;
278 parton = fOutput->GetParton(0);
279
ab01dff2 280 fhPartonEta2->Fill( parton->Eta() );
281 fhPartonPhi2->Fill( parton->Phi() );
185da5d3 282
283 //hJetEtDiff->Fill( jet->Energy() - parton->Energy() );
ab01dff2 284 fhEtaDiff2->Fill( jethighest->Eta() - parton->Eta() );
285 fhPhiDiff2->Fill( jethighest->Phi() - parton->Phi() );
286 fhEtaPhiSpread2->Fill(jethighest->Eta()-parton->Eta(),jethighest->Phi() - parton->Phi());
287 fhJetEtSecond2->Fill(jetsecond->Energy());
288 fhJetEtRatio2->Fill(jetsecond->Energy()/jethighest->Energy());
289 fhEtaPhiDist2->Fill( TMath::Sqrt((jethighest->Eta() - jetsecond->Eta())*(jethighest->Eta() - jetsecond->Eta())
185da5d3 290 + (jethighest->Phi() - jetsecond->Phi())*(jethighest->Phi() - jetsecond->Phi()) ));
291 /*
292 Float_t *pt,*phi,*eta;
293 Int_t *pdg;
294 pt = new Float_t[parton->GetNTracks()];
295 eta = new Float_t[parton->GetNTracks()];
296 phi = new Float_t[parton->GetNTracks()];
297 pdg = new Int_t[parton->GetNTracks()];*/
298
299
300
301 Float_t pt[2000];
302 Float_t eta[2000];
303 Float_t phi[2000];
304 Int_t pdg[2000];
305
306 parton->GetTrackList(pt,eta,phi,pdg);
307 for(Int_t iT=0; iT< parton->GetNTracks() ; iT++ )
308 {
309 if ( (eta[iT]-parton->Eta())*(eta[iT]-parton->Eta())+
310 (phi[iT]-parton->Phi())*(phi[iT]-parton->Phi()) >fConeRadius * fConeRadius ) continue;
311 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
312 Double_t rt = 2.0*atan(exp(-parton->Eta()));
313 Double_t ctt = cos(tt);
314 Double_t crt = cos(rt);
315 Double_t stt = sin(tt);
316 Double_t srt = sin(rt);
317 Double_t ctp = cos(phi[iT]);
318 Double_t crp = cos(parton->Phi());
319 Double_t stp = sin(phi[iT]);
320 Double_t srp = sin(parton->Phi());
321 Double_t alpha = acos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
322 Double_t correctp = pt[iT]/stt;
ab01dff2 323 fhPartonPL2->Fill( correctp*cos(alpha));
185da5d3 324 if ( (parton->Eta()-eta[iT])*(parton->Eta()-eta[iT]) +
325 (parton->Phi()-phi[iT])*(parton->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 326 fhPartonJT2->Fill( correctp*sin(alpha));
185da5d3 327 if (fNominalEnergy == 0.0) {
ab01dff2 328 fhPartonFragmFcn2->Fill( correctp*sin(tt)/parton->Energy() );
185da5d3 329 }else
330 {
ab01dff2 331 fhPartonFragmFcn2->Fill(correctp*sin(tt)/fNominalEnergy);
185da5d3 332 }
333 }// loop over tracks
334
335/*
336 pt = new Float_t[jet->NTracks()];
337 eta = new Float_t[jet->NTracks()];
338 phi = new Float_t[jet->NTracks()];
339 pdg = new Int_t[jet->NTracks()];*/
340 jethighest->TrackList(pt,eta,phi,pdg);
341 for(Int_t iT=0; iT< jethighest->NTracks() ; iT++ )
342 {
343 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
344 Double_t rt = 2.0*atan(exp(-jethighest->Eta()));
345 Double_t ctt = cos(tt);
346 Double_t crt = cos(rt);
347 Double_t stt = sin(tt);
348 Double_t srt = sin(rt);
349 Double_t ctp = cos(phi[iT]);
350 Double_t crp = cos(jethighest->Phi());
351 Double_t stp = sin(phi[iT]);
352 Double_t srp = sin(jethighest->Phi());
353 Double_t alpha = acos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
354 Double_t correctp = pt[iT]/stt;
ab01dff2 355 fhJetPL2->Fill( correctp*cos(alpha));
185da5d3 356 if ( (jethighest->Eta()-eta[iT])*(jethighest->Eta()-eta[iT]) +
357 (jethighest->Phi()-phi[iT])*(jethighest->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 358 fhJetJT2->Fill( correctp*sin(alpha));
185da5d3 359 if (fNominalEnergy==0.0){
ab01dff2 360 fhFragmFcn2->Fill( correctp*sin(tt)/parton->Energy() );
185da5d3 361 } else
362 {
ab01dff2 363 fhFragmFcn2->Fill( correctp*sin(tt)/fNominalEnergy );
185da5d3 364 }
365 }// loop over tracks
366 }
367
368 if (fOutput->GetNJets()==1)
369 {
370
371//========================= CASE 1 ===========================
372 Int_t nPartons = fOutput->GetNPartons();
f7d5860b 373 if (fOutput->GetNJets()!=1) return;
374 AliEMCALParton* parton;
375 AliEMCALJet* jet;
376 jet = fOutput->GetJet(0);
ab01dff2 377 fhJetEt->Fill(jet->Energy());
378 fhJetEta->Fill(jet->Eta() );
379 fhJetPhi->Fill(jet->Phi() );
f7d5860b 380 if (nPartons ==0) return;
381 parton = fOutput->GetParton(0);
185da5d3 382
ab01dff2 383 fhPartonEta->Fill( parton->Eta() );
384 fhPartonPhi->Fill( parton->Phi() );
f7d5860b 385
386 //hJetEtDiff->Fill( jet->Energy() - parton->Energy() );
ab01dff2 387 fhEtaDiff->Fill( jet->Eta() - parton->Eta() );
388 fhPhiDiff->Fill( jet->Phi() - parton->Phi() );
389 fhEtaPhiSpread->Fill(jet->Eta()-parton->Eta(),jet->Phi() - parton->Phi());
f7d5860b 390 /*
391 Float_t *pt,*phi,*eta;
392 Int_t *pdg;
393 pt = new Float_t[parton->GetNTracks()];
394 eta = new Float_t[parton->GetNTracks()];
395 phi = new Float_t[parton->GetNTracks()];
396 pdg = new Int_t[parton->GetNTracks()];*/
397
398
399
400 Float_t pt[2000];
401 Float_t eta[2000];
402 Float_t phi[2000];
403 Int_t pdg[2000];
404
405 parton->GetTrackList(pt,eta,phi,pdg);
406 for(Int_t iT=0; iT< parton->GetNTracks() ; iT++ )
407 {
408 if ( (eta[iT]-parton->Eta())*(eta[iT]-parton->Eta())+
409 (phi[iT]-parton->Phi())*(phi[iT]-parton->Phi()) >fConeRadius * fConeRadius ) continue;
410 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
411 Double_t rt = 2.0*atan(exp(-parton->Eta()));
412 Double_t ctt = cos(tt);
413 Double_t crt = cos(rt);
414 Double_t stt = sin(tt);
415 Double_t srt = sin(rt);
416 Double_t ctp = cos(phi[iT]);
417 Double_t crp = cos(parton->Phi());
418 Double_t stp = sin(phi[iT]);
419 Double_t srp = sin(parton->Phi());
420 Double_t alpha = acos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
421 Double_t correctp = pt[iT]/stt;
ab01dff2 422 fhPartonPL->Fill( correctp*cos(alpha));
f7d5860b 423 if ( (parton->Eta()-eta[iT])*(parton->Eta()-eta[iT]) +
424 (parton->Phi()-phi[iT])*(parton->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 425 fhPartonJT->Fill( correctp*sin(alpha));
f7d5860b 426 if (fNominalEnergy == 0.0) {
ab01dff2 427 fhPartonFragmFcn->Fill( correctp*sin(tt)/parton->Energy() );
f7d5860b 428 }else
429 {
ab01dff2 430 fhPartonFragmFcn->Fill(correctp*sin(tt)/fNominalEnergy);
f7d5860b 431 }
432 }// loop over tracks
433
434/*
435 pt = new Float_t[jet->NTracks()];
436 eta = new Float_t[jet->NTracks()];
437 phi = new Float_t[jet->NTracks()];
438 pdg = new Int_t[jet->NTracks()];*/
439 jet->TrackList(pt,eta,phi,pdg);
440 for(Int_t iT=0; iT< jet->NTracks() ; iT++ )
441 {
442 Double_t tt = 2.0*atan(exp(-eta[iT])); // These names are short to make the equation manageable
443 Double_t rt = 2.0*atan(exp(-jet->Eta()));
444 Double_t ctt = cos(tt);
445 Double_t crt = cos(rt);
446 Double_t stt = sin(tt);
447 Double_t srt = sin(rt);
448 Double_t ctp = cos(phi[iT]);
449 Double_t crp = cos(jet->Phi());
450 Double_t stp = sin(phi[iT]);
451 Double_t srp = sin(jet->Phi());
452 Double_t alpha = acos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);
453 Double_t correctp = pt[iT]/stt;
ab01dff2 454 fhJetPL->Fill( correctp*cos(alpha));
f7d5860b 455 if ( (jet->Eta()-eta[iT])*(jet->Eta()-eta[iT]) +
456 (jet->Phi()-phi[iT])*(jet->Phi()-phi[iT]) < 0.2*0.2 )
ab01dff2 457 fhJetJT->Fill( correctp*sin(alpha));
f7d5860b 458 if (fNominalEnergy==0.0){
ab01dff2 459 fhFragmFcn->Fill( correctp*sin(tt)/parton->Energy() );
f7d5860b 460 } else
461 {
ab01dff2 462 fhFragmFcn->Fill( correctp*sin(tt)/fNominalEnergy );
f7d5860b 463 }
464 }// loop over tracks
185da5d3 465 }
f7d5860b 466}
467
468