Argghh... small typo... corrected
[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   fRecPointDataPtr(0),\r
44   fDigitIndexPtr(0),\r
45   fEmcClusteringThreshold(0),\r
46   fEmcMinEnergyThreshold(0),\r
47   fEmcTimeGate(0),\r
48   fDigitsInCluster(0),\r
49   fDigitsPointerArray(0),\r
50   fDigitContainerPtr(0),\r
51   fMaxDigitIndexDiff(0),\r
52   fNDigits(0)\r
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
60 \r
61 }//end\r
62 \r
63 \r
64 //BALLE how do you set the right detector?\r
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
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
92 Int_t \r
93 AliHLTCaloClusterizer::ClusterizeEvent(Int_t nDigits, UInt_t availableSize, UInt_t& totSize)\r
94 {\r
95   //see header file for documentation\r
96   Int_t nRecPoints = 0;\r
97 \r
98   fNDigits = nDigits;\r
99 \r
100   UInt_t maxRecPointSize = sizeof(AliHLTCaloRecPointDataStruct) + (sizeof(AliHLTCaloDigitDataStruct) << 7); //Reasonable estimate... \r
101 \r
102   //Clusterization starts\r
103   for(Int_t i = 0; i < nDigits; i++)\r
104     { \r
105       fDigitsInCluster = 0;\r
106       //      printf("ENERGY: %f\n", fDigitsPointerArray[i]->fEnergy);\r
107       if(fDigitsPointerArray[i]->fEnergy < fEmcClusteringThreshold)\r
108         {\r
109           continue;\r
110         }\r
111       if(availableSize < (totSize + maxRecPointSize)) \r
112         {\r
113           return -1; //Might get out of buffer, exiting\r
114         }\r
115       //            printf("cluster candidate!\n");\r
116       // First digit is placed at the fDigits member variable in the recpoint\r
117       fDigitIndexPtr = &(fRecPointDataPtr->fDigits);\r
118 \r
119       fRecPointDataPtr->fAmp = 0;\r
120       fRecPointDataPtr->fModule = fDigitsPointerArray[i]->fModule;\r
121 \r
122       // Assigning the digit to this rec point\r
123       fRecPointDataPtr->fDigits = i;\r
124 \r
125       // Incrementing the pointer to be ready for new entry\r
126       fDigitIndexPtr++;\r
127 \r
128       fRecPointDataPtr->fAmp += fDigitsPointerArray[i]->fEnergy;\r
129       fDigitsPointerArray[i]->fEnergy = 0;\r
130       fDigitsInCluster++;\r
131       nRecPoints++;\r
132 \r
133       // Scanning for the neighbours\r
134       ScanForNeighbourDigits(i, fRecPointDataPtr);\r
135 \r
136       totSize += sizeof(AliHLTCaloRecPointDataStruct) + (fDigitsInCluster-1)*sizeof(AliHLTCaloDigitDataStruct);   \r
137       fRecPointDataPtr->fMultiplicity = fDigitsInCluster;     \r
138       //      printf("Rec point energy: %f\n", fRecPointDataPtr->fAmp);\r
139       fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fDigitIndexPtr);\r
140 \r
141     }//end of clusterization\r
142 \r
143    return nRecPoints;\r
144 }\r
145 \r
146 void\r
147 AliHLTCaloClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTCaloRecPointDataStruct* recPoint)\r
148 {\r
149   //see header file for documentation\r
150   Int_t max = TMath::Min(fNDigits, (Int_t)fMaxDigitIndexDiff+index);\r
151   Int_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff));\r
152 \r
153   max = fNDigits;\r
154   min = 0;\r
155   for(Int_t j = min; j < max; j++)\r
156     {\r
157       if(fDigitsPointerArray[j]->fEnergy > fEmcMinEnergyThreshold)\r
158         {\r
159           if(j != index)\r
160             {\r
161               if(AreNeighbours(fDigitsPointerArray[index],\r
162                                fDigitsPointerArray[j]))\r
163                 {\r
164                   // Assigning value to digit ptr\r
165                   *fDigitIndexPtr = j;\r
166                   // Incrementing digit pointer to be ready for new entry\r
167                   fDigitIndexPtr++;\r
168 \r
169                   recPoint->fAmp += fDigitsPointerArray[j]->fEnergy;\r
170                   fDigitsPointerArray[j]->fEnergy = 0;        \r
171                   fDigitsInCluster++;\r
172                   ScanForNeighbourDigits(j, recPoint);\r
173                 }\r
174             }\r
175         }\r
176     }\r
177   return;\r
178 }\r
179 \r
180 Int_t \r
181 AliHLTCaloClusterizer::AreNeighbours(AliHLTCaloDigitDataStruct* digit1, \r
182                                             AliHLTCaloDigitDataStruct* digit2)\r
183 {\r
184   //see header file for documentation\r
185   if ( (digit1->fModule == digit2->fModule) /*&& (coord1[1]==coord2[1])*/ ) // inside the same PHOS module\r
186     { \r
187 //       Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );  \r
188 //       Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX ); \r
189 //       if (( coldiff <= 1   &&  rowdiff == 0 ) || ( coldiff == 0 &&  rowdiff <= 1 ))\r
190 //      {\r
191 //        cout << "Are neighbours: digit (E = "  << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r
192 //          " is neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r
193 \r
194 //        if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r
195 //          {\r
196 //            return 1; \r
197 //          }\r
198 //      }\r
199 \r
200       Float_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );  \r
201       Float_t coldiff = TMath::Abs( digit1->fX - digit2->fX ); \r
202       if (( coldiff <= 2.4   &&  rowdiff < 0.4 ) || ( coldiff < 0.4 &&  rowdiff <= 2.4 ))\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       else\r
213         {\r
214           //  cout << "Not neighbours: digit (E = "  << digit1->fEnergy << ") with x = " << digit1->fX << " and z = " << digit1->fZ << \r
215           //  " is not neighbour with digit (E = " << digit2->fEnergy << ") with x = " << digit2->fX << " and z = " << digit2->fZ << endl;\r
216         }\r
217     }\r
218   return 0;\r
219 }\r