]>
Commit | Line | Data |
---|---|---|
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 |
39 | ClassImp(AliHLTCaloClusterizer);\r | |
40 | \r | |
41 | AliHLTCaloClusterizer::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 | |
7cc065a2 | 62 | fSortDigits(false),\r |
36b00077 | 63 | fIsEMCAL(false),\r |
7cc065a2 | 64 | fBuffer(0)\r |
36b00077 | 65 | \r |
ef44ec64 | 66 | {\r |
31b89da4 | 67 | //See header file for documentation\r |
68 | //fEmcClusteringThreshold = 0.2;\r | |
69 | //fEmcMinEnergyThreshold = 0.03;\r | |
70 | \r | |
71 | fEmcClusteringThreshold = 0.1;\r | |
72 | fEmcMinEnergyThreshold = 0.01;\r | |
73 | fEmcTimeGate = 1.e-6 ;\r | |
74 | \r | |
75 | fMaxDigitIndexDiff = 2*fCaloConstants->GetNZROWSMOD();\r | |
76 | \r | |
77 | \r | |
78 | fArraySize = 10;\r | |
79 | fRecPointArray = new AliHLTCaloRecPointDataStruct*[fArraySize];\r | |
80 | \r | |
81 | fAvailableSize = sizeof(AliHLTCaloRecPointDataStruct) * 20;\r | |
7cc065a2 | 82 | fBuffer = new UChar_t[fAvailableSize]; //FR\r |
dd71d443 | 83 | \r |
7cc065a2 | 84 | fFirstRecPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fBuffer);\r |
31b89da4 | 85 | fRecPointDataPtr = fFirstRecPointPtr;\r |
c375e15d | 86 | \r |
c375e15d | 87 | }//end\r |
88 | \r | |
31b89da4 | 89 | AliHLTCaloClusterizer::~AliHLTCaloClusterizer()\r |
ef44ec64 | 90 | {\r |
31b89da4 | 91 | //See header file for documentation\r |
dd71d443 | 92 | delete [] fBuffer; //FR\r |
93 | fBuffer = NULL; //FR\r | |
94 | fAvailableSize = 0; //FR\r | |
95 | \r | |
9991396c | 96 | delete [] fRecPointArray;\r |
ef44ec64 | 97 | }\r |
98 | \r | |
31b89da4 | 99 | void\r |
ef44ec64 | 100 | AliHLTCaloClusterizer::SetRecPointDataPtr(AliHLTCaloRecPointDataStruct* recPointDataPtr)\r |
101 | {\r | |
31b89da4 | 102 | // See header file for documentation\r |
103 | fRecPointDataPtr = recPointDataPtr;\r | |
ef44ec64 | 104 | }\r |
105 | \r | |
31b89da4 | 106 | Int_t\r |
7c80a370 | 107 | AliHLTCaloClusterizer::ClusterizeEvent(Int_t nDigits)\r |
ef44ec64 | 108 | {\r |
31b89da4 | 109 | //see header file for documentation\r |
110 | Int_t nRecPoints = 0;\r | |
111 | fNRecPoints = 0;\r | |
112 | fUsedSize = 0;\r | |
113 | fNDigits = nDigits;\r | |
114 | fRecPointDataPtr = fFirstRecPointPtr;\r | |
115 | \r | |
116 | // Sort our digits\r | |
117 | SortDigits();\r | |
118 | \r | |
119 | //Clusterization starts\r | |
120 | for (Int_t i = 0; i < nDigits; i++)\r | |
121 | {\r | |
122 | fDigitsInCluster = 0;\r | |
123 | \r | |
f92543f5 | 124 | HLTDebug("Digit with energy: %f", fDigitsPointerArray[i]->fEnergy);\r |
125 | \r | |
126 | if (fDigitsPointerArray[i]->fEnergy < fEmcClusteringThreshold && fSortedByEnergy)\r | |
31b89da4 | 127 | {\r |
f92543f5 | 128 | // Since we have sorted by energy the next digit will have even lower energy, so we return \r |
129 | return fNRecPoints;\r | |
130 | }\r | |
131 | \r | |
132 | if(fDigitsPointerArray[i]->fAssociatedCluster != -1)\r | |
133 | {\r | |
134 | // The digit is added to a previous cluster, continue\r | |
135 | continue;\r | |
136 | }\r | |
31b89da4 | 137 | \r |
f92543f5 | 138 | CheckArray();\r |
31b89da4 | 139 | CheckBuffer();\r |
140 | \r | |
141 | // First digit is placed at the fDigits member variable in the recpoint\r | |
142 | fDigitIndexPtr = &(fRecPointDataPtr->fDigits);\r | |
143 | \r | |
144 | fRecPointDataPtr->fAmp = 0;\r | |
145 | fRecPointDataPtr->fModule = fDigitsPointerArray[i]->fModule;\r | |
afcce849 | 146 | \r |
31b89da4 | 147 | // Assigning the digit to this rec point\r |
148 | fRecPointDataPtr->fDigits = i;\r | |
149 | fUsedSize += sizeof(AliHLTCaloRecPointDataStruct);\r | |
150 | \r | |
151 | // Incrementing the pointer to be ready for new entry\r | |
152 | fDigitIndexPtr++;\r | |
153 | \r | |
154 | fRecPointDataPtr->fAmp += fDigitsPointerArray[i]->fEnergy;\r | |
f92543f5 | 155 | \r |
156 | \r | |
157 | //fDigitsPointerArray[i]->fEnergy = 0;\r | |
158 | fDigitsPointerArray[i]->fAssociatedCluster = fNRecPoints;\r | |
159 | \r | |
160 | \r | |
161 | fDigitsInCluster++;\r | |
31b89da4 | 162 | nRecPoints++;\r |
163 | \r | |
164 | // Scanning for the neighbours\r | |
165 | if (ScanForNeighbourDigits(i, fRecPointDataPtr) != 0)\r | |
166 | {\r | |
167 | return -1;\r | |
168 | }\r | |
169 | \r | |
170 | //fUsedSize += sizeof(AliHLTCaloRecPointDataStruct) + (fDigitsInCluster-1)*sizeof(AliHLTCaloDigitDataStruct);\r | |
171 | \r | |
172 | fRecPointDataPtr->fMultiplicity = fDigitsInCluster;\r | |
173 | fRecPointArray[fNRecPoints] = fRecPointDataPtr;\r | |
174 | \r | |
175 | fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fDigitIndexPtr);\r | |
176 | \r | |
177 | fNRecPoints++;\r | |
178 | \r | |
ef44ec64 | 179 | }//end of clusterization\r |
70c86cde | 180 | \r |
181 | return nRecPoints;\r | |
ef44ec64 | 182 | }\r |
183 | \r | |
9bd6f29b | 184 | Int_t\r |
ef44ec64 | 185 | AliHLTCaloClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTCaloRecPointDataStruct* recPoint)\r |
186 | {\r | |
31b89da4 | 187 | //see header file for documentation\r |
ef44ec64 | 188 | \r |
31b89da4 | 189 | // The following cuts can be used if we sort by posisiton. Not tested, but it should be fine...\r |
190 | Int_t max = TMath::Min(fNDigits, (Int_t)fMaxDigitIndexDiff+index);\r | |
191 | Int_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff));\r | |
192 | \r | |
193 | // All digits for now\r | |
194 | max = fNDigits;\r | |
195 | min = 0;\r | |
196 | \r | |
197 | for (Int_t j = min; j < max; j++)\r | |
ef44ec64 | 198 | {\r |
f92543f5 | 199 | if (fDigitsPointerArray[j]->fAssociatedCluster == -1 && fDigitsPointerArray[j]->fEnergy > fEmcMinEnergyThreshold)\r |
31b89da4 | 200 | {\r |
201 | if (j != index)\r | |
202 | {\r | |
203 | if (AreNeighbours(fDigitsPointerArray[index],\r | |
204 | fDigitsPointerArray[j]))\r | |
205 | {\r | |
206 | // Check that the buffer is large enough for adding a digit (can be heavily improved wrt performance)\r | |
207 | CheckBuffer();\r | |
208 | \r | |
209 | // Assigning index to digit\r | |
210 | *fDigitIndexPtr = j;\r | |
211 | fUsedSize += sizeof(Int_t);\r | |
212 | \r | |
213 | // Incrementing digit pointer to be ready for new entry\r | |
214 | fDigitIndexPtr++;\r | |
215 | \r | |
216 | // Adding the digit energy to the rec point\r | |
217 | fRecPointDataPtr->fAmp += fDigitsPointerArray[j]->fEnergy;\r | |
218 | \r | |
f92543f5 | 219 | // Setting energy to 0\r |
220 | //fDigitsPointerArray[j]->fEnergy = 0;\r | |
221 | \r | |
222 | // Setting the associated cluster \r | |
223 | fDigitsPointerArray[j]->fAssociatedCluster = fNRecPoints;\r | |
224 | \r | |
225 | HLTDebug("Added digit with index: %d, energy: %f, to associated cluster: %d", fDigitsPointerArray[j]->fID, fDigitsPointerArray[j]->fEnergy, fDigitsPointerArray[j]->fAssociatedCluster);\r | |
226 | \r | |
31b89da4 | 227 | fDigitsInCluster++;\r |
228 | \r | |
229 | // Scan for neighbours of this digit\r | |
230 | ScanForNeighbourDigits(j, recPoint);\r | |
231 | }\r | |
232 | }\r | |
233 | }\r | |
ef44ec64 | 234 | }\r |
31b89da4 | 235 | return 0;\r |
ef44ec64 | 236 | }\r |
237 | \r | |
31b89da4 | 238 | Int_t\r |
239 | AliHLTCaloClusterizer::AreNeighbours(AliHLTCaloDigitDataStruct* digit1,\r | |
240 | AliHLTCaloDigitDataStruct* digit2)\r | |
ef44ec64 | 241 | {\r |
31b89da4 | 242 | //see header file for documentation\r |
36b00077 | 243 | if ( (digit1->fModule == digit2->fModule) || AreEdgeCells(digit1, digit2))\r |
31b89da4 | 244 | {\r |
245 | Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );\r | |
246 | Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX );\r | |
247 | \r | |
36b00077 | 248 | // Common edge defines neighbour\r |
249 | //if (( coldiff <= 1 && rowdiff == 0 ) || ( coldiff == 0 && rowdiff <= 1 ))\r | |
250 | // Common edge and corner defines neighbour\r | |
251 | if (( coldiff <= 1 && rowdiff <= 1 ))\r | |
31b89da4 | 252 | {\r |
253 | // Check also for time\r | |
254 | if (TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r | |
255 | {\r | |
256 | return 1;\r | |
257 | }\r | |
258 | }\r | |
ef44ec64 | 259 | }\r |
31b89da4 | 260 | return 0;\r |
ef44ec64 | 261 | }\r |
7c80a370 | 262 | \r |
263 | \r | |
264 | \r | |
265 | Int_t AliHLTCaloClusterizer::CheckArray()\r | |
266 | {\r | |
31b89da4 | 267 | // See header file for class documentation\r |
268 | if (fArraySize == fNRecPoints)\r | |
269 | {\r | |
270 | fArraySize *= 2;\r | |
271 | AliHLTCaloRecPointDataStruct **tmp = new AliHLTCaloRecPointDataStruct*[fArraySize];\r | |
272 | memcpy(tmp, fRecPointArray, fArraySize/2 * sizeof(AliHLTCaloRecPointDataStruct*));\r | |
273 | delete [] fRecPointArray;\r | |
274 | fRecPointArray = tmp;\r | |
275 | }\r | |
276 | return 0;\r | |
7c80a370 | 277 | }\r |
278 | \r | |
279 | Int_t AliHLTCaloClusterizer::CheckBuffer()\r | |
280 | {\r | |
31b89da4 | 281 | // See header file for class documentation\r |
282 | if ((fAvailableSize - fUsedSize) < (Int_t)sizeof(AliHLTCaloRecPointDataStruct))\r | |
283 | {\r | |
284 | Int_t recPointOffset = reinterpret_cast<UChar_t*>(fRecPointDataPtr) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr);\r | |
285 | Int_t digitIndexOffset = reinterpret_cast<UChar_t*>(fDigitIndexPtr) - reinterpret_cast<UChar_t*>(fRecPointDataPtr);\r | |
286 | UChar_t *tmp = new UChar_t[fAvailableSize*2];\r | |
dd71d443 | 287 | \r |
288 | if (tmp == NULL)\r | |
289 | {\r | |
290 | HLTError("Pointer error");\r | |
291 | return(-1);\r | |
292 | }\r | |
293 | \r | |
31b89da4 | 294 | memcpy(tmp, fFirstRecPointPtr, fUsedSize);\r |
295 | fAvailableSize *= 2;\r | |
296 | for (Int_t n = 0; n < fNRecPoints; n++)\r | |
297 | {\r | |
298 | fRecPointArray[n] = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(reinterpret_cast<UChar_t*>(fRecPointArray[n]) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr) + reinterpret_cast<UChar_t*>(tmp));\r | |
299 | }\r | |
dd71d443 | 300 | delete [] fBuffer; //FR\r |
301 | fBuffer = tmp; //FR\r | |
302 | \r | |
31b89da4 | 303 | fFirstRecPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp);\r |
304 | fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp + recPointOffset);\r | |
305 | fDigitIndexPtr = reinterpret_cast<Int_t*>(reinterpret_cast<UChar_t*>(fRecPointDataPtr) + digitIndexOffset);\r | |
306 | //fUsedSize = 0;\r | |
307 | }\r | |
308 | return 0;\r | |
7c80a370 | 309 | }\r |
310 | \r | |
31b89da4 | 311 | void AliHLTCaloClusterizer::SetSortDigitsByPosition()\r |
7c80a370 | 312 | {\r |
31b89da4 | 313 | // Sort the digit pointers by position\r |
314 | fCompareFunction = &CompareDigitsByPosition;\r | |
315 | fSortDigits = true;\r | |
316 | fSortedByPosition = true;\r | |
7c80a370 | 317 | }\r |
318 | \r | |
31b89da4 | 319 | void AliHLTCaloClusterizer::SetSortDigitsByEnergy()\r |
320 | {\r | |
321 | // See header file for class documentation\r | |
322 | fCompareFunction = &CompareDigitsByEnergy;\r | |
323 | fSortDigits = true;\r | |
324 | fSortedByEnergy = true;\r | |
325 | }\r | |
326 | \r | |
327 | void AliHLTCaloClusterizer::SortDigits()\r | |
328 | {\r | |
329 | // See header file for class documentation\r | |
330 | if (fSortDigits) qsort(fDigitsPointerArray, fNDigits, sizeof(AliHLTCaloDigitDataStruct*), fCompareFunction);\r | |
331 | }\r | |
332 | \r | |
333 | Int_t\r | |
334 | AliHLTCaloClusterizer::CompareDigitsByPosition(const void *dig0, const void *dig1)\r | |
335 | {\r | |
336 | // See header file for documentation\r | |
337 | return (*((AliHLTCaloDigitDataStruct**)(dig0)))->fID - (*((AliHLTCaloDigitDataStruct**)(dig1)))->fID;\r | |
338 | }\r | |
339 | \r | |
340 | Int_t\r | |
341 | AliHLTCaloClusterizer::CompareDigitsByEnergy(const void *dig0, const void *dig1)\r | |
7c80a370 | 342 | {\r |
31b89da4 | 343 | // See header file for documentation\r |
02f2218d | 344 | if ( ((*((AliHLTCaloDigitDataStruct**)(dig1)))->fEnergy - (*((AliHLTCaloDigitDataStruct**)(dig0)))->fEnergy) < 0) return -1;\r |
345 | return 1;\r | |
2a24cbbe | 346 | }\r |
36b00077 | 347 | \r |
348 | void AliHLTCaloClusterizer::SetDetector(TString det)\r | |
349 | {\r | |
350 | if(det.CompareTo("EMCAL"))\r | |
351 | {\r | |
352 | fIsEMCAL = true;\r | |
353 | }\r | |
354 | else\r | |
355 | {\r | |
356 | fIsEMCAL = false;\r | |
357 | }\r | |
358 | }\r | |
359 | \r | |
360 | Bool_t AliHLTCaloClusterizer::AreEdgeCells(AliHLTCaloDigitDataStruct *digit0, AliHLTCaloDigitDataStruct *digit1)\r | |
361 | {\r | |
362 | if(fIsEMCAL)\r | |
363 | {\r | |
364 | Int_t modDiff = digit0->fModule - digit1->fModule;\r | |
365 | if(TMath::Abs(modDiff) > 1) return kFALSE;\r | |
366 | if(digit0->fModule > digit1->fModule && digit1->fModule%2 == 0) \r | |
367 | {\r | |
368 | if(digit0->fZ == 0 && digit1->fZ == (fCaloConstants->GetNZROWSMOD()-1))\r | |
369 | return kTRUE;\r | |
370 | }\r | |
371 | if(digit1->fModule > digit0->fModule && digit0->fModule%2 == 0) \r | |
372 | {\r | |
373 | if(digit1->fZ == 0 && digit0->fZ == (fCaloConstants->GetNZROWSMOD()-1))\r | |
374 | return kTRUE;\r | |
375 | }\r | |
376 | }\r | |
377 | \r | |
378 | return false;\r | |
379 | \r | |
380 | }\r |