-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Authors: Øystein Djuvsland <oysteind@ift.uib.no> *
- * *
- * 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 *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-
-
-/** @file AliHLTPHOSClusterizer.cxx
- @author Øystein Djuvsland
- @date
- @brief A temporary clusterizer for PHOS
-*/
-
+//Insert copyright
-//#include "AliHLTPHOSDefinitions.h"
#include "AliHLTPHOSClusterizer.h"
-//#include "AliHLTPHOSCommonDefs.h"
-#include "TVector3.h"
+#include "AliHLTPHOSBase.h"
#include "TMath.h"
-#include "AliHLTPHOSRcuCellEnergyDataStruct.h"
-#include "AliHLTPHOSRecPointListDataStruct.h"
-#include "AliHLTPHOSValidCellDataStruct.h"
+#include "AliHLTPHOSRecPointContainerStruct.h"
#include "AliHLTPHOSRecPointDataStruct.h"
-#include "AliHLTPHOSClusterDataStruct.h"
-//#include "AliHLTPHOSConstants.h"
-
-//using namespace PhosHLTConst;
-
-
-
-ClassImp(AliHLTPHOSClusterizer)
-
-/**
-* Main constructor
-**/
-AliHLTPHOSClusterizer::AliHLTPHOSClusterizer():AliHLTPHOSBase(),fPHOSModule(0), fThreshold(0), fClusterThreshold(0),
- fHighGainFactor(0.005), fLowGainFactor(0.08),
- fArraySize(3), fMultiplicity(fArraySize*fArraySize)
-//AliHLTPHOSClusterizer::AliHLTPHOSClusterizer():AliHLTPHOSProcessor(),fPHOSModule(0), fThreshold(0), fClusterThreshold(0),
-// fHighGainFactor(0.005), fLowGainFactor(0.08),
-// fArraySize(3), fMultiplicity(fArraySize*fArraySize)
+#include "AliHLTPHOSDigitDataStruct.h"
+#include "AliHLTPHOSDigitContainerDataStruct.h"
+#include "TClonesArray.h"
+#include "AliPHOSGeometry.h"
+#include "AliPHOSDigit.h"
+#include "AliPHOSRecoParamEmc.h"
+
+ClassImp(AliHLTPHOSClusterizer);
+
+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)
{
//See header file for documentation
+ fPHOSGeometry = AliPHOSGeometry::GetInstance("noCPV");
+ fEmcClusteringThreshold = 0.2;
+ fEmcMinEnergyThreshold = 0.03;
+ fEmcTimeGate = 1.e-8 ;
+ fLogWeight = 4.5;
+
}//end
-// PTH AliHLTPHOSClusterizer::AliHLTPHOSClusterizer(const AliHLTPHOSClusterizer &):AliHLTPHOSBase(),fPHOSModule(0), fThreshold(0), fClusterThreshold(0),
-// fHighGainFactor(0.005), fLowGainFactor(0.08),
-// fArraySize(3), fMultiplicity(fArraySize*fArraySize)
-//{
+AliHLTPHOSClusterizer::AliHLTPHOSClusterizer(const AliHLTPHOSClusterizer &)
+{
//Copy constructor, not implemented
-//}//end
+}//end
-AliHLTPHOSClusterizer:: ~AliHLTPHOSClusterizer()
+AliHLTPHOSClusterizer::~AliHLTPHOSClusterizer()
{
//Destructor
}
-/**
-* 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::SetRecoParameters(AliPHOSRecoParamEmc* params)
+{
+ fEmcClusteringThreshold = params->GetClusteringThreshold();
+ fEmcMinEnergyThreshold = params->GetMinE();
+ fLogWeight = params->GetLogWeight();
+}
-int
-AliHLTPHOSClusterizer::BuildCellEnergyArray(AliHLTPHOSRcuCellEnergyDataStruct* cellData,
- AliHLTPHOSRecPointListDataStruct* recPointList)
+void
+AliHLTPHOSClusterizer::SetOfflineMode(AliPHOSGetter* getter)
{
- //BuildCellEnergyArray
+ fRecPointContainerPtr = new AliHLTPHOSRecPointContainerStruct();
+ fDigitContainerPtr = new AliHLTPHOSDigitContainerDataStruct();
+ fGetterPtr = getter;
+ fDigitArrayPtr = fGetterPtr->Digits();
+ fEmcRecPointsPtr = fGetterPtr->EmcRecPoints();
+ fOnlineMode = false;
+}
- 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;
+Int_t
+AliHLTPHOSClusterizer::GetEvent(Int_t i)
+{
+ Int_t coord[4];
- for(Int_t cell = 0; cell < cellData->fCnt; cell++)
+ fGetterPtr->Event(i, "D");
+ for(Int_t j = 0; j < fDigitArrayPtr->GetEntries(); j++)
{
-
- z = (cellData->fValidData[cell]).fZ;
- x = (cellData->fValidData[cell]).fX;
- gain = (cellData->fValidData[cell]).fGain;
-
- zMod = z + (cellData->fRcuZ)*N_ZROWS_RCU;
- xMod = x + (cellData->fRcuX)*N_XCOLUMNS_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++;
- }
- }
+ fDigitPtr = (AliPHOSDigit*)fDigitArrayPtr->At(j);
+ fPHOSGeometry->AbsToRelNumbering(fDigitPtr->GetId(), coord);
+ fDigitContainerPtr->fDigitDataStruct[j].fX = coord[3];
+ fDigitContainerPtr->fDigitDataStruct[j].fZ = coord[2];
+ fDigitContainerPtr->fDigitDataStruct[j].fModule = coord[0];
+ fDigitContainerPtr->fDigitDataStruct[j].fEnergy = fDigitPtr->GetEnergy();
+ fDigitContainerPtr->fDigitDataStruct[j].fTime = fDigitPtr->GetTime();
+ }
+ fDigitContainerPtr->fNDigits = fDigitArrayPtr->GetEntriesFast();
+ return 0;
+}
- }
+Int_t
+AliHLTPHOSClusterizer::GetNEvents()
+{
+ if(fOnlineMode)
+ {
+ printf("Number of events not available in online mode!\n");
+ return -1;
+ }
+ return fGetterPtr->MaxEvent();
+}
- 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
-AliHLTPHOSClusterizer::CreateRecPointStructArray(AliHLTPHOSRecPointDataStruct* recPointStructArrayPtr,
- AliHLTPHOSRecPointListDataStruct* list,
- int nPoints)
+Int_t
+AliHLTPHOSClusterizer::ClusterizeEvent()
{
- //CreateRecPointStructArray
-
- Int_t flag = 0;
- Int_t edgeFlagRows = 0;
- Int_t edgeFlagCols = 0;
- Int_t k = 0;
+ // Steering method to construct the clusters stored in a list of Reconstructed Points
+ // A cluster is defined as a list of neighbour digits
Int_t nRecPoints = 0;
- Int_t z = 0;
- Int_t x = 0;
- Int_t i = 0;
- Int_t j = 0;
- Int_t a = 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_ZROWS_MOD/2*/ || (z+fArraySize/2) >= N_ZROWS_MOD/*/2*/)
- {
- edgeFlagRows = 1;
- continue;
- }
- if((x-fArraySize/2) < 0/*= - N_COLUMNS_MOD/2*/ || (x+fArraySize/2) >= N_XCOLUMNS_MOD)
- {
- edgeFlagCols = 1;
- continue;
- }
-
+ UInt_t i = 0;
- if(!flag) energiesList = new Float_t[fMultiplicity];
- flag = 0;
- k = 0;
-
- for(i = -fArraySize/2; i <= fArraySize/2; i++)
- {
- if(flag) break;
- for(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++;
- }
- }
-
+ AliHLTPHOSRecPointDataStruct *recPoint = 0;
+
+ //Clusterization starts
+ for(i = 0; i < fDigitContainerPtr->fNDigits; i++)
+ {
+
+ fDigitsInCluster = 0;
- if(!flag && k)
+ if(fDigitContainerPtr->fDigitDataStruct[i].fEnergy < fEmcClusteringThreshold)
{
- for(a = 0; a < k; a++)
- {
- (recPointStructArrayPtr[nRecPoints].fEnergiesListPtr)[a] = energiesList[a];
- }
- 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++;
+ continue;
}
- }
- if(energiesList)
- {
- delete [] energiesList;
- energiesList = NULL;
- }
- return nRecPoints;
+ recPoint = &(fRecPointContainerPtr->fRecPointArray[nRecPoints]);
+ recPoint->fAmp = 0;
+ //recPoint->
+ recPoint->fDigitsList[fDigitsInCluster] = fDigitContainerPtr->fDigitDataStruct[i];
+ recPoint->fAmp += fDigitContainerPtr->fDigitDataStruct[i].fEnergy;
+ fDigitContainerPtr->fDigitDataStruct[i].fEnergy = 0;
+ fDigitsInCluster++;
+ nRecPoints++;
+ ScanForNeighbourDigits(i, recPoint);
+
+ recPoint->fMultiplicity = fDigitsInCluster;
+ }//end of clusterization
+ fRecPointContainerPtr->fNRecPoints = nRecPoints;
-}//end CreateRecPointStructArray
+ return nRecPoints;
+}
+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
-AliHLTPHOSClusterizer::CalculateCenterOfGravity(AliHLTPHOSRecPointDataStruct* recPointPtr)
{
- //CalculateCenterOfGravity
- //Copied from offline code, and modified
-
- 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++)
+
+ for(UInt_t j = 0; j < fDigitContainerPtr->fNDigits; 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]) ) ) ;
- //printf("log in cog: %f\n", log( fEnergyArray[x+i][z+j] / fEnergyArray[x][z]));
- xt += xi * w ;
- zt += zj * w ;
- wtot += w ;
+ if(fDigitContainerPtr->fDigitDataStruct[j].fEnergy > fEmcMinEnergyThreshold)
+ {
+ switch(AreNeighbours(&(fDigitContainerPtr->fDigitDataStruct[index]),
+ &(fDigitContainerPtr->fDigitDataStruct[j])))
+ {
+ case 0:
+ break;
+ case 1:
+ recPoint->fDigitsList[fDigitsInCluster] = fDigitContainerPtr->fDigitDataStruct[j];
+ recPoint->fAmp += fDigitContainerPtr->fDigitDataStruct[j].fEnergy;
+ fDigitContainerPtr->fDigitDataStruct[j].fEnergy = 0;
+ fDigitsInCluster++;
+ ScanForNeighbourDigits(j, recPoint);
+ break;
+ case 2:
+ break;
+ }
}
- }
- xt /= wtot ;
- zt /= wtot ;
-
- //printf("wtot cog: %f\n", wtot);
+ }
+ return;
+}
- recPointPtr->fX = xt;
- recPointPtr->fZ = zt;
-
- return 0;
-}//end CalculateCenterOfGravity
-
-int
-AliHLTPHOSClusterizer::CalculateMoments(AliHLTPHOSRecPointDataStruct* recPointPtr, Bool_t axisOnly)
+Int_t
+AliHLTPHOSClusterizer::AreNeighbours(AliHLTPHOSDigitDataStruct* digit1,
+ AliHLTPHOSDigitDataStruct* digit2)
{
- //Copied from offline code, and modified
+ //Int_t coord1[4];
+ //Int_t coord2[4];
- // Calculate the shower moments in the eigen reference system
- // M2x, M2z, M3x, M4z
- // Calculate the angle between the shower position vector and the eigen vector
+ // fPHOSGeometry->AbsToRelNumbering(digit1->fID, coord1);
+ //fPHOSGeometry->AbsToRelNumbering(digit2->fID, coord2);
- Double_t wtot = 0. ;
- Double_t x = 0.;
- Double_t z = 0.;
- Double_t dxx = 0.;
- Double_t dzz = 0.;
- Double_t dxz = 0.;
- Double_t lambda0=0, lambda1=0;
- Float_t logWeight = 4.5;
-
- //Int_t iDigit;
-
- // 1) Find covariance matrix elements:
- // || dxx dxz ||
- // || dxz dzz ||
-
- Int_t xi;
- Int_t zj;
- Double_t w;
-
- Int_t xc = recPointPtr->fCoordinatesPtr[0];
- Int_t zc = recPointPtr->fCoordinatesPtr[1];
-
- for(Int_t i = -fArraySize/2; i <= fArraySize/2; i++)
- {
- xi = xc + i;
- for(Int_t j = -fArraySize/2; j <= fArraySize/2; j++)
- {
- zj = zc + j;
- if (fEnergyArray[xi][zj] > 0)
- {
- // w = TMath::Max(0.,logWeight + log( fEnergyArray[xi][zj]/fEnergyArray[xc][zc] ) ) ;
+ if ( (digit1->fModule == digit2->fModule) /*&& (coord1[1]==coord2[1])*/ ) // inside the same PHOS module
+ {
- //printf("log in mom: %f\n", TMath::Log( fEnergyArray[xi][zj] / fEnergyArray[xc][zc]));
- x += w * xi ;
- z += w * zj ;
- dxx += w * xi * xi ;
- dzz += w * zj * zj ;
- dxz += w * xi * zj ;
- wtot += w ;
- }
- }
- }
- //printf("wtot: %f\n", wtot);
- if (wtot>0)
- {
- x /= wtot ;
- z /= wtot ;
- dxx /= wtot ;
- dzz /= wtot ;
- dxz /= wtot ;
- dxx -= x * x ;
- dzz -= z * z ;
- dxz -= x * z ;
-
- // 2) Find covariance matrix eigen values lambda0 and lambda1
-
- lambda0 = 0.5 * (dxx + dzz) + TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz ) ;
- lambda1 = 0.5 * (dxx + dzz) - TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz ) ;
- }
-
- recPointPtr->fM2x = lambda0;
- recPointPtr->fM2z = lambda1;
- // printf("lambda0 = %f -- lambda1 = %f\n", lambda0, lambda1);
- // 3) Find covariance matrix eigen vectors e0 and e1
- if(!axisOnly)
- {
- TVector2 e0,e1;
- if (dxz != 0)
- e0.Set(1.,(lambda0-dxx)/dxz);
- else
- e0.Set(0.,1.);
-
- e0 = e0.Unit();
- e1.Set(-e0.Y(),e0.X());
-
- // 4) Rotate cluster tensor from (x,z) to (e0,e1) system
- // and calculate moments M3x and M4z
-
- Float_t cosPhi = e0.X();
- Float_t sinPhi = e0.Y();
-
- Float_t xiPHOS ;
- Float_t zjPHOS ;
- Double_t dx3, dz3, dz4;
- wtot = 0.;
- x = 0.;
- z = 0.;
- dxx = 0.;
- dzz = 0.;
- dxz = 0.;
- dx3 = 0.;
- dz3 = 0.;
- dz4 = 0.;
- for(Int_t i = -fArraySize/2; i <= fArraySize/2; i++)
+ Int_t rowdiff = TMath::Abs( digit1->fZ - digit2->fZ );
+ Int_t coldiff = TMath::Abs( digit1->fX - digit2->fX );
+
+ if (( coldiff <= 1 ) && ( rowdiff <= 1 ))
{
- for(Int_t j = -fArraySize/2; j <= fArraySize/2; j++)
- {
- xi = xc + i;
- zj = zc + j;
- xiPHOS = (Float_t)xi*cosPhi + (Float_t)zj*sinPhi;
- zjPHOS = (Float_t)zj*cosPhi - (Float_t)xi*sinPhi;
- // xi = (Float_t)xi*cosPhi + (Float_t)zj*sinPhi;
- // zj = (Float_t)zj*cosPhi - (Float_t)xi*sinPhi;
- if (fEnergyArray[xi][zj] > 0)
- {
- w = TMath::Max(0.,logWeight+TMath::Log(fEnergyArray[xi][zj]/fEnergyArray[xc][zc] ) ) ;
- // printf("log in mom: %f\n", TMath::Log( fEnergyArray[xi][zj] / fEnergyArray[xc][zc]));
-
- x += w * xi ;
- z += w * zj ;
- dxx += w * xi * xi ;
- dzz += w * zj * zj ;
- dxz += w * xi * zj ;
- dx3 += w * xi * xi * xi;
- dz3 += w * zj * zj * zj ;
- dz4 += w * zj * zj * zj * zj ;
- wtot += w ;
- }
- }
+ if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)
+ return 1;
}
- if (wtot>0)
- {
- x /= wtot ;
- z /= wtot ;
- dxx /= wtot ;
- dzz /= wtot ;
- dxz /= wtot ;
- dx3 /= wtot ;
- dz3 /= wtot ;
- dz4 /= wtot ;
- dx3 += -3*dxx*x + 2*x*x*x;
- dz4 += -4*dz3*z + 6*dzz*z*z -3*z*z*z*z;
- dxx -= x * x ;
- dzz -= z * z ;
- dxz -= x * z ;
- }
-
- // 5) Find an angle between cluster center vector and eigen vector e0
-
- Float_t phi = 0;
- if ( (x*x+z*z) > 0 )
- phi = TMath::ACos ((x*e0.X() + z*e0.Y()) / sqrt(x*x + z*z));
-
- recPointPtr->fM3x = dx3;
- recPointPtr->fM4z = dz4;
- recPointPtr->fPhixe = phi;
- // printf("%f\n",phi);
}
return 0;
}
-
-/**
-* 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
-AliHLTPHOSClusterizer::ClusterizeStruct(AliHLTPHOSRecPointDataStruct* recPointPtr, AliHLTPHOSClusterDataStruct* clusterStructPtr)
+void
+AliHLTPHOSClusterizer::CalculateCenterOfGravity()
{
- //ClusterizeStruct
-
- 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
+ // Calculates the center of gravity in the local PHOS-module coordinates
+ Float_t wtot = 0.;
+
+ Int_t relid[4];
+ Float_t x = 0.;
+ Float_t z = 0.;
+ Float_t xi = 0.;
+ Float_t zi = 0.;
+ AliHLTPHOSRecPointDataStruct *recPoint = 0;
+ AliHLTPHOSDigitDataStruct *digit = 0;
-/**
-* Resets the cell energy array
-**/
-int
-AliHLTPHOSClusterizer::ResetCellEnergyArray()
+ AliPHOSGeometry * phosgeom = AliPHOSGeometry::GetInstance() ;
-{
- //ResetCellEnergyArray
+ Int_t iDigit = 0;
+ Int_t iRecPoint = 0;
- for(Int_t x = 0; x < N_ZROWS_MOD; x++)
+ for(iRecPoint=0; iRecPoint<fRecPointContainerPtr->fNRecPoints; iRecPoint++)
{
- for(Int_t z = 0; z < N_XCOLUMNS_MOD; z++)
+ recPoint = &(fRecPointContainerPtr->fRecPointArray[iRecPoint]);
+ for(iDigit = 0; iDigit < recPoint->fMultiplicity; iDigit++)
{
- fEnergyArray[x][z] = 0;
+ digit = &(recPoint->fDigitsList[iDigit]);
+
+ //fPHOSGeometry->AbsToRelNumbering(digit->fID, relid) ;
+ // fPHOSGeometry->RelPosInModule(relid, xi, zi);
+ xi = digit->fX;
+ zi = digit->fZ;
+
+ if (recPoint->fAmp > 0 && digit->fEnergy > 0)
+ {
+ Float_t w = TMath::Max( 0., fLogWeight + TMath::Log( digit->fEnergy / recPoint->fAmp ) ) ;
+ x += xi * w ;
+ z += zi * w ;
+ wtot += w ;
+ }
+ }
+
+ if (wtot>0)
+ {
+ recPoint->fX = x/wtot ;
+ recPoint->fZ = z/wtot ;
+ }
+ else
+ {
+ recPoint->fAmp = 0;
}
}
+
+}
- return 0;
-
-}//end ResetCellEnergyArray
#ifndef ALIHLTPHOSCLUSTERIZER_H
#define ALIHLTPHOSCLUSTERIZER_H
+#include "AliHLTPHOSBase.h"
+#include "AliHLTPHOSRecPointContainerStruct.h"
+#include "AliHLTPHOSDigitContainerDataStruct.h"
+#include "AliPHOSGeometry.h"
+#include "TClonesArray.h"
+#include "AliPHOSDigit.h"
+#include "AliPHOSGetter.h"
+#include "AliPHOSRecoParamEmc.h"
-
-//#include "AliHLTPHOSProcessor.h"
-#include "AliHLTPHOSBase.h"
-
-
-//#include "AliHLTPHOSCommonDefs.h"
-//#include "AliHLTPHOSConstants.h"
-
-//using namespace PhosHLTConst;
-
-struct AliHLTPHOSClusterDataStruct;
-struct AliHLTPHOSRecPointDataStruct;
-struct AliHLTPHOSValidCellDataStruct;
-struct AliHLTPHOSRecPointListDataStruct;
-struct AliHLTPHOSRcuCellEnergyDataStruct;
-
-//class AliHLTPHOSClusterizer: public AliHLTPHOSProcessor
-
-class AliHLTPHOSClusterizer: public AliHLTPHOSBase
+class AliHLTPHOSClusterizer : public AliHLTPHOSBase
{
- public:
+public:
+
+ AliHLTPHOSClusterizer();
- AliHLTPHOSClusterizer();
virtual ~AliHLTPHOSClusterizer();
- // AliHLTPHOSClusterizer(const AliHLTPHOSClusterizer &);
- //AliHLTPHOSClusterizer & operator = (const AliHLTPHOSClusterizer &) {return *this;}
-
- void SetThreshold(float threshold) {fThreshold = threshold;}
- void SetClusterThreshold(float clusterThreshold) {fClusterThreshold = clusterThreshold;}
- void SetHighGainFactor(float highGain) {fHighGainFactor = highGain;}
- void SetLowGainFactor(float lowGain) {fLowGainFactor = lowGain;}
- void SetArraySize(int size)
- {
- fArraySize = size;
- fMultiplicity = fArraySize * fArraySize;
- }
- float GetThreshold() {return fThreshold;}
- float GetClusterThreshold() {return fClusterThreshold;}
- float GetHighGainFactor() {return fHighGainFactor;}
- float GetLowGainFactor() {return fLowGainFactor;}
- float GetArraySize() {return fArraySize;}
- float GetMultiplicity() {return fMultiplicity;}
+ AliHLTPHOSClusterizer(const AliHLTPHOSClusterizer &);
+ AliHLTPHOSClusterizer & operator = (const AliHLTPHOSClusterizer &) {return *this;}
+
+ void SetRecPointContainer(AliHLTPHOSRecPointContainerStruct *RecPointContainerPtr)
+ { fRecPointContainerPtr = RecPointContainerPtr; }
+
+ void SetRecoParameters(AliPHOSRecoParamEmc*);
+
+ void SetEmcClusteringThreshold(Float_t threshold) { fEmcClusteringThreshold = threshold; }
+ void SetEmcMinEnergyThreshold(Float_t threshold) { fEmcMinEnergyThreshold = threshold; }
+ void SetEmcTimeGate(Float_t gate) { fEmcTimeGate = gate; }
+ void SetLogWeight(Float_t weight) { fLogWeight = weight; }
+
+ void SetOfflineMode(AliPHOSGetter*);
+
+ virtual Int_t ClusterizeEvent();
+ virtual Int_t GetEvent(Int_t);
- int BuildCellEnergyArray(AliHLTPHOSRcuCellEnergyDataStruct *structPtr, AliHLTPHOSRecPointListDataStruct* recPointList);
- int CreateRecPointStructArray(AliHLTPHOSRecPointDataStruct* rectStructsPtr, AliHLTPHOSRecPointListDataStruct* list, int nPoints);
- int CalculateCenterOfGravity(AliHLTPHOSRecPointDataStruct* recPointPtr);
- int CalculateMoments(AliHLTPHOSRecPointDataStruct* recPointPtr, Bool_t axisOnly);
- int ClusterizeStruct(AliHLTPHOSRecPointDataStruct* recArrayPtr, AliHLTPHOSClusterDataStruct* clusterArrayPtr);
- int ResetCellEnergyArray();
+ Int_t GetNEvents();
+ virtual void ScanForNeighbourDigits(Int_t, AliHLTPHOSRecPointDataStruct*);
+ virtual Int_t AreNeighbours(AliHLTPHOSDigitDataStruct*, AliHLTPHOSDigitDataStruct*);
+ virtual void CalculateCenterOfGravity();
+
+private:
- private:
+ Float_t fEmcClusteringThreshold;
+ Float_t fEmcMinEnergyThreshold;
+ Float_t fEmcTimeGate;
+ Float_t fLogWeight;
+ Int_t fDigitsInCluster;
- AliHLTUInt8_t fPHOSModule; /**<Number of the PHOSModule*/
- Float_t fEnergyArray[N_XCOLUMNS_MOD][N_ZROWS_MOD]; /**<2D array of cell energies*/
- Float_t fThreshold; /**<Energy threshold*/
- Float_t fClusterThreshold; /**<Cluster threshold*/
- Float_t fHighGainFactor; /**<High gain factor*/
- Float_t fLowGainFactor; /**<Low gain factor*/
- Int_t fArraySize; /**<Size of the array which the energies are summed*/
- Int_t fMultiplicity; /**<Number of crystals the energies are summed for*/
+ Bool_t fOnlineMode;
+
+ TClonesArray *fDigitArrayPtr;
+ TObjArray *fEmcRecPointsPtr;
+ AliPHOSDigit *fDigitPtr;
+ AliHLTPHOSDigitContainerDataStruct *fDigitContainerPtr;
+ AliHLTPHOSRecPointContainerStruct *fRecPointContainerPtr;
+ AliPHOSGeometry *fPHOSGeometry;
+
+ AliPHOSGetter *fGetterPtr;
+
ClassDef(AliHLTPHOSClusterizer, 1);
};
#include "AliHLTPHOSRecPointDataStruct.h"
#include "AliHLTPHOSClusterDataStruct.h"
#include "AliHLTPHOSRecPointListDataStruct.h"
+#include "AliHLTPHOSDigitContainerDataStruct.h"
using namespace std;
-const AliHLTComponentDataType AliHLTPHOSClusterizerComponent::fgkInputDataTypes[]={kAliHLTVoidDataType,{0,"",""}};
+const AliHLTComponentDataType AliHLTPHOSClusterizerComponent::fgkInputDataTypes[]=
+ {
+ kAliHLTVoidDataType,{0,"",""}
+ };
AliHLTPHOSClusterizerComponent gAliHLTPHOSClusterizerComponent;
-//AliHLTPHOSClusterizerComponent::AliHLTPHOSClusterizerComponent(): AliHLTPHOSBase(), AliHLTProcessor(), fClusterizerPtr(0), fOutPtr(0),
+//AliHLTPHOSClusterizerComponent::AliHLTPHOSClusterizerComponent(): AliHLTPHOSBase(), AliHLTProcessor(), fClusterizerPtr(0), fOutPtr(0),
// fRecPointStructArrayPtr(0), fRecPointListPtr(0)
-AliHLTPHOSClusterizerComponent::AliHLTPHOSClusterizerComponent(): AliHLTPHOSProcessor(), fClusterizerPtr(0), fOutPtr(0),
- fRecPointStructArrayPtr(0), fRecPointListPtr(0)
+AliHLTPHOSClusterizerComponent::AliHLTPHOSClusterizerComponent(): AliHLTPHOSProcessor(), fClusterizerPtr(0), fOutPtr(0),
+ fRecPointStructArrayPtr(0), fRecPointListPtr(0)
{
//Constructor
}
{
//Destructor
- if(fClusterizerPtr)
+ if (fClusterizerPtr)
{
delete fClusterizerPtr;
fClusterizerPtr = 0;
}
-
- if(fRecPointListPtr)
+
+ if (fRecPointListPtr)
{
delete fRecPointListPtr;
fRecPointListPtr = 0;
}
- if(fRecPointStructArrayPtr)
+ if (fRecPointStructArrayPtr)
{
- for(int i = 0; i < 1000; i++)
- {
- fRecPointStructArrayPtr[i].Del();
- }
+ for (int i = 0; i < 1000; i++)
+ {
+ // fRecPointStructArrayPtr[i].Del();
+ }
delete fRecPointStructArrayPtr;
fRecPointStructArrayPtr = 0;
}
-
+
}
/*
-int
+int
AliHLTPHOSClusterizerComponent::AliHLTPHOSClusterizerComponent::Deinit()
{
////////// PTH WARNING you should Define a class AliHLTPHOSModuleProcessor
}
*/
-// PTH AliHLTPHOSClusterizerComponent::AliHLTPHOSClusterizerComponent(const AliHLTPHOSClusterizerComponent &):AliHLTProcessor(),
-// fClusterizerPtr(0),
+// PTH AliHLTPHOSClusterizerComponent::AliHLTPHOSClusterizerComponent(const AliHLTPHOSClusterizerComponent &):AliHLTProcessor(),
+// fClusterizerPtr(0),
// fOutPtr(0),
// fRecPointStructArrayPtr(0),
// fRecPointListPtr(0)
//{
- //Copy constructor, not implemented
+//Copy constructor, not implemented
//}
int
{
//Deinitialization
- if(fClusterizerPtr)
+ if (fClusterizerPtr)
{
delete fClusterizerPtr;
fClusterizerPtr = 0;
}
-
- if(fRecPointListPtr)
+
+ if (fRecPointListPtr)
{
delete fRecPointListPtr;
fRecPointListPtr = 0;
}
-
- for(int i = 0; i < 1000; i++)
+
+ for (int i = 0; i < 1000; i++)
{
- fRecPointStructArrayPtr[i].Del();
+ // fRecPointStructArrayPtr[i].Del();
}
-
- if(fRecPointStructArrayPtr)
+
+ if (fRecPointStructArrayPtr)
{
- for(int i = 0; i < 1000; i++)
- {
- fRecPointStructArrayPtr[i].Del();
- }
+ for (int i = 0; i < 1000; i++)
+ {
+ // fRecPointStructArrayPtr[i].Del();
+ }
delete fRecPointStructArrayPtr;
fRecPointStructArrayPtr = 0;
}
-
-
-
-
return 0;
}
}
-const Char_t*
+const Char_t*
AliHLTPHOSClusterizerComponent::GetComponentID()
{
return "AliHltPhosClusterizer";
{
//Get datatypes for input
const AliHLTComponentDataType* pType=fgkInputDataTypes;
- while (pType->fID!=0) {
- list.push_back(*pType);
- pType++;
- }
+ while (pType->fID!=0)
+ {
+ list.push_back(*pType);
+ pType++;
+ }
}
-AliHLTComponentDataType
+AliHLTComponentDataType
AliHLTPHOSClusterizerComponent::GetOutputDataType()
{
// return AliHLTPHOSPhysicsDefinitions::fgkAliHLTClusterDataType;
inputMultiplier = 0.2;
}
-int
+int
AliHLTPHOSClusterizerComponent::DoEvent(const AliHLTComponentEventData& evtData, const AliHLTComponentBlockData* blocks,
- AliHLTComponentTriggerData& /*trigData*/, AliHLTUInt8_t* outputPtr, AliHLTUInt32_t& size,
- std::vector<AliHLTComponentBlockData>& outputBlocks)
+ AliHLTComponentTriggerData& /*trigData*/, AliHLTUInt8_t* outputPtr, AliHLTUInt32_t& size,
+
std::vector<AliHLTComponentBlockData>& outputBlocks)
{
//Do event
-
+
UInt_t tSize = 0;
- UInt_t offset = 0;
+ UInt_t offset = 0;
UInt_t mysize = 0;
Int_t nRecPoints = 0;
+ Int_t nDigits = 0;
Int_t index = 0;
+ Int_t j =0;
AliHLTUInt8_t* outBPtr;
outBPtr = outputPtr;
- const AliHLTComponentBlockData* iter = 0;
- unsigned long ndx;
+ const AliHLTComponentBlockData* iter = 0;
+ unsigned long ndx;
- for( ndx = 0; ndx < evtData.fBlockCnt; ndx++ )
+ AliHLTPHOSDigitContainerDataStruct *digitContainerPtr = 0;
+
+ //AliHLTPHOSRecPointContainerStruct *recPointContainerPtr = (AliHLTPHOSRecPointContainerStruct*)outBPtr;
+ fClusterizerPtr->SetRecPointContainer((AliHLTPHOSRecPointContainerStruct*)outBPtr);
+ for ( ndx = 0; ndx < evtData.fBlockCnt; ndx++ )
{
iter = blocks+ndx;
-
- if(iter->fDataType != AliHLTPHOSDefinitions::fgkCellEnergyDataType)
- {
- cout << "Warning: data type is not fgkCellEnergyDataType " << endl;
- continue;
- }
- index = fClusterizerPtr->BuildCellEnergyArray( reinterpret_cast<AliHLTPHOSRcuCellEnergyDataStruct*>(iter->fPtr),
- fRecPointListPtr);
-
+ digitContainerPtr = reinterpret_cast<AliHLTPHOSDigitContainerDataStruct*>(iter->fPtr);
+ if (iter->fDataType != AliHLTPHOSDefinitions::fgkAliHLTDigitDataType)
+ {
+ // cout << "Warning: data type is not fgkAliHLTDigitDataType " << endl;
+ continue;
+ }
+ for (Int_t i = 0; i < digitContainerPtr->fNDigits; i++)
+ {
+ if(fNoCrazyness && digitContainerPtr->fDigitDataStruct[i].fCrazyness)
+ continue;
+
+ fAllDigitsPtr->fDigitDataStruct[j+nDigits].fX = digitContainerPtr->fDigitDataStruct[i].fX;
+ fAllDigitsPtr->fDigitDataStruct[j+nDigits].fZ = digitContainerPtr->fDigitDataStruct[i].fZ;
+ fAllDigitsPtr->fDigitDataStruct[j+nDigits].fAmplitude = digitContainerPtr->fDigitDataStruct[i].fAmplitude;
+ fAllDigitsPtr->fDigitDataStruct[j+nDigits].fTime = digitContainerPtr->fDigitDataStruct[i].fTime;
+ // fAllDigitsPtr->fDigitDataStruct[i+nDigits].fCrazyness = digitContainerPtr->fDigitDataStruct[i].fCrazyness;
+ j++;
+ }
+ nDigits++;
}
-
- nRecPoints = fClusterizerPtr->CreateRecPointStructArray(fRecPointStructArrayPtr, fRecPointListPtr, index);
-
- cout << "Number of clusters found: " << nRecPoints << endl;
-
- for(Int_t i = 0; i < nRecPoints; i++)
+
+// nRecPoints = fClusterizerPtr->CreateRecPointStructArray(fRecPointStructArrayPtr, fRecPointListPtr, index);
+
+ fOutPtr = (AliHLTPHOSRecPointContainerStruct*)outBPtr;
+ nRecPoints = fClusterizerPtr->ClusterizeEvent();
+ //nRecPoints = fClusterizerPtr->CalculateCenterOfGravity(&fRecPointStructArrayPtr[i]);
+ cout << "Number of clusters found: " << nRecPoints << " extracted from " << nDigits << " digits" << endl;
+
+ mysize = 0;
+ offset = tSize;
+
+ // fClusterizerPtr->CalculateMoments(&fRecPointStructArrayPtr[i], 0);
+ // fClusterizerPtr->ClusterizeStruct(&fRecPointStructArrayPtr[i], fOutPtr);
+
+ mysize += sizeof(AliHLTPHOSClusterDataStruct);
+
+ AliHLTComponentBlockData bd;
+ FillBlockData( bd );
+ bd.fOffset = offset;
+ bd.fSize = mysize;
+ // PTH bd.fDataType = AliHLTPHOSPhysicsDefinitions::fgkAliHLTClusterDataType;
+ bd.fDataType = AliHLTPHOSDefinitions::fgkAliHLTClusterDataType;
+ bd.fSpecification = 0xFFFFFFFF;
+ outputBlocks.push_back( bd );
+
+ tSize += mysize;
+ outBPtr += mysize;
+
+ if ( tSize > size )
{
- mysize = 0;
- offset = tSize;
-
- fOutPtr = (AliHLTPHOSClusterDataStruct*)outBPtr;
- fClusterizerPtr->CalculateCenterOfGravity(&fRecPointStructArrayPtr[i]);
- // fClusterizerPtr->CalculateMoments(&fRecPointStructArrayPtr[i], 0);
- // fClusterizerPtr->ClusterizeStruct(&fRecPointStructArrayPtr[i], fOutPtr);
-
- mysize += sizeof(AliHLTPHOSClusterDataStruct);
-
- AliHLTComponentBlockData bd;
- FillBlockData( bd );
- bd.fOffset = offset;
- bd.fSize = mysize;
- // PTH bd.fDataType = AliHLTPHOSPhysicsDefinitions::fgkAliHLTClusterDataType;
- bd.fDataType = AliHLTPHOSDefinitions::fgkAliHLTClusterDataType;
- bd.fSpecification = 0xFFFFFFFF;
- outputBlocks.push_back( bd );
-
- tSize += mysize;
- outBPtr += mysize;
-
- if( tSize > size )
- {
- Logging( kHLTLogFatal, "HLT::AliHLTPHOSClusterizerComponent::DoEvent", "Too much data",
- "Data written over allowed buffer. Amount written: %lu, allowed amount: %lu."
- , tSize, size );
- return EMSGSIZE;
- }
+ Logging( kHLTLogFatal, "HLT::AliHLTPHOSClusterizerComponent::DoEvent", "Too much data",
+ "Data written over allowed buffer. Amount written: %lu, allowed amount: %lu."
+ , tSize, size );
+ return EMSGSIZE;
}
+
size = tSize;
- fClusterizerPtr->ResetCellEnergyArray();
+// fClusterizerPtr->ResetCellEnergyArray();
return 0;
AliHLTPHOSClusterizerComponent::DoInit(int argc, const char** argv )
{
//Do initialization
+ fAllDigitsPtr = new AliHLTPHOSDigitContainerDataStruct();
fClusterizerPtr = new AliHLTPHOSClusterizer();
- for(int i = 0; i < argc; i++)
+ //fClusterizerPtr->SetNoCrazyness(true);
+ //
+ for (int i = 0; i < argc; i++)
{
+ /*
if(!strcmp("-threshold", argv[i]))
- fClusterizerPtr->SetThreshold(atof(argv[i+1]));
+ fClusterizerPtr->SetThreshold(atof(argv[i+1]));
if(!strcmp("-clusterthreshold", argv[i]))
- fClusterizerPtr->SetClusterThreshold(atof(argv[i+1]));
+ fClusterizerPtr->SetClusterThreshold(atof(argv[i+1]));
if(!strcmp("-highgain", argv[i]))
- fClusterizerPtr->SetHighGainFactor(atof(argv[i+1]));
+ fClusterizerPtr->SetHighGainFactor(atof(argv[i+1]));
if(!strcmp("-lowgain", argv[i]))
- fClusterizerPtr->SetLowGainFactor(atof(argv[i+1]));
+ fClusterizerPtr->SetLowGainFactor(atof(argv[i+1]));
if(!strcmp("-arraysize", argv[i]))
- fClusterizerPtr->SetArraySize(atoi(argv[i+1]));
+ fClusterizerPtr->SetArraySize(atoi(argv[i+1]));*/
}
- fClusterizerPtr->ResetCellEnergyArray();
+ // fClusterizerPtr->ResetCellEnergyArray();
fRecPointListPtr = new AliHLTPHOSRecPointListDataStruct[N_ZROWS_MOD*N_XCOLUMNS_MOD];
fRecPointStructArrayPtr = new AliHLTPHOSRecPointDataStruct[1000];
- for(int i = 0; i < 1000; i++)
+ for (int i = 0; i < 1000; i++)
{
fRecPointStructArrayPtr[i].fMultiplicity = atoi(argv[4])* atoi(argv[4]);
- fRecPointStructArrayPtr[i].New();
+ // fRecPointStructArrayPtr[i].New();
}
+ /*
printf("Clusterizer component started with:\n");
printf(" Cell threshold: %f\n", fClusterizerPtr->GetThreshold());
printf(" Cluster threshold: %f\n", fClusterizerPtr->GetClusterThreshold());
printf(" High gain factor: %f\n", fClusterizerPtr->GetHighGainFactor());
printf(" Low gain factor: %f\n", fClusterizerPtr->GetLowGainFactor());
printf(" Cluster array size: %d\n\n", fClusterizerPtr->GetArraySize());
-
+ */
return 0;
}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff