Events by default are selected.
[u/mrichter/AliRoot.git] / HLT / CALO / AliHLTCaloClusterizer.cxx
CommitLineData
ef44ec64 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
ef44ec64 32#include "AliHLTLogging.h"\r
33#include "TMath.h"\r
ef44ec64 34#include "AliHLTCaloRecPointDataStruct.h"\r
35#include "AliHLTCaloDigitDataStruct.h"\r
36#include "AliHLTCaloDigitContainerDataStruct.h"\r
4f4b7ba4 37#include "AliHLTCaloConstantsHandler.h"\r
ef44ec64 38\r
39ClassImp(AliHLTCaloClusterizer);\r
40\r
41AliHLTCaloClusterizer::AliHLTCaloClusterizer(TString det):\r
4f4b7ba4 42 AliHLTCaloConstantsHandler(det),\r
7c80a370 43 fRecPointArray(0),\r
ef44ec64 44 fRecPointDataPtr(0),\r
7c80a370 45 fFirstRecPointPtr(0),\r
46 fArraySize(0),\r
47 fAvailableSize(0),\r
48 fUsedSize(0),\r
49 fNRecPoints(0),\r
ad44d760 50 fDigitIndexPtr(0),\r
ef44ec64 51 fEmcClusteringThreshold(0),\r
52 fEmcMinEnergyThreshold(0),\r
53 fEmcTimeGate(0),\r
54 fDigitsInCluster(0),\r
f92dcf21 55 fDigitsPointerArray(0),\r
ef44ec64 56 fDigitContainerPtr(0),\r
ad44d760 57 fMaxDigitIndexDiff(0),\r
58 fNDigits(0)\r
ef44ec64 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
7c80a370 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
c375e15d 74\r
c375e15d 75}//end\r
76\r
ef44ec64 77AliHLTCaloClusterizer::~AliHLTCaloClusterizer() \r
78{\r
79 //See header file for documentation\r
80}\r
81\r
82void \r
83AliHLTCaloClusterizer::SetRecPointDataPtr(AliHLTCaloRecPointDataStruct* recPointDataPtr)\r
84{\r
85 // See header file for documentation\r
86 fRecPointDataPtr = recPointDataPtr;\r
87}\r
88\r
ef44ec64 89Int_t \r
7c80a370 90AliHLTCaloClusterizer::ClusterizeEvent(Int_t nDigits)\r
ef44ec64 91{\r
92 //see header file for documentation\r
93 Int_t nRecPoints = 0;\r
7c80a370 94 fNRecPoints = 0;\r
95 fUsedSize = 0;\r
ad44d760 96 fNDigits = nDigits;\r
7c80a370 97 fRecPointDataPtr = fFirstRecPointPtr;\r
ef44ec64 98 //Clusterization starts\r
f92dcf21 99 for(Int_t i = 0; i < nDigits; i++)\r
ef44ec64 100 { \r
101 fDigitsInCluster = 0;\r
98baf84d 102\r
ad44d760 103 if(fDigitsPointerArray[i]->fEnergy < fEmcClusteringThreshold)\r
ef44ec64 104 {\r
105 continue;\r
106 }\r
7c80a370 107 CheckArray();\r
108 CheckBuffer();\r
98baf84d 109\r
ef44ec64 110 // First digit is placed at the fDigits member variable in the recpoint\r
ad44d760 111 fDigitIndexPtr = &(fRecPointDataPtr->fDigits);\r
b4479a87 112 fUsedSize += sizeof(AliHLTCaloRecPointDataStruct);\r
113 \r
ef44ec64 114 fRecPointDataPtr->fAmp = 0;\r
ad44d760 115 fRecPointDataPtr->fModule = fDigitsPointerArray[i]->fModule;\r
ef44ec64 116\r
ad44d760 117 // Assigning the digit to this rec point\r
118 fRecPointDataPtr->fDigits = i;\r
7c80a370 119 fUsedSize += sizeof(AliHLTCaloRecPointDataStruct);\r
120 \r
ef44ec64 121 // Incrementing the pointer to be ready for new entry\r
ad44d760 122 fDigitIndexPtr++;\r
ef44ec64 123\r
ad44d760 124 fRecPointDataPtr->fAmp += fDigitsPointerArray[i]->fEnergy;\r
125 fDigitsPointerArray[i]->fEnergy = 0;\r
ef44ec64 126 fDigitsInCluster++;\r
127 nRecPoints++;\r
128\r
129 // Scanning for the neighbours\r
9bd6f29b 130 if(ScanForNeighbourDigits(i, fRecPointDataPtr) != 0)\r
131 {\r
132 return -1;\r
133 }\r
ef44ec64 134\r
7c80a370 135 //fUsedSize += sizeof(AliHLTCaloRecPointDataStruct) + (fDigitsInCluster-1)*sizeof(AliHLTCaloDigitDataStruct); \r
136 \r
ef44ec64 137 fRecPointDataPtr->fMultiplicity = fDigitsInCluster; \r
7c80a370 138 fRecPointArray[fNRecPoints] = fRecPointDataPtr; \r
139 \r
ad44d760 140 fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fDigitIndexPtr);\r
b4479a87 141\r
142 fNRecPoints++;\r
7c80a370 143 \r
ef44ec64 144 }//end of clusterization\r
b4479a87 145\r
ef44ec64 146 return nRecPoints;\r
147}\r
148\r
9bd6f29b 149Int_t\r
ef44ec64 150AliHLTCaloClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTCaloRecPointDataStruct* recPoint)\r
151{\r
152 //see header file for documentation\r
ad44d760 153 Int_t max = TMath::Min(fNDigits, (Int_t)fMaxDigitIndexDiff+index);\r
ef44ec64 154 Int_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff));\r
155\r
ad44d760 156 max = fNDigits;\r
ef44ec64 157 min = 0;\r
158 for(Int_t j = min; j < max; j++)\r
159 {\r
ad44d760 160 if(fDigitsPointerArray[j]->fEnergy > fEmcMinEnergyThreshold)\r
ef44ec64 161 {\r
162 if(j != index)\r
163 {\r
ad44d760 164 if(AreNeighbours(fDigitsPointerArray[index],\r
165 fDigitsPointerArray[j]))\r
ef44ec64 166 {\r
7c80a370 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
9bd6f29b 180 // Assigning index to digit\r
ad44d760 181 *fDigitIndexPtr = j;\r
7c80a370 182 fUsedSize += sizeof(Int_t);\r
183 \r
ef44ec64 184 // Incrementing digit pointer to be ready for new entry\r
ad44d760 185 fDigitIndexPtr++;\r
ef44ec64 186\r
ad44d760 187 recPoint->fAmp += fDigitsPointerArray[j]->fEnergy;\r
188 fDigitsPointerArray[j]->fEnergy = 0; \r
ef44ec64 189 fDigitsInCluster++;\r
190 ScanForNeighbourDigits(j, recPoint);\r
191 }\r
192 }\r
193 }\r
194 }\r
9bd6f29b 195 return 0;\r
ef44ec64 196}\r
197\r
198Int_t \r
199AliHLTCaloClusterizer::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
7c80a370 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
ef44ec64 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
7c80a370 230 */\r
ef44ec64 231 else\r
232 {\r
7c80a370 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
ef44ec64 235 }\r
236 }\r
237 return 0;\r
238}\r
7c80a370 239\r
240\r
241\r
242Int_t AliHLTCaloClusterizer::CheckArray()\r
243{\r
7c80a370 244 if(fArraySize == fNRecPoints)\r
245 {\r
7c80a370 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
255Int_t AliHLTCaloClusterizer::CheckBuffer()\r
256{\r
257 // See header file for class documentation \r
98baf84d 258 if((fAvailableSize - fUsedSize) < sizeof(AliHLTCaloRecPointDataStruct))\r
7c80a370 259 {\r
7c80a370 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
276Int_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
7c80a370 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
7c80a370 302 for(Int_t j = 0; j < multiplicity; j++)\r
303 {\r
304 //AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[j];\r
98baf84d 305// AliHLTCaloDigitDataStruct *digit = digits[*digitIndexPtr];\r
7c80a370 306 digitIndexPtr++;\r
307 //recPoint = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(digitIndexPtr);\r
308 }\r
309 }\r
310 \r
98baf84d 311 return 0;\r
7c80a370 312}\r
313\r
314Int_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
7c80a370 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
7c80a370 336 for(Int_t j = 0; j < multiplicity; j++)\r
337 {\r
338 //AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[j];\r
339 AliHLTCaloDigitDataStruct digit = digits[*digitIndexPtr];\r
7c80a370 340 digitIndexPtr++;\r
341 //recPoint = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(digitIndexPtr);\r
342 }\r
343 }\r
98baf84d 344 return 0;\r
7c80a370 345 \r
346}