#include "AliHLTPHOSClusterizer.h"
#include "AliHLTPHOSBase.h"
+#include "AliHLTLogging.h"
#include "TMath.h"
#include "AliHLTPHOSRecPointContainerStruct.h"
#include "AliHLTPHOSRecPointDataStruct.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 <iostream>
+
+using namespace std;
ClassImp(AliHLTPHOSClusterizer);
fDigitContainerPtr(0),
fRecPointContainerPtr(0),
fPHOSGeometry(0),
- fGetterPtr(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
//See header file for documentation
}
+void
+AliHLTPHOSClusterizer::SetRecPointContainer(AliHLTPHOSRecPointContainerStruct* recPointContainerPtr)
+ {
+ fRecPointContainerPtr = recPointContainerPtr;
+ fRecPointContainerPtr->fNRecPoints = 0;
+ }
void
-AliHLTPHOSClusterizer::SetRecoParameters(AliPHOSRecoParamEmc* params)
+AliHLTPHOSClusterizer::SetRecoParameters(AliPHOSRecoParam* params)
{
//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();
-}
-
-void
-AliHLTPHOSClusterizer::SetOfflineMode(AliPHOSGetter* getter)
-{
- //see header file for documentation
- fRecPointContainerPtr = new AliHLTPHOSRecPointContainerStruct();
- fDigitContainerPtr = new AliHLTPHOSDigitContainerDataStruct();
- fGetterPtr = getter;
- fDigitArrayPtr = fGetterPtr->Digits();
- fEmcRecPointsPtr = fGetterPtr->EmcRecPoints();
- fOnlineMode = false;
-}
-
-Int_t
-AliHLTPHOSClusterizer::GetEvent(Int_t i)
-{
- //see header file for documentation
- Int_t coord[4];
-
- fGetterPtr->Event(i, "D");
- for(Int_t j = 0; j < fDigitArrayPtr->GetEntries(); j++)
- {
- 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()
-{
- //see header file for documentation
- if(fOnlineMode)
- {
- printf("Number of events not available in online mode!\n");
- return -1;
- }
- return fGetterPtr->MaxEvent();
-}
-
+#endif
+}
Int_t
AliHLTPHOSClusterizer::ClusterizeEvent()
//see header file for documentation
Int_t nRecPoints = 0;
UInt_t i = 0;
-
- //printf("Starting clusterisation of event.\n");
-
AliHLTPHOSRecPointDataStruct *recPoint = 0;
-
- //printf("Number of digits in event: %d\n", fDigitContainerPtr->fNDigits);
-
//Clusterization starts
for(i = 0; i < fDigitContainerPtr->fNDigits; i++)
{
-
fDigitsInCluster = 0;
- //printf("Digit energy: %f\n", fDigitContainerPtr->fDigitDataStruct[i].fEnergy);
+
if(fDigitContainerPtr->fDigitDataStruct[i].fEnergy < fEmcClusteringThreshold)
{
continue;
}
- //printf("Got rec point above clustering threshold!\n");
recPoint = &(fRecPointContainerPtr->fRecPointArray[nRecPoints]);
recPoint->fAmp = 0;
- //recPoint->
+ 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)
+ {
+ // HLTWarning("Too many rec points in event. Aborting clusterisation");
+ break;
+ }
ScanForNeighbourDigits(i, recPoint);
-
recPoint->fMultiplicity = fDigitsInCluster;
-
+
}//end of clusterization
fRecPointContainerPtr->fNRecPoints = nRecPoints;
-
return nRecPoints;
}
{
//see header file for documentation
-
- for(UInt_t j = 0; j < fDigitContainerPtr->fNDigits; j++)
+ 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++)
{
if(fDigitContainerPtr->fDigitDataStruct[j].fEnergy > fEmcMinEnergyThreshold)
{
- switch(AreNeighbours(&(fDigitContainerPtr->fDigitDataStruct[index]),
+ if(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;
+ fDigitContainerPtr->fDigitDataStruct[j].fEnergy = 0;
fDigitsInCluster++;
ScanForNeighbourDigits(j, recPoint);
- break;
- case 2:
- break;
}
}
}
return;
}
-
Int_t
AliHLTPHOSClusterizer::AreNeighbours(AliHLTPHOSDigitDataStruct* digit1,
AliHLTPHOSDigitDataStruct* digit2)
{
//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 );
-
- //cout << "x1: " << digit1->fX << " x2: " << digit2->fX << " z1: " << digit1->fZ << " z2: " << digit2->fZ <<
- //" t1: " << digit1->fTime << " t2: " << digit1->fTime << " --> ";
-
- if (( coldiff <= 1 ) && ( rowdiff <= 1 ))
+ if (( coldiff <= 1 && rowdiff < 1 ) || ( coldiff < 1 && rowdiff <= 1 ))
{
if(TMath::Abs(digit1->fTime - digit2->fTime ) < fEmcTimeGate)
{
- //cout << "OK\n";
return 1;
}
- // cout << "Time difference\n";
}
- //else
- //{
- // cout << "Space difference\n";
- //}
}
return 0;
}
-
-void
-AliHLTPHOSClusterizer::CalculateCenterOfGravity()
-{
- //see header file for documentation
-
- 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;
-
- //AliPHOSGeometry * phosgeom = AliPHOSGeometry::GetInstance() ;
-
- UInt_t iDigit = 0;
- UInt_t iRecPoint = 0;
-
- for(iRecPoint=0; iRecPoint<fRecPointContainerPtr->fNRecPoints; iRecPoint++)
- {
- recPoint = &(fRecPointContainerPtr->fRecPointArray[iRecPoint]);
- for(iDigit = 0; iDigit < recPoint->fMultiplicity; iDigit++)
- {
- 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;
- }
- }
-
-}
-
-