reverting r45444 to disentangle modules and make porting possible
[u/mrichter/AliRoot.git] / HLT / CALO / AliHLTCaloClusterizer.cxx
CommitLineData
ef44ec64 1// $Id$\r
2\r
3/**************************************************************************\r
31b89da4 4 * This file is property of and copyright by the ALICE HLT Project *\r
ef44ec64 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
31b89da4 14 * about the suitability of this software for any purpose. It is *\r
ef44ec64 15 * provided "as is" without express or implied warranty. *\r
16 **************************************************************************/\r
17\r
31b89da4 18/**\r
ef44ec64 19 * @file AliHLTCaloClusterizer.cxx\r
20 * @author Oystein Djuvsland\r
31b89da4 21 * @date\r
22 * @brief Clusterizer for PHOS HLT\r
ef44ec64 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
31b89da4 42 AliHLTCaloConstantsHandler(det),\r
43 fCompareFunction(CompareDigitsByPosition),\r
44 fRecPointArray(0),\r
45 fRecPointDataPtr(0),\r
46 fFirstRecPointPtr(0),\r
47 fArraySize(0),\r
48 fAvailableSize(0),\r
49 fUsedSize(0),\r
50 fNRecPoints(0),\r
51 fDigitIndexPtr(0),\r
52 fEmcClusteringThreshold(0),\r
53 fEmcMinEnergyThreshold(0),\r
54 fEmcTimeGate(0),\r
55 fDigitsInCluster(0),\r
56 fDigitsPointerArray(0),\r
57 fDigitContainerPtr(0),\r
58 fMaxDigitIndexDiff(0),\r
59 fNDigits(0),\r
60 fSortedByPosition(false),\r
61 fSortedByEnergy(false),\r
62 fSortDigits(false)\r
ef44ec64 63{\r
31b89da4 64 //See header file for documentation\r
65 //fEmcClusteringThreshold = 0.2;\r
66 //fEmcMinEnergyThreshold = 0.03;\r
67\r
68 fEmcClusteringThreshold = 0.1;\r
69 fEmcMinEnergyThreshold = 0.01;\r
70 fEmcTimeGate = 1.e-6 ;\r
71\r
72 fMaxDigitIndexDiff = 2*fCaloConstants->GetNZROWSMOD();\r
73\r
74\r
75 fArraySize = 10;\r
76 fRecPointArray = new AliHLTCaloRecPointDataStruct*[fArraySize];\r
77\r
78 fAvailableSize = sizeof(AliHLTCaloRecPointDataStruct) * 20;\r
79 fFirstRecPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(new UChar_t[fAvailableSize]);\r
80 fRecPointDataPtr = fFirstRecPointPtr;\r
c375e15d 81\r
c375e15d 82}//end\r
83\r
31b89da4 84AliHLTCaloClusterizer::~AliHLTCaloClusterizer()\r
ef44ec64 85{\r
31b89da4 86 //See header file for documentation\r
ef44ec64 87}\r
88\r
31b89da4 89void\r
ef44ec64 90AliHLTCaloClusterizer::SetRecPointDataPtr(AliHLTCaloRecPointDataStruct* recPointDataPtr)\r
91{\r
31b89da4 92 // See header file for documentation\r
93 fRecPointDataPtr = recPointDataPtr;\r
ef44ec64 94}\r
95\r
31b89da4 96Int_t\r
7c80a370 97AliHLTCaloClusterizer::ClusterizeEvent(Int_t nDigits)\r
ef44ec64 98{\r
31b89da4 99 //see header file for documentation\r
100 Int_t nRecPoints = 0;\r
101 fNRecPoints = 0;\r
102 fUsedSize = 0;\r
103 fNDigits = nDigits;\r
104 fRecPointDataPtr = fFirstRecPointPtr;\r
105\r
106 // Sort our digits\r
107 SortDigits();\r
108\r
109 //Clusterization starts\r
110 for (Int_t i = 0; i < nDigits; i++)\r
111 {\r
112 fDigitsInCluster = 0;\r
113\r
f92543f5 114 HLTDebug("Digit with energy: %f", fDigitsPointerArray[i]->fEnergy);\r
115 \r
116 if (fDigitsPointerArray[i]->fEnergy < fEmcClusteringThreshold && fSortedByEnergy)\r
31b89da4 117 {\r
f92543f5 118 // Since we have sorted by energy the next digit will have even lower energy, so we return \r
119 return fNRecPoints;\r
120 }\r
121\r
122 if(fDigitsPointerArray[i]->fAssociatedCluster != -1)\r
123 {\r
124 // The digit is added to a previous cluster, continue\r
125 continue;\r
126 }\r
31b89da4 127\r
f92543f5 128 CheckArray();\r
31b89da4 129 CheckBuffer();\r
130\r
131 // First digit is placed at the fDigits member variable in the recpoint\r
132 fDigitIndexPtr = &(fRecPointDataPtr->fDigits);\r
133\r
134 fRecPointDataPtr->fAmp = 0;\r
135 fRecPointDataPtr->fModule = fDigitsPointerArray[i]->fModule;\r
136\r
137 // Assigning the digit to this rec point\r
138 fRecPointDataPtr->fDigits = i;\r
139 fUsedSize += sizeof(AliHLTCaloRecPointDataStruct);\r
140\r
141 // Incrementing the pointer to be ready for new entry\r
142 fDigitIndexPtr++;\r
143\r
144 fRecPointDataPtr->fAmp += fDigitsPointerArray[i]->fEnergy;\r
f92543f5 145 \r
146 \r
147 //fDigitsPointerArray[i]->fEnergy = 0;\r
148 fDigitsPointerArray[i]->fAssociatedCluster = fNRecPoints;\r
149 \r
150 \r
151 fDigitsInCluster++;\r
31b89da4 152 nRecPoints++;\r
153\r
154 // Scanning for the neighbours\r
155 if (ScanForNeighbourDigits(i, fRecPointDataPtr) != 0)\r
156 {\r
157 return -1;\r
158 }\r
159\r
160 //fUsedSize += sizeof(AliHLTCaloRecPointDataStruct) + (fDigitsInCluster-1)*sizeof(AliHLTCaloDigitDataStruct);\r
161\r
162 fRecPointDataPtr->fMultiplicity = fDigitsInCluster;\r
163 fRecPointArray[fNRecPoints] = fRecPointDataPtr;\r
164\r
165 fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fDigitIndexPtr);\r
166\r
167 fNRecPoints++;\r
168\r
ef44ec64 169 }//end of clusterization\r
70c86cde 170\r
171 return nRecPoints;\r
ef44ec64 172}\r
173\r
9bd6f29b 174Int_t\r
ef44ec64 175AliHLTCaloClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTCaloRecPointDataStruct* recPoint)\r
176{\r
31b89da4 177 //see header file for documentation\r
ef44ec64 178\r
31b89da4 179 // The following cuts can be used if we sort by posisiton. Not tested, but it should be fine...\r
180 Int_t max = TMath::Min(fNDigits, (Int_t)fMaxDigitIndexDiff+index);\r
181 Int_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff));\r
182\r
183 // All digits for now\r
184 max = fNDigits;\r
185 min = 0;\r
186\r
187 for (Int_t j = min; j < max; j++)\r
ef44ec64 188 {\r
f92543f5 189 if (fDigitsPointerArray[j]->fAssociatedCluster == -1 && fDigitsPointerArray[j]->fEnergy > fEmcMinEnergyThreshold)\r
31b89da4 190 {\r
191 if (j != index)\r
192 {\r
193 if (AreNeighbours(fDigitsPointerArray[index],\r
194 fDigitsPointerArray[j]))\r
195 {\r
196 // Check that the buffer is large enough for adding a digit (can be heavily improved wrt performance)\r
197 CheckBuffer();\r
198\r
199 // Assigning index to digit\r
200 *fDigitIndexPtr = j;\r
201 fUsedSize += sizeof(Int_t);\r
202\r
203 // Incrementing digit pointer to be ready for new entry\r
204 fDigitIndexPtr++;\r
205\r
206 // Adding the digit energy to the rec point\r
207 fRecPointDataPtr->fAmp += fDigitsPointerArray[j]->fEnergy;\r
208\r
f92543f5 209 // Setting energy to 0\r
210 //fDigitsPointerArray[j]->fEnergy = 0;\r
211 \r
212 // Setting the associated cluster \r
213 fDigitsPointerArray[j]->fAssociatedCluster = fNRecPoints;\r
214 \r
215 HLTDebug("Added digit with index: %d, energy: %f, to associated cluster: %d", fDigitsPointerArray[j]->fID, fDigitsPointerArray[j]->fEnergy, fDigitsPointerArray[j]->fAssociatedCluster);\r
216 \r
31b89da4 217 fDigitsInCluster++;\r
218\r
219 // Scan for neighbours of this digit\r
220 ScanForNeighbourDigits(j, recPoint);\r
221 }\r
222 }\r
223 }\r
ef44ec64 224 }\r
31b89da4 225 return 0;\r
ef44ec64 226}\r
227\r
31b89da4 228Int_t\r
229AliHLTCaloClusterizer::AreNeighbours(AliHLTCaloDigitDataStruct* digit1,\r
230 AliHLTCaloDigitDataStruct* digit2)\r
ef44ec64 231{\r
31b89da4 232 //see header file for documentation\r
233 if ( (digit1->fModule == digit2->fModule) /*&& (coord1[1]==coord2[1])*/ ) // inside the same PHOS module\r
234 {\r
235 Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );\r
236 Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX );\r
237\r
238 // As in the offline code we define neighbours as cells that share an edge, a corner is not enough\r
c5688f76 239 // if (( coldiff <= 1 && rowdiff == 0 ) || ( coldiff == 0 && rowdiff <= 1 ))\r
240 if (( coldiff <= 1) || ( rowdiff <= 1 ))\r
31b89da4 241 {\r
242 // Check also for time\r
243 if (TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r
244 {\r
245 return 1;\r
246 }\r
247 }\r
ef44ec64 248 }\r
31b89da4 249 return 0;\r
ef44ec64 250}\r
7c80a370 251\r
252\r
253\r
254Int_t AliHLTCaloClusterizer::CheckArray()\r
255{\r
31b89da4 256 // See header file for class documentation\r
257 if (fArraySize == fNRecPoints)\r
258 {\r
259 fArraySize *= 2;\r
260 AliHLTCaloRecPointDataStruct **tmp = new AliHLTCaloRecPointDataStruct*[fArraySize];\r
261 memcpy(tmp, fRecPointArray, fArraySize/2 * sizeof(AliHLTCaloRecPointDataStruct*));\r
262 delete [] fRecPointArray;\r
263 fRecPointArray = tmp;\r
264 }\r
265 return 0;\r
7c80a370 266}\r
267\r
268Int_t AliHLTCaloClusterizer::CheckBuffer()\r
269{\r
31b89da4 270 // See header file for class documentation\r
271 if ((fAvailableSize - fUsedSize) < (Int_t)sizeof(AliHLTCaloRecPointDataStruct))\r
272 {\r
273 Int_t recPointOffset = reinterpret_cast<UChar_t*>(fRecPointDataPtr) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr);\r
274 Int_t digitIndexOffset = reinterpret_cast<UChar_t*>(fDigitIndexPtr) - reinterpret_cast<UChar_t*>(fRecPointDataPtr);\r
275 UChar_t *tmp = new UChar_t[fAvailableSize*2];\r
276\r
277 memcpy(tmp, fFirstRecPointPtr, fUsedSize);\r
278 fAvailableSize *= 2;\r
279 for (Int_t n = 0; n < fNRecPoints; n++)\r
280 {\r
281 fRecPointArray[n] = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(reinterpret_cast<UChar_t*>(fRecPointArray[n]) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr) + reinterpret_cast<UChar_t*>(tmp));\r
282 }\r
283 delete [] fFirstRecPointPtr;\r
284 fFirstRecPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp);\r
285 fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp + recPointOffset);\r
286 fDigitIndexPtr = reinterpret_cast<Int_t*>(reinterpret_cast<UChar_t*>(fRecPointDataPtr) + digitIndexOffset);\r
287 //fUsedSize = 0;\r
288 }\r
289 return 0;\r
7c80a370 290}\r
291\r
31b89da4 292void AliHLTCaloClusterizer::SetSortDigitsByPosition()\r
7c80a370 293{\r
31b89da4 294 // Sort the digit pointers by position\r
295 fCompareFunction = &CompareDigitsByPosition;\r
296 fSortDigits = true;\r
297 fSortedByPosition = true;\r
7c80a370 298}\r
299\r
31b89da4 300void AliHLTCaloClusterizer::SetSortDigitsByEnergy()\r
301{\r
302 // See header file for class documentation\r
303 fCompareFunction = &CompareDigitsByEnergy;\r
304 fSortDigits = true;\r
305 fSortedByEnergy = true;\r
306}\r
307\r
308void AliHLTCaloClusterizer::SortDigits()\r
309{\r
310 // See header file for class documentation\r
311 if (fSortDigits) qsort(fDigitsPointerArray, fNDigits, sizeof(AliHLTCaloDigitDataStruct*), fCompareFunction);\r
312}\r
313\r
314Int_t\r
315AliHLTCaloClusterizer::CompareDigitsByPosition(const void *dig0, const void *dig1)\r
316{\r
317 // See header file for documentation\r
318 return (*((AliHLTCaloDigitDataStruct**)(dig0)))->fID - (*((AliHLTCaloDigitDataStruct**)(dig1)))->fID;\r
319}\r
320\r
321Int_t\r
322AliHLTCaloClusterizer::CompareDigitsByEnergy(const void *dig0, const void *dig1)\r
7c80a370 323{\r
31b89da4 324 // See header file for documentation\r
80d83a65 325 return (*((AliHLTCaloDigitDataStruct**)(dig1)))->fEnergy - (*((AliHLTCaloDigitDataStruct**)(dig0)))->fEnergy;\r
2a24cbbe 326}\r