]> git.uio.no Git - u/mrichter/AliRoot.git/blame - HLT/CALO/AliHLTCaloClusterAnalyser.cxx
added DigitPublisher class and relative modifications to the Geometry
[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),
4b759e43 60#ifndef HAVE_NOT_ALIVCLUSTER // backward compatibility for r42844
c8fe2783 61 fClusterType(AliVCluster::kPHOSNeutral),
4b759e43 62#else
63 fClusterType(AliESDCaloCluster::kPHOSCluster),
64#endif
38a3b2ad 65 fRecoParamsPtr(0),
66 fCutOnSingleCellClusters(false),
67 fSingleCellEnergyCut(0.5)
ea367f6a 68{
69 //See header file for documentation
ea367f6a 70}
71
72AliHLTCaloClusterAnalyser::~AliHLTCaloClusterAnalyser()
73{
98baf84d 74 // See header file for class documentation
ea367f6a 75}
76
77void
7c80a370 78AliHLTCaloClusterAnalyser::SetCaloClusterData(AliHLTCaloClusterDataStruct *caloClusterDataPtr)
ea367f6a 79{
80 //see header file for documentation
81 fCaloClusterDataPtr = caloClusterDataPtr;
82}
7c80a370 83
ea367f6a 84void
7c80a370 85AliHLTCaloClusterAnalyser::SetRecPointArray(AliHLTCaloRecPointDataStruct **recPointDataPtr, Int_t nRecPoints)
ea367f6a 86{
7c80a370 87 fRecPointArray = recPointDataPtr;
88 fNRecPoints = nRecPoints;
89}
90
91void
f92543f5 92AliHLTCaloClusterAnalyser::SetDigitDataArray(AliHLTCaloDigitDataStruct **digits)
7c80a370 93{
94// AliHLTCaloClusterizer cl("PHOS");
95 // cl.CheckDigits(fRecPointArray, digits, fNRecPoints);
96 fDigitDataArray = digits;
97 //cl.CheckDigits(fRecPointArray, fDigitDataArray, fNRecPoints);
ea367f6a 98}
99
100Int_t
101AliHLTCaloClusterAnalyser::CalculateCenterOfGravity()
102{
103 //see header file for documentation
104 Float_t wtot = 0.;
105 Float_t x = 0.;
106 Float_t z = 0.;
107 Float_t xi = 0.;
108 Float_t zi = 0.;
109
110 AliHLTCaloDigitDataStruct *digit = 0;
ea367f6a 111
112 UInt_t iDigit = 0;
113
114 for(Int_t iRecPoint=0; iRecPoint < fNRecPoints; iRecPoint++)
115 {
7c80a370 116 AliHLTCaloRecPointDataStruct *recPoint = fRecPointArray[iRecPoint];
ad44d760 117 // digit = &(recPoint->fDigits);
98baf84d 118
7c80a370 119 Int_t *digitIndexPtr = &(recPoint->fDigits);
256f1bf0 120 wtot = 0;
121 x = 0;
122 z = 0;
123
c1e4a18c 124/* Float_t maxAmp = 0;
256f1bf0 125 Int_t maxX = 0;
c1e4a18c 126 Int_t maxZ = 0;*/
127 for(iDigit = 0; iDigit < recPoint->fMultiplicity; iDigit++)
ea367f6a 128 {
7c80a370 129
f92543f5 130 digit = fDigitDataArray[*digitIndexPtr];
98baf84d 131
1be6beff 132 xi = digit->fX;
133 zi = digit->fZ;
134
135 //xi = digit->fX+0.5;
136 //zi = digit->fZ+0.5;
98baf84d 137
ea367f6a 138 if (recPoint->fAmp > 0 && digit->fEnergy > 0)
139 {
140 Float_t w = TMath::Max( 0., fLogWeight + TMath::Log( digit->fEnergy / recPoint->fAmp ) ) ;
141 x += xi * w ;
142 z += zi * w ;
143 wtot += w ;
256f1bf0 144/* if(digit->fEnergy > maxAmp)
145 {
146 maxAmp = digit->fEnergy;
c1e4a18c 147 maxX = digit->fX;// + 0.5;
148 maxZ = digit->fZ;// + 0.5;
256f1bf0 149 }*/
ea367f6a 150 }
7c80a370 151 digitIndexPtr++;
ea367f6a 152 }
98baf84d 153
ea367f6a 154 if (wtot>0)
155 {
156 recPoint->fX = x/wtot ;
157 recPoint->fZ = z/wtot ;
158 }
159 else
160 {
1be6beff 161 recPoint->fX = -9999;
162 recPoint->fZ =-9999;
163 // no good crashes depth with FP exception
164 //recPoint->fAmp = 0;
ea367f6a 165 }
256f1bf0 166// 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 167 }
168 return 0;
169}
170
171
172Int_t
173AliHLTCaloClusterAnalyser::CalculateRecPointMoments()
cd9e2797 174{
ea367f6a 175 //See header file for documentation
176 return 0;
177}
178
179Int_t
180AliHLTCaloClusterAnalyser::CalculateClusterMoments(AliHLTCaloRecPointDataStruct */*recPointPtr*/, AliHLTCaloClusterDataStruct* /*clusterPtr*/)
181{
182 //See header file for documentation
183 return 0;
184}
185
186
187Int_t
188AliHLTCaloClusterAnalyser::DeconvoluteClusters()
189{
190 //See header file for documentation
191 return 0;
192}
193
194Int_t
7c80a370 195AliHLTCaloClusterAnalyser::CreateClusters(Int_t nRecPoints, UInt_t availableSize, UInt_t& totSize)
ea367f6a 196{
197 //See header file for documentation
198
c1e4a18c 199 if(fRecoParamsPtr)
200 {
201 fLogWeight = fRecoParamsPtr->GetLogWeight();
202 }
73d6f579 203
7c80a370 204 fNRecPoints = nRecPoints;
205
206 if(fGeometry == 0)
207 {
208 HLTError("No geometry object is initialised, creation of clusters stopped");
73cdb797 209 return -1;
7c80a370 210 }
211
212 CalculateCenterOfGravity();
213
ad44d760 214 // AliHLTCaloDigitDataStruct* digitPtr = &(recPointPtr->fDigits);
215 AliHLTCaloDigitDataStruct* digitPtr = 0;
216
73d6f579 217 AliHLTCaloClusterDataStruct* caloClusterPtr = 0;
7fdb2519 218
219 // Int_t id = -1;
ea367f6a 220 TVector3 globalPos;
221
222 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
223 {
73d6f579 224 if((availableSize - totSize) < sizeof(AliHLTCaloClusterDataStruct))
225 {
9a0cbaba 226 HLTError("Out of buffer: available size is: %d, total size used: %d", availableSize, totSize);
73d6f579 227 return -ENOBUFS;
228 }
38a3b2ad 229
230 AliHLTCaloRecPointDataStruct *recPointPtr = fRecPointArray[i];
231
232 if(fCutOnSingleCellClusters && recPointPtr->fAmp > fSingleCellEnergyCut && recPointPtr->fMultiplicity == 1) continue;
233
7fdb2519 234 totSize += sizeof(AliHLTCaloClusterDataStruct);
73d6f579 235
236 caloClusterPtr = fCaloClusterDataPtr;
34eb3256 237 caloClusterPtr->fChi2 = 0;
238 caloClusterPtr->fClusterType = kUndef;
239 caloClusterPtr->fDispersion = 0;
240 caloClusterPtr->fDistanceToBadChannel = 0;
241 caloClusterPtr->fDistToBadChannel = 0;
242 caloClusterPtr->fEmcCpvDistance = 0;
243 caloClusterPtr->fEnergy = 0;
244 caloClusterPtr->fFitQuality = 0;
245 caloClusterPtr->fID = 0;
246 caloClusterPtr->fM02 = 0;
247 caloClusterPtr->fM20 = 0;
248 caloClusterPtr->fNCells = 0;
249 caloClusterPtr->fNExMax = 0;
250 caloClusterPtr->fTOF = 0;
fe1a0255 251 caloClusterPtr->fTrackDx = -999;
252 caloClusterPtr->fTrackDz = -999;
7c80a370 253
254 AliHLTCaloGlobalCoordinate globalCoord;
03a0ff8a 255
256 // 0 = assume photon
257 fGeometry->GetGlobalCoordinates(*recPointPtr, globalCoord, 0);
ea367f6a 258
353c33d3 259 caloClusterPtr->fModule = recPointPtr->fModule;
a6db7bcd 260 caloClusterPtr->fGlobalPos[0] = globalCoord.fX;
cd9e2797 261 caloClusterPtr->fGlobalPos[1] = globalCoord.fY;
a6db7bcd 262 caloClusterPtr->fGlobalPos[2] = globalCoord.fZ;
ea367f6a 263
5259027d 264 HLTDebug("Cluster local position: x = %f, z = %f, module = %d", recPointPtr->fX, recPointPtr->fZ, recPointPtr->fModule);
265 HLTDebug("Cluster global position: x = %f, y = %f, z = %f", globalCoord.fX, globalCoord.fY, globalCoord.fZ);
ee72b4e5 266
7fdb2519 267 //caloClusterPtr->fNCells = 0;//recPointPtr->fMultiplicity;
268 caloClusterPtr->fNCells = recPointPtr->fMultiplicity;
cd9e2797 269
a6db7bcd 270 caloClusterPtr->fClusterType = fClusterType;
7fdb2519 271// Int_t tmpSize = 0;//totSize + (caloClusterPtr->fNCells-1)*(sizeof(Short_t) + sizeof(Float_t));
a6db7bcd 272
2a24cbbe 273 //TODO remove hardcoded 10;
274 memset(caloClusterPtr->fTracksMatched, 0xff, sizeof(Int_t)*10);
275
9a0cbaba 276 //Int_t tmpSize = totSize + (caloClusterPtr->fNCells-1)*(sizeof(Short_t) + sizeof(Float_t));
cefad3b2 277 UInt_t tmpSize = (caloClusterPtr->fNCells-1)*sizeof(AliHLTCaloCellDataStruct);
cd9e2797 278
73d6f579 279 if((availableSize - totSize) < tmpSize)
280 {
9a0cbaba 281 HLTError("Out of buffer, available size is: %d, total size used: %d, extra size needed: %d", availableSize, totSize, tmpSize);
73d6f579 282 return -ENOBUFS;
283 }
c1e4a18c 284
cd9e2797 285 Int_t *digitIndexPtr = &(recPointPtr->fDigits);
286 Int_t id = 0;
287
7fdb2519 288 AliHLTCaloCellDataStruct *cellPtr = &(caloClusterPtr->fCaloCells);
289
ea367f6a 290 for(UInt_t j = 0; j < caloClusterPtr->fNCells; j++)
291 {
f92543f5 292 digitPtr = fDigitDataArray[*digitIndexPtr];
936dbe15 293 fGeometry->GetCellAbsId(recPointPtr->fModule, digitPtr->fX, digitPtr->fZ, id);
7fdb2519 294
295 cellPtr->fCellsAbsId= id;
296 cellPtr->fCellsAmpFraction = digitPtr->fEnergy/recPointPtr->fAmp;
cd9e2797 297 //printf("Cell ID pointer: %x\n", cellIDPtr);
298 //printf("Cell Amp Pointer: %x\n", cellAmpFracPtr);
9a0cbaba 299 //printf("Cell pos: x = %d, z = %d\n", digitPtr->fX, digitPtr->fZ);
34eb3256 300// printf("Cell ID: %d\n", cellPtr->fCellsAbsId);
9a0cbaba 301 //printf("Cell Amp: %f, pointer: %x\n", *cellAmpFracPtr, cellAmpFracPtr);
7fdb2519 302 cellPtr++;
cd9e2797 303 digitIndexPtr++;
304
ea367f6a 305 }
306
73d6f579 307 totSize += tmpSize;
308
c1e4a18c 309 if(fRecoParamsPtr)
310 {
311 caloClusterPtr->fEnergy = fRecoParamsPtr->GetCorrectedEnergy(recPointPtr->fAmp);
312 }
313 else
314 {
315 caloClusterPtr->fEnergy = recPointPtr->fAmp;
316 }
34eb3256 317
afcce849 318 // Set the time of the maximum digit as cluster time
319 caloClusterPtr->fTOF = recPointPtr->fTime;
320
321 HLTDebug("Cluster global position: x = %f, y = %f, z = %f, energy: %f, time: %f, number of cells: %d, cluster pointer: %x", globalCoord.fX, globalCoord.fY, globalCoord.fZ, caloClusterPtr->fEnergy, caloClusterPtr->fTOF, caloClusterPtr->fNCells, caloClusterPtr);
7fdb2519 322
ea367f6a 323 if(fDoClusterFit)
324 {
325 FitCluster(recPointPtr);
326 }
327 else
328 {
329 caloClusterPtr->fDispersion = 0;
330 caloClusterPtr->fFitQuality = 0;
331 caloClusterPtr->fM20 = 0;
332 caloClusterPtr->fM02 = 0;
333
334 }
335 if(fHaveCPVInfo)
336 {
337 caloClusterPtr->fEmcCpvDistance = GetCPVDistance(recPointPtr);
338 }
339 else
340 {
341 caloClusterPtr->fEmcCpvDistance = -1;
342 }
343 if(fDoPID)
344 {
345 DoParticleIdentification(caloClusterPtr);
346 }
347 else
348 {
349 for(Int_t k = 0; k < AliPID::kSPECIESN; k++)
350 {
351 caloClusterPtr->fPID[k] = 0;
352 }
353 }
354 if(fHaveDistanceToBadChannel)
355 {
356 caloClusterPtr->fDistanceToBadChannel = GetDistanceToBadChannel(caloClusterPtr);
357 }
358 else
359 {
360 caloClusterPtr->fDistanceToBadChannel = -1;
361 }
362
73d6f579 363 caloClusterPtr->fClusterType = fClusterType;
cd9e2797 364 //totSize += sizeof(AliHLTCaloClusterDataStruct) + (caloClusterPtr->fNCells)*(sizeof(Short_t) +sizeof(Float_t)-1);
73d6f579 365 //totSize += sizeof(AliHLTCaloClusterDataStruct) + (caloClusterPtr->fNCells-1)*(sizeof(Short_t) + sizeof(Float_t));
cd9e2797 366 //printf("CaloClusterPtr: %x, energy: %f\n", caloClusterPtr, caloClusterPtr->fEnergy);
367
ea367f6a 368 // caloClusterPtr = reinterpret_cast<AliHLTCaloClusterDataStruct*>(cellAmpFracPtr);
cd9e2797 369 //caloClusterPtr = reinterpret_cast<AliHLTCaloClusterDataStruct*>(cellIDPtr);
7fdb2519 370 fCaloClusterDataPtr = reinterpret_cast<AliHLTCaloClusterDataStruct*>(cellPtr);
ea367f6a 371 recPointPtr = reinterpret_cast<AliHLTCaloRecPointDataStruct*>(digitPtr);
73d6f579 372 //digitPtr = &(recPointPtr->fDigits);
ea367f6a 373 }
98baf84d 374
375return fNRecPoints;
ea367f6a 376
377}
378