HLT/PHOS/AliHLTPHOSClusterizer.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Authors: Øystein Djuvsland <oysteind@ift.uib.no>                       *
   5  *                                                                        *
   6  * Permission to use, copy, modify and distribute this software and its   *
   7  * documentation strictly for non-commercial purposes is hereby granted   *
   8  * without fee, provided that the above copyright notice appears in all   *
   9  * copies and that both the copyright notice and this permission notice   *
  10  * appear in the supporting documentation. The authors make no claims     *
  11  * about the suitability of this software for any purpose. It is          *
  12  * provided "as is" without express or implied warranty.                  *
  13  **************************************************************************/
  14
  15
  16 /** @file   AliHLTPHOSClusterizer.cxx
  17     @author Øystein Djuvsland
  18     @date
  19     @brief  A temporary clusterizer for PHOS
  20 */
  21
  22
  23
  24 #include "AliHLTPHOSPhysicsDefinitions.h"
  25 #include "AliHLTPHOSClusterizer.h"
  26 #include "AliHLTPHOSCommonDefs.h"
  27 #include "TVector3.h"
  28 #include "TMath.h"
  29 #include "AliHLTPHOSRcuCellEnergyDataStruct.h"
  30 #include "AliHLTPHOSRecPointListDataStruct.h"
  31 #include "AliHLTPHOSValidCellDataStruct.h"
  32 #include "AliHLTPHOSRecPointDataStruct.h"
  33 #include "AliHLTPHOSClusterDataStruct.h"
  34
  35
  36 ClassImp(AliHLTPHOSClusterizer);
  37
  38 /**
  39 * Main constructor
  40 **/
  41 AliHLTPHOSClusterizer::AliHLTPHOSClusterizer():fPHOSModule(-1), fThreshold(0), fClusterThreshold(0),
  42                                                fHighGainFactor(0.005), fLowGainFactor(0.08),
  43                                                fArraySize(3), fMultiplicity(fArraySize*fArraySize)
  44 {
  45   //Main constructor
  46
  47 }//end
  48
  49 AliHLTPHOSClusterizer::AliHLTPHOSClusterizer(const AliHLTPHOSClusterizer &):fPHOSModule(-1), fThreshold(0), fClusterThreshold(0),
  50                                                                             fHighGainFactor(0.005), fLowGainFactor(0.08),
  51                                                                             fArraySize(3), fMultiplicity(fArraySize*fArraySize)
  52 {
  53   //Copy constructor, not implemented
  54 }//end
  55
  56 AliHLTPHOSClusterizer:: ~AliHLTPHOSClusterizer()
  57 {
  58   //Destructor
  59 }
  60
  61 /**
  62 * Building a 2D array of cell energies of the PHOS detector
  63 * @param cellData object containing the cell energies from one event
  64 * @param recPointList list to be filled with coordinates of local maxima in the detector
  65 **/
  66
  67 Int_t
  68 AliHLTPHOSClusterizer::BuildCellEnergyArray(AliHLTPHOSRcuCellEnergyDataStruct* cellData,
  69                                             AliHLTPHOSRecPointListDataStruct* recPointList)
  70 {
  71   //Build the cell energy array of the detector
  72
  73   Int_t x = 0;
  74   Int_t z = 0;
  75   Int_t gain = 0;
  76   Int_t xMod = 0;
  77   Int_t zMod = 0;
  78   Int_t index = 0;
  79   Double_t energyCount = 0;
  80
  81   for(Int_t cell = 0; cell < cellData->fCnt; cell++)
  82     {
  83
  84       z = (cellData->fValidData[cell]).fZ;
  85       x = (cellData->fValidData[cell]).fX;
  86       gain = (cellData->fValidData[cell]).fGain;
  87
  88       zMod = z + (cellData->fRcuZ)*N_ROWS_RCU;
  89       xMod = x + (cellData->fRcuX)*N_COLUMNS_RCU;
  90
  91       energyCount = (cellData->fValidData[cell]).fEnergy;
  92
  93       if(gain == 0 && energyCount < 1023)
  94         {
  95           fEnergyArray[xMod][zMod] = fHighGainFactor * energyCount;
  96
  97           if(fEnergyArray[xMod][zMod] > fClusterThreshold)
  98             {
  99               recPointList[index].fZ = zMod;
 100               recPointList[index].fX = xMod;
 101
 102               for(Int_t j = 0; j < index; j++)
 103                 {
 104                   if(recPointList[j].fZ == zMod && recPointList[j].fX == xMod)
 105                     recPointList[j].fZ = -1;
 106                 }
 107               index++;
 108             }
 109
 110         }
 111       else if(gain == 1 && fEnergyArray[xMod][zMod] == 0)
 112         {
 113           fEnergyArray[xMod][zMod] = fLowGainFactor * energyCount;
 114            if(fEnergyArray[xMod][zMod] > fClusterThreshold)
 115             {
 116               recPointList[index].fZ = zMod;
 117               recPointList[index].fX = xMod;
 118               recPointList[index].fModule = cellData->fModuleID;
 119               index++;
 120             }
 121         }
 122
 123     }
 124
 125   fPHOSModule = cellData->fModuleID;
 126
 127   return index;
 128 }//end BuildCellEnergyArray
 129
 130
 131
 132 /**
 133 * Creating an array of rec points
 134 * @param recPointStructArrayPtr array to store the rec points
 135 * @param list list of rec points
 136 * @param nPoints number of points
 137 **/
 138 Int_t
 139 AliHLTPHOSClusterizer::CreateRecPointStructArray(AliHLTPHOSRecPointDataStruct* recPointStructArrayPtr,
 140                                                  AliHLTPHOSRecPointListDataStruct* list,
 141                                                  Int_t nPoints)
 142
 143 {
 144   //Create the rec point struct array
 145
 146   Int_t flag = 0;
 147   Int_t edgeFlagRows = 0;
 148   Int_t edgeFlagCols = 0;
 149   Int_t k = 0;
 150   Int_t nRecPoints = 0;
 151   Int_t z = 0;
 152   Int_t x = 0;
 153
 154   Float_t* energiesList = NULL;
 155
 156   for(Int_t point = 0; point < nPoints; point++)
 157     {
 158       z = list[point].fZ;
 159       x = list[point].fX;
 160       if(z == -1) continue;
 161
 162       if((z-fArraySize/2) < 0/*= - N_ROWS_MOD/2*/ || (z+fArraySize/2) >= N_ROWS_MOD/*/2*/)
 163         {
 164           edgeFlagRows = 1;
 165           continue;
 166         }
 167       if((x-fArraySize/2) < 0/*= - N_COLUMNS_MOD/2*/ || (x+fArraySize/2) >= N_COLUMNS_MOD)
 168         {
 169           edgeFlagCols = 1;
 170           continue;
 171         }
 172
 173
 174       if(!flag) energiesList = new Float_t[fMultiplicity];
 175       flag = 0;
 176       k = 0;
 177
 178       for(Int_t i = -fArraySize/2; i <= fArraySize/2; i++)
 179         {
 180           if(flag) break;
 181           for(Int_t j = -fArraySize/2; j <= fArraySize/2; j++)
 182             {
 183
 184               if(fEnergyArray[x+i][z+j] > fEnergyArray[x][z] && abs(i) < 2 && abs(j) < 2)
 185                 {
 186                   flag = 1;
 187                   break;
 188                 }
 189               energiesList[k] = fEnergyArray[x+i][z+j];
 190               k++;
 191             }
 192         }
 193
 194
 195       if(!flag && k)
 196         {
 197           recPointStructArrayPtr[nRecPoints].fEnergiesListPtr = energiesList;
 198           recPointStructArrayPtr[nRecPoints].fCoordinatesPtr[0] = x;
 199           recPointStructArrayPtr[nRecPoints].fCoordinatesPtr[1] = z;
 200           recPointStructArrayPtr[nRecPoints].fX = x;
 201           recPointStructArrayPtr[nRecPoints].fZ = z;
 202           recPointStructArrayPtr[nRecPoints].fMultiplicity = fMultiplicity;
 203           recPointStructArrayPtr[nRecPoints].fPHOSModule = list[point].fModule;
 204           nRecPoints++;
 205         }
 206     }
 207
 208
 209
 210   energiesList = NULL;
 211
 212   return nRecPoints;
 213
 214 }//end CreateRecPointStructArray
 215
 216
 217 /**
 218 * Calculating the center of gravity of a rec point
 219 * Not working well at this point!
 220 * @param recPointPtr pointer to the rec point
 221 **/
 222 Int_t
 223 AliHLTPHOSClusterizer::CalculateCenterOfGravity(AliHLTPHOSRecPointDataStruct* recPointPtr)
 224 {
 225   //Calculate the center of gravity
 226
 227   Float_t xt = 0;
 228   Float_t zt = 0;
 229   Float_t xi = 0;
 230   Float_t zj = 0;
 231   Float_t w = 0;
 232   Float_t w0 = 4.5;
 233   Float_t wtot = 0;
 234
 235   Int_t x = recPointPtr->fCoordinatesPtr[0];
 236   Int_t z = recPointPtr->fCoordinatesPtr[1];
 237
 238   for(Int_t i = -fArraySize/2; i <= fArraySize/2; i++)
 239     {
 240       for(Int_t j = -fArraySize/2; j <= fArraySize/2; j++)
 241         {
 242           xi = x + i;
 243           zj = z + j;
 244           w = TMath::Max( (Float_t)0., (Float_t)(w0 + log( fEnergyArray[x+i][z+j] / fEnergyArray[x][z]) ) ) ;
 245           xt += xi * w ;
 246           zt += zj * w ;
 247           wtot += w ;
 248
 249         }
 250     }
 251   xt /= wtot ;
 252   zt /= wtot ;
 253
 254   recPointPtr->fX = xt;
 255   recPointPtr->fZ = zt;
 256
 257   return 0;
 258
 259 }//end CalculateCenterOfGravity
 260
 261
 262
 263 /**
 264 * Create a simpler data struct for a rec point, containing
 265 * only coordinates, energy and timing
 266 * @param recPointPtr pointer to the rec point
 267 * @clusterStructPtr pointer to the cluster struct
 268 **/
 269 Int_t
 270 AliHLTPHOSClusterizer::ClusterizeStruct(AliHLTPHOSRecPointDataStruct* recPointPtr, AliHLTPHOSClusterDataStruct* clusterStructPtr)
 271 {
 272   //Simplify the rec points
 273
 274   Float_t clusterEnergy = 0;
 275   Float_t* energiesListPtr = recPointPtr->fEnergiesListPtr;
 276
 277   for(Int_t i = 0; i < recPointPtr->fMultiplicity; i++)
 278     {
 279       clusterEnergy += energiesListPtr[i];
 280     }
 281
 282   clusterStructPtr->fLocalPositionPtr[0] = recPointPtr->fX;
 283   clusterStructPtr->fLocalPositionPtr[1] = recPointPtr->fZ;
 284   clusterStructPtr->fClusterEnergy = clusterEnergy;
 285   clusterStructPtr->fPHOSModule = recPointPtr->fPHOSModule;
 286
 287   return 0;
 288
 289 }//end ClusterizeStruct
 290
 291
 292
 293
 294 /**
 295 * Resets the cell energy array
 296 **/
 297 Int_t
 298 AliHLTPHOSClusterizer::ResetCellEnergyArray()
 299
 300 {
 301   //Reset the cell energy array
 302
 303   for(Int_t x = 0; x < N_ROWS_MOD; x++)
 304     {
 305       for(Int_t z = 0; z < N_COLUMNS_MOD; z++)
 306         {
 307           fEnergyArray[x][z] = 0;
 308         }
 309     }
 310
 311   return 0;
 312
 313 }//end ResetCellEnergyArray
 314