X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=HLT%2FPHOS%2FAliHLTPHOSClusterizer.cxx;h=3cf7a23ba9f8b47feae01ab4764a3dc3072a1305;hb=1af38adf7128b4501530551fa272929028f63e49;hp=377f58eec576d49defeebed89fb024689918b191;hpb=9be2600f6614107c11a1d269dbac57242e5854f4;p=u%2Fmrichter%2FAliRoot.git diff --git a/HLT/PHOS/AliHLTPHOSClusterizer.cxx b/HLT/PHOS/AliHLTPHOSClusterizer.cxx index 377f58eec57..3cf7a23ba9f 100644 --- a/HLT/PHOS/AliHLTPHOSClusterizer.cxx +++ b/HLT/PHOS/AliHLTPHOSClusterizer.cxx @@ -1,314 +1,185 @@ /************************************************************************** - * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * This file is property of and copyright by the ALICE HLT Project * + * All rights reserved. * * * - * Authors: Øystein Djuvsland * + * Primary Authors: Oystein Djuvsland * * * * Permission to use, copy, modify and distribute this software and its * * documentation strictly for non-commercial purposes is hereby granted * * without fee, provided that the above copyright notice appears in all * * copies and that both the copyright notice and this permission notice * * appear in the supporting documentation. The authors make no claims * - * about the suitability of this software for any purpose. It is * + * about the suitability of this software for any purpose. It is * * provided "as is" without express or implied warranty. * **************************************************************************/ +/** + * @file AliHLTPHOSClusterizer.cxx + * @author Oystein Djuvsland + * @date + * @brief Clusterizer for PHOS HLT + */ -/** @file AliHLTPHOSClusterizer.cxx - @author Øystein Djuvsland - @date - @brief A temporary clusterizer for PHOS -*/ +// see header file for class documentation +// or +// refer to README to build package +// or +// visit http://web.ift.uib.no/~kjeks/doc/alice-hlt - - -#include "AliHLTPHOSPhysicsDefinitions.h" #include "AliHLTPHOSClusterizer.h" -#include "AliHLTPHOSCommonDefs.h" -#include "TVector3.h" +#include "AliHLTPHOSBase.h" +#include "AliHLTLogging.h" #include "TMath.h" -#include "AliHLTPHOSRcuCellEnergyDataStruct.h" -#include "AliHLTPHOSRecPointListDataStruct.h" -#include "AliHLTPHOSValidCellDataStruct.h" +#include "AliHLTPHOSRecPointContainerStruct.h" #include "AliHLTPHOSRecPointDataStruct.h" -#include "AliHLTPHOSClusterDataStruct.h" - +#include "AliHLTPHOSDigitDataStruct.h" +#include "AliHLTPHOSDigitContainerDataStruct.h" +#include "TClonesArray.h" +#include "AliPHOSGeometry.h" +#include "AliPHOSDigit.h" +#ifndef HAVE_NOT_PHOSRECOPARAMEMC // set from configure if EMC functionality not available in AliPHOSRecoParam +#include "AliPHOSRecoParam.h" +#else +#include "AliPHOSRecoParamEmc.h" +#endif +#include + +using namespace std; ClassImp(AliHLTPHOSClusterizer); -/** -* Main constructor -**/ -AliHLTPHOSClusterizer::AliHLTPHOSClusterizer():fPHOSModule(-1), fThreshold(0), fClusterThreshold(0), - fHighGainFactor(0.005), fLowGainFactor(0.08), - fArraySize(3), fMultiplicity(fArraySize*fArraySize) -{ - //Main constructor - +AliHLTPHOSClusterizer::AliHLTPHOSClusterizer(): + AliHLTPHOSBase(), + fEmcClusteringThreshold(0), + fEmcMinEnergyThreshold(0), + fEmcTimeGate(0), + fLogWeight(0), + fDigitsInCluster(0), + fOnlineMode(true), + fDigitArrayPtr(0), + fEmcRecPointsPtr(0), + fDigitPtr(0), + fDigitContainerPtr(0), + fRecPointContainerPtr(0), + fPHOSGeometry(0), + fGetterPtr(0), + fMaxDigitIndexDiff(2*N_ZROWS_MOD) + { + //See header file for documentation + fPHOSGeometry = AliPHOSGeometry::GetInstance("noCPV"); + fEmcClusteringThreshold = 0.2; + fEmcMinEnergyThreshold = 0.03; + fEmcTimeGate = 1.e-6 ; + fLogWeight = 4.5; }//end -AliHLTPHOSClusterizer::AliHLTPHOSClusterizer(const AliHLTPHOSClusterizer &):fPHOSModule(-1), fThreshold(0), fClusterThreshold(0), - fHighGainFactor(0.005), fLowGainFactor(0.08), - fArraySize(3), fMultiplicity(fArraySize*fArraySize) -{ - //Copy constructor, not implemented -}//end -AliHLTPHOSClusterizer:: ~AliHLTPHOSClusterizer() +AliHLTPHOSClusterizer::~AliHLTPHOSClusterizer() { - //Destructor + //See header file for documentation } -/** -* Building a 2D array of cell energies of the PHOS detector -* @param cellData object containing the cell energies from one event -* @param recPointList list to be filled with coordinates of local maxima in the detector -**/ +void +AliHLTPHOSClusterizer::SetRecPointContainer(AliHLTPHOSRecPointContainerStruct* recPointContainerPtr) + { + fRecPointContainerPtr = recPointContainerPtr; + fRecPointContainerPtr->fNRecPoints = 0; + } -Int_t -AliHLTPHOSClusterizer::BuildCellEnergyArray(AliHLTPHOSRcuCellEnergyDataStruct* cellData, - AliHLTPHOSRecPointListDataStruct* recPointList) +void +AliHLTPHOSClusterizer::SetRecoParameters(AliPHOSRecoParam* params) { - //Build the cell energy array of the detector - - Int_t x = 0; - Int_t z = 0; - Int_t gain = 0; - Int_t xMod = 0; - Int_t zMod = 0; - Int_t index = 0; - Double_t energyCount = 0; - - for(Int_t cell = 0; cell < cellData->fCnt; cell++) - { - - z = (cellData->fValidData[cell]).fZ; - x = (cellData->fValidData[cell]).fX; - gain = (cellData->fValidData[cell]).fGain; - - zMod = z + (cellData->fRcuZ)*N_ROWS_RCU; - xMod = x + (cellData->fRcuX)*N_COLUMNS_RCU; - - energyCount = (cellData->fValidData[cell]).fEnergy; - - if(gain == 0 && energyCount < 1023) - { - fEnergyArray[xMod][zMod] = fHighGainFactor * energyCount; - - if(fEnergyArray[xMod][zMod] > fClusterThreshold) - { - recPointList[index].fZ = zMod; - recPointList[index].fX = xMod; - - for(Int_t j = 0; j < index; j++) - { - if(recPointList[j].fZ == zMod && recPointList[j].fX == xMod) - recPointList[j].fZ = -1; - } - index++; - } - - } - else if(gain == 1 && fEnergyArray[xMod][zMod] == 0) - { - fEnergyArray[xMod][zMod] = fLowGainFactor * energyCount; - if(fEnergyArray[xMod][zMod] > fClusterThreshold) - { - recPointList[index].fZ = zMod; - recPointList[index].fX = xMod; - recPointList[index].fModule = cellData->fModuleID; - index++; - } - } + //see header file for documentation +#ifndef HAVE_NOT_PHOSRECOPARAMEMC // set from configure if EMC functionality not available in AliPHOSRecoParam + // the new AliPHOSRecoParam functions, available from revision + // fEmcClusteringThreshold = params->GetEMCClusteringThreshold(); + // fEmcMinEnergyThreshold = params->GetEMCMinE(); + // fLogWeight = params->GetEMCLogWeight(); +#else + fEmcClusteringThreshold = params->GetClusteringThreshold(); + fEmcMinEnergyThreshold = params->GetMinE(); + fLogWeight = params->GetLogWeight(); +#endif +} - } - - fPHOSModule = cellData->fModuleID; - - return index; -}//end BuildCellEnergyArray - - - -/** -* Creating an array of rec points -* @param recPointStructArrayPtr array to store the rec points -* @param list list of rec points -* @param nPoints number of points -**/ Int_t -AliHLTPHOSClusterizer::CreateRecPointStructArray(AliHLTPHOSRecPointDataStruct* recPointStructArrayPtr, - AliHLTPHOSRecPointListDataStruct* list, - Int_t nPoints) - +AliHLTPHOSClusterizer::ClusterizeEvent() { - //Create the rec point struct array - - Int_t flag = 0; - Int_t edgeFlagRows = 0; - Int_t edgeFlagCols = 0; - Int_t k = 0; + //see header file for documentation Int_t nRecPoints = 0; - Int_t z = 0; - Int_t x = 0; - - Float_t* energiesList = NULL; - - for(Int_t point = 0; point < nPoints; point++) - { - z = list[point].fZ; - x = list[point].fX; - if(z == -1) continue; - - if((z-fArraySize/2) < 0/*= - N_ROWS_MOD/2*/ || (z+fArraySize/2) >= N_ROWS_MOD/*/2*/) + UInt_t i = 0; + AliHLTPHOSRecPointDataStruct *recPoint = 0; + //Clusterization starts + for(i = 0; i < fDigitContainerPtr->fNDigits; i++) + { + fDigitsInCluster = 0; + + if(fDigitContainerPtr->fDigitDataStruct[i].fEnergy < fEmcClusteringThreshold) { - edgeFlagRows = 1; continue; } - if((x-fArraySize/2) < 0/*= - N_COLUMNS_MOD/2*/ || (x+fArraySize/2) >= N_COLUMNS_MOD) + recPoint = &(fRecPointContainerPtr->fRecPointArray[nRecPoints]); + recPoint->fAmp = 0; + recPoint->fModule = fDigitContainerPtr->fDigitDataStruct[i].fModule; + recPoint->fDigitsList[fDigitsInCluster] = fDigitContainerPtr->fDigitDataStruct[i]; + recPoint->fAmp += fDigitContainerPtr->fDigitDataStruct[i].fEnergy; + fDigitContainerPtr->fDigitDataStruct[i].fEnergy = 0; + fDigitsInCluster++; + nRecPoints++; + if(nRecPoints == 100) { - edgeFlagCols = 1; - continue; + // HLTWarning("Too many rec points in event. Aborting clusterisation"); + break; } - - - if(!flag) energiesList = new Float_t[fMultiplicity]; - flag = 0; - k = 0; - - for(Int_t i = -fArraySize/2; i <= fArraySize/2; i++) - { - if(flag) break; - for(Int_t j = -fArraySize/2; j <= fArraySize/2; j++) - { - - if(fEnergyArray[x+i][z+j] > fEnergyArray[x][z] && abs(i) < 2 && abs(j) < 2) - { - flag = 1; - break; - } - energiesList[k] = fEnergyArray[x+i][z+j]; - k++; - } - } - - - if(!flag && k) - { - recPointStructArrayPtr[nRecPoints].fEnergiesListPtr = energiesList; - recPointStructArrayPtr[nRecPoints].fCoordinatesPtr[0] = x; - recPointStructArrayPtr[nRecPoints].fCoordinatesPtr[1] = z; - recPointStructArrayPtr[nRecPoints].fX = x; - recPointStructArrayPtr[nRecPoints].fZ = z; - recPointStructArrayPtr[nRecPoints].fMultiplicity = fMultiplicity; - recPointStructArrayPtr[nRecPoints].fPHOSModule = list[point].fModule; - nRecPoints++; - } - } - - - - energiesList = NULL; - + ScanForNeighbourDigits(i, recPoint); + recPoint->fMultiplicity = fDigitsInCluster; + + }//end of clusterization + fRecPointContainerPtr->fNRecPoints = nRecPoints; return nRecPoints; - -}//end CreateRecPointStructArray +} +void +AliHLTPHOSClusterizer::ScanForNeighbourDigits(Int_t index, AliHLTPHOSRecPointDataStruct* recPoint) -/** -* Calculating the center of gravity of a rec point -* Not working well at this point! -* @param recPointPtr pointer to the rec point -**/ -Int_t -AliHLTPHOSClusterizer::CalculateCenterOfGravity(AliHLTPHOSRecPointDataStruct* recPointPtr) { - //Calculate the center of gravity - - Float_t xt = 0; - Float_t zt = 0; - Float_t xi = 0; - Float_t zj = 0; - Float_t w = 0; - Float_t w0 = 4.5; - Float_t wtot = 0; - - Int_t x = recPointPtr->fCoordinatesPtr[0]; - Int_t z = recPointPtr->fCoordinatesPtr[1]; - - for(Int_t i = -fArraySize/2; i <= fArraySize/2; i++) + //see header file for documentation + UInt_t max = TMath::Min((Int_t)fDigitContainerPtr->fNDigits, (Int_t)fMaxDigitIndexDiff+index); + UInt_t min = TMath::Max(0, (Int_t)(index - (Int_t)fMaxDigitIndexDiff)); + for(UInt_t j = min; j < max; j++) { - for(Int_t j = -fArraySize/2; j <= fArraySize/2; j++) - { - xi = x + i; - zj = z + j; - w = TMath::Max( (Float_t)0., (Float_t)(w0 + log( fEnergyArray[x+i][z+j] / fEnergyArray[x][z]) ) ) ; - xt += xi * w ; - zt += zj * w ; - wtot += w ; - + if(fDigitContainerPtr->fDigitDataStruct[j].fEnergy > fEmcMinEnergyThreshold) + { + if(AreNeighbours(&(fDigitContainerPtr->fDigitDataStruct[index]), + &(fDigitContainerPtr->fDigitDataStruct[j]))) + { + recPoint->fDigitsList[fDigitsInCluster] = fDigitContainerPtr->fDigitDataStruct[j]; + recPoint->fAmp += fDigitContainerPtr->fDigitDataStruct[j].fEnergy; + fDigitContainerPtr->fDigitDataStruct[j].fEnergy = 0; + fDigitsInCluster++; + ScanForNeighbourDigits(j, recPoint); + } } - } - xt /= wtot ; - zt /= wtot ; - - recPointPtr->fX = xt; - recPointPtr->fZ = zt; - - return 0; - -}//end CalculateCenterOfGravity - - - -/** -* Create a simpler data struct for a rec point, containing -* only coordinates, energy and timing -* @param recPointPtr pointer to the rec point -* @param clusterStructPtr pointer to the cluster struct -**/ -Int_t -AliHLTPHOSClusterizer::ClusterizeStruct(AliHLTPHOSRecPointDataStruct* recPointPtr, AliHLTPHOSClusterDataStruct* clusterStructPtr) -{ - //Simplify the rec points - - Float_t clusterEnergy = 0; - Float_t* energiesListPtr = recPointPtr->fEnergiesListPtr; - - for(Int_t i = 0; i < recPointPtr->fMultiplicity; i++) - { - clusterEnergy += energiesListPtr[i]; - } - - clusterStructPtr->fLocalPositionPtr[0] = recPointPtr->fX; - clusterStructPtr->fLocalPositionPtr[1] = recPointPtr->fZ; - clusterStructPtr->fClusterEnergy = clusterEnergy; - clusterStructPtr->fPHOSModule = recPointPtr->fPHOSModule; - - return 0; - -}//end ClusterizeStruct - - - + } + return; +} -/** -* Resets the cell energy array -**/ Int_t -AliHLTPHOSClusterizer::ResetCellEnergyArray() - +AliHLTPHOSClusterizer::AreNeighbours(AliHLTPHOSDigitDataStruct* digit1, + AliHLTPHOSDigitDataStruct* digit2) { - //Reset the cell energy array - - for(Int_t x = 0; x < N_ROWS_MOD; x++) - { - for(Int_t z = 0; z < N_COLUMNS_MOD; z++) + //see header file for documentation + if ( (digit1->fModule == digit2->fModule) /*&& (coord1[1]==coord2[1])*/ ) // inside the same PHOS module + { + Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ ); + Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX ); + if (( coldiff <= 1 && rowdiff < 1 ) || ( coldiff < 1 && rowdiff <= 1 )) { - fEnergyArray[x][z] = 0; + if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate) + { + return 1; + } } } - return 0; - -}//end ResetCellEnergyArray - +}