new analysis classes, not yet in build system
[u/mrichter/AliRoot.git] / HLT / JET / analysis / AliHLTJETAnalysisJets.cxx
1 //-*- Mode: C++ -*-
2 // $Id: AliHLTJETAnalysisJets.cxx  $
3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE HLT Project        * 
5 //* ALICE Experiment at CERN, All rights reserved.                         *
6 //*                                                                        *
7 //* Primary Authors: Jochen Thaeder <thaeder@kip.uni-heidelberg.de>        *
8 //*                  for The ALICE HLT Project.                            *
9 //*                                                                        *
10 //* Permission to use, copy, modify and distribute this software and its   *
11 //* documentation strictly for non-commercial purposes is hereby granted   *
12 //* without fee, provided that the above copyright notice appears in all   *
13 //* copies and that both the copyright notice and this permission notice   *
14 //* appear in the supporting documentation. The authors make no claims     *
15 //* about the suitability of this software for any purpose. It is          *
16 //* provided "as is" without express or implied warranty.                  *
17 //**************************************************************************
18
19 /** @file   AliHLTJETAnalysisJets.cxx
20     @author Jochen Thaeder
21     @date   
22     @brief  Container holding analysis objects
23 */
24
25 // see header file for class documentation
26 // or
27 // refer to README to build package
28 // or
29 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
30
31 #if __GNUC__>= 3
32    using namespace std;
33 #endif
34
35 #include "TH2F.h"
36
37 #include "AliHLTJETAnalysisJets.h"
38
39 /** ROOT macro for the implementation of ROOT specific class methods */
40 ClassImp(AliHLTJETAnalysisJets)
41
42 /*
43  * ---------------------------------------------------------------------------------
44  *                            Constructor / Destructor
45  * ---------------------------------------------------------------------------------
46  */
47
48 //##################################################################################
49 AliHLTJETAnalysisJets::AliHLTJETAnalysisJets() :
50   fJets(NULL ),
51   fJetsMC(NULL ),
52   fMatchedJets(NULL),
53   fMatchedJetsMC(NULL),
54   fDeltaEt(NULL),
55   fDeltaEta(NULL),
56   fDeltaPhi(NULL),
57   fDeltaEtaDeltaPhi(NULL),
58   fSpectraEt(NULL),
59   fSpectraEta(NULL),
60   fSpectraPhi(NULL),
61   fCorrelationsJetEt(NULL),
62   fResolutionsJetEt(NULL),
63   fResolutionsDiJetEt(NULL) {
64   // see header file for class documentation
65   // or
66   // refer to README to build package
67   // or
68   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
69   
70 }
71
72 //##################################################################################
73 AliHLTJETAnalysisJets::~AliHLTJETAnalysisJets() {
74   // see header file for class documentation
75
76   if ( fJetsMC )
77     delete fJetsMC;
78   fJetsMC = NULL;
79
80   if ( fMatchedJets )
81     delete fMatchedJets;
82   fMatchedJets = NULL;
83
84   if ( fMatchedJetsMC )
85     delete fMatchedJetsMC;
86   fMatchedJetsMC = NULL;
87
88   if ( fDeltaEt ) {
89     fDeltaEt->Clear();
90     delete fDeltaEt;
91   }
92   fDeltaEt = NULL;
93
94   if ( fDeltaEta ) {
95     fDeltaEta->Clear();
96     delete fDeltaEta;
97   }
98   fDeltaEta = NULL;
99
100   if ( fDeltaPhi ) {
101     fDeltaPhi->Clear();
102     delete fDeltaPhi;
103   }
104   fDeltaPhi = NULL;
105   
106   if ( fDeltaEtaDeltaPhi ) {
107     fDeltaEtaDeltaPhi->Clear();
108     delete fDeltaEtaDeltaPhi;
109   }
110   fDeltaEtaDeltaPhi = NULL;
111
112   if ( fSpectraEt ) {
113     fSpectraEt->Clear();
114     delete fSpectraEt;
115   }
116   fSpectraEt = NULL;
117
118   if ( fSpectraEta ) {
119     fSpectraEta->Clear();
120     delete fSpectraEta;
121   }
122   fSpectraEta = NULL;
123
124   if ( fSpectraPhi ) {
125     fSpectraPhi->Clear();
126     delete fSpectraPhi;
127   }
128   fSpectraPhi = NULL;
129
130   if ( fCorrelationsJetEt ) {
131     fCorrelationsJetEt->Clear();
132     delete fCorrelationsJetEt;   
133   }
134   fCorrelationsJetEt = NULL;
135
136   if ( fResolutionsJetEt ) {
137     fResolutionsJetEt->Clear();
138     delete fResolutionsJetEt;   
139   }
140   fResolutionsJetEt = NULL;
141
142   if ( fResolutionsDiJetEt ) {
143     fResolutionsDiJetEt->Clear();
144     delete fResolutionsDiJetEt;   
145   }
146   fResolutionsDiJetEt = NULL;
147
148 }
149
150 /*
151  * ---------------------------------------------------------------------------------
152  *                                   Initialize / Reset
153  * ---------------------------------------------------------------------------------
154  */
155
156 //##################################################################################
157 Int_t AliHLTJETAnalysisJets::Initialize() {
158   // see header file for class documentation  
159
160   Int_t iResult = 0;
161
162   // -- Setup match arrays
163   fMatchedJets = new TArrayI(100);
164   fMatchedJetsMC = new TArrayI(100);
165
166   // -- Setup Delta histograms  
167   SetupDeltaHistograms();
168
169   // -- Setup Spectra histograms
170   SetupSpectraHistograms();
171
172   // -- Setup Matched histograms
173   SetupMatchedHistograms();
174
175   return iResult;
176 }
177
178 //##################################################################################
179 void AliHLTJETAnalysisJets::ResetEvent() {
180   // see header file for class documentation  
181
182   if ( fJetsMC )
183     delete fJetsMC;
184   fJetsMC = NULL;
185
186   fJets = NULL;
187
188   return;
189 }
190
191 /*
192  * ---------------------------------------------------------------------------------
193  *                                Setter - public
194  * ---------------------------------------------------------------------------------
195  */
196
197 //##################################################################################
198 void AliHLTJETAnalysisJets::SetHLTMC( AliHLTMCEvent* mcEvent ) {
199   // see header file for class documentation
200
201   // -- New MC jets 
202   if ( fJetsMC )
203     delete fJetsMC;
204   fJetsMC = new AliHLTJETJets();
205   
206   AliAODJet* jet = NULL;
207
208   while ( (jet = mcEvent->NextGenJet()) )
209     fJetsMC->AddJet(jet);
210
211   // -- Sort jets
212   fJetsMC->Sort();
213
214   return;
215 }
216
217 //##################################################################################
218 void AliHLTJETAnalysisJets::SetMC( AliMCEvent* /*mcEvent*/ ) {
219   // see header file for class documentation
220
221   HLTFatal("No implemented!");
222   
223   return;
224 }
225
226 //##################################################################################
227 void AliHLTJETAnalysisJets::SetJets( AliHLTJETJets* jets ) {
228   // see header file for class documentation
229
230   fJets = jets;
231
232   // -- Sort jets
233   fJets->Sort();
234
235   return;
236 }
237
238 /*
239  * ---------------------------------------------------------------------------------
240  *                                 Getter - public
241  * ---------------------------------------------------------------------------------
242  */
243
244 //##################################################################################
245 TH1* AliHLTJETAnalysisJets::GetHistogram( Int_t histIdx, Int_t plotIdx ) {
246   // see header file for class documentation
247   
248   if ( histIdx >= AliHLTJETAnalysisBase::kHistMax ) {
249     HLTError("Histogram index out of bound : %d - max = %d", histIdx, AliHLTJETAnalysisBase::kHistMax );
250     return NULL; 
251   }
252   
253   TClonesArray *hist = NULL;
254   Int_t maxIdx = -1;
255   
256   switch (histIdx ) {
257   case AliHLTJETAnalysisBase::kHistDeltaEt : 
258     hist = fDeltaEt;
259     maxIdx = AliHLTJETAnalysisBase::kDeltaMax;
260     break;
261   case AliHLTJETAnalysisBase::kHistDeltaEta : 
262     hist = fDeltaEta;
263     maxIdx = AliHLTJETAnalysisBase::kDeltaMax;
264     break;
265   case AliHLTJETAnalysisBase::kHistDeltaPhi : 
266     hist = fDeltaPhi;
267     maxIdx = AliHLTJETAnalysisBase::kDeltaMax;
268     break;
269   case AliHLTJETAnalysisBase::kHistDeltaEtaDeltaPhi : 
270     hist = fDeltaEtaDeltaPhi;
271     maxIdx = AliHLTJETAnalysisBase::kDeltaMax;
272     break;
273   case AliHLTJETAnalysisBase::kHistSpectraEt :
274     hist = fSpectraEt;
275     maxIdx = AliHLTJETAnalysisBase::kSpectraMax;
276     break;
277   case AliHLTJETAnalysisBase::kHistSpectraEta :
278     hist = fSpectraEta;
279     maxIdx = AliHLTJETAnalysisBase::kSpectraMax;
280     break;
281   case AliHLTJETAnalysisBase::kHistSpectraPhi :
282     hist = fSpectraPhi;
283     maxIdx = AliHLTJETAnalysisBase::kSpectraMax;
284     break;
285   case AliHLTJETAnalysisBase::kHistCorrelationsJetEt :
286     hist = fCorrelationsJetEt;
287     maxIdx = AliHLTJETAnalysisBase::kPlotMax;
288     break;
289   case AliHLTJETAnalysisBase::kHistResolutionsJetEt :
290     hist = fResolutionsJetEt;
291     maxIdx = AliHLTJETAnalysisBase::kPlotMax;
292     break;
293   case AliHLTJETAnalysisBase::kHistResolutionsDiJetEt :
294     hist = fResolutionsDiJetEt;
295     maxIdx = AliHLTJETAnalysisBase::kPlotMax;
296     break;
297   }
298  
299   // -- Check Boundaries
300   if ( plotIdx >= maxIdx ) {
301     HLTError("Index out of bound : %d - max = %d", plotIdx, maxIdx );
302     return NULL; 
303   }
304   else if ( maxIdx == -1 ) {
305     HLTError("Histogram index not found." );
306     return NULL; 
307   }
308
309   // -- Retrieve histogram
310   return reinterpret_cast<TH1*>((*hist)[plotIdx]);  
311 }
312
313 /*
314  * ---------------------------------------------------------------------------------
315  *                              Analysis - public
316  * ---------------------------------------------------------------------------------
317  */
318
319 //##################################################################################
320 Int_t AliHLTJETAnalysisJets::Analyze() {
321   // see header file for class documentation
322   
323   Int_t iResult = 0;
324
325   if ( !fJets ) {
326     HLTError("No input jets set.");
327     iResult = -1;
328   }
329   if ( !fJetsMC ) {
330     HLTError("No input MC jets set.");
331     iResult = -1;
332   }
333
334   if ( !iResult) {    
335     
336     // -- Fill unmatched jets into histograms
337     FillBasicSpectraHistograms();
338     FillUnmatchedDeltaHistograms();
339     
340     // -- Match jets
341     MatchJets();
342    
343     // -- Fill matched jets into histograms
344     FillMatchedDeltaHistograms();  
345     FillMatchedSpectraHistograms();  
346     FillMatchedHistograms();  
347   }
348
349   return iResult;
350 }
351
352 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
353 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
354 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
355 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
356 ///                                                                                                              ///
357 //////                                               PRIVATE                                                  //////
358 ///                                                                                                              ///
359 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
360 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
361 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
362 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
363
364 /*
365  * ---------------------------------------------------------------------------------
366  *                             Setup / Reset - private
367  * ---------------------------------------------------------------------------------
368  */
369
370 //##################################################################################
371 void AliHLTJETAnalysisJets::SetupDeltaHistograms() {
372   // see header file for class documentation
373   
374   // ---------------------------------------------------
375   // -- Difference in reconstruction
376   //    All / Matched jets
377   // ---------------------------------------------------
378
379   fDeltaEt                 = new TClonesArray( "TH1F", AliHLTJETAnalysisBase::kDeltaMax );
380   fDeltaEta                = new TClonesArray( "TH1F", AliHLTJETAnalysisBase::kDeltaMax );
381   fDeltaPhi                = new TClonesArray( "TH1F", AliHLTJETAnalysisBase::kDeltaMax );
382   fDeltaEtaDeltaPhi        = new TClonesArray( "TH2F", AliHLTJETAnalysisBase::kDeltaMax );
383
384   for ( Int_t idx = 0; idx < AliHLTJETAnalysisBase::kDeltaMax; ++idx ) {
385
386     // -- Delta Et -------------------------------------------------
387     new ((*fDeltaEt)[idx]) TH1F(Form("delta E_{t} : %s", 
388                                      AliHLTJETAnalysisBase::fgkDeltaType[idx] ), 
389                                 Form("#DeltaE_{t} : %s;#DeltaE_{t};dN/d#DeltaE_{t}", 
390                                      AliHLTJETAnalysisBase::fgkDeltaType[idx] ),
391                                 100, -200., 200.);
392     SetupHist(reinterpret_cast<TH1F*>((*fDeltaEt)[idx]));
393         
394     // -- Delta Eta ------------------------------------------------
395     new ((*fDeltaEta)[idx]) TH1F(Form("delta Eta : %s", 
396                                       AliHLTJETAnalysisBase::fgkDeltaType[idx] ), 
397                                  Form("#Delta#eta : %s;#Delta#eta;dN/d#Delta#eta", 
398                                       AliHLTJETAnalysisBase::fgkDeltaType[idx] ),
399                                  100, -1.2, 1.2);
400     SetupHist(reinterpret_cast<TH1F*>((*fDeltaEta)[idx]));
401     
402     // -- Delta Phi ------------------------------------------------
403     new ((*fDeltaPhi )[idx]) TH1F(Form("delta Phi : %s", 
404                                        AliHLTJETAnalysisBase::fgkDeltaType[idx] ), 
405                                   Form("#Delta#phi : %s;#Delta #phi;dN/d#Delta#phi", 
406                                        AliHLTJETAnalysisBase::fgkDeltaType[idx] ), 
407                                   100, -7., 7.);
408     SetupHist(reinterpret_cast<TH1F*>((*fDeltaPhi)[idx]));
409     
410     // -- Delta Eta Delta Phi --------------------------------------
411     new ((*fDeltaEtaDeltaPhi) [idx] ) TH2F(Form("delta Eta delta Phi : %s", 
412                                                 AliHLTJETAnalysisBase::fgkDeltaType[idx] ), 
413                                            Form("#Delta#eta #Delta#phi : %s;#Delta#eta;#Delta#phi", 
414                                                 AliHLTJETAnalysisBase::fgkDeltaType[idx] ), 
415                                            100, -1.2, 1.2, 100, -7., 7.);
416     SetupHist(reinterpret_cast<TH2F*>((*fDeltaEtaDeltaPhi)[idx]));
417
418   } //  for ( Int_t idx = 0; idx < AliHLTJETAnalysisBase::kDeltaMax; ++idx ) {
419
420   return;
421 }
422
423 //##################################################################################
424 void AliHLTJETAnalysisJets::SetupSpectraHistograms() {
425   // see header file for class documentation
426     
427   // ---------------------------------------------------
428   // -- Jet spectra
429   // ---------------------------------------------------
430
431   fSpectraEt  = new TClonesArray( "TH1F", AliHLTJETAnalysisBase::kSpectraMax );
432   fSpectraEta = new TClonesArray( "TH1F", AliHLTJETAnalysisBase::kSpectraMax );
433   fSpectraPhi = new TClonesArray( "TH1F", AliHLTJETAnalysisBase::kSpectraMax );
434
435   for ( Int_t idx = 0; idx < AliHLTJETAnalysisBase::kSpectraMax; ++idx ) {
436     
437     // -- Spectra Et -----------------------------------------------
438     new ((*fSpectraEt)[idx]) TH1F(Form("E_{t} : %s", 
439                                        AliHLTJETAnalysisBase::fgkSpectraType[idx] ), 
440                                   Form("E_{t} : %s;E_{t} (GeV/c);dN/dE_{t}", 
441                                        AliHLTJETAnalysisBase::fgkSpectraType[idx] ),
442                                   100, 0., 200.);
443     SetupHist(reinterpret_cast<TH1F*>((*fSpectraEt)[idx]));
444     
445     // -- Spectra Eta ----------------------------------------------
446     new ((*fSpectraEta)[idx]) TH1F(Form("#eta : %s", 
447                                         AliHLTJETAnalysisBase::fgkSpectraType[idx] ), 
448                                    Form("#eta : %s;#eta;dN/d#eta",
449                                         AliHLTJETAnalysisBase::fgkSpectraType[idx] ),
450                                    80, -0.9, 0.9);
451     SetupHist(reinterpret_cast<TH1F*>((*fSpectraEta)[idx]));
452     
453     // -- Spectra Phi ----------------------------------------------
454     new ((*fSpectraPhi)[idx]) TH1F(Form("#phi : %s", 
455                                         AliHLTJETAnalysisBase::fgkSpectraType[idx] ), 
456                                    Form("#phi : %s;#phi;dN/d#phi",
457                                         AliHLTJETAnalysisBase::fgkSpectraType[idx] ),
458                                    50, 0., 7.);
459     SetupHist(reinterpret_cast<TH1F*>((*fSpectraPhi)[idx]));
460     
461   } // for ( Int_t idx = 0; idx < AliHLTJETAnalysisBase::kSpectraMax; ++idx ) {
462
463   return;
464 }
465
466 //##################################################################################
467 void AliHLTJETAnalysisJets::SetupMatchedHistograms() {
468   // see header file for class documentation
469   
470   // ---------------------------------------------------
471   // -- Correlations
472   // -- Resolutions
473   // ---------------------------------------------------
474   
475   fCorrelationsJetEt  = new TClonesArray( "TH2F", AliHLTJETAnalysisBase::kPlotMax );
476   fResolutionsJetEt   = new TClonesArray( "TH2F", AliHLTJETAnalysisBase::kPlotMax );
477   fResolutionsDiJetEt = new TClonesArray( "TH2F", AliHLTJETAnalysisBase::kPlotMax );
478
479   for ( Int_t idx = 0; idx < AliHLTJETAnalysisBase::kPlotMax; ++idx ) {
480
481     // -- Correlations ---------------------------------------------
482     new ((*fCorrelationsJetEt)[idx]) TH2F(Form("Correlations : %s", 
483                                                AliHLTJETAnalysisBase::fgkPlotType[idx] ), 
484                                           Form("Correlations : %s; E_{t,Pythia} (GeV/c); E_{t,Rec} (GeV/c)", 
485                                                AliHLTJETAnalysisBase::fgkPlotType[idx] ),
486                                           100, 0., 200.,100, 0., 200.);
487     SetupHist(reinterpret_cast<TH1F*>((*fCorrelationsJetEt)[idx]));
488     
489     // -- Jet Resolutions ------------------------------------------
490     new ((*fResolutionsJetEt)[idx]) TH2F(Form("Resolutions : %s", 
491                                               AliHLTJETAnalysisBase::fgkPlotType[idx] ), 
492                                          Form("Resolutions : %s; E_{t,Pythia} (GeV/c); f", 
493                                               AliHLTJETAnalysisBase::fgkPlotType[idx] ),
494                                          200, 0., 200.,200, -2., 2.);
495                                         
496     SetupHist(reinterpret_cast<TH1F*>((*fResolutionsJetEt)[idx]));
497
498     // -- Di-Jet Resolutions ---------------------------------------
499     new ((*fResolutionsDiJetEt)[idx]) TH2F(Form("Di Jet Resolutions : %s", 
500                                                 AliHLTJETAnalysisBase::fgkPlotType[idx] ), 
501                                            Form("Di Jet Resolutions : %s; E_{t,Pythia} (GeV/c); f", 
502                                                 AliHLTJETAnalysisBase::fgkPlotType[idx] ),
503                                            100, -200., 200.,100, -200., 200.);
504                                         
505     SetupHist(reinterpret_cast<TH1F*>((*fResolutionsDiJetEt)[idx]));
506
507   } // for ( Int_t idx = 0; idx < kJetPlotMax; ++idx ) {
508
509   return;
510 }
511
512 /*
513  * ---------------------------------------------------------------------------------
514  *                             Analysis - private
515  * ---------------------------------------------------------------------------------
516  */
517
518 //##################################################################################
519 Int_t AliHLTJETAnalysisJets::MatchJets() {
520   // see header file for class documentation
521
522   Int_t iResult = 0;
523
524   // -- Reset match arrays
525   Int_t maxNJets;
526   if ( fJets->GetNAODJets() > fJetsMC->GetNAODJets() )
527     maxNJets = fJets->GetNAODJets();
528   else
529     maxNJets = fJetsMC->GetNAODJets();
530
531   for ( Int_t idx = 0; idx < maxNJets; ++idx ) {
532     (*fMatchedJets)[idx] = -1;
533     (*fMatchedJetsMC)[idx] = -1;
534   }
535
536   // -- Match Jets - pythia fixed
537   for ( Int_t jetIterMC = 0; jetIterMC < fJetsMC->GetNAODJets(); ++jetIterMC ) {
538
539     Float_t minDistance2 = 100.;
540     Int_t idxClosest = -1;
541
542     for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
543
544       // -- very Jet matched only once
545       if ( (*fMatchedJets)[jetIter] != -1 )
546         continue;
547                           
548       Float_t distance2 = GetDistance2( fJetsMC->GetJet(jetIterMC), 
549                                         fJets->GetJet(jetIter) );
550       
551       // -- Get Closest
552       if ( distance2 < minDistance2 ) {
553         minDistance2 = distance2;
554         idxClosest = jetIter;
555       }
556
557     } //     for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
558  
559     // -- Match found
560     if  ( minDistance2 < 10. ) {
561       (*fMatchedJetsMC)[jetIterMC] = idxClosest;
562       (*fMatchedJets)[idxClosest] = jetIterMC;
563     }
564
565   } //  for ( Int_t jetIterMC = 0; jetIterMC < fJetsMC->GetNAODJets(); ++jetIterMC ) {
566
567   return iResult;
568 }
569
570 /*
571  * ---------------------------------------------------------------------------------
572  *                                Fill - private
573  * ---------------------------------------------------------------------------------
574  */
575
576 //##################################################################################
577 void AliHLTJETAnalysisJets::FillBasicSpectraHistograms() {
578   // see header file for class documentation
579
580   // -- Fill basic MC spectras
581   // ---------------------------
582   for ( Int_t jetIter = 0; jetIter < fJetsMC->GetNAODJets(); ++jetIter ) {
583     FillHist(fSpectraEt, AliHLTJETAnalysisBase::kSpectraPythiaAll, fJetsMC->GetJet(jetIter)->Pt());
584     FillHist(fSpectraEta, AliHLTJETAnalysisBase::kSpectraPythiaAll, fJetsMC->GetJet(jetIter)->Eta());
585     FillHist(fSpectraPhi, AliHLTJETAnalysisBase::kSpectraPythiaAll, fJetsMC->GetJet(jetIter)->Phi());
586   } // for ( Int_t jetIter = 0; jetIter < fJetsMC->GetNAODJets(); ++jetIter ) {
587
588   // -- Fill basic Reco spectras
589   // -----------------------------
590   for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
591     FillHist(fSpectraEt, AliHLTJETAnalysisBase::kSpectraRecoAll, fJets->GetJet(jetIter)->Pt());
592     FillHist(fSpectraEta, AliHLTJETAnalysisBase::kSpectraRecoAll, fJets->GetJet(jetIter)->Eta());
593     FillHist(fSpectraPhi, AliHLTJETAnalysisBase::kSpectraRecoAll, fJets->GetJet(jetIter)->Phi());
594   } // for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
595
596   // -- Fill basic Reco leading spectras
597   // -------------------------------------
598   if ( fJets->GetNAODJets() > 0 ) {
599     FillHist(fSpectraEt, AliHLTJETAnalysisBase::kSpectraRecoLeadAll, fJets->GetJet(0)->Pt());
600     FillHist(fSpectraEta, AliHLTJETAnalysisBase::kSpectraRecoLeadAll, fJets->GetJet(0)->Eta());
601     FillHist(fSpectraPhi, AliHLTJETAnalysisBase::kSpectraRecoLeadAll, fJets->GetJet(0)->Phi());
602   } 
603
604   return;
605 }
606
607 //##################################################################################
608 void AliHLTJETAnalysisJets::FillUnmatchedDeltaHistograms() {
609   // see header file for class documentation
610
611   for ( Int_t jetIterMC = 0; jetIterMC < fJetsMC->GetNAODJets(); ++jetIterMC ) {
612     AliAODJet *jetMC = fJetsMC->GetJet(jetIterMC);
613     
614     for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
615       AliAODJet *jet = fJets->GetJet(jetIter);
616
617       FillHist(fDeltaEt, AliHLTJETAnalysisBase::kDeltaAll, jetMC->Pt()-jet->Pt());
618       FillHist(fDeltaEta, AliHLTJETAnalysisBase::kDeltaAll, jetMC->Eta()-jet->Eta());
619       FillHist(fDeltaPhi, AliHLTJETAnalysisBase::kDeltaAll, jetMC->Phi()-jet->Phi());
620       FillHist(fDeltaEtaDeltaPhi, AliHLTJETAnalysisBase::kDeltaAll,
621                jetMC->Eta()-jet->Eta(), jetMC->Phi()-jet->Phi());
622       
623     } // for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
624   } // for ( Int_t jetIterMC = 0; jetIterMC < fJetsMC->GetNAODJets(); ++jetIterMC ) {
625
626   // -- Leading Jets
627   // -----------------
628   if ( fJets->GetNAODJets() > 0 && fJetsMC->GetNAODJets() > 0 ) {
629     AliAODJet *jetMC = fJetsMC->GetJet(0);
630     AliAODJet *jet = fJets->GetJet(0);
631
632     FillHist(fDeltaEt, AliHLTJETAnalysisBase::kDeltaLead, jetMC->Pt()-jet->Pt());
633     FillHist(fDeltaEta, AliHLTJETAnalysisBase::kDeltaLead, jetMC->Eta()-jet->Eta());
634     FillHist(fDeltaPhi, AliHLTJETAnalysisBase::kDeltaLead, jetMC->Phi()-jet->Phi());
635     FillHist(fDeltaEtaDeltaPhi, AliHLTJETAnalysisBase::kDeltaLead,
636              jetMC->Eta()-jet->Eta(), jetMC->Phi()-jet->Phi());
637   }
638
639   return;
640 }
641
642 //##################################################################################
643 void AliHLTJETAnalysisJets::FillMatchedDeltaHistograms() {
644   // see header file for class documentation
645
646   for ( Int_t jetIterMC = 0; jetIterMC < fJetsMC->GetNAODJets(); ++jetIterMC ) {
647
648     // -- Not matched jet
649     if ( (*fMatchedJetsMC)[jetIterMC] == -1 )
650       continue;
651     
652     AliAODJet *jetMC = fJetsMC->GetJet(jetIterMC);
653     AliAODJet *jet = fJets->GetJet((*fMatchedJetsMC)[jetIterMC]);
654
655     FillHist(fDeltaEt, AliHLTJETAnalysisBase::kDeltaMatchedAll, jetMC->Pt()-jet->Pt());
656     FillHist(fDeltaEta, AliHLTJETAnalysisBase::kDeltaMatchedAll, jetMC->Eta()-jet->Eta());
657     FillHist(fDeltaPhi, AliHLTJETAnalysisBase::kDeltaMatchedAll, jetMC->Phi()-jet->Phi());
658     FillHist(fDeltaEtaDeltaPhi, AliHLTJETAnalysisBase::kDeltaMatchedAll,
659              jetMC->Eta()-jet->Eta(), jetMC->Phi()-jet->Phi());
660     
661   } // for ( Int_t jetIter = 0; jetIter < fJetsMC->GetNAODJets(); ++jetIter ) {
662   
663   // -- Leading Jets
664   // -----------------
665   if ((*fMatchedJetsMC)[0] != -1 ) {
666
667     AliAODJet *jetMC = fJetsMC->GetJet(0);
668     AliAODJet *jet = fJets->GetJet((*fMatchedJetsMC)[0]);
669     
670     FillHist(fDeltaEt, AliHLTJETAnalysisBase::kDeltaMatchedLead, jetMC->Pt()-jet->Pt());
671     FillHist(fDeltaEta, AliHLTJETAnalysisBase::kDeltaMatchedLead, jetMC->Eta()-jet->Eta());
672     FillHist(fDeltaPhi, AliHLTJETAnalysisBase::kDeltaMatchedLead, jetMC->Phi()-jet->Phi());
673     FillHist(fDeltaEtaDeltaPhi, AliHLTJETAnalysisBase::kDeltaMatchedLead,
674              jetMC->Eta()-jet->Eta(), jetMC->Phi()-jet->Phi());
675   }
676   
677   return;
678 }
679  
680 //##################################################################################
681 void AliHLTJETAnalysisJets::FillMatchedSpectraHistograms() {
682   // see header file for class documentation
683
684   Int_t idx = 0;
685
686   // -- Fill matched MC spectras
687   // ---------------------------
688   for ( Int_t jetIter = 0; jetIter < fJetsMC->GetNAODJets(); ++jetIter ) {
689     if ( (*fMatchedJetsMC)[jetIter] == -1 )
690       idx = AliHLTJETAnalysisBase::kSpectraPythiaUnmatched;
691     else
692       idx = AliHLTJETAnalysisBase::kSpectraPythiaMatched;
693
694     FillHist(fSpectraEt, idx, fJetsMC->GetJet(jetIter)->Pt());
695     FillHist(fSpectraEta, idx, fJetsMC->GetJet(jetIter)->Eta());
696     FillHist(fSpectraPhi, idx, fJetsMC->GetJet(jetIter)->Phi());
697   } // for ( Int_t jetIter = 0; jetIter < fJetsMC->GetNAODJets(); ++jetIter ) {
698
699   // -- Fill matched Reco spectras
700   // -----------------------------
701   for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
702     if ( (*fMatchedJets)[jetIter] == -1 )
703       idx = AliHLTJETAnalysisBase::kSpectraRecoUnmatched;
704     else
705       idx = AliHLTJETAnalysisBase::kSpectraRecoMatched;
706     
707     FillHist(fSpectraEt, idx, fJets->GetJet(jetIter)->Pt());
708     FillHist(fSpectraEta, idx, fJets->GetJet(jetIter)->Eta());
709     FillHist(fSpectraPhi, idx, fJets->GetJet(jetIter)->Phi());
710   } // for ( Int_t jetIter = 0; jetIter < fJets->GetNAODJets(); ++jetIter ) {
711
712   // -- Fill matched Reco leading spectras
713   // -------------------------------------
714   if ( fJets->GetNAODJets() > 0 ) {
715     if ( (*fMatchedJets)[0] == -1 )
716       idx = AliHLTJETAnalysisBase::kSpectraRecoLeadUnmatched;
717     else
718       idx = AliHLTJETAnalysisBase::kSpectraRecoLeadMatched;
719
720     FillHist(fSpectraEt, idx, fJets->GetJet(0)->Pt());
721     FillHist(fSpectraEta, idx, fJets->GetJet(0)->Eta());
722     FillHist(fSpectraPhi, idx, fJets->GetJet(0)->Phi());
723   } 
724   
725   return;
726 }
727
728 //##################################################################################
729 void AliHLTJETAnalysisJets::FillMatchedHistograms() {
730   // see header file for class documentation
731
732   for ( Int_t jetIterMC = 0; jetIterMC < fJetsMC->GetNAODJets(); ++jetIterMC ) {
733
734     // -- Not matched jet
735     if ( (*fMatchedJetsMC)[jetIterMC] == -1 )
736       continue;
737     
738     AliAODJet *jetMC = fJetsMC->GetJet(jetIterMC);
739     AliAODJet *jet = fJets->GetJet((*fMatchedJetsMC)[jetIterMC]);
740
741     // -- Correlations    
742     FillHist(fCorrelationsJetEt, AliHLTJETAnalysisBase::kPlotAll, jetMC->Pt(), jet->Pt());
743
744     // -- Resolutions
745     FillHist(fResolutionsJetEt, AliHLTJETAnalysisBase::kPlotAll, jetMC->Pt(),
746              (jetMC->Pt()-jet->Pt())/(jetMC->Pt()+jet->Pt()));
747
748   } // for ( Int_t jetIter = 0; jetIter < fJetsMC->GetNAODJets(); ++jetIter ) {
749   
750   // -- Leading Jets
751   // -----------------
752   if ((*fMatchedJetsMC)[0] != -1 ) {
753
754     AliAODJet *jetMC = fJetsMC->GetJet(0);
755     AliAODJet *jet = fJets->GetJet((*fMatchedJetsMC)[0]);
756     
757     // -- Correlations    
758     FillHist(fCorrelationsJetEt, AliHLTJETAnalysisBase::kPlotLead, jetMC->Pt(), jet->Pt());
759     
760     // -- Resolutions
761     FillHist(fResolutionsJetEt, AliHLTJETAnalysisBase::kPlotLead, jetMC->Pt(),
762              (jetMC->Pt()-jet->Pt())/(jetMC->Pt()+jet->Pt()));
763   }
764
765   return;
766 }
767
768 /*
769  * ---------------------------------------------------------------------------------
770  *                               Helper - private
771  * ---------------------------------------------------------------------------------
772  */
773
774 //##################################################################################
775 void AliHLTJETAnalysisJets::SetupHist( TH1* hist ) {
776   // see header file for class documentation
777   
778   hist->GetXaxis()->CenterTitle();
779   hist->GetYaxis()->CenterTitle();
780   hist->SetMarkerStyle(kFullCircle);
781   //  hist->SetStats(kFALSE);
782
783   return;
784 }
785
786 //##################################################################################
787 void AliHLTJETAnalysisJets::FillHist( TClonesArray* array, Int_t idx, Float_t valueX ) {
788     // see header file for class documentation
789   
790   reinterpret_cast<TH1F*>((*array)[idx])->Fill(valueX);
791   
792   return;
793 }
794
795 //##################################################################################
796 void AliHLTJETAnalysisJets::FillHist( TClonesArray* array, Int_t idx, Float_t valueX, Float_t valueY ) {
797   // see header file for class documentation
798   
799   reinterpret_cast<TH2F*>((*array)[idx])->Fill(valueX, valueY);
800   
801   return;
802 }
803
804 //################################################################################## 
805 Float_t AliHLTJETAnalysisJets::GetDistance2( AliAODJet *jet1, AliAODJet *jet2 ) {
806   // see header file for class documentation
807   
808   return ( (jet1->Eta()-jet2->Eta())*(jet1->Eta()-jet2->Eta()) ) + 
809     ( (jet1->Phi()-jet2->Phi())*(jet1->Phi()-jet2->Phi()) );
810 }