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