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         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
  63 {\r
  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
  81 \r
  82 }//end\r
  83 \r
  84 AliHLTCaloClusterizer::~AliHLTCaloClusterizer()\r
  85 {\r
  86     //See header file for documentation\r
  87 }\r
  88 \r
  89 void\r
  90 AliHLTCaloClusterizer::SetRecPointDataPtr(AliHLTCaloRecPointDataStruct* recPointDataPtr)\r
  91 {\r
  92     // See header file for documentation\r
  93     fRecPointDataPtr = recPointDataPtr;\r
  94 }\r
  95 \r
  96 Int_t\r
  97 AliHLTCaloClusterizer::ClusterizeEvent(Int_t nDigits)\r
  98 {\r
  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
 114         if (fDigitsPointerArray[i]->fEnergy < fEmcClusteringThreshold)\r
 115         {\r
 116             continue;\r
 117         }\r
 118 \r
 119         CheckArray();\r
 120         CheckBuffer();\r
 121 \r
 122         // First digit is placed at the fDigits member variable in the recpoint\r
 123         fDigitIndexPtr = &(fRecPointDataPtr->fDigits);\r
 124 \r
 125         fRecPointDataPtr->fAmp = 0;\r
 126         fRecPointDataPtr->fModule = fDigitsPointerArray[i]->fModule;\r
 127 \r
 128         // Assigning the digit to this rec point\r
 129         fRecPointDataPtr->fDigits = i;\r
 130         fUsedSize += sizeof(AliHLTCaloRecPointDataStruct);\r
 131 \r
 132         // Incrementing the pointer to be ready for new entry\r
 133         fDigitIndexPtr++;\r
 134 \r
 135         fRecPointDataPtr->fAmp += fDigitsPointerArray[i]->fEnergy;\r
 136         fDigitsPointerArray[i]->fEnergy = 0;\r
 137         fDigitsInCluster++;\r
 138         nRecPoints++;\r
 139 \r
 140         // Scanning for the neighbours\r
 141         if (ScanForNeighbourDigits(i, fRecPointDataPtr) != 0)\r
 142         {\r
 143             return -1;\r
 144         }\r
 145 \r
 146         //fUsedSize += sizeof(AliHLTCaloRecPointDataStruct) + (fDigitsInCluster-1)*sizeof(AliHLTCaloDigitDataStruct);\r
 147 \r
 148         fRecPointDataPtr->fMultiplicity = fDigitsInCluster;\r
 149         fRecPointArray[fNRecPoints] = fRecPointDataPtr;\r
 150 \r
 151         fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(fDigitIndexPtr);\r
 152 \r
 153         fNRecPoints++;\r
 154 \r
 155     }//end of clusterization\r
 156 \r
 157     return nRecPoints;\r
 158 }\r
 159 \r
 160 Int_t\r
 161 AliHLTCaloClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTCaloRecPointDataStruct* recPoint)\r
 162 {\r
 163     //see header file for documentation\r
 164 \r
 165     // The following cuts can be used if we sort by posisiton. Not tested, but it should be fine...\r
 166     Int_t max = TMath::Min(fNDigits, (Int_t)fMaxDigitIndexDiff+index);\r
 167     Int_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff));\r
 168 \r
 169     // All digits for now\r
 170     max = fNDigits;\r
 171     min = 0;\r
 172 \r
 173     for (Int_t j = min; j < max; j++)\r
 174     {\r
 175         if (fDigitsPointerArray[j]->fEnergy > fEmcMinEnergyThreshold)\r
 176         {\r
 177             if (j != index)\r
 178             {\r
 179                 if (AreNeighbours(fDigitsPointerArray[index],\r
 180                                   fDigitsPointerArray[j]))\r
 181                 {\r
 182                     // Check that the buffer is large enough for adding a digit (can be heavily improved wrt performance)\r
 183                     CheckBuffer();\r
 184 \r
 185                     // Assigning index to digit\r
 186                     *fDigitIndexPtr = j;\r
 187                     fUsedSize += sizeof(Int_t);\r
 188 \r
 189                     // Incrementing digit pointer to be ready for new entry\r
 190                     fDigitIndexPtr++;\r
 191 \r
 192                     // Adding the digit energy to the rec point\r
 193                     fRecPointDataPtr->fAmp += fDigitsPointerArray[j]->fEnergy;\r
 194 \r
 195                     // Setting the digit energy to 0 (not included in clusterisation anymore)\r
 196                     fDigitsPointerArray[j]->fEnergy = 0;\r
 197 \r
 198                     fDigitsInCluster++;\r
 199 \r
 200                     // Scan for neighbours of this digit\r
 201                     ScanForNeighbourDigits(j, recPoint);\r
 202                 }\r
 203             }\r
 204         }\r
 205     }\r
 206     return 0;\r
 207 }\r
 208 \r
 209 Int_t\r
 210 AliHLTCaloClusterizer::AreNeighbours(AliHLTCaloDigitDataStruct* digit1,\r
 211                                      AliHLTCaloDigitDataStruct* digit2)\r
 212 {\r
 213     //see header file for documentation\r
 214     if ( (digit1->fModule == digit2->fModule) /*&& (coord1[1]==coord2[1])*/ ) // inside the same PHOS module\r
 215     {\r
 216         Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );\r
 217         Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX );\r
 218 \r
 219         // As in the offline code we define neighbours as cells that share an edge, a corner is not  enough\r
 220         if (( coldiff <= 1   &&  rowdiff == 0 ) || ( coldiff == 0 &&  rowdiff <= 1 ))\r
 221         {\r
 222             // Check also for time\r
 223             if (TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)\r
 224             {\r
 225                 return 1;\r
 226             }\r
 227         }\r
 228     }\r
 229     return 0;\r
 230 }\r
 231 \r
 232 \r
 233 \r
 234 Int_t AliHLTCaloClusterizer::CheckArray()\r
 235 {\r
 236     // See header file for class documentation\r
 237     if (fArraySize == fNRecPoints)\r
 238     {\r
 239         fArraySize *= 2;\r
 240         AliHLTCaloRecPointDataStruct **tmp = new AliHLTCaloRecPointDataStruct*[fArraySize];\r
 241         memcpy(tmp, fRecPointArray, fArraySize/2 * sizeof(AliHLTCaloRecPointDataStruct*));\r
 242         delete [] fRecPointArray;\r
 243         fRecPointArray = tmp;\r
 244     }\r
 245     return 0;\r
 246 }\r
 247 \r
 248 Int_t AliHLTCaloClusterizer::CheckBuffer()\r
 249 {\r
 250     // See header file for class documentation\r
 251     if ((fAvailableSize - fUsedSize) < (Int_t)sizeof(AliHLTCaloRecPointDataStruct))\r
 252     {\r
 253         Int_t recPointOffset = reinterpret_cast<UChar_t*>(fRecPointDataPtr) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr);\r
 254         Int_t digitIndexOffset = reinterpret_cast<UChar_t*>(fDigitIndexPtr) - reinterpret_cast<UChar_t*>(fRecPointDataPtr);\r
 255         UChar_t *tmp = new UChar_t[fAvailableSize*2];\r
 256 \r
 257         memcpy(tmp, fFirstRecPointPtr, fUsedSize);\r
 258         fAvailableSize *= 2;\r
 259         for (Int_t n = 0; n < fNRecPoints; n++)\r
 260         {\r
 261             fRecPointArray[n] = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(reinterpret_cast<UChar_t*>(fRecPointArray[n]) - reinterpret_cast<UChar_t*>(fFirstRecPointPtr) + reinterpret_cast<UChar_t*>(tmp));\r
 262         }\r
 263         delete [] fFirstRecPointPtr;\r
 264         fFirstRecPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp);\r
 265         fRecPointDataPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(tmp + recPointOffset);\r
 266         fDigitIndexPtr = reinterpret_cast<Int_t*>(reinterpret_cast<UChar_t*>(fRecPointDataPtr) + digitIndexOffset);\r
 267         //fUsedSize = 0;\r
 268     }\r
 269     return 0;\r
 270 }\r
 271 \r
 272 void AliHLTCaloClusterizer::SetSortDigitsByPosition()\r
 273 {\r
 274     // Sort the digit pointers by position\r
 275     fCompareFunction = &CompareDigitsByPosition;\r
 276     fSortDigits = true;\r
 277     fSortedByPosition = true;\r
 278 }\r
 279 \r
 280 void AliHLTCaloClusterizer::SetSortDigitsByEnergy()\r
 281 {\r
 282     // See header file for class documentation\r
 283     fCompareFunction = &CompareDigitsByEnergy;\r
 284     fSortDigits = true;\r
 285     fSortedByEnergy = true;\r
 286 }\r
 287 \r
 288 void AliHLTCaloClusterizer::SortDigits()\r
 289 {\r
 290     // See header file for class documentation\r
 291     if (fSortDigits) qsort(fDigitsPointerArray, fNDigits, sizeof(AliHLTCaloDigitDataStruct*), fCompareFunction);\r
 292 }\r
 293 \r
 294 Int_t\r
 295 AliHLTCaloClusterizer::CompareDigitsByPosition(const void *dig0, const void *dig1)\r
 296 {\r
 297     // See header file for documentation\r
 298     return (*((AliHLTCaloDigitDataStruct**)(dig0)))->fID - (*((AliHLTCaloDigitDataStruct**)(dig1)))->fID;\r
 299 }\r
 300 \r
 301 Int_t\r
 302 AliHLTCaloClusterizer::CompareDigitsByEnergy(const void *dig0, const void *dig1)\r
 303 {\r
 304     // See header file for documentation\r
 305     return (*((AliHLTCaloDigitDataStruct**)(dig1)))->fEnergy - (*((AliHLTCaloDigitDataStruct**)(dig0)))->fEnergy;\r
 306 }\r