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 |
39 | ClassImp(AliHLTCaloClusterizer);\r |
40 | \r |
41 | AliHLTCaloClusterizer::AliHLTCaloClusterizer(TString det):\r |
4f4b7ba4 |
42 | AliHLTCaloConstantsHandler(det),\r |
ef44ec64 |
43 | fRecPointDataPtr(0),\r |
ad44d760 |
44 | fDigitIndexPtr(0),\r |
ef44ec64 |
45 | fEmcClusteringThreshold(0),\r |
46 | fEmcMinEnergyThreshold(0),\r |
47 | fEmcTimeGate(0),\r |
48 | fDigitsInCluster(0),\r |
f92dcf21 |
49 | fDigitsPointerArray(0),\r |
ef44ec64 |
50 | fDigitContainerPtr(0),\r |
ad44d760 |
51 | fMaxDigitIndexDiff(0),\r |
52 | fNDigits(0)\r |
ef44ec64 |
53 | {\r |
54 | //See header file for documentation\r |
55 | fEmcClusteringThreshold = 0.2;\r |
56 | fEmcMinEnergyThreshold = 0.03;\r |
57 | fEmcTimeGate = 1.e-6 ;\r |
58 | \r |
59 | fMaxDigitIndexDiff = 2*fCaloConstants->GetNZROWSMOD();\r |
c375e15d |
60 | \r |
c375e15d |
61 | }//end\r |
62 | \r |
4f4b7ba4 |
63 | \r |
64 | //BALLE how do you set the right detector?\r |
ad44d760 |
65 | // AliHLTCaloClusterizer::AliHLTCaloClusterizer(const AliHLTCaloClusterizer &) :\r |
66 | // AliHLTCaloConstantsHandler("BALLE"),\r |
67 | // fRecPointDataPtr(0),\r |
68 | // fDigitDataPtr(0),\r |
69 | // fEmcClusteringThreshold(0),\r |
70 | // fEmcMinEnergyThreshold(0),\r |
71 | // fEmcTimeGate(0),\r |
72 | // fDigitsInCluster(0),\r |
73 | // fDigitContainerPtr(0),\r |
74 | // fMaxDigitIndexDiff(0)\r |
75 | // {\r |
76 | // // dummy copy constructor\r |
77 | // }//end\r |
ef44ec64 |
78 | \r |
79 | \r |
80 | AliHLTCaloClusterizer::~AliHLTCaloClusterizer() \r |
81 | {\r |
82 | //See header file for documentation\r |
83 | }\r |
84 | \r |
85 | void \r |
86 | AliHLTCaloClusterizer::SetRecPointDataPtr(AliHLTCaloRecPointDataStruct* recPointDataPtr)\r |
87 | {\r |
88 | // See header file for documentation\r |
89 | fRecPointDataPtr = recPointDataPtr;\r |
90 | }\r |
91 | \r |
ef44ec64 |
92 | Int_t \r |
ad44d760 |
93 | AliHLTCaloClusterizer::ClusterizeEvent(Int_t nDigits, UInt_t availableSize, UInt_t& totSize)\r |
ef44ec64 |
94 | {\r |
95 | //see header file for documentation\r |
96 | Int_t nRecPoints = 0;\r |
9bd6f29b |
97 | \r |
98 | fAvailableSize = availableSize;\r |
99 | \r |
ad44d760 |
100 | fNDigits = nDigits;\r |
101 | \r |
ef44ec64 |
102 | UInt_t maxRecPointSize = sizeof(AliHLTCaloRecPointDataStruct) + (sizeof(AliHLTCaloDigitDataStruct) << 7); //Reasonable estimate... \r |
103 | \r |
104 | //Clusterization starts\r |
f92dcf21 |
105 | for(Int_t i = 0; i < nDigits; i++)\r |
ef44ec64 |
106 | { \r |
107 | fDigitsInCluster = 0;\r |
ad44d760 |
108 | // printf("ENERGY: %f\n", fDigitsPointerArray[i]->fEnergy);\r |
109 | if(fDigitsPointerArray[i]->fEnergy < fEmcClusteringThreshold)\r |
ef44ec64 |
110 | {\r |
111 | continue;\r |
112 | }\r |
9bd6f29b |
113 | \r |
114 | if(fAvailableSize < (sizeof(AliHLTCaloRecPointDataStruct)))\r |
ef44ec64 |
115 | {\r |
9bd6f29b |
116 | HLTError("Out of buffer, stopping clusterisation");\r |
117 | return -1; \r |
ef44ec64 |
118 | }\r |
9bd6f29b |
119 | \r |
ad44d760 |
120 | // printf("cluster candidate!\n");\r |
ef44ec64 |
121 | // First digit is placed at the fDigits member variable in the recpoint\r |
ad44d760 |
122 | fDigitIndexPtr = &(fRecPointDataPtr->fDigits);\r |
ef44ec64 |
123 | \r |
124 | fRecPointDataPtr->fAmp = 0;\r |
ad44d760 |
125 | fRecPointDataPtr->fModule = fDigitsPointerArray[i]->fModule;\r |
ef44ec64 |
126 | \r |
ad44d760 |
127 | // Assigning the digit to this rec point\r |
128 | fRecPointDataPtr->fDigits = i;\r |
ef44ec64 |
129 | \r |
130 | // Incrementing the pointer to be ready for new entry\r |
ad44d760 |
131 | fDigitIndexPtr++;\r |
ef44ec64 |
132 | \r |
ad44d760 |
133 | fRecPointDataPtr->fAmp += fDigitsPointerArray[i]->fEnergy;\r |
134 | fDigitsPointerArray[i]->fEnergy = 0;\r |
ef44ec64 |
135 | fDigitsInCluster++;\r |
136 | nRecPoints++;\r |
137 | \r |
138 | // Scanning for the neighbours\r |
9bd6f29b |
139 | if(ScanForNeighbourDigits(i, fRecPointDataPtr) != 0)\r |
140 | {\r |
141 | return -1;\r |
142 | }\r |
ef44ec64 |
143 | \r |
144 | totSize += sizeof(AliHLTCaloRecPointDataStruct) + (fDigitsInCluster-1)*sizeof(AliHLTCaloDigitDataStruct); \r |
145 | fRecPointDataPtr->fMultiplicity = fDigitsInCluster; \r |
ad44d760 |
146 | // printf("Rec point energy: %f\n", fRecPointDataPtr->fAmp);\r |
147 | fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fDigitIndexPtr);\r |
ef44ec64 |
148 | \r |
ef44ec64 |
149 | }//end of clusterization\r |
150 | \r |
151 | return nRecPoints;\r |
152 | }\r |
153 | \r |
9bd6f29b |
154 | Int_t\r |
ef44ec64 |
155 | AliHLTCaloClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTCaloRecPointDataStruct* recPoint)\r |
156 | {\r |
157 | //see header file for documentation\r |
ad44d760 |
158 | Int_t max = TMath::Min(fNDigits, (Int_t)fMaxDigitIndexDiff+index);\r |
ef44ec64 |
159 | Int_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff));\r |
160 | \r |
ad44d760 |
161 | max = fNDigits;\r |
ef44ec64 |
162 | min = 0;\r |
163 | for(Int_t j = min; j < max; j++)\r |
164 | {\r |
ad44d760 |
165 | if(fDigitsPointerArray[j]->fEnergy > fEmcMinEnergyThreshold)\r |
ef44ec64 |
166 | {\r |
167 | if(j != index)\r |
168 | {\r |
ad44d760 |
169 | if(AreNeighbours(fDigitsPointerArray[index],\r |
170 | fDigitsPointerArray[j]))\r |
ef44ec64 |
171 | {\r |
9bd6f29b |
172 | if(fAvailableSize < (sizeof(Int_t)))\r |
173 | {\r |
174 | HLTError("Out of buffer, stopping clusterisation");\r |
175 | return -1; \r |
176 | } \r |
177 | // Assigning index to digit\r |
ad44d760 |
178 | *fDigitIndexPtr = j;\r |
ef44ec64 |
179 | // Incrementing digit pointer to be ready for new entry\r |
ad44d760 |
180 | fDigitIndexPtr++;\r |
ef44ec64 |
181 | \r |
ad44d760 |
182 | recPoint->fAmp += fDigitsPointerArray[j]->fEnergy;\r |
183 | fDigitsPointerArray[j]->fEnergy = 0; \r |
ef44ec64 |
184 | fDigitsInCluster++;\r |
185 | ScanForNeighbourDigits(j, recPoint);\r |
186 | }\r |
187 | }\r |
188 | }\r |
189 | }\r |
9bd6f29b |
190 | return 0;\r |
ef44ec64 |
191 | }\r |
192 | \r |
193 | Int_t \r |
194 | AliHLTCaloClusterizer::AreNeighbours(AliHLTCaloDigitDataStruct* digit1, \r |
195 | AliHLTCaloDigitDataStruct* digit2)\r |
196 | {\r |
197 | //see header file for documentation\r |
198 | if ( (digit1->fModule == digit2->fModule) /*&& (coord1[1]==coord2[1])*/ ) // inside the same PHOS module\r |
199 | { \r |
200 | // Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ ); \r |
201 | // Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX ); \r |
202 | // if (( coldiff <= 1 && rowdiff == 0 ) || ( coldiff == 0 && rowdiff <= 1 ))\r |
203 | // {\r |
204 | // cout << "Are neighbours: digit (E = " << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r |
205 | // " is neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r |
206 | \r |
207 | // if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r |
208 | // {\r |
209 | // return 1; \r |
210 | // }\r |
211 | // }\r |
212 | \r |
213 | Float_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ ); \r |
214 | Float_t coldiff = TMath::Abs( digit1->fX - digit2->fX ); \r |
215 | if (( coldiff <= 2.4 && rowdiff < 0.4 ) || ( coldiff < 0.4 && rowdiff <= 2.4 ))\r |
216 | {\r |
217 | // cout << "Are neighbours: digit (E = " << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r |
218 | // " is neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r |
219 | \r |
220 | if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r |
221 | {\r |
222 | return 1; \r |
223 | }\r |
224 | }\r |
225 | else\r |
226 | {\r |
227 | // cout << "Not neighbours: digit (E = " << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r |
228 | // " is not neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r |
229 | }\r |
230 | }\r |
231 | return 0;\r |
232 | }\r |