e21c56114890f7b1c1e65cef023ecaed74cd936d
[u/mrichter/AliRoot.git] / HLT / CALO / AliHLTCaloClusterizer.cxx
1 // $Id$\r
2 \r
3 /**************************************************************************\r
4  * This file is property of and copyright by the ALICE HLT Project        * \r
5  * All rights reserved.                                                   *\r
6  *                                                                        *\r
7  * Primary Authors: Oystein Djuvsland                                     *\r
8  *                                                                        *\r
9  * Permission to use, copy, modify and distribute this software and its   *\r
10  * documentation strictly for non-commercial purposes is hereby granted   *\r
11  * without fee, provided that the above copyright notice appears in all   *\r
12  * copies and that both the copyright notice and this permission notice   *\r
13  * appear in the supporting documentation. The authors make no claims     *\r
14  * about the suitability of this software for any purpose. It is          * \r
15  * provided "as is" without express or implied warranty.                  *\r
16  **************************************************************************/\r
17 \r
18 /** \r
19  * @file   AliHLTCaloClusterizer.cxx\r
20  * @author Oystein Djuvsland\r
21  * @date \r
22  * @brief  Clusterizer for PHOS HLT \r
23  */\r
24 \r
25 // see header file for class documentation\r
26 // or\r
27 // refer to README to build package\r
28 // or\r
29 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt\r
30 \r
31 #include "AliHLTCaloClusterizer.h"\r
32 #include "AliHLTLogging.h"\r
33 #include "TMath.h"\r
34 #include "AliHLTCaloRecPointDataStruct.h"\r
35 #include "AliHLTCaloDigitDataStruct.h"\r
36 #include "AliHLTCaloDigitContainerDataStruct.h"\r
37 #include "AliHLTCaloConstantsHandler.h"\r
38 \r
39 ClassImp(AliHLTCaloClusterizer);\r
40 \r
41 AliHLTCaloClusterizer::AliHLTCaloClusterizer(TString det):\r
42   AliHLTCaloConstantsHandler(det),\r
43   fRecPointArray(0),\r
44   fRecPointDataPtr(0),\r
45   fFirstRecPointPtr(0),\r
46   fArraySize(0),\r
47   fAvailableSize(0),\r
48   fUsedSize(0),\r
49   fNRecPoints(0),\r
50   fDigitIndexPtr(0),\r
51   fEmcClusteringThreshold(0),\r
52   fEmcMinEnergyThreshold(0),\r
53   fEmcTimeGate(0),\r
54   fDigitsInCluster(0),\r
55   fDigitsPointerArray(0),\r
56   fDigitContainerPtr(0),\r
57   fMaxDigitIndexDiff(0),\r
58   fNDigits(0)\r
59 {\r
60   //See header file for documentation\r
61   fEmcClusteringThreshold = 0.2;\r
62   fEmcMinEnergyThreshold = 0.03;\r
63   fEmcTimeGate = 1.e-6 ;\r
64   \r
65   fMaxDigitIndexDiff = 2*fCaloConstants->GetNZROWSMOD();\r
66   \r
67   \r
68   fArraySize = 10;\r
69   fRecPointArray = new AliHLTCaloRecPointDataStruct*[fArraySize];\r
70   \r
71   fAvailableSize = sizeof(AliHLTCaloRecPointDataStruct) * 20;\r
72   fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(new UChar_t[fAvailableSize]);\r
73   fFirstRecPointPtr = fRecPointDataPtr;  \r
74 \r
75 }//end\r
76 \r
77 AliHLTCaloClusterizer::~AliHLTCaloClusterizer()  \r
78 {\r
79   //See header file for documentation\r
80 }\r
81 \r
82 void \r
83 AliHLTCaloClusterizer::SetRecPointDataPtr(AliHLTCaloRecPointDataStruct* recPointDataPtr)\r
84 {\r
85   // See header file for documentation\r
86   fRecPointDataPtr = recPointDataPtr;\r
87 }\r
88 \r
89 Int_t \r
90 AliHLTCaloClusterizer::ClusterizeEvent(Int_t nDigits)\r
91 {\r
92   //see header file for documentation\r
93   Int_t nRecPoints = 0;\r
94   fNRecPoints = 0;\r
95   fUsedSize = 0;\r
96   fNDigits = nDigits;\r
97   fRecPointDataPtr = fFirstRecPointPtr;\r
98   //Clusterization starts\r
99   for(Int_t i = 0; i < nDigits; i++)\r
100     { \r
101       fDigitsInCluster = 0;\r
102 \r
103       if(fDigitsPointerArray[i]->fEnergy < fEmcClusteringThreshold)\r
104         {\r
105           continue;\r
106         }\r
107       CheckArray();\r
108       CheckBuffer();\r
109 \r
110       // First digit is placed at the fDigits member variable in the recpoint\r
111       fDigitIndexPtr = &(fRecPointDataPtr->fDigits);\r
112       fUsedSize += sizeof(AliHLTCaloRecPointDataStruct);\r
113       \r
114       fRecPointDataPtr->fAmp = 0;\r
115       fRecPointDataPtr->fModule = fDigitsPointerArray[i]->fModule;\r
116 \r
117       // Assigning the digit to this rec point\r
118       fRecPointDataPtr->fDigits = i;\r
119       fUsedSize += sizeof(AliHLTCaloRecPointDataStruct);\r
120       \r
121       // Incrementing the pointer to be ready for new entry\r
122       fDigitIndexPtr++;\r
123 \r
124       fRecPointDataPtr->fAmp += fDigitsPointerArray[i]->fEnergy;\r
125       fDigitsPointerArray[i]->fEnergy = 0;\r
126       fDigitsInCluster++;\r
127       nRecPoints++;\r
128 \r
129       // Scanning for the neighbours\r
130       if(ScanForNeighbourDigits(i, fRecPointDataPtr) != 0)\r
131       {\r
132          return -1;\r
133       }\r
134 \r
135       //fUsedSize += sizeof(AliHLTCaloRecPointDataStruct) + (fDigitsInCluster-1)*sizeof(AliHLTCaloDigitDataStruct);   \r
136       \r
137       fRecPointDataPtr->fMultiplicity = fDigitsInCluster;     \r
138       fRecPointArray[fNRecPoints] = fRecPointDataPtr; \r
139       \r
140       fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fDigitIndexPtr);\r
141 \r
142       fNRecPoints++;\r
143       \r
144     }//end of clusterization\r
145 \r
146    return nRecPoints;\r
147 }\r
148 \r
149 Int_t\r
150 AliHLTCaloClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTCaloRecPointDataStruct* recPoint)\r
151 {\r
152   //see header file for documentation\r
153   Int_t max = TMath::Min(fNDigits, (Int_t)fMaxDigitIndexDiff+index);\r
154   Int_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff));\r
155 \r
156   max = fNDigits;\r
157   min = 0;\r
158   for(Int_t j = min; j < max; j++)\r
159     {\r
160       if(fDigitsPointerArray[j]->fEnergy > fEmcMinEnergyThreshold)\r
161         {\r
162           if(j != index)\r
163             {\r
164               if(AreNeighbours(fDigitsPointerArray[index],\r
165                                fDigitsPointerArray[j]))\r
166                 {\r
167 //                if((fAvailableSize - fUsedSize) < sizeof(Int_t))\r
168 //                   {\r
169 //                      UChar_t *tmp = new UChar_t[fAvailableSize*2];\r
170 //                      memcpy(tmp, fRecPointDataPtr, fAvailableSize);\r
171 //                      for(Int_t n = 0; n < fNRecPoints; n++)\r
172 //                      {\r
173 //                         fRecPointArray[n] = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(reinterpret_cast<UChar_t*>(fRecPointArray[n]) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr) + reinterpret_cast<UChar_t*>(tmp));\r
174 //                      }\r
175 //                      fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp);\r
176 //                      fFirstRecPointPtr = fRecPointDataPtr;\r
177 //                      fUsedSize = 0;\r
178 //                   }  \r
179                   CheckBuffer();\r
180                   // Assigning index to digit\r
181                   *fDigitIndexPtr = j;\r
182                   fUsedSize += sizeof(Int_t);\r
183                   \r
184                   // Incrementing digit pointer to be ready for new entry\r
185                   fDigitIndexPtr++;\r
186 \r
187                   recPoint->fAmp += fDigitsPointerArray[j]->fEnergy;\r
188                   fDigitsPointerArray[j]->fEnergy = 0;        \r
189                   fDigitsInCluster++;\r
190                   ScanForNeighbourDigits(j, recPoint);\r
191                 }\r
192             }\r
193         }\r
194     }\r
195   return 0;\r
196 }\r
197 \r
198 Int_t \r
199 AliHLTCaloClusterizer::AreNeighbours(AliHLTCaloDigitDataStruct* digit1, \r
200                                             AliHLTCaloDigitDataStruct* digit2)\r
201 {\r
202   //see header file for documentation\r
203   if ( (digit1->fModule == digit2->fModule) /*&& (coord1[1]==coord2[1])*/ ) // inside the same PHOS module\r
204     { \r
205       Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );  \r
206       Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX ); \r
207       if (( coldiff <= 1   &&  rowdiff == 0 ) || ( coldiff == 0 &&  rowdiff <= 1 ))\r
208         {\r
209 //        cout << "Are neighbours: digit (E = "  << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r
210 //          " is neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r
211 \r
212           if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r
213             {\r
214               return 1; \r
215             }\r
216         }\r
217 \r
218    /*   Float_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );  \r
219       Float_t coldiff = TMath::Abs( digit1->fX - digit2->fX ); \r
220       if (( coldiff <= 2.4   &&  rowdiff < 0.4 ) || ( coldiff < 0.4 &&  rowdiff <= 2.4 ))\r
221         {\r
222           //      cout << "Are neighbours: digit (E = "  << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r
223           //        " is neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r
224 \r
225           if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r
226             {\r
227               return 1; \r
228             }\r
229         }\r
230    */\r
231       else\r
232         {\r
233 //          cout << "Not neighbours: digit (E = "  << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r
234 //          " is not neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r
235         }\r
236     }\r
237   return 0;\r
238 }\r
239 \r
240 \r
241 \r
242 Int_t AliHLTCaloClusterizer::CheckArray()\r
243 {\r
244       if(fArraySize == fNRecPoints)\r
245         {\r
246            fArraySize *= 2;\r
247            AliHLTCaloRecPointDataStruct **tmp = new AliHLTCaloRecPointDataStruct*[fArraySize];\r
248            memcpy(tmp, fRecPointArray, fArraySize/2 * sizeof(AliHLTCaloRecPointDataStruct*));\r
249            delete fRecPointArray;\r
250            fRecPointArray = tmp;\r
251         }\r
252    return 0;\r
253 }\r
254 \r
255 Int_t AliHLTCaloClusterizer::CheckBuffer()\r
256 {\r
257    // See header file for class documentation \r
258         if((fAvailableSize - fUsedSize) < sizeof(AliHLTCaloRecPointDataStruct))\r
259         {\r
260             Int_t recPointOffset = reinterpret_cast<UChar_t*>(fRecPointDataPtr) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr);\r
261             fAvailableSize *= 2;\r
262             UChar_t *tmp = new UChar_t[fAvailableSize];\r
263             memcpy(tmp, fRecPointDataPtr, fAvailableSize/2);\r
264             for(Int_t n = 0; n < fNRecPoints; n++)\r
265             {\r
266                fRecPointArray[n] = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(reinterpret_cast<UChar_t*>(fRecPointArray[n]) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr) + reinterpret_cast<UChar_t*>(tmp));\r
267             }\r
268             delete fRecPointDataPtr;\r
269             fFirstRecPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp);\r
270             fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp + recPointOffset);\r
271             fUsedSize = 0;\r
272         }\r
273    return 0;\r
274 }\r
275 \r
276 Int_t AliHLTCaloClusterizer::CheckDigits(AliHLTCaloRecPointDataStruct** recArray, AliHLTCaloDigitDataStruct** digitArray, Int_t nRP)\r
277 {\r
278   AliHLTCaloRecPointDataStruct **recpoints = recArray;\r
279   AliHLTCaloDigitDataStruct **digits = digitArray;\r
280   Int_t nRecPoints = nRP;\r
281   \r
282   if(recArray == 0)\r
283   {\r
284      recpoints = fRecPointArray;\r
285   }\r
286   if(digitArray == 0)\r
287   {\r
288      digits = fDigitsPointerArray;\r
289   }\r
290   if(nRP == 0)\r
291   {\r
292      nRecPoints = fNRecPoints;\r
293   }\r
294   for(Int_t i = 0; i < nRecPoints; i++)\r
295   {\r
296           \r
297      AliHLTCaloRecPointDataStruct *recPoint = recpoints[i];\r
298 \r
299      //AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[0];\r
300      Int_t multiplicity = recPoint->fMultiplicity;\r
301      Int_t *digitIndexPtr = &(recPoint->fDigits);\r
302      for(Int_t j = 0; j < multiplicity; j++)\r
303      {\r
304         //AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[j];\r
305 //      AliHLTCaloDigitDataStruct *digit = digits[*digitIndexPtr];\r
306         digitIndexPtr++;\r
307         //recPoint = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(digitIndexPtr);\r
308      }\r
309   }\r
310      \r
311      return 0;\r
312 }\r
313 \r
314 Int_t AliHLTCaloClusterizer::CheckDigits(AliHLTCaloRecPointDataStruct** recArray, AliHLTCaloDigitDataStruct* digitArray, Int_t nRP)\r
315 {\r
316   AliHLTCaloRecPointDataStruct **recpoints = recArray;\r
317   AliHLTCaloDigitDataStruct *digits = digitArray;\r
318   Int_t nRecPoints = nRP;\r
319   \r
320   if(recArray == 0)\r
321   {\r
322      recpoints = fRecPointArray;\r
323   }\r
324     if(nRP == 0)\r
325   {\r
326      nRecPoints = fNRecPoints;\r
327   }\r
328   for(Int_t i = 0; i < nRecPoints; i++)\r
329   {\r
330           \r
331      AliHLTCaloRecPointDataStruct *recPoint = recpoints[i];\r
332 \r
333      //AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[0];\r
334      Int_t multiplicity = recPoint->fMultiplicity;\r
335      Int_t *digitIndexPtr = &(recPoint->fDigits);\r
336      for(Int_t j = 0; j < multiplicity; j++)\r
337      {\r
338         //AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[j];\r
339         AliHLTCaloDigitDataStruct digit = digits[*digitIndexPtr];\r
340         digitIndexPtr++;\r
341         //recPoint = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(digitIndexPtr);\r
342      }\r
343   }\r
344    return 0;\r
345      \r
346 }