Increase default size of ObjArray to avoid crash in PbPb data
[u/mrichter/AliRoot.git] / PWG / CaloTrackCorrBase / AliCalorimeterUtils.cxx
CommitLineData
765d44e7 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15/* $Id: $ */
16
17//_________________________________________________________________________
18// Class utility for Calorimeter specific selection methods ///
19//
20//
21//
22//-- Author: Gustavo Conesa (LPSC-Grenoble)
23//////////////////////////////////////////////////////////////////////////////
24
25
26// --- ROOT system ---
27#include "TGeoManager.h"
28
c5693f62 29// --- ANALYSIS system ---
765d44e7 30#include "AliCalorimeterUtils.h"
4b892846 31#include "AliESDEvent.h"
46a3cde6 32#include "AliMCEvent.h"
33#include "AliStack.h"
765d44e7 34#include "AliAODPWG4Particle.h"
c8fe2783 35#include "AliVCluster.h"
36#include "AliVCaloCells.h"
37#include "AliMixedEvent.h"
765d44e7 38
c5693f62 39// --- Detector ---
40#include "AliEMCALGeometry.h"
41#include "AliPHOSGeoUtils.h"
42
765d44e7 43ClassImp(AliCalorimeterUtils)
44
45
cb5780f4 46//____________________________________________
765d44e7 47 AliCalorimeterUtils::AliCalorimeterUtils() :
48 TObject(), fDebug(0),
7e4e67ab 49 fEMCALGeoName("EMCAL_COMPLETEV1"),fPHOSGeoName("PHOSgeo"),
9e8998b1 50 fEMCALGeo(0x0), fPHOSGeo(0x0),
51 fEMCALGeoMatrixSet(kFALSE), fPHOSGeoMatrixSet(kFALSE),
52 fLoadEMCALMatrices(kFALSE), fLoadPHOSMatrices(kFALSE),
53 fRemoveBadChannels(kFALSE), fPHOSBadChannelMap(0x0),
a5fb4114 54 fNCellsFromPHOSBorder(0),
9e8998b1 55 fNMaskCellColumns(0), fMaskCellColumns(0x0),
56 fRecalibration(kFALSE), fPHOSRecalibrationFactors(),
247abff4 57 fEMCALRecoUtils(new AliEMCALRecoUtils),
9e8998b1 58 fRecalculatePosition(kFALSE), fCorrectELinearity(kFALSE),
59 fRecalculateMatching(kFALSE),
60 fCutR(20), fCutZ(20),
61 fCutEta(20), fCutPhi(20)
765d44e7 62{
63 //Ctor
64
65 //Initialize parameters
66 InitParameters();
3b13c34c 67 for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
68 for(Int_t i = 0; i < 5 ; i++) fPHOSMatrix[i] = 0 ;
7cc7d3f8 69
765d44e7 70}
765d44e7 71
cb5780f4 72//_________________________________________
73AliCalorimeterUtils::~AliCalorimeterUtils()
74{
765d44e7 75 //Dtor
76
4df35693 77 //if(fPHOSGeo) delete fPHOSGeo ;
765d44e7 78 if(fEMCALGeo) delete fEMCALGeo ;
7cc7d3f8 79
765d44e7 80 if(fPHOSBadChannelMap) {
81 fPHOSBadChannelMap->Clear();
82 delete fPHOSBadChannelMap;
83 }
84
09e819c9 85 if(fPHOSRecalibrationFactors) {
78219bac 86 fPHOSRecalibrationFactors->Clear();
87 delete fPHOSRecalibrationFactors;
09e819c9 88 }
89
a5fb4114 90 if(fEMCALRecoUtils) delete fEMCALRecoUtils ;
91 if(fNMaskCellColumns) delete [] fMaskCellColumns;
9584c261 92
765d44e7 93}
94
cb5780f4 95//_____________________________________________________________________________________
96Bool_t AliCalorimeterUtils::CheckCellFiducialRegion(AliVCluster* cluster,
97 AliVCaloCells* cells,
98 AliVEvent * event, Int_t iev) const
99{
7cc7d3f8 100
765d44e7 101 // Given the list of AbsId of the cluster, get the maximum cell and
102 // check if there are fNCellsFromBorder from the calorimeter border
103
7cc7d3f8 104 //If the distance to the border is 0 or negative just exit accept all clusters
247abff4 105 if(cells->GetType()==AliVCaloCells::kEMCALCell && fEMCALRecoUtils->GetNumberOfCellsFromEMCALBorder() <= 0 ) return kTRUE;
c8fe2783 106 if(cells->GetType()==AliVCaloCells::kPHOSCell && fNCellsFromPHOSBorder <= 0 ) return kTRUE;
7cc7d3f8 107
c8fe2783 108 Int_t absIdMax = -1;
765d44e7 109 Float_t ampMax = -1;
c8fe2783 110
111 AliMixedEvent * mixEvent = dynamic_cast<AliMixedEvent*> (event);
112 Int_t nMixedEvents = 0 ;
113 Int_t * cellsCumul = NULL ;
114 Int_t numberOfCells = 0 ;
115 if (mixEvent){
116 nMixedEvents = mixEvent->GetNumberOfEvents() ;
117 if (cells->GetType()==AliVCaloCells::kEMCALCell) {
118 cellsCumul = mixEvent->GetEMCALCellsCumul() ;
119 numberOfCells = mixEvent->GetNumberOfEMCALCells() ;
120 }
7cc7d3f8 121
c8fe2783 122 else if (cells->GetType()==AliVCaloCells::kPHOSCell) {
123 cellsCumul = mixEvent->GetPHOSCellsCumul() ;
124 numberOfCells = mixEvent->GetNumberOfPHOSCells() ;
125 }
898c9d44 126
127 if(cellsCumul){
128
129 Int_t startCell = cellsCumul[iev] ;
130 Int_t endCell = (iev+1 < nMixedEvents)?cellsCumul[iev+1]:numberOfCells;
c8fe2783 131 //Find cells with maximum amplitude
898c9d44 132 for(Int_t i = 0; i < cluster->GetNCells() ; i++){
133 Int_t absId = cluster->GetCellAbsId(i) ;
134 for (Int_t j = startCell; j < endCell ; j++) {
135 Short_t cellNumber;
136 Double_t amp ;
137 Double_t time;
138 cells->GetCell(j, cellNumber, amp, time) ;
139 if (absId == cellNumber) {
140 if(amp > ampMax){
141 ampMax = amp;
142 absIdMax = absId;
143 }
144 }
c8fe2783 145 }
898c9d44 146 }//loop on cluster cells
147 }// cells cumul available
148 else {
149 printf("AliCalorimeterUtils::CheckCellFiducialRegion() - CellsCumul is NULL!!!\n");
150 abort();
c8fe2783 151 }
898c9d44 152 } else {//Normal SE Events
c8fe2783 153 for(Int_t i = 0; i < cluster->GetNCells() ; i++){
154 Int_t absId = cluster->GetCellAbsId(i) ;
155 Float_t amp = cells->GetCellAmplitude(absId);
156 if(amp > ampMax){
157 ampMax = amp;
158 absIdMax = absId;
159 }
160 }
161 }
765d44e7 162
163 if(fDebug > 1)
280e6766 164 printf("AliCalorimeterUtils::CheckCellFiducialRegion() - Cluster Max AbsId %d, Cell Energy %2.2f, Cluster Energy %2.2f\n",
7cc7d3f8 165 absIdMax, ampMax, cluster->E());
765d44e7 166
167 if(absIdMax==-1) return kFALSE;
168
169 //Check if the cell is close to the borders:
170 Bool_t okrow = kFALSE;
171 Bool_t okcol = kFALSE;
7cc7d3f8 172
c8fe2783 173 if(cells->GetType()==AliVCaloCells::kEMCALCell){
7cc7d3f8 174
765d44e7 175 Int_t iTower = -1, iIphi = -1, iIeta = -1, iphi = -1, ieta = -1, iSM = -1;
176 fEMCALGeo->GetCellIndex(absIdMax,iSM,iTower,iIphi,iIeta);
177 fEMCALGeo->GetCellPhiEtaIndexInSModule(iSM,iTower,iIphi, iIeta,iphi,ieta);
280e6766 178 if(iSM < 0 || iphi < 0 || ieta < 0 ) {
179 Fatal("CheckCellFidutialRegion","Negative value for super module: %d, or cell ieta: %d, or cell iphi: %d, check EMCAL geometry name\n",iSM,ieta,iphi);
180 }
7cc7d3f8 181
765d44e7 182 //Check rows/phi
247abff4 183 Int_t nborder = fEMCALRecoUtils->GetNumberOfCellsFromEMCALBorder();
765d44e7 184 if(iSM < 10){
247abff4 185 if(iphi >= nborder && iphi < 24-nborder) okrow =kTRUE;
7cc7d3f8 186 }
765d44e7 187 else{
247abff4 188 if(iphi >= nborder && iphi < 12-nborder) okrow =kTRUE;
765d44e7 189 }
190
247abff4 191 //Check columns/eta
192 if(!fEMCALRecoUtils->IsEMCALNoBorderAtEta0()){
193 if(ieta > nborder && ieta < 48-nborder) okcol =kTRUE;
765d44e7 194 }
195 else{
196 if(iSM%2==0){
247abff4 197 if(ieta >= nborder) okcol = kTRUE;
765d44e7 198 }
199 else {
247abff4 200 if(ieta < 48-nborder) okcol = kTRUE;
765d44e7 201 }
202 }//eta 0 not checked
203 if(fDebug > 1)
204 {
280e6766 205 printf("AliCalorimeterUtils::CheckCellFiducialRegion() - EMCAL Cluster in %d cells fiducial volume: ieta %d, iphi %d, SM %d ?",
7cc7d3f8 206 nborder, ieta, iphi, iSM);
765d44e7 207 if (okcol && okrow ) printf(" YES \n");
208 else printf(" NO: column ok? %d, row ok? %d \n",okcol,okrow);
209 }
210 }//EMCAL
c8fe2783 211 else if(cells->GetType()==AliVCaloCells::kPHOSCell){
765d44e7 212 Int_t relId[4];
213 Int_t irow = -1, icol = -1;
214 fPHOSGeo->AbsToRelNumbering(absIdMax,relId);
215 irow = relId[2];
216 icol = relId[3];
217 //imod = relId[0]-1;
218 if(irow >= fNCellsFromPHOSBorder && irow < 64-fNCellsFromPHOSBorder) okrow =kTRUE;
219 if(icol >= fNCellsFromPHOSBorder && icol < 56-fNCellsFromPHOSBorder) okcol =kTRUE;
220 if(fDebug > 1)
221 {
280e6766 222 printf("AliCalorimeterUtils::CheckCellFiducialRegion() - PHOS Cluster in %d cells fiducial volume: icol %d, irow %d, Module %d?",
7cc7d3f8 223 fNCellsFromPHOSBorder, icol, irow, relId[0]-1);
765d44e7 224 if (okcol && okrow ) printf(" YES \n");
225 else printf(" NO: column ok? %d, row ok? %d \n",okcol,okrow);
226 }
227 }//PHOS
228
229 if (okcol && okrow) return kTRUE;
230 else return kFALSE;
231
232}
233
765d44e7 234//_________________________________________________________________________________________________________
cb5780f4 235Bool_t AliCalorimeterUtils::ClusterContainsBadChannel(TString calorimeter,UShort_t* cellList, Int_t nCells)
236{
765d44e7 237 // Check that in the cluster cells, there is no bad channel of those stored
238 // in fEMCALBadChannelMap or fPHOSBadChannelMap
239
240 if (!fRemoveBadChannels) return kFALSE;
36037088 241 //printf("fEMCALBadChannelMap %p, fPHOSBadChannelMap %p \n",fEMCALBadChannelMap,fPHOSBadChannelMap);
247abff4 242 if(calorimeter == "EMCAL" && !fEMCALRecoUtils->GetEMCALChannelStatusMap(0)) return kFALSE;
78219bac 243 if(calorimeter == "PHOS" && !fPHOSBadChannelMap) return kFALSE;
7cc7d3f8 244
765d44e7 245 Int_t icol = -1;
246 Int_t irow = -1;
247 Int_t imod = -1;
248 for(Int_t iCell = 0; iCell<nCells; iCell++){
7cc7d3f8 249
765d44e7 250 //Get the column and row
251 if(calorimeter == "EMCAL"){
247abff4 252 return fEMCALRecoUtils->ClusterContainsBadChannel((AliEMCALGeometry*)fEMCALGeo,cellList,nCells);
765d44e7 253 }
254 else if(calorimeter=="PHOS"){
255 Int_t relId[4];
256 fPHOSGeo->AbsToRelNumbering(cellList[iCell],relId);
257 irow = relId[2];
258 icol = relId[3];
259 imod = relId[0]-1;
260 if(fPHOSBadChannelMap->GetEntries() <= imod)continue;
c5693f62 261 //printf("PHOS bad channels imod %d, icol %d, irow %d\n",imod, irow, icol);
262 if(GetPHOSChannelStatus(imod, irow, icol)) return kTRUE;
765d44e7 263 }
264 else return kFALSE;
265
266 }// cell cluster loop
7cc7d3f8 267
765d44e7 268 return kFALSE;
7cc7d3f8 269
765d44e7 270}
271
cb5780f4 272//_______________________________________________________________
273void AliCalorimeterUtils::CorrectClusterEnergy(AliVCluster *clus)
274{
9584c261 275 // Correct cluster energy non linearity
13cd2872 276
9584c261 277 clus->SetE(fEMCALRecoUtils->CorrectClusterEnergyLinearity(clus));
13cd2872 278
279}
280
c5693f62 281//________________________________________________________________________________________
282Int_t AliCalorimeterUtils::GetMaxEnergyCell(AliVCaloCells* cells, const AliVCluster* clu,
283 Float_t & clusterFraction) const
cb5780f4 284{
13cd2872 285
286 //For a given CaloCluster gets the absId of the cell
287 //with maximum energy deposit.
288
289 if( !clu || !cells ){
290 AliInfo("Cluster or cells pointer is null!");
291 return -1;
292 }
293
294 Double_t eMax =-1.;
295 Double_t eTot = 0.;
296 Double_t eCell =-1.;
297 Float_t fraction = 1.;
298 Float_t recalFactor = 1.;
299 Int_t cellAbsId =-1 , absId =-1 ;
300 Int_t iSupMod =-1 , ieta =-1 , iphi = -1, iRCU = -1;
301
302 TString calo = "EMCAL";
303 if(clu->IsPHOS()) calo = "PHOS";
304
305 for (Int_t iDig=0; iDig< clu->GetNCells(); iDig++) {
306
307 cellAbsId = clu->GetCellAbsId(iDig);
308
309 fraction = clu->GetCellAmplitudeFraction(iDig);
310 if(fraction < 1e-4) fraction = 1.; // in case unfolding is off
311
b08dd6d5 312 iSupMod = GetModuleNumberCellIndexes(cellAbsId, calo, ieta, iphi, iRCU);
13cd2872 313
314 if(IsRecalibrationOn()) {
315 if(calo=="EMCAL") recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
c5693f62 316 else recalFactor = GetPHOSChannelRecalibrationFactor (iSupMod,iphi,ieta);
13cd2872 317 }
318
319 eCell = cells->GetCellAmplitude(cellAbsId)*fraction*recalFactor;
320
321 if(eCell > eMax) {
322 eMax = eCell;
323 absId = cellAbsId;
324 }
325
326 eTot+=eCell;
327
328 }// cell loop
329
330 clusterFraction = (eTot-eMax)/eTot; //Do not use cluster energy in case it was corrected for non linearity.
331
332 return absId;
333
9584c261 334}
335
cb5780f4 336//_____________________________________________________________________________________________________
765d44e7 337Int_t AliCalorimeterUtils::GetModuleNumber(AliAODPWG4Particle * particle, AliVEvent * inputEvent) const
338{
339 //Get the EMCAL/PHOS module number that corresponds to this particle
340
341 Int_t absId = -1;
342 if(particle->GetDetector()=="EMCAL"){
343 fEMCALGeo->GetAbsCellIdFromEtaPhi(particle->Eta(),particle->Phi(), absId);
344 if(fDebug > 2)
345 printf("AliCalorimeterUtils::GetModuleNumber(PWG4AOD) - EMCAL: cluster eta %f, phi %f, absid %d, SuperModule %d\n",
7cc7d3f8 346 particle->Eta(), particle->Phi()*TMath::RadToDeg(),absId, fEMCALGeo->GetSuperModuleNumber(absId));
765d44e7 347 return fEMCALGeo->GetSuperModuleNumber(absId) ;
348 }//EMCAL
349 else if(particle->GetDetector()=="PHOS"){
46a3cde6 350 // In case we use the MC reader, the input are TParticles,
351 // in this case use the corresponing method in PHOS Geometry to get the particle.
352 if(strcmp(inputEvent->ClassName(), "AliMCEvent") == 0 )
353 {
354 Int_t mod =-1;
355 Double_t z = 0., x=0.;
356 TParticle* primary = 0x0;
357 AliStack * stack = ((AliMCEvent*)inputEvent)->Stack();
358 if(stack) {
359 primary = stack->Particle(particle->GetCaloLabel(0));
360 }
361 else {
280e6766 362 Fatal("GetModuleNumber(PWG4AOD)", "Stack not available, stop!");
46a3cde6 363 }
898c9d44 364
365 if(primary){
366 fPHOSGeo->ImpactOnEmc(primary,mod,z,x) ;
367 }
368 else{
369 Fatal("GetModuleNumber(PWG4AOD)", "Primary not available, stop!");
370 }
46a3cde6 371 return mod;
372 }
373 // Input are ESDs or AODs, get the PHOS module number like this.
374 else{
2206e918 375 //FIXME
376 //AliVCluster *cluster = inputEvent->GetCaloCluster(particle->GetCaloLabel(0));
377 //return GetModuleNumber(cluster);
378 //MEFIX
379 return -1;
46a3cde6 380 }
765d44e7 381 }//PHOS
382
383 return -1;
384}
385
cb5780f4 386//_____________________________________________________________________
c8fe2783 387Int_t AliCalorimeterUtils::GetModuleNumber(AliVCluster * cluster) const
765d44e7 388{
280e6766 389 //Get the EMCAL/PHOS module number that corresponds to this cluster
765d44e7 390 TLorentzVector lv;
391 Double_t v[]={0.,0.,0.}; //not necessary to pass the real vertex.
9cbbc28b 392 if(!cluster){
393 if(fDebug > 1) printf("AliCalorimeterUtils::GetModuleNumber() - NUL Cluster, please check!!!");
394 return -1;
395 }
765d44e7 396 cluster->GetMomentum(lv,v);
397 Float_t phi = lv.Phi();
398 if(phi < 0) phi+=TMath::TwoPi();
399 Int_t absId = -1;
c8fe2783 400 if(cluster->IsEMCAL()){
765d44e7 401 fEMCALGeo->GetAbsCellIdFromEtaPhi(lv.Eta(),phi, absId);
402 if(fDebug > 2)
280e6766 403 printf("AliCalorimeterUtils::GetModuleNumber() - EMCAL: cluster eta %f, phi %f, absid %d, SuperModule %d\n",
7cc7d3f8 404 lv.Eta(), phi*TMath::RadToDeg(),absId, fEMCALGeo->GetSuperModuleNumber(absId));
765d44e7 405 return fEMCALGeo->GetSuperModuleNumber(absId) ;
406 }//EMCAL
c8fe2783 407 else if(cluster->IsPHOS()) {
765d44e7 408 Int_t relId[4];
409 if ( cluster->GetNCells() > 0) {
410 absId = cluster->GetCellAbsId(0);
411 if(fDebug > 2)
280e6766 412 printf("AliCalorimeterUtils::GetModuleNumber() - PHOS: cluster eta %f, phi %f, e %f, absId %d\n",
7cc7d3f8 413 lv.Eta(), phi*TMath::RadToDeg(), lv.E(), absId);
765d44e7 414 }
415 else return -1;
416
417 if ( absId >= 0) {
418 fPHOSGeo->AbsToRelNumbering(absId,relId);
419 if(fDebug > 2)
280e6766 420 printf("AliCalorimeterUtils::GetModuleNumber() - PHOS: Module %d\n",relId[0]-1);
765d44e7 421 return relId[0]-1;
422 }
423 else return -1;
424 }//PHOS
425
426 return -1;
427}
428
cb5780f4 429//___________________________________________________________________________________________________
430Int_t AliCalorimeterUtils::GetModuleNumberCellIndexes(const Int_t absId, const TString calo,
431 Int_t & icol, Int_t & irow, Int_t & iRCU) const
765d44e7 432{
433 //Get the EMCAL/PHOS module, columns, row and RCU number that corresponds to this absId
434 Int_t imod = -1;
435 if ( absId >= 0) {
436 if(calo=="EMCAL"){
437 Int_t iTower = -1, iIphi = -1, iIeta = -1;
438 fEMCALGeo->GetCellIndex(absId,imod,iTower,iIphi,iIeta);
439 fEMCALGeo->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol);
280e6766 440 if(imod < 0 || irow < 0 || icol < 0 ) {
441 Fatal("GetModuleNumberCellIndexes()","Negative value for super module: %d, or cell icol: %d, or cell irow: %d, check EMCAL geometry name\n",imod,icol,irow);
442 }
443
765d44e7 444 //RCU0
445 if (0<=irow&&irow<8) iRCU=0; // first cable row
446 else if (8<=irow&&irow<16 && 0<=icol&&icol<24) iRCU=0; // first half;
447 //second cable row
448 //RCU1
449 else if(8<=irow&&irow<16 && 24<=icol&&icol<48) iRCU=1; // second half;
450 //second cable row
451 else if(16<=irow&&irow<24) iRCU=1; // third cable row
452
453 if (imod%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same
454 if (iRCU<0) {
280e6766 455 Fatal("GetModuleNumberCellIndexes()","Wrong EMCAL RCU number = %d\n", iRCU);
765d44e7 456 }
457
458 return imod ;
459 }//EMCAL
460 else{//PHOS
461 Int_t relId[4];
462 fPHOSGeo->AbsToRelNumbering(absId,relId);
463 irow = relId[2];
464 icol = relId[3];
465 imod = relId[0]-1;
466 iRCU= (Int_t)(relId[2]-1)/16 ;
467 //Int_t iBranch= (Int_t)(relid[3]-1)/28 ; //0 to 1
468 if (iRCU >= 4) {
280e6766 469 Fatal("GetModuleNumberCellIndexes()","Wrong PHOS RCU number = %d\n", iRCU);
765d44e7 470 }
471 return imod;
472 }//PHOS
473 }
474
475 return -1;
476}
477
cb5780f4 478//________________________________________
765d44e7 479void AliCalorimeterUtils::InitParameters()
480{
481 //Initialize the parameters of the analysis.
7cc7d3f8 482 fEMCALGeoName = "EMCAL_COMPLETEV1";
483 fPHOSGeoName = "PHOSgeo";
484 fEMCALGeoMatrixSet = kFALSE;
485 fPHOSGeoMatrixSet = kFALSE;
486 fRemoveBadChannels = kFALSE;
487 fNCellsFromPHOSBorder = 0;
a5fb4114 488
7cc7d3f8 489 // fMaskCellColumns = new Int_t[fNMaskCellColumns];
490 // fMaskCellColumns[0] = 6 ; fMaskCellColumns[1] = 7 ; fMaskCellColumns[2] = 8 ;
491 // fMaskCellColumns[3] = 35; fMaskCellColumns[4] = 36; fMaskCellColumns[5] = 37;
492 // fMaskCellColumns[6] = 12+AliEMCALGeoParams::fgkEMCALCols; fMaskCellColumns[7] = 13+AliEMCALGeoParams::fgkEMCALCols;
493 // fMaskCellColumns[8] = 40+AliEMCALGeoParams::fgkEMCALCols; fMaskCellColumns[9] = 41+AliEMCALGeoParams::fgkEMCALCols;
494 // fMaskCellColumns[10]= 42+AliEMCALGeoParams::fgkEMCALCols;
a5fb4114 495
765d44e7 496}
497
765d44e7 498
cb5780f4 499//_____________________________________________________
500void AliCalorimeterUtils::InitPHOSBadChannelStatusMap()
501{
765d44e7 502 //Init PHOS bad channels map
503 if(fDebug > 0 )printf("AliCalorimeterUtils::InitPHOSBadChannelStatusMap()\n");
504 //In order to avoid rewriting the same histograms
505 Bool_t oldStatus = TH1::AddDirectoryStatus();
506 TH1::AddDirectory(kFALSE);
7cc7d3f8 507
78219bac 508 fPHOSBadChannelMap = new TObjArray(5);
c5693f62 509 for (int i = 0; i < 5; i++)fPHOSBadChannelMap->Add(new TH2I(Form("PHOS_BadMap_mod%d",i),Form("PHOS_BadMap_mod%d",i), 64, 0, 64, 56, 0, 56));
7cc7d3f8 510
765d44e7 511 fPHOSBadChannelMap->SetOwner(kTRUE);
512 fPHOSBadChannelMap->Compress();
513
514 //In order to avoid rewriting the same histograms
515 TH1::AddDirectory(oldStatus);
516}
517
cb5780f4 518//______________________________________________________
519void AliCalorimeterUtils::InitPHOSRecalibrationFactors()
520{
09e819c9 521 //Init EMCAL recalibration factors
522 if(fDebug > 0 )printf("AliCalorimeterUtils::InitPHOSRecalibrationFactors()\n");
523 //In order to avoid rewriting the same histograms
524 Bool_t oldStatus = TH1::AddDirectoryStatus();
525 TH1::AddDirectory(kFALSE);
7cc7d3f8 526
78219bac 527 fPHOSRecalibrationFactors = new TObjArray(5);
c5693f62 528 for (int i = 0; i < 5; i++)fPHOSRecalibrationFactors->Add(new TH2F(Form("PHOSRecalFactors_Mod%d",i),Form("PHOSRecalFactors_Mod%d",i), 64, 0, 64, 56, 0, 56));
09e819c9 529 //Init the histograms with 1
530 for (Int_t m = 0; m < 5; m++) {
531 for (Int_t i = 0; i < 56; i++) {
532 for (Int_t j = 0; j < 64; j++) {
c5693f62 533 SetPHOSChannelRecalibrationFactor(m,j,i,1.);
09e819c9 534 }
535 }
536 }
537 fPHOSRecalibrationFactors->SetOwner(kTRUE);
538 fPHOSRecalibrationFactors->Compress();
539
540 //In order to avoid rewriting the same histograms
541 TH1::AddDirectory(oldStatus);
542}
543
544
cb5780f4 545//___________________________________________
765d44e7 546void AliCalorimeterUtils::InitEMCALGeometry()
547{
548 //Initialize EMCAL geometry if it did not exist previously
549 if (!fEMCALGeo){
a38a48f2 550 fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName);
765d44e7 551 if(fDebug > 0){
552 printf("AliCalorimeterUtils::InitEMCALGeometry()");
553 if (!gGeoManager) printf(" - Careful!, gGeoManager not loaded, load misalign matrices");
554 printf("\n");
555 }
556 }
557}
558
cb5780f4 559//__________________________________________
765d44e7 560void AliCalorimeterUtils::InitPHOSGeometry()
561{
562 //Initialize PHOS geometry if it did not exist previously
563 if (!fPHOSGeo){
564 fPHOSGeo = new AliPHOSGeoUtils(fPHOSGeoName);
565 if(fDebug > 0){
566 printf("AliCalorimeterUtils::InitPHOSGeometry()");
567 if (!gGeoManager) printf(" - Careful!, gGeoManager not loaded, load misalign matrices");
568 printf("\n");
569 }
570 }
571}
572
cb5780f4 573//_________________________________________________________
765d44e7 574void AliCalorimeterUtils::Print(const Option_t * opt) const
575{
7cc7d3f8 576
765d44e7 577 //Print some relevant parameters set for the analysis
578 if(! opt)
579 return;
7cc7d3f8 580
765d44e7 581 printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
582 printf("Remove Clusters with bad channels? %d\n",fRemoveBadChannels);
583 printf("Remove Clusters with max cell at less than %d cells from EMCAL border and %d cells from PHOS border\n",
7cc7d3f8 584 fEMCALRecoUtils->GetNumberOfCellsFromEMCALBorder(), fNCellsFromPHOSBorder);
247abff4 585 if(fEMCALRecoUtils->IsEMCALNoBorderAtEta0()) printf("Do not remove EMCAL clusters at Eta = 0\n");
09e819c9 586 printf("Recalibrate Clusters? %d\n",fRecalibration);
9584c261 587 printf("Recalculate Clusters Position? %d\n",fRecalculatePosition);
588 printf("Recalculate Clusters Energy? %d\n",fCorrectELinearity);
f2ccb5b8 589 printf("Matching criteria: dR < %2.2f[cm], dZ < %2.2f[cm]\n",fCutR,fCutZ);
7cc7d3f8 590
765d44e7 591 printf(" \n") ;
592}
593
cb5780f4 594//__________________________________________________________________
595Bool_t AliCalorimeterUtils::MaskFrameCluster(const Int_t iSM,
596 const Int_t ieta) const
597{
a5fb4114 598 //Check if cell is in one of the regions where we have significant amount
599 //of material in front. Only EMCAL
600
601 Int_t icol = ieta;
602 if(iSM%2) icol+=48; // Impair SM, shift index [0-47] to [48-96]
603
604 if (fNMaskCellColumns && fMaskCellColumns) {
605 for (Int_t imask = 0; imask < fNMaskCellColumns; imask++) {
606 if(icol==fMaskCellColumns[imask]) return kTRUE;
607 }
608 }
609
610 return kFALSE;
611
612}
613
614
cb5780f4 615//__________________________________________________________________________
616Float_t AliCalorimeterUtils::RecalibrateClusterEnergy(AliVCluster * cluster,
617 AliVCaloCells * cells)
618{
09e819c9 619 // Recalibrate the cluster energy, considering the recalibration map and the energy of the cells that compose the cluster.
7cc7d3f8 620
898c9d44 621 //Initialize some used variables
622 Float_t energy = 0;
623 Int_t absId = -1;
624 Int_t icol = -1, irow = -1, iRCU = -1, module=1;
625 Float_t factor = 1, frac = 0;
626
627 if(cells) {
7cc7d3f8 628
629 //Get the cluster number of cells and list of absId, check what kind of cluster do we have.
630 UShort_t * index = cluster->GetCellsAbsId() ;
631 Double_t * fraction = cluster->GetCellsAmplitudeFraction() ;
632 Int_t ncells = cluster->GetNCells();
633 TString calo = "EMCAL";
634 if(cluster->IsPHOS()) calo = "PHOS";
635
636 //Loop on the cells, get the cell amplitude and recalibration factor, multiply and and to the new energy
637 for(Int_t icell = 0; icell < ncells; icell++){
638 absId = index[icell];
639 frac = fraction[icell];
640 if(frac < 1e-3) frac = 1; //in case of EMCAL, this is set as 0, not used.
641 module = GetModuleNumberCellIndexes(absId,calo,icol,irow,iRCU);
c5693f62 642 if(cluster->IsPHOS()) factor = GetPHOSChannelRecalibrationFactor (module,irow,icol);
7cc7d3f8 643 else factor = GetEMCALChannelRecalibrationFactor(module,icol,irow);
644 if(fDebug>2)
645 printf("AliCalorimeterUtils::RecalibrateClusterEnergy() - recalibrate cell: %s, module %d, col %d, row %d, cell fraction %f, recalibration factor %f, cell energy %f\n",
646 calo.Data(),module,icol,irow,frac,factor,cells->GetCellAmplitude(absId));
647
648 energy += cells->GetCellAmplitude(absId)*factor*frac;
649 }
650
651 if(fDebug>1)
652 printf("AliCalorimeterUtils::RecalibrateClusterEnergy() - Energy before %f, after %f\n",cluster->E(),energy);
898c9d44 653
654 }// cells available
655 else{
656 Fatal("RecalibrateClusterEnergy()","Cells pointer does not exist!");
657 }
658
09e819c9 659 return energy;
09e819c9 660}
661
cb5780f4 662//________________________________________________________________________________
765d44e7 663void AliCalorimeterUtils::SetGeometryTransformationMatrices(AliVEvent* inputEvent)
664{
665 //Set the calorimeters transformation matrices
7cc7d3f8 666
765d44e7 667 //Get the EMCAL transformation geometry matrices from ESD
3b13c34c 668 if(!fEMCALGeoMatrixSet && fEMCALGeo){
669 if(fLoadEMCALMatrices){
c5693f62 670 printf("AliCalorimeterUtils::SetGeometryTransformationMatrices() - Load user defined EMCAL geometry matrices\n");
3b13c34c 671 for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
672 if(fEMCALMatrix[mod]){
673 if(fDebug > 1)
674 fEMCALMatrix[mod]->Print();
675 fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod) ;
676 }
677 }//SM loop
678 fEMCALGeoMatrixSet = kTRUE;//At least one, so good
679
680 }//Load matrices
681 else if (!gGeoManager) {
7cc7d3f8 682
3b13c34c 683 if(fDebug > 1)
684 printf(" AliCalorimeterUtils::SetGeometryTransformationMatrices() - Load EMCAL misalignment matrices. \n");
685 if(!strcmp(inputEvent->GetName(),"AliESDEvent")) {
686 for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
7cc7d3f8 687 //printf("Load ESD matrix %d, %p\n",mod,((AliESDEvent*)inputEvent)->GetEMCALMatrix(mod));
4b892846 688 if(((AliESDEvent*)inputEvent)->GetEMCALMatrix(mod)) {
689 fEMCALGeo->SetMisalMatrix(((AliESDEvent*)inputEvent)->GetEMCALMatrix(mod),mod) ;
3b13c34c 690 }
691 }// loop over super modules
c5693f62 692
3b13c34c 693 fEMCALGeoMatrixSet = kTRUE;//At least one, so good
694
695 }//ESD as input
696 else {
697 if(fDebug > 1)
698 printf("AliCalorimeterUtils::SetGeometryTransformationMatrices() - Setting of EMCAL transformation matrixes for AODs not implemented yet. \n Import geometry.root file\n");
699 }//AOD as input
700 }//Get matrix from data
04131edb 701 else if(gGeoManager){
702 fEMCALGeoMatrixSet = kTRUE;
703 }
765d44e7 704 }//EMCAL geo && no geoManager
c5693f62 705
765d44e7 706 //Get the PHOS transformation geometry matrices from ESD
3b13c34c 707 if(!fPHOSGeoMatrixSet && fPHOSGeo){
c5693f62 708
3b13c34c 709 if(fLoadPHOSMatrices){
c5693f62 710 printf("AliCalorimeterUtils::SetGeometryTransformationMatrices() - Load user defined PHOS geometry matrices\n");
711 for(Int_t mod = 0 ; mod < 5 ; mod++){
3b13c34c 712 if(fPHOSMatrix[mod]){
c5693f62 713 if(fDebug > 1 )
3b13c34c 714 fPHOSMatrix[mod]->Print();
dbf54f1e 715 fPHOSGeo->SetMisalMatrix(fPHOSMatrix[mod],mod) ;
3b13c34c 716 }
717 }//SM loop
718 fPHOSGeoMatrixSet = kTRUE;//At least one, so good
c5693f62 719
3b13c34c 720 }//Load matrices
721 else if (!gGeoManager) {
722 if(fDebug > 1)
723 printf(" AliCalorimeterUtils::SetGeometryTransformationMatrices() - Load PHOS misalignment matrices. \n");
765d44e7 724 if(!strcmp(inputEvent->GetName(),"AliESDEvent")) {
dbf54f1e 725 for(Int_t mod = 0; mod < 5; mod++){
4b892846 726 if( ((AliESDEvent*)inputEvent)->GetPHOSMatrix(mod)) {
765d44e7 727 //printf("PHOS: mod %d, matrix %p\n",mod, ((AliESDEvent*)inputEvent)->GetPHOSMatrix(mod));
4b892846 728 fPHOSGeo->SetMisalMatrix( ((AliESDEvent*)inputEvent)->GetPHOSMatrix(mod),mod) ;
765d44e7 729 }
3b13c34c 730 }// loop over modules
731 fPHOSGeoMatrixSet = kTRUE; //At least one so good
765d44e7 732 }//ESD as input
733 else {
734 if(fDebug > 1)
735 printf("AliCalorimeterUtils::SetGeometryTransformationMatrices() - Setting of EMCAL transformation matrixes for AODs not implemented yet. \n Import geometry.root file\n");
3b13c34c 736 }//AOD as input
737 }// get matrix from data
04131edb 738 else if(gGeoManager){
739 fPHOSGeoMatrixSet = kTRUE;
740 }
765d44e7 741 }//PHOS geo and geoManager was not set
c5693f62 742
765d44e7 743}
744
cb5780f4 745//__________________________________________________________________________________________
746void AliCalorimeterUtils::RecalculateClusterPosition(AliVCaloCells* cells, AliVCluster* clu)
747{
9584c261 748
749 //Recalculate EMCAL cluster position
750
19db8f8c 751 fEMCALRecoUtils->RecalculateClusterPosition((AliEMCALGeometry*)fEMCALGeo, cells,clu);
7cc7d3f8 752
9584c261 753}
cb5780f4 754
755//________________________________________________________________________________
756void AliCalorimeterUtils::RecalculateClusterTrackMatching(AliVEvent * event,
757 TObjArray* clusterArray)
758{
759 //Recalculate track matching
760
761 if (fRecalculateMatching) {
762 fEMCALRecoUtils->FindMatches(event,clusterArray,fEMCALGeo) ;
763 //AliESDEvent* esdevent = dynamic_cast<AliESDEvent*> (event);
764 //if(esdevent){
765 // fEMCALRecoUtils->SetClusterMatchedToTrack(esdevent) ;
766 // fEMCALRecoUtils->SetTracksMatchedToCluster(esdevent) ;
767 //}
768 }
769}