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