]> git.uio.no Git - u/mrichter/AliRoot.git/blame - HLT/CALO/AliHLTCaloClusterAnalyser.cxx
Refactoring of handling of constants.
[u/mrichter/AliRoot.git] / HLT / CALO / AliHLTCaloClusterAnalyser.cxx
CommitLineData
ea367f6a 1// $Id: AliHLTCaloClusterAnalyser.cxx 35107 2009-09-30 01:45:06Z phille $
2
3/**************************************************************************
4 * This file is property of and copyright by the ALICE HLT Project *
5 * All rights reserved. *
6 * *
7 * Primary Authors: Oystein Djuvsland *
8 * *
9 * Permission to use, copy, modify and distribute this software and its *
10 * documentation strictly for non-commercial purposes is hereby granted *
11 * without fee, provided that the above copyright notice appears in all *
12 * copies and that both the copyright notice and this permission notice *
13 * appear in the supporting documentation. The authors make no claims *
14 * about the suitability of this software for any purpose. It is *
15 * provided "as is" without express or implied warranty. *
16 **************************************************************************/
17
18/**
19 * @file AliHLTCaloClusterAnalyser.cxx
20 * @author Oystein Djuvsland
21 * @date
22 * @brief Cluster analyser for Calo HLT
23 */
24
25// see header file for class documentation
26// or
27// refer to README to build package
28// or
29// visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
30
31#include "AliHLTCaloClusterAnalyser.h"
32#include "AliHLTCaloRecPointHeaderStruct.h"
33#include "AliHLTCaloRecPointDataStruct.h"
34#include "AliHLTCaloClusterDataStruct.h"
35//#include "AliHLTCaloPhysicsAnalyzer.h"
7c80a370 36#include "AliHLTCaloGeometry.h"
ea367f6a 37#include "AliESDCaloCluster.h"
38#include "TMath.h"
39#include "TVector3.h"
73d6f579 40#include "TH1F.h"
41#include "TFile.h"
c1e4a18c 42#include "AliHLTCaloRecoParamHandler.h"
ea367f6a 43
44ClassImp(AliHLTCaloClusterAnalyser);
45
7c80a370 46AliHLTCaloClusterAnalyser::AliHLTCaloClusterAnalyser() :
ea367f6a 47 // AliHLTCaloBase(),
7c80a370 48 fLogWeight(4.5),
49 fRecPointArray(0),
50 fDigitDataArray(0),
ea367f6a 51 fNRecPoints(0),
52 fCaloClusterDataPtr(0),
53 fCaloClusterHeaderPtr(0),
ea367f6a 54 //fAnalyzerPtr(0),
55 fDoClusterFit(false),
56 fHaveCPVInfo(false),
57 fDoPID(false),
7c80a370 58 fHaveDistanceToBadChannel(false),
73d6f579 59 fGeometry(0),
c1e4a18c 60 fClusterType(AliESDCaloCluster::kPHOSCluster),
38a3b2ad 61 fRecoParamsPtr(0),
62 fCutOnSingleCellClusters(false),
63 fSingleCellEnergyCut(0.5)
ea367f6a 64{
65 //See header file for documentation
ea367f6a 66}
67
68AliHLTCaloClusterAnalyser::~AliHLTCaloClusterAnalyser()
69{
98baf84d 70 // See header file for class documentation
ea367f6a 71}
72
73void
7c80a370 74AliHLTCaloClusterAnalyser::SetCaloClusterData(AliHLTCaloClusterDataStruct *caloClusterDataPtr)
ea367f6a 75{
76 //see header file for documentation
77 fCaloClusterDataPtr = caloClusterDataPtr;
78}
7c80a370 79
ea367f6a 80void
7c80a370 81AliHLTCaloClusterAnalyser::SetRecPointArray(AliHLTCaloRecPointDataStruct **recPointDataPtr, Int_t nRecPoints)
ea367f6a 82{
7c80a370 83 fRecPointArray = recPointDataPtr;
84 fNRecPoints = nRecPoints;
85}
86
87void
88AliHLTCaloClusterAnalyser::SetDigitDataArray(AliHLTCaloDigitDataStruct *digits)
89{
90// AliHLTCaloClusterizer cl("PHOS");
91 // cl.CheckDigits(fRecPointArray, digits, fNRecPoints);
92 fDigitDataArray = digits;
93 //cl.CheckDigits(fRecPointArray, fDigitDataArray, fNRecPoints);
ea367f6a 94}
95
96Int_t
97AliHLTCaloClusterAnalyser::CalculateCenterOfGravity()
98{
99 //see header file for documentation
100 Float_t wtot = 0.;
101 Float_t x = 0.;
102 Float_t z = 0.;
103 Float_t xi = 0.;
104 Float_t zi = 0.;
105
106 AliHLTCaloDigitDataStruct *digit = 0;
ea367f6a 107
108 UInt_t iDigit = 0;
109
110 for(Int_t iRecPoint=0; iRecPoint < fNRecPoints; iRecPoint++)
111 {
7c80a370 112 AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[iRecPoint];
ad44d760 113 // digit = &(recPoint->fDigits);
98baf84d 114
7c80a370 115 Int_t *digitIndexPtr = &(recPoint->fDigits);
256f1bf0 116 wtot = 0;
117 x = 0;
118 z = 0;
119
c1e4a18c 120/* Float_t maxAmp = 0;
256f1bf0 121 Int_t maxX = 0;
c1e4a18c 122 Int_t maxZ = 0;*/
123 for(iDigit = 0; iDigit < recPoint->fMultiplicity; iDigit++)
ea367f6a 124 {
7c80a370 125
126 digit = &(fDigitDataArray[*digitIndexPtr]);
98baf84d 127
a6db7bcd 128 xi = digit->fX+0.5;
129 zi = digit->fZ+0.5;
98baf84d 130
ea367f6a 131 if (recPoint->fAmp > 0 && digit->fEnergy > 0)
132 {
133 Float_t w = TMath::Max( 0., fLogWeight + TMath::Log( digit->fEnergy / recPoint->fAmp ) ) ;
134 x += xi * w ;
135 z += zi * w ;
136 wtot += w ;
256f1bf0 137/* if(digit->fEnergy > maxAmp)
138 {
139 maxAmp = digit->fEnergy;
c1e4a18c 140 maxX = digit->fX;// + 0.5;
141 maxZ = digit->fZ;// + 0.5;
256f1bf0 142 }*/
ea367f6a 143 }
7c80a370 144 digitIndexPtr++;
ea367f6a 145 }
98baf84d 146
ea367f6a 147 if (wtot>0)
148 {
149 recPoint->fX = x/wtot ;
150 recPoint->fZ = z/wtot ;
151 }
152 else
153 {
154 recPoint->fAmp = 0;
155 }
256f1bf0 156// printf("Max digit: E = %f, x = %d, z= %d, cluster: E = %f, x = %f, z = %f\n" , maxAmp, maxX, maxZ, recPoint->fAmp, recPoint->fX, recPoint->fZ);
ea367f6a 157 }
158 return 0;
159}
160
161
162Int_t
163AliHLTCaloClusterAnalyser::CalculateRecPointMoments()
cd9e2797 164{
ea367f6a 165 //See header file for documentation
166 return 0;
167}
168
169Int_t
170AliHLTCaloClusterAnalyser::CalculateClusterMoments(AliHLTCaloRecPointDataStruct */*recPointPtr*/, AliHLTCaloClusterDataStruct* /*clusterPtr*/)
171{
172 //See header file for documentation
173 return 0;
174}
175
176
177Int_t
178AliHLTCaloClusterAnalyser::DeconvoluteClusters()
179{
180 //See header file for documentation
181 return 0;
182}
183
184Int_t
7c80a370 185AliHLTCaloClusterAnalyser::CreateClusters(Int_t nRecPoints, UInt_t availableSize, UInt_t& totSize)
ea367f6a 186{
187 //See header file for documentation
188
c1e4a18c 189 if(fRecoParamsPtr)
190 {
191 fLogWeight = fRecoParamsPtr->GetLogWeight();
192 }
73d6f579 193
7c80a370 194 fNRecPoints = nRecPoints;
195
196 if(fGeometry == 0)
197 {
198 HLTError("No geometry object is initialised, creation of clusters stopped");
73cdb797 199 return -1;
7c80a370 200 }
201
202 CalculateCenterOfGravity();
203
ad44d760 204 // AliHLTCaloDigitDataStruct* digitPtr = &(recPointPtr->fDigits);
205 AliHLTCaloDigitDataStruct* digitPtr = 0;
206
73d6f579 207 AliHLTCaloClusterDataStruct* caloClusterPtr = 0;
7fdb2519 208
209 // Int_t id = -1;
ea367f6a 210 TVector3 globalPos;
211
212 for(Int_t i = 0; i < fNRecPoints; i++) //TODO needs fix when we start unfolding (number of clusters not necessarily same as number of recpoints gotten from the clusterizer
213 {
73d6f579 214 if((availableSize - totSize) < sizeof(AliHLTCaloClusterDataStruct))
215 {
9a0cbaba 216 HLTError("Out of buffer: available size is: %d, total size used: %d", availableSize, totSize);
73d6f579 217 return -ENOBUFS;
218 }
38a3b2ad 219
220 AliHLTCaloRecPointDataStruct *recPointPtr = fRecPointArray[i];
221
222 if(fCutOnSingleCellClusters && recPointPtr->fAmp > fSingleCellEnergyCut && recPointPtr->fMultiplicity == 1) continue;
223
7fdb2519 224 totSize += sizeof(AliHLTCaloClusterDataStruct);
73d6f579 225
226 caloClusterPtr = fCaloClusterDataPtr;
34eb3256 227 caloClusterPtr->fChi2 = 0;
228 caloClusterPtr->fClusterType = kUndef;
229 caloClusterPtr->fDispersion = 0;
230 caloClusterPtr->fDistanceToBadChannel = 0;
231 caloClusterPtr->fDistToBadChannel = 0;
232 caloClusterPtr->fEmcCpvDistance = 0;
233 caloClusterPtr->fEnergy = 0;
234 caloClusterPtr->fFitQuality = 0;
235 caloClusterPtr->fID = 0;
236 caloClusterPtr->fM02 = 0;
237 caloClusterPtr->fM20 = 0;
238 caloClusterPtr->fNCells = 0;
239 caloClusterPtr->fNExMax = 0;
240 caloClusterPtr->fTOF = 0;
241 caloClusterPtr->fTrackDx = 0;
242 caloClusterPtr->fTrackDz = 0;
7c80a370 243
244 AliHLTCaloGlobalCoordinate globalCoord;
245 fGeometry->GetGlobalCoordinates(*recPointPtr, globalCoord);
ea367f6a 246
353c33d3 247 caloClusterPtr->fModule = recPointPtr->fModule;
a6db7bcd 248 caloClusterPtr->fGlobalPos[0] = globalCoord.fX;
cd9e2797 249 caloClusterPtr->fGlobalPos[1] = globalCoord.fY;
a6db7bcd 250 caloClusterPtr->fGlobalPos[2] = globalCoord.fZ;
ea367f6a 251
ee72b4e5 252 HLTDebug("Cluster local position: x = %f, z = %f, module = %d", recPointPtr->fX, recPointPtr->fZ, recPointPtr->fModule);
253 HLTDebug("Cluster global position: x = %f, y = %f, z = %f", globalCoord.fX, globalCoord.fY, globalCoord.fZ);
254
7fdb2519 255 //caloClusterPtr->fNCells = 0;//recPointPtr->fMultiplicity;
256 caloClusterPtr->fNCells = recPointPtr->fMultiplicity;
cd9e2797 257
a6db7bcd 258 caloClusterPtr->fClusterType = fClusterType;
7fdb2519 259// Int_t tmpSize = 0;//totSize + (caloClusterPtr->fNCells-1)*(sizeof(Short_t) + sizeof(Float_t));
a6db7bcd 260
2a24cbbe 261 //TODO remove hardcoded 10;
262 memset(caloClusterPtr->fTracksMatched, 0xff, sizeof(Int_t)*10);
263
9a0cbaba 264 //Int_t tmpSize = totSize + (caloClusterPtr->fNCells-1)*(sizeof(Short_t) + sizeof(Float_t));
cefad3b2 265 UInt_t tmpSize = (caloClusterPtr->fNCells-1)*sizeof(AliHLTCaloCellDataStruct);
cd9e2797 266
73d6f579 267 if((availableSize - totSize) < tmpSize)
268 {
9a0cbaba 269 HLTError("Out of buffer, available size is: %d, total size used: %d, extra size needed: %d", availableSize, totSize, tmpSize);
73d6f579 270 return -ENOBUFS;
271 }
c1e4a18c 272
cd9e2797 273 Int_t *digitIndexPtr = &(recPointPtr->fDigits);
274 Int_t id = 0;
275
7fdb2519 276 AliHLTCaloCellDataStruct *cellPtr = &(caloClusterPtr->fCaloCells);
277
ea367f6a 278 for(UInt_t j = 0; j < caloClusterPtr->fNCells; j++)
279 {
cd9e2797 280 digitPtr = &(fDigitDataArray[*digitIndexPtr]);
936dbe15 281 fGeometry->GetCellAbsId(recPointPtr->fModule, digitPtr->fX, digitPtr->fZ, id);
7fdb2519 282
283 cellPtr->fCellsAbsId= id;
284 cellPtr->fCellsAmpFraction = digitPtr->fEnergy/recPointPtr->fAmp;
cd9e2797 285 //printf("Cell ID pointer: %x\n", cellIDPtr);
286 //printf("Cell Amp Pointer: %x\n", cellAmpFracPtr);
9a0cbaba 287 //printf("Cell pos: x = %d, z = %d\n", digitPtr->fX, digitPtr->fZ);
34eb3256 288// printf("Cell ID: %d\n", cellPtr->fCellsAbsId);
9a0cbaba 289 //printf("Cell Amp: %f, pointer: %x\n", *cellAmpFracPtr, cellAmpFracPtr);
7fdb2519 290 cellPtr++;
cd9e2797 291 digitIndexPtr++;
292
ea367f6a 293 }
294
73d6f579 295 totSize += tmpSize;
296
c1e4a18c 297 if(fRecoParamsPtr)
298 {
299 caloClusterPtr->fEnergy = fRecoParamsPtr->GetCorrectedEnergy(recPointPtr->fAmp);
300 }
301 else
302 {
303 caloClusterPtr->fEnergy = recPointPtr->fAmp;
304 }
34eb3256 305
7fdb2519 306 HLTDebug("Cluster global position: x = %f, y = %f, z = %f, energy: %f, number of cells: %d, cluster pointer: %x", globalCoord.fX, globalCoord.fY, globalCoord.fZ, caloClusterPtr->fEnergy, caloClusterPtr->fNCells, caloClusterPtr);
307
ea367f6a 308 if(fDoClusterFit)
309 {
310 FitCluster(recPointPtr);
311 }
312 else
313 {
314 caloClusterPtr->fDispersion = 0;
315 caloClusterPtr->fFitQuality = 0;
316 caloClusterPtr->fM20 = 0;
317 caloClusterPtr->fM02 = 0;
318
319 }
320 if(fHaveCPVInfo)
321 {
322 caloClusterPtr->fEmcCpvDistance = GetCPVDistance(recPointPtr);
323 }
324 else
325 {
326 caloClusterPtr->fEmcCpvDistance = -1;
327 }
328 if(fDoPID)
329 {
330 DoParticleIdentification(caloClusterPtr);
331 }
332 else
333 {
334 for(Int_t k = 0; k < AliPID::kSPECIESN; k++)
335 {
336 caloClusterPtr->fPID[k] = 0;
337 }
338 }
339 if(fHaveDistanceToBadChannel)
340 {
341 caloClusterPtr->fDistanceToBadChannel = GetDistanceToBadChannel(caloClusterPtr);
342 }
343 else
344 {
345 caloClusterPtr->fDistanceToBadChannel = -1;
346 }
347
73d6f579 348 caloClusterPtr->fClusterType = fClusterType;
cd9e2797 349 //totSize += sizeof(AliHLTCaloClusterDataStruct) + (caloClusterPtr->fNCells)*(sizeof(Short_t) +sizeof(Float_t)-1);
73d6f579 350 //totSize += sizeof(AliHLTCaloClusterDataStruct) + (caloClusterPtr->fNCells-1)*(sizeof(Short_t) + sizeof(Float_t));
cd9e2797 351 //printf("CaloClusterPtr: %x, energy: %f\n", caloClusterPtr, caloClusterPtr->fEnergy);
352
ea367f6a 353 // caloClusterPtr = reinterpret_cast<AliHLTCaloClusterDataStruct*>(cellAmpFracPtr);
cd9e2797 354 //caloClusterPtr = reinterpret_cast<AliHLTCaloClusterDataStruct*>(cellIDPtr);
7fdb2519 355 fCaloClusterDataPtr = reinterpret_cast<AliHLTCaloClusterDataStruct*>(cellPtr);
ea367f6a 356 recPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(digitPtr);
73d6f579 357 //digitPtr = &(recPointPtr->fDigits);
ea367f6a 358 }
98baf84d 359
360return fNRecPoints;
ea367f6a 361
362}
363