Split: fixed more module incpaths
[u/mrichter/AliRoot.git] / HLT / JET / cone / AliHLTJETConeGrid.cxx
1 //-*- Mode: C++ -*-
2 // $Id: AliHLTJETConeGrid.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   AliHLTJETConeGrid.cxx
20     @author Jochen Thaeder
21     @date   
22     @brief  Eta-Phi grid of the cone finder
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 #include "AliHLTJETConeGrid.h"
32 #include "AliHLTJETConeEtaPhiCell.h"
33
34 using namespace std;
35
36 /** ROOT macro for the implementation of ROOT specific class methods */
37 ClassImp(AliHLTJETConeGrid)
38
39 /*
40  * ---------------------------------------------------------------------------------
41  *                            Constructor / Destructor
42  * ---------------------------------------------------------------------------------
43  */
44   
45 // #################################################################################
46 AliHLTJETConeGrid::AliHLTJETConeGrid()
47   : 
48   fGrid(NULL),
49   fEtaMin(-0.9),
50   fEtaMax(0.9),
51   fPhiMin(0.0),
52   fPhiMax(6.3),
53   fEtaRange(1.8),
54   fPhiRange(6.3),
55   fEtaBinning(0.05),
56   fPhiBinning(0.05),
57   fEtaNGridBins(-1),
58   fPhiNGridBins(-1),
59   fNBins(-1),
60   fEtaNRBins(-1),
61   fPhiNRBins(-1),
62   fEtaIdxCurrent(0),
63   fEtaIdxMin(0),
64   fEtaIdxMax(0),
65   fPhiIdxCurrent(0),
66   fPhiIdxMin(0),
67   fPhiIdxMax(0),
68   fConeRadius(0.0) {
69   // see header file for class documentation
70   // or
71   // refer to README to build package
72   // or
73   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
74
75 }
76
77 // #################################################################################
78 AliHLTJETConeGrid::~AliHLTJETConeGrid() {
79   // see header file for class documentation
80  
81   if ( fGrid ) {
82     fGrid->Clear("C");
83     delete fGrid;
84   }
85   fGrid = NULL;
86  
87 }
88
89 /*
90  * ---------------------------------------------------------------------------------
91  *                                   Setup / Reset
92  * ---------------------------------------------------------------------------------
93  */
94
95 // #################################################################################
96 Int_t AliHLTJETConeGrid::Initialize() {
97   // see header file for class documentation
98
99   Int_t iResult = 0;
100
101   // -- Set total N bins in eta and phi
102   fEtaNGridBins = TMath::CeilNint( fEtaRange / fEtaBinning );
103   fPhiNGridBins = TMath::CeilNint( fPhiRange / fPhiBinning );
104
105   // -- Set total number of bins
106   fNBins = fEtaNGridBins * fPhiNGridBins;
107
108   // -- Set cells around R in eta and phi
109   fEtaNRBins = TMath::CeilNint( fConeRadius / fEtaBinning );
110   fPhiNRBins = TMath::CeilNint( fConeRadius / fPhiBinning );
111
112   HLTInfo(" -= Grid =- ");
113   HLTInfo(" NGridBins (%d,%d)", fEtaNGridBins, fPhiNGridBins);
114   HLTInfo(" NRBins    (%d,%d)", fEtaNRBins   , fPhiNRBins);
115   HLTInfo(" NBins      %d", fNBins );
116
117   fGrid = new TClonesArray("AliHLTJETConeEtaPhiCell", fNBins );
118   
119   if ( ! fGrid ) {
120     HLTError( "Error: Setup search grid with size %d .", fNBins );
121     iResult = 1;
122   }
123
124   return iResult;
125 }
126
127 // #################################################################################
128 void AliHLTJETConeGrid::Reset() { 
129   // see header file for class documentation
130
131   if ( fGrid )
132     fGrid->Clear("C");
133
134   return;
135 }
136
137 /*
138  * ---------------------------------------------------------------------------------
139  *                                    Process
140  * ---------------------------------------------------------------------------------
141  */
142
143 //##################################################################################
144 Int_t AliHLTJETConeGrid::FillTrack( TParticle* particle, const Float_t* aEtaPhi, Int_t* aGridIdx ) {
145   // see header file for class documentation
146
147   Int_t iResult = 0;
148
149   // ---------------------------
150   // -- Get Cell Indices
151   // ---------------------------
152   iResult = GetCellIndex( aEtaPhi, aGridIdx );
153   if ( iResult < 0 ) {
154     HLTError("Error getting cell index");
155     return iResult;
156   }
157
158   // ---------------------------
159   // -- Fill track in primary region
160   // ---------------------------
161   
162   // -- Create new cell and add track to cell
163   if (! fGrid->UncheckedAt( aGridIdx[kIdxPrimary] ) ) {
164     new( (*fGrid ) [aGridIdx[kIdxPrimary]] ) AliHLTJETConeEtaPhiCell( aGridIdx[kIdxEtaPrimary], 
165                                                                       aGridIdx[kIdxPhiPrimary], 
166                                                                       particle );
167   }
168   // -- Add track to existing cell
169   else {
170     (reinterpret_cast<AliHLTJETConeEtaPhiCell*> ((*fGrid)[aGridIdx[kIdxPrimary]]))->AddTrack(particle);
171   }
172
173   // ---------------------------
174   // -- Fill track in outter region
175   // ---------------------------
176   
177   // -- if it has to be filled
178   if ( iResult == 1 ) {
179
180     // -- Create new cell and add track to cell
181     if (! fGrid->UncheckedAt( aGridIdx[kIdxOutter] ) ) {
182       new( (*fGrid) [aGridIdx[kIdxOutter]] ) AliHLTJETConeEtaPhiCell( aGridIdx[kIdxEtaPrimary], 
183                                                                       aGridIdx[kIdxPhiOutter], 
184                                                                       particle );
185     }
186     // -- Add track to existing cell
187     else {
188       (reinterpret_cast<AliHLTJETConeEtaPhiCell*> ((*fGrid)[aGridIdx[kIdxOutter]]))->AddTrack(particle);
189     }
190   }
191
192   return 0;
193 }
194
195 //##################################################################################
196 Int_t AliHLTJETConeGrid::FillTrack( AliESDtrack *esdTrack, const Float_t* aEtaPhi, Int_t* aGridIdx ) {
197   // see header file for class documentation
198
199   Int_t iResult = 0;
200  
201   // ---------------------------
202   // -- Get Cell Indices
203   // ---------------------------
204
205   iResult = GetCellIndex( aEtaPhi, aGridIdx );
206   if ( iResult < 0 ) {
207     return iResult;
208   }
209     
210   // ---------------------------
211   // -- Fill track in primary region
212   // ---------------------------
213   
214   // -- Create new cell and add track to cell
215   if (! fGrid->UncheckedAt( aGridIdx[kIdxPrimary] ) ) {
216     new( (*fGrid ) [aGridIdx[kIdxPrimary]] ) AliHLTJETConeEtaPhiCell( aGridIdx[kIdxEtaPrimary], 
217                                                                       aGridIdx[kIdxPhiPrimary], 
218                                                                       esdTrack );
219   }
220   // -- Add track to existing cell
221   else {
222     (reinterpret_cast<AliHLTJETConeEtaPhiCell*> ((*fGrid)[aGridIdx[kIdxPrimary]]))->AddTrack(esdTrack);
223   }
224    
225   // ---------------------------
226   // -- Fill track in outter region
227   // ---------------------------
228   
229   // -- if it has to be filled
230   if ( iResult == 1 ) {
231     
232     // -- Create new cell and add track to cell
233     if (! fGrid->UncheckedAt( aGridIdx[kIdxOutter] ) ) {
234       new( (*fGrid) [aGridIdx[kIdxOutter]] ) AliHLTJETConeEtaPhiCell( aGridIdx[kIdxEtaPrimary], 
235                                                                       aGridIdx[kIdxPhiOutter], 
236                                                                       esdTrack );
237     }
238     // -- Add track to existing cell
239     else {
240       (reinterpret_cast<AliHLTJETConeEtaPhiCell*> ((*fGrid)[aGridIdx[kIdxOutter]]))->AddTrack(esdTrack);
241     }
242   }
243   
244   return 0;
245 }
246
247 /*
248  * ---------------------------------------------------------------------------------
249  *                             Helper - public
250  * ---------------------------------------------------------------------------------
251  */
252
253 // #################################################################################
254 Int_t AliHLTJETConeGrid::NextCell() {
255   // see header file for class documentation
256
257   Int_t cellIdx = 0;
258
259   ++fPhiIdxCurrent;
260   
261   if ( fPhiIdxCurrent > fPhiIdxMax ) {
262     fPhiIdxCurrent = fPhiIdxMin;
263
264     ++fEtaIdxCurrent;
265     
266     if ( fEtaIdxCurrent > fEtaIdxMax )
267       cellIdx = -1;
268   }
269
270   if ( cellIdx != -1 ) 
271     cellIdx = fEtaIdxCurrent + ( fPhiIdxCurrent * fEtaNGridBins );
272
273   if ( cellIdx > fNBins ) {
274     HLTError("Idx out of bound (%d,%d) -> %d", fEtaIdxCurrent, fPhiIdxCurrent, cellIdx );
275     HLTError("MAX %d,%d - %d", fEtaNGridBins, fPhiNGridBins, fNBins );
276
277     cellIdx = -1;
278   }
279
280   return cellIdx;
281 }
282
283 // #################################################################################
284 void AliHLTJETConeGrid::SetCellIter( const Int_t etaIdx, const Int_t phiIdx ) {
285   // see header file for class documentation
286
287   fEtaIdxMax = etaIdx + fEtaNRBins;
288   fEtaIdxMin = etaIdx - fEtaNRBins;
289   fEtaIdxCurrent = fEtaIdxMin;
290
291   fPhiIdxMax = phiIdx + fPhiNRBins;
292   fPhiIdxMin = phiIdx - fPhiNRBins;
293   fPhiIdxCurrent = fPhiIdxMin - 1;
294
295   return;
296 }
297
298 /*
299  * ---------------------------------------------------------------------------------
300  *                             Helper - private
301  * ---------------------------------------------------------------------------------
302  */
303
304 //##################################################################################
305 Int_t AliHLTJETConeGrid::GetCellIndex( const Float_t* aEtaPhi, Int_t* aGridIdx ) {
306   // see header file for class documentation
307
308   Int_t iResult = 0;
309
310   // -- Prime is relative to (0,0) in grid which is (-fEtaMax,-fConeRadius)
311   Float_t etaPrime = fEtaMax + aEtaPhi[kIdxEta];
312   Float_t phiPrime = aEtaPhi[kIdxPhi] + fConeRadius;
313
314   HLTDebug("eta    : %f - phi :%f", aEtaPhi[kIdxEta], aEtaPhi[kIdxPhi] );
315   HLTDebug("eta\'   : %f - phi\':%f", etaPrime, phiPrime );
316
317   // -- Index in 2D (idxEta,idxPhi) 
318   aGridIdx[kIdxEtaPrimary] = TMath::FloorNint( etaPrime / fEtaBinning );
319   aGridIdx[kIdxPhiPrimary] = TMath::FloorNint( phiPrime / fPhiBinning );
320  
321   // -- Index in 1D 
322   aGridIdx[kIdxPrimary] = aGridIdx[kIdxEtaPrimary] + ( aGridIdx[kIdxPhiPrimary] * fEtaNGridBins );
323   
324   // -- Boundery Check 2D
325   if ( (aGridIdx[kIdxEtaPrimary] < 0) || 
326        (aGridIdx[kIdxPhiPrimary] < 0) || 
327        (aGridIdx[kIdxEtaPrimary] >= fEtaNGridBins) || 
328        (aGridIdx[kIdxPhiPrimary] >= fPhiNGridBins) ) {
329     HLTError ( "Index out of range: idxEta %d - max %d, idxPhi %d - max %d", 
330                aGridIdx[kIdxEtaPrimary], fEtaNGridBins, 
331                aGridIdx[kIdxPhiPrimary], fPhiNGridBins);
332     iResult = -1;
333   }
334
335   // -- Boundery Check 1D
336   if ( (aGridIdx[kIdxPrimary] < 0) ||
337        (aGridIdx[kIdxPrimary] >= fNBins) ) {
338     HLTError( "Index out of range: 1D idx %d - max %d", 
339               aGridIdx[kIdxPrimary], fNBins );
340     iResult = -2;
341   }
342   
343   HLTDebug( "idxEta %d - max %d, idxPhi %d - max %d", 
344             aGridIdx[kIdxEtaPrimary], fEtaNGridBins, 
345             aGridIdx[kIdxPhiPrimary], fPhiNGridBins);
346   HLTDebug( "1D idx %d - max %d", aGridIdx[kIdxPrimary], fNBins );
347   
348   if ( iResult ) 
349     return iResult;
350
351   // -- check if to map from border region to outter region
352   Float_t phiOutterPrime = 0.0;
353   
354   // -- upper border
355   if ( aEtaPhi[kIdxPhi] > ( fPhiMax - fConeRadius ) ) {
356     phiOutterPrime = phiPrime - fPhiMax;
357     iResult = 1;  
358
359     HLTDebug("eta    : %f - phiOutter :%f .", aEtaPhi[kIdxEta], aEtaPhi[kIdxPhi] - fPhiMax );
360     HLTDebug("eta\'   : %f - phiOutter\':%f .", etaPrime, phiOutterPrime );
361   }
362   // -- lower border
363   else if ( aEtaPhi[kIdxPhi] < fConeRadius ) {
364     phiOutterPrime = phiPrime + fPhiMax;
365     iResult = 1;
366
367     HLTDebug("eta    : %f - phiOutter :%f .", aEtaPhi[kIdxEta], aEtaPhi[kIdxPhi] + fPhiMax );
368     HLTDebug("eta\'   : %f - phiOutter\':%f .", etaPrime, phiOutterPrime );
369   }
370   
371   // -- if outter phi present
372   if ( iResult == 1 ) {
373
374     // -- Index in 2D (idxEta,idxPhiOutter) 
375     aGridIdx[kIdxPhiOutter] = TMath::FloorNint( phiOutterPrime / fPhiBinning );
376
377     // -- Index in 1D
378     aGridIdx[kIdxOutter] = aGridIdx[kIdxEtaPrimary] + ( aGridIdx[kIdxPhiOutter] * fEtaNGridBins );
379   
380     // -- Boundery Check 2D
381     if ( aGridIdx[kIdxPhiOutter] < 0 || aGridIdx[kIdxPhiOutter] >= fPhiNGridBins ) {
382       HLTError( "Index out of range ( Outter Phi ): idxPhiOutter %d - max %d", 
383                 aGridIdx[kIdxPhiOutter], fPhiNGridBins);
384       iResult = -3;
385     }
386     
387     // -- Boundery Check 1D
388     if ( aGridIdx[kIdxOutter] >= fNBins ) {
389       HLTError( "Index out of range ( Outter Phi ): 1D idx %d - max %d", 
390                 aGridIdx[kIdxOutter], fNBins );
391       iResult = -4;
392     }
393
394     HLTDebug( "idxEta %d - max %d, idxPhiOutter %d - max %d", 
395              aGridIdx[kIdxEtaPrimary], fEtaNGridBins, aGridIdx[kIdxPhiOutter], fPhiNGridBins);
396     HLTDebug( "1D idx %d - max %d", aGridIdx[kIdxOutter], fNBins );
397   }
398
399   return iResult;
400 }