correct filling and initializaton of transformation matrices for ESD analysis (Dmitri)
[u/mrichter/AliRoot.git] / PWG4 / PartCorrBase / 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
29//---- ANALYSIS system ----
30#include "AliCalorimeterUtils.h"
31#include "AliAODEvent.h"
32#include "AliESDEvent.h"
33#include "AliAODPWG4Particle.h"
34#include "AliAODCaloCluster.h"
35#include "AliESDCaloCluster.h"
36#include "AliAODCaloCells.h"
37#include "AliESDCaloCells.h"
38
39
40ClassImp(AliCalorimeterUtils)
41
42
43//____________________________________________________________________________
44 AliCalorimeterUtils::AliCalorimeterUtils() :
45 TObject(), fDebug(0),
46 fEMCALGeoName("EMCAL_COMPLETE"),fPHOSGeoName("PHOSgeo"),
47 fEMCALGeo(0x0), fPHOSGeo(0x0),
48 fEMCALGeoMatrixSet(kFALSE), fPHOSGeoMatrixSet(kFALSE),
49 fRemoveBadChannels(kFALSE),
78219bac 50 fEMCALBadChannelMap(0x0),fPHOSBadChannelMap(0x0),
765d44e7 51 fNCellsFromEMCALBorder(0), fNCellsFromPHOSBorder(0),
09e819c9 52 fNoEMCALBorderAtEta0(kFALSE),fRecalibration(kFALSE),
78219bac 53 fEMCALRecalibrationFactors(), fPHOSRecalibrationFactors()
765d44e7 54{
55 //Ctor
56
57 //Initialize parameters
58 InitParameters();
59}
78219bac 60/*
765d44e7 61//____________________________________________________________________________
62AliCalorimeterUtils::AliCalorimeterUtils(const AliCalorimeterUtils & calo) :
63 TObject(calo), fDebug(calo.fDebug),
64 fEMCALGeoName(calo.fEMCALGeoName), fPHOSGeoName(calo.fPHOSGeoName),
65 fEMCALGeo(new AliEMCALGeoUtils(*calo.fEMCALGeo)),
66 fPHOSGeo(new AliPHOSGeoUtils(*calo.fPHOSGeo)),
67 fEMCALGeoMatrixSet(calo.fEMCALGeoMatrixSet),
68 fPHOSGeoMatrixSet(calo.fPHOSGeoMatrixSet),
69 fRemoveBadChannels(calo.fRemoveBadChannels),
70 fEMCALBadChannelMap(new TObjArray(*calo.fEMCALBadChannelMap)),
71 fPHOSBadChannelMap(new TObjArray(*calo.fPHOSBadChannelMap)),
72 fNCellsFromEMCALBorder(calo.fNCellsFromEMCALBorder),
73 fNCellsFromPHOSBorder(calo.fNCellsFromPHOSBorder),
09e819c9 74 fNoEMCALBorderAtEta0(calo.fNoEMCALBorderAtEta0),
75 fRecalibration(calo.fRecalibration),
76 fEMCALRecalibrationFactors(new TObjArray(*calo.fEMCALRecalibrationFactors)),
77 fPHOSRecalibrationFactors(new TObjArray(*calo.fEMCALRecalibrationFactors))
765d44e7 78{
79 // cpy ctor
80}
78219bac 81*/
765d44e7 82
83//_________________________________
84AliCalorimeterUtils::~AliCalorimeterUtils() {
85 //Dtor
86
4df35693 87 //if(fPHOSGeo) delete fPHOSGeo ;
765d44e7 88 if(fEMCALGeo) delete fEMCALGeo ;
89
90 if(fEMCALBadChannelMap) {
91 fEMCALBadChannelMap->Clear();
92 delete fEMCALBadChannelMap;
93 }
94 if(fPHOSBadChannelMap) {
95 fPHOSBadChannelMap->Clear();
96 delete fPHOSBadChannelMap;
97 }
98
09e819c9 99 if(fEMCALRecalibrationFactors) {
78219bac 100 fEMCALRecalibrationFactors->Clear();
101 delete fEMCALRecalibrationFactors;
09e819c9 102 }
103 if(fPHOSRecalibrationFactors) {
78219bac 104 fPHOSRecalibrationFactors->Clear();
105 delete fPHOSRecalibrationFactors;
09e819c9 106 }
107
765d44e7 108}
109
110//_______________________________________________________________
111Bool_t AliCalorimeterUtils::CheckCellFiducialRegion(AliAODCaloCluster* cluster, AliAODCaloCells* cells) const {
112 // Given the list of AbsId of the cluster, get the maximum cell and
113 // check if there are fNCellsFromBorder from the calorimeter border
114
115 //If the distance to the border is 0 or negative just exit accept all clusters
116 if(cells->GetType()==AliAODCaloCells::kEMCAL && fNCellsFromEMCALBorder <= 0 ) return kTRUE;
117 if(cells->GetType()==AliAODCaloCells::kPHOS && fNCellsFromPHOSBorder <= 0 ) return kTRUE;
118
119 //Find cells with maximum amplitude
120 Int_t absIdMax = -1;
121 Float_t ampMax = -1;
122 for(Int_t i = 0; i < cluster->GetNCells() ; i++){
123 Int_t absId = cluster->GetCellAbsId(i) ;
124 Float_t amp = cells->GetCellAmplitude(absId);
125 if(amp > ampMax){
126 ampMax = amp;
127 absIdMax = absId;
128 }
129 }
130
131 if(fDebug > 1)
132 printf("AliCalorimeterUtils::CheckCellFiducialRegion(AOD) - Cluster Max AbsId %d, Cell Energy %2.2f, Cluster Energy %2.2f\n",
133 absIdMax, ampMax, cluster->E());
134
135 if(absIdMax==-1) return kFALSE;
136
137 //Check if the cell is close to the borders:
138 Bool_t okrow = kFALSE;
139 Bool_t okcol = kFALSE;
140
141 if(cells->GetType()==AliAODCaloCells::kEMCAL){
142
143 Int_t iTower = -1, iIphi = -1, iIeta = -1, iphi = -1, ieta = -1, iSM = -1;
144 fEMCALGeo->GetCellIndex(absIdMax,iSM,iTower,iIphi,iIeta);
145 fEMCALGeo->GetCellPhiEtaIndexInSModule(iSM,iTower,iIphi, iIeta,iphi,ieta);
146
147 //Check rows/phi
148 if(iSM < 10){
149 if(iphi >= fNCellsFromEMCALBorder && iphi < 24-fNCellsFromEMCALBorder) okrow =kTRUE;
150 }
151 else{
152 if(iphi >= fNCellsFromEMCALBorder && iphi < 12-fNCellsFromEMCALBorder) okrow =kTRUE;
153 }
154
155 //Check collumns/eta
156 if(!fNoEMCALBorderAtEta0){
157 if(ieta > fNCellsFromEMCALBorder && ieta < 48-fNCellsFromEMCALBorder) okcol =kTRUE;
158 }
159 else{
160 if(iSM%2==0){
161 if(ieta >= fNCellsFromEMCALBorder) okcol = kTRUE;
162 }
163 else {
164 if(ieta < 48-fNCellsFromEMCALBorder) okcol = kTRUE;
165 }
166 }//eta 0 not checked
167 if(fDebug > 1)
168 {
169 printf("AliCalorimeterUtils::CheckCellFiducialRegion(AOD) - EMCAL Cluster in %d cells fiducial volume: ieta %d, iphi %d, SM %d ?",
170 fNCellsFromEMCALBorder, ieta, iphi, iSM);
171 if (okcol && okrow ) printf(" YES \n");
172 else printf(" NO: column ok? %d, row ok? %d \n",okcol,okrow);
173 }
174 }//EMCAL
175 else if(cells->GetType()==AliAODCaloCells::kPHOS){
176 Int_t relId[4];
177 Int_t irow = -1, icol = -1;
178 fPHOSGeo->AbsToRelNumbering(absIdMax,relId);
179 irow = relId[2];
180 icol = relId[3];
181 //imod = relId[0]-1;
182 if(irow >= fNCellsFromPHOSBorder && irow < 64-fNCellsFromPHOSBorder) okrow =kTRUE;
183 if(icol >= fNCellsFromPHOSBorder && icol < 56-fNCellsFromPHOSBorder) okcol =kTRUE;
184 if(fDebug > 1)
185 {
186 printf("AliCalorimeterUtils::CheckCellFiducialRegion(AOD) - PHOS Cluster in %d cells fiducial volume: icol %d, irow %d, Module %d?",
187 fNCellsFromPHOSBorder, icol, irow, relId[0]-1);
188 if (okcol && okrow ) printf(" YES \n");
189 else printf(" NO: column ok? %d, row ok? %d \n",okcol,okrow);
190 }
191 }//PHOS
192
193 if (okcol && okrow) return kTRUE;
194 else return kFALSE;
195
196}
197
198//_______________________________________________________________
199Bool_t AliCalorimeterUtils::CheckCellFiducialRegion(AliESDCaloCluster* cluster, AliESDCaloCells* cells) const {
200 // Given the list of AbsId of the cluster, get the maximum cell and
201 // check if there are fNCellsFromBorder from the calorimeter border
202
203 //If the distance to the border is 0 or negative just exit accept all clusters
204 if(cells->GetType()==AliESDCaloCells::kEMCALCell && fNCellsFromEMCALBorder <= 0 ) return kTRUE;
205 if(cells->GetType()==AliESDCaloCells::kPHOSCell && fNCellsFromPHOSBorder <= 0 ) return kTRUE;
206
207 //Find cell with maximum amplitude
208 Int_t absIdMax = -1;
209 Float_t ampMax = -1;
210 for(Int_t i = 0; i < cluster->GetNCells() ; i++){
211 Int_t absId = cluster->GetCellAbsId(i) ;
212 Float_t amp = cells->GetCellAmplitude(absId);
213 if(amp > ampMax){
214 ampMax = amp;
215 absIdMax = absId;
216 }
217 }
218
219 if(fDebug > 1)
220 printf("AliCalorimeterUtils::CheckCellFiducialRegion(ESD) - Cluster Max AbsId %d, Cell Energy %2.2f, Cluster Energy %2.2f\n",
221 absIdMax, ampMax, cluster->E());
222 if(absIdMax==-1) return kFALSE;
223
224 //Check if the cell is close to the borders:
225 Bool_t okrow = kFALSE;
226 Bool_t okcol = kFALSE;
227
228 if(cells->GetType()==AliESDCaloCells::kEMCALCell){
229
230 Int_t iTower = -1, iIphi = -1, iIeta = -1, iphi = -1, ieta = -1, iSM = -1;
231 fEMCALGeo->GetCellIndex(absIdMax,iSM,iTower,iIphi,iIeta);
232 fEMCALGeo->GetCellPhiEtaIndexInSModule(iSM,iTower,iIphi, iIeta,iphi,ieta);
233
234 //Check rows/phi
235 if(iSM < 10){
236 if(iphi >= fNCellsFromEMCALBorder && iphi < 24-fNCellsFromEMCALBorder) okrow =kTRUE;
237 }
238 else{
239 if(iphi >= fNCellsFromEMCALBorder && iphi < 12-fNCellsFromEMCALBorder) okrow =kTRUE;
240 }
241
242 //Check collumns/eta
243 if(!fNoEMCALBorderAtEta0){
244 if(ieta > fNCellsFromEMCALBorder && ieta < 48-fNCellsFromEMCALBorder) okcol =kTRUE;
245 }
246 else{
247 if(iSM%2==0){
248 if(ieta >= fNCellsFromEMCALBorder) okcol = kTRUE;
249 }
250 else {
251 if(ieta < 48-fNCellsFromEMCALBorder) okcol = kTRUE;
252 }
253 }//eta 0 not checked
254 if(fDebug > 1)
255 {
256 printf("AliCalorimeterUtils::CheckCellFiducialRegion(ESD) - EMCAL Cluster in %d cells fiducial volume: ieta %d, iphi %d, SM %d ?",
257 fNCellsFromEMCALBorder, ieta, iphi, iSM);
258 if (okcol && okrow ) printf(" YES \n");
259 else printf(" NO: column ok? %d, row ok? %d \n",okcol,okrow);
260 }
261 }//EMCAL
262 else if(cells->GetType()==AliESDCaloCells::kPHOSCell){
263 Int_t relId[4];
264 Int_t irow = -1, icol = -1;
265 fPHOSGeo->AbsToRelNumbering(absIdMax,relId);
266 irow = relId[2];
267 icol = relId[3];
268 //imod = relId[0]-1;
269 if(irow >= fNCellsFromPHOSBorder && irow < 64-fNCellsFromPHOSBorder) okrow =kTRUE;
270 if(icol >= fNCellsFromPHOSBorder && icol < 56-fNCellsFromPHOSBorder) okcol =kTRUE;
271 if(fDebug > 1)
272 {
273 printf("AliCalorimeterUtils::CheckCellFiducialRegion(ESD) - PHOS Cluster in %d cells fiducial volume: icol %d, irow %d, Module %d ?",
274 fNCellsFromPHOSBorder, icol, irow,relId[0]-1);
275 if (okcol && okrow ) printf(" YES \n");
276 else printf(" NO: column ok? %d, row ok? %d \n",okcol,okrow);
277 }
278 }//PHOS
279
280 if (okcol && okrow) return kTRUE;
281 else return kFALSE;
282
283}
284
285
286//_________________________________________________________________________________________________________
287Bool_t AliCalorimeterUtils::ClusterContainsBadChannel(TString calorimeter,UShort_t* cellList, Int_t nCells){
288 // Check that in the cluster cells, there is no bad channel of those stored
289 // in fEMCALBadChannelMap or fPHOSBadChannelMap
290
291 if (!fRemoveBadChannels) return kFALSE;
292
78219bac 293 if(calorimeter == "EMCAL" && !fEMCALBadChannelMap) return kFALSE;
294 if(calorimeter == "PHOS" && !fPHOSBadChannelMap) return kFALSE;
765d44e7 295
296 Int_t icol = -1;
297 Int_t irow = -1;
298 Int_t imod = -1;
299 for(Int_t iCell = 0; iCell<nCells; iCell++){
300
301 //Get the column and row
302 if(calorimeter == "EMCAL"){
303 Int_t iTower = -1, iIphi = -1, iIeta = -1;
304 fEMCALGeo->GetCellIndex(cellList[iCell],imod,iTower,iIphi,iIeta);
305 if(fEMCALBadChannelMap->GetEntries() <= imod) continue;
306 fEMCALGeo->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol);
307 if(GetEMCALChannelStatus(imod, icol, irow))return kTRUE;
308 }
309 else if(calorimeter=="PHOS"){
310 Int_t relId[4];
311 fPHOSGeo->AbsToRelNumbering(cellList[iCell],relId);
312 irow = relId[2];
313 icol = relId[3];
314 imod = relId[0]-1;
315 if(fPHOSBadChannelMap->GetEntries() <= imod)continue;
316 if(GetPHOSChannelStatus(imod, icol, irow)) return kTRUE;
317 }
318 else return kFALSE;
319
320 }// cell cluster loop
321
322 return kFALSE;
323
324}
325
326//____________________________________________________________________________________________________________________________________________________
327Int_t AliCalorimeterUtils::GetModuleNumber(AliAODPWG4Particle * particle, AliVEvent * inputEvent) const
328{
329 //Get the EMCAL/PHOS module number that corresponds to this particle
330
331 Int_t absId = -1;
332 if(particle->GetDetector()=="EMCAL"){
333 fEMCALGeo->GetAbsCellIdFromEtaPhi(particle->Eta(),particle->Phi(), absId);
334 if(fDebug > 2)
335 printf("AliCalorimeterUtils::GetModuleNumber(PWG4AOD) - EMCAL: cluster eta %f, phi %f, absid %d, SuperModule %d\n",
336 particle->Eta(), particle->Phi()*TMath::RadToDeg(),absId, fEMCALGeo->GetSuperModuleNumber(absId));
337 return fEMCALGeo->GetSuperModuleNumber(absId) ;
338 }//EMCAL
339 else if(particle->GetDetector()=="PHOS"){
340 if(!strcmp(inputEvent->GetName(),"AliESDEvent")) {
341 AliESDCaloCluster *cluster = ((AliESDEvent*)inputEvent)->GetCaloCluster(particle->GetCaloLabel(0));
342 return GetModuleNumber(cluster);
343 }//ESDs
344 else{
345 AliAODCaloCluster *cluster = ((AliAODEvent*)inputEvent)->GetCaloCluster(particle->GetCaloLabel(0));
346 return GetModuleNumber(cluster);
347 }//AODs
348 }//PHOS
349
350 return -1;
351}
352
353//____________________________________________________________________________________________________________________________________________________
354Int_t AliCalorimeterUtils::GetModuleNumber(AliAODCaloCluster * cluster) const
355{
356 //Get the EMCAL/PHOS module number that corresponds to this cluster, input are AODs
357 TLorentzVector lv;
358 Double_t v[]={0.,0.,0.}; //not necessary to pass the real vertex.
359 cluster->GetMomentum(lv,v);
360 Float_t phi = lv.Phi();
361 if(phi < 0) phi+=TMath::TwoPi();
362 Int_t absId = -1;
363 if(cluster->IsEMCALCluster()){
364 fEMCALGeo->GetAbsCellIdFromEtaPhi(lv.Eta(),phi, absId);
365 if(fDebug > 2)
366 printf("AliCalorimeterUtils::GetModuleNumber(ESD) - EMCAL: cluster eta %f, phi %f, absid %d, SuperModule %d\n",
367 lv.Eta(), phi*TMath::RadToDeg(),absId, fEMCALGeo->GetSuperModuleNumber(absId));
368 return fEMCALGeo->GetSuperModuleNumber(absId) ;
369 }//EMCAL
370 else if(cluster->IsPHOSCluster()) {
371 Int_t relId[4];
372 if ( cluster->GetNCells() > 0) {
373 absId = cluster->GetCellAbsId(0);
374 if(fDebug > 2)
375 printf("AliCalorimeterUtils::GetModuleNumber(AOD) - PHOS: cluster eta %f, phi %f, e %f, absId %d\n",
376 lv.Eta(), phi*TMath::RadToDeg(), lv.E(), absId);
377 }
378 else return -1;
379
380 if ( absId >= 0) {
381 fPHOSGeo->AbsToRelNumbering(absId,relId);
382 if(fDebug > 2)
383 printf("AliCalorimeterUtils::GetModuleNumber(AOD) - PHOS: Module %d\n",relId[0]-1);
384 return relId[0]-1;
385 }
386 else return -1;
387 }//PHOS
388
389 return -1;
390}
391
392//____________________________________________________________________________________________________________________________________________________
393Int_t AliCalorimeterUtils::GetModuleNumber(AliESDCaloCluster * cluster) const
394{
395 //Get the EMCAL/PHOS module number that corresponds to this cluster, input are ESDs
396 TLorentzVector lv;
397 Double_t v[]={0.,0.,0.}; //not necessary to pass the real vertex.
398 cluster->GetMomentum(lv,v);
399 Float_t phi = lv.Phi();
400 if(phi < 0) phi+=TMath::TwoPi();
401 Int_t absId = -1;
402 if(cluster->IsEMCAL()){
403 fEMCALGeo->GetAbsCellIdFromEtaPhi(lv.Eta(),phi, absId);
404 if(fDebug > 2)
405 printf("AliCalorimeterUtils::GetModuleNumber(ESD) - EMCAL: cluster eta %f, phi %f, absid %d, SuperModule %d\n",
406 lv.Eta(), phi*TMath::RadToDeg(),absId, fEMCALGeo->GetSuperModuleNumber(absId));
407 return fEMCALGeo->GetSuperModuleNumber(absId) ;
408 }//EMCAL
409 else if(cluster->IsPHOS()){
410 Int_t relId[4];
411 if ( cluster->GetNCells() > 0) {
412 absId = cluster->GetCellAbsId(0);
413 if(fDebug > 2)
414 printf("AliCalorimeterUtils::GetModuleNumber(ESD) - PHOS: cluster eta %f, phi %f, e %f, absId %d\n",
415 lv.Eta(), phi*TMath::RadToDeg(), lv.E(), absId);
416 }
417 else return -1;
418
419 if ( absId >= 0) {
420 fPHOSGeo->AbsToRelNumbering(absId,relId);
421 if(fDebug > 2)
422 printf("AliCalorimeterUtils::GetModuleNumber(ESD) - PHOS: Module %d\n",relId[0]-1);
423 return relId[0]-1;
424 }
425 else return -1;
426 }//PHOS
427
428 return -1;
429}
430
431
432//_____________________________________________________________________________________________________________
433Int_t AliCalorimeterUtils::GetModuleNumberCellIndexes(const Int_t absId, const TString calo, Int_t & icol, Int_t & irow, Int_t & iRCU) const
434{
435 //Get the EMCAL/PHOS module, columns, row and RCU number that corresponds to this absId
436 Int_t imod = -1;
437 if ( absId >= 0) {
438 if(calo=="EMCAL"){
439 Int_t iTower = -1, iIphi = -1, iIeta = -1;
440 fEMCALGeo->GetCellIndex(absId,imod,iTower,iIphi,iIeta);
441 fEMCALGeo->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol);
442
443 //RCU0
444 if (0<=irow&&irow<8) iRCU=0; // first cable row
445 else if (8<=irow&&irow<16 && 0<=icol&&icol<24) iRCU=0; // first half;
446 //second cable row
447 //RCU1
448 else if(8<=irow&&irow<16 && 24<=icol&&icol<48) iRCU=1; // second half;
449 //second cable row
450 else if(16<=irow&&irow<24) iRCU=1; // third cable row
451
452 if (imod%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same
453 if (iRCU<0) {
454 printf("AliCalorimeterUtils::GetModuleNumberCellIndexes() - Wrong EMCAL RCU number = %d\n", iRCU);
455 abort();
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) {
469 printf("AliCalorimeterUtils::GetModuleNumberCellIndexes() - Wrong PHOS RCU number = %d\n", iRCU);
470 abort();
471 }
472 return imod;
473 }//PHOS
474 }
475
476 return -1;
477}
478
479//_______________________________________________________________
480//void AliCalorimeterUtils::Init()
481//{
482// //Init reader. Method to be called in AliAnaPartCorrMaker
483//
484// fEMCALBadChannelMap->SetName(Form("EMCALBadMap_%s",fTaskName.Data()));
485// fPHOSBadChannelMap->SetName(Form("PHOSBadMap_%s",fTaskName.Data()));
486//}
487
488//_______________________________________________________________
489void AliCalorimeterUtils::InitParameters()
490{
491 //Initialize the parameters of the analysis.
492 fEMCALGeoName = "EMCAL_COMPLETE";
493 fPHOSGeoName = "PHOSgeo";
494
495 if(gGeoManager) {// geoManager was set
496 if(fDebug > 2)printf("AliCalorimeterUtils::InitParameters() - Geometry manager available\n");
497 fEMCALGeoMatrixSet = kTRUE;
498 fPHOSGeoMatrixSet = kTRUE;
499 }
500 else{
501 fEMCALGeoMatrixSet = kFALSE;
502 fPHOSGeoMatrixSet = kFALSE;
503 }
504
505 fRemoveBadChannels = kFALSE;
506
507 fNCellsFromEMCALBorder = 0;
508 fNCellsFromPHOSBorder = 0;
509 fNoEMCALBorderAtEta0 = kFALSE;
510}
511
512//________________________________________________________________
513void AliCalorimeterUtils::InitEMCALBadChannelStatusMap(){
514 //Init EMCAL bad channels map
515 if(fDebug > 0 )printf("AliCalorimeterUtils::InitEMCALBadChannelStatusMap()\n");
516 //In order to avoid rewriting the same histograms
517 Bool_t oldStatus = TH1::AddDirectoryStatus();
518 TH1::AddDirectory(kFALSE);
78219bac 519
520 fEMCALBadChannelMap = new TObjArray(12);
521 //TH2F * hTemp = new TH2I("EMCALBadChannelMap","EMCAL SuperModule bad channel map", 48, 0, 48, 24, 0, 24);
522 for (int i = 0; i < 12; i++) {
523 fEMCALBadChannelMap->Add(new TH2I(Form("EMCALBadChannelMap_Mod%d",i),Form("EMCALBadChannelMap_Mod%d",i), 48, 0, 48, 24, 0, 24));
524 //fEMCALBadChannelMap->Add((TH2I*)hTemp->Clone(Form("EMCALBadChannelMap_Mod%d",i)));
525 }
526
527 //delete hTemp;
765d44e7 528
529 fEMCALBadChannelMap->SetOwner(kTRUE);
530 fEMCALBadChannelMap->Compress();
78219bac 531
765d44e7 532 //In order to avoid rewriting the same histograms
533 TH1::AddDirectory(oldStatus);
534}
535
536//________________________________________________________________
537void AliCalorimeterUtils::InitPHOSBadChannelStatusMap(){
538 //Init PHOS bad channels map
539 if(fDebug > 0 )printf("AliCalorimeterUtils::InitPHOSBadChannelStatusMap()\n");
540 //In order to avoid rewriting the same histograms
541 Bool_t oldStatus = TH1::AddDirectoryStatus();
542 TH1::AddDirectory(kFALSE);
78219bac 543
544 fPHOSBadChannelMap = new TObjArray(5);
765d44e7 545 for (int i = 0; i < 5; i++)fPHOSBadChannelMap->Add(new TH2I(Form("PHOSBadChannelMap_Mod%d",i),Form("PHOSBadChannelMap_Mod%d",i), 56, 0, 56, 64, 0, 64));
546
547 fPHOSBadChannelMap->SetOwner(kTRUE);
548 fPHOSBadChannelMap->Compress();
549
550 //In order to avoid rewriting the same histograms
551 TH1::AddDirectory(oldStatus);
552}
553
554//________________________________________________________________
09e819c9 555void AliCalorimeterUtils::InitEMCALRecalibrationFactors(){
556 //Init EMCAL recalibration factors
557 if(fDebug > 0 )printf("AliCalorimeterUtils::InitEMCALRecalibrationFactors()\n");
558 //In order to avoid rewriting the same histograms
559 Bool_t oldStatus = TH1::AddDirectoryStatus();
560 TH1::AddDirectory(kFALSE);
78219bac 561
562 fEMCALRecalibrationFactors = new TObjArray(12);
09e819c9 563 for (int i = 0; i < 12; i++) fEMCALRecalibrationFactors->Add(new TH2F(Form("EMCALRecalFactors_SM%d",i),Form("EMCALRecalFactors_SM%d",i), 48, 0, 48, 24, 0, 24));
564 //Init the histograms with 1
565 for (Int_t sm = 0; sm < 12; sm++) {
566 for (Int_t i = 0; i < 48; i++) {
567 for (Int_t j = 0; j < 24; j++) {
568 SetEMCALChannelRecalibrationFactor(sm,i,j,1.);
569 }
570 }
571 }
572 fEMCALRecalibrationFactors->SetOwner(kTRUE);
573 fEMCALRecalibrationFactors->Compress();
574
575 //In order to avoid rewriting the same histograms
576 TH1::AddDirectory(oldStatus);
577}
578
579//________________________________________________________________
580void AliCalorimeterUtils::InitPHOSRecalibrationFactors(){
581 //Init EMCAL recalibration factors
582 if(fDebug > 0 )printf("AliCalorimeterUtils::InitPHOSRecalibrationFactors()\n");
583 //In order to avoid rewriting the same histograms
584 Bool_t oldStatus = TH1::AddDirectoryStatus();
585 TH1::AddDirectory(kFALSE);
78219bac 586
587 fPHOSRecalibrationFactors = new TObjArray(5);
09e819c9 588 for (int i = 0; i < 5; i++)fPHOSRecalibrationFactors->Add(new TH2F(Form("PHOSRecalFactors_Mod%d",i),Form("PHOSRecalFactors_Mod%d",i), 56, 0, 56, 64, 0, 64));
589 //Init the histograms with 1
590 for (Int_t m = 0; m < 5; m++) {
591 for (Int_t i = 0; i < 56; i++) {
592 for (Int_t j = 0; j < 64; j++) {
593 SetPHOSChannelRecalibrationFactor(m,i,j,1.);
594 }
595 }
596 }
597 fPHOSRecalibrationFactors->SetOwner(kTRUE);
598 fPHOSRecalibrationFactors->Compress();
599
600 //In order to avoid rewriting the same histograms
601 TH1::AddDirectory(oldStatus);
602}
603
604
605//________________________________________________________________
765d44e7 606void AliCalorimeterUtils::InitEMCALGeometry()
607{
608 //Initialize EMCAL geometry if it did not exist previously
609 if (!fEMCALGeo){
610 fEMCALGeo = new AliEMCALGeoUtils(fEMCALGeoName);
611 if(fDebug > 0){
612 printf("AliCalorimeterUtils::InitEMCALGeometry()");
613 if (!gGeoManager) printf(" - Careful!, gGeoManager not loaded, load misalign matrices");
614 printf("\n");
615 }
616 }
617}
618
619//________________________________________________________________
620void AliCalorimeterUtils::InitPHOSGeometry()
621{
622 //Initialize PHOS geometry if it did not exist previously
623 if (!fPHOSGeo){
624 fPHOSGeo = new AliPHOSGeoUtils(fPHOSGeoName);
625 if(fDebug > 0){
626 printf("AliCalorimeterUtils::InitPHOSGeometry()");
627 if (!gGeoManager) printf(" - Careful!, gGeoManager not loaded, load misalign matrices");
628 printf("\n");
629 }
630 }
631}
632
633//________________________________________________________________
634void AliCalorimeterUtils::Print(const Option_t * opt) const
635{
636
637 //Print some relevant parameters set for the analysis
638 if(! opt)
639 return;
640
641 printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
642 printf("Remove Clusters with bad channels? %d\n",fRemoveBadChannels);
643 printf("Remove Clusters with max cell at less than %d cells from EMCAL border and %d cells from PHOS border\n",
644 fNCellsFromEMCALBorder, fNCellsFromPHOSBorder);
645 if(fNoEMCALBorderAtEta0) printf("Do not remove EMCAL clusters at Eta = 0\n");
09e819c9 646 printf("Recalibrate Clusters? %d\n",fRecalibration);
765d44e7 647
648 printf(" \n") ;
649}
650
651//________________________________________________________________
09e819c9 652Float_t AliCalorimeterUtils::RecalibrateClusterEnergy(AliESDCaloCluster * cluster, AliESDCaloCells * cells){
653 // Recalibrate the cluster energy, considering the recalibration map and the energy of the cells that compose the cluster.
654 // ESD case
655
656 if(!cells) {
657 printf("AliCalorimeterUtils::RecalibrateClusterEnergy(ESD) - Cells pointer does not exist, stop!");
658 abort();
659 }
660 //Get the cluster number of cells and list of absId, check what kind of cluster do we have.
661 UShort_t * index = cluster->GetCellsAbsId() ;
662 Double_t * fraction = cluster->GetCellsAmplitudeFraction() ;
663 Int_t ncells = cluster->GetNCells();
664 TString calo = "EMCAL";
665 if(cluster->IsPHOS()) calo = "PHOS";
666 //Initialize some used variables
667 Float_t energy = 0;
668 Int_t absId = -1;
669 Int_t icol = -1, irow = -1, iRCU = -1, module=1;
670 Float_t factor = 1, frac = 0;
671
672 //Loop on the cells, get the cell amplitude and recalibration factor, multiply and and to the new energy
673 for(Int_t icell = 0; icell < ncells; icell++){
674 absId = index[icell];
675 frac = fraction[icell];
676 if(frac < 1e-3) frac = 1; //in case of EMCAL, this is set as 0, not used.
677 module = GetModuleNumberCellIndexes(absId,calo,icol,irow,iRCU);
678 if(cluster->IsPHOS()) factor = GetPHOSChannelRecalibrationFactor (module,icol,irow);
679 else factor = GetEMCALChannelRecalibrationFactor(module,icol,irow);
680 if(fDebug>2)
681 printf("AliCalorimeterUtils::RecalibrateClusterEnergy(ESD) - recalibrate cell: %s, module %d, col %d, row %d, cell fraction %f, recalibration factor %f, cell energy %f\n",
682 calo.Data(),module,icol,irow,frac,factor,cells->GetCellAmplitude(absId));
683
684 energy += cells->GetCellAmplitude(absId)*factor*frac;
685 }
686
687 if(fDebug>1)
688 printf("AliCalorimeterUtils::RecalibrateClusterEnergy(ESD) - Energy before %f, after %f\n",cluster->E(),energy);
689
690 return energy;
691
692}
693
694//________________________________________________________________
695Float_t AliCalorimeterUtils::RecalibrateClusterEnergy(AliAODCaloCluster * cluster, AliAODCaloCells * cells){
696 // Recalibrate the cluster energy, considering the recalibration map and the energy of the cells that compose the cluster.
697 // AOD case
698
699 if(!cells) {
700 printf("AliCalorimeterUtils::RecalibrateClusterEnergy(AOD) - Cells pointer does not exist, stop!");
701 abort();
702 }
703
704 //Get the cluster number of cells and list of absId, check what kind of cluster do we have.
705 UShort_t * index = cluster->GetCellsAbsId() ;
706 Double_t * fraction = cluster->GetCellsAmplitudeFraction() ;
707 Int_t ncells = cluster->GetNCells();
708 TString calo = "EMCAL";
709 if(cluster->IsPHOSCluster()) calo = "PHOS";
710
711 //Initialize some used variables
712 Float_t energy = 0;
713 Int_t absId = -1;
714 Int_t icol = -1, irow = -1, iRCU = -1, module=1;
715 Float_t factor = 1, frac = 0;
716
717 //Loop on the cells, get the cell amplitude and recalibration factor, multiply and and to the new energy
718 for(Int_t icell = 0; icell < ncells; icell++){
719 absId = index[icell];
720 frac = fraction[icell];
721 if(frac < 1e-3) frac = 1; //in case of EMCAL, this is set as 0, not used.
722 module = GetModuleNumberCellIndexes(absId,calo,icol,irow,iRCU);
723 if(cluster->IsPHOSCluster()) factor = GetPHOSChannelRecalibrationFactor (module,icol,irow);
724 else factor = GetEMCALChannelRecalibrationFactor(module,icol,irow);
725 if(fDebug>2)
726 printf("AliCalorimeterUtils::RecalibrateClusterEnergy(ESD) - recalibrate cell: %s, module %d, col %d, row %d, cell fraction %f,recalibration factor %f, cell energy %f\n",
727 calo.Data(),module,icol,irow,frac,factor,cells->GetCellAmplitude(absId));
728
729 energy += cells->GetCellAmplitude(absId)*factor*frac;
730 }
731
732 if(fDebug>1)
733 printf("AliCalorimeterUtils::RecalibrateClusterEnergy(ESD) - Energy before %f, after %f\n",cluster->E(),energy);
734
735 return energy;
736
737}
738
739
740//________________________________________________________________
765d44e7 741void AliCalorimeterUtils::SetGeometryTransformationMatrices(AliVEvent* inputEvent)
742{
743 //Set the calorimeters transformation matrices
744
745 //Get the EMCAL transformation geometry matrices from ESD
746 if (!gGeoManager && fEMCALGeo) {//&& !fEMCALGeoMatrixSet) { FIXME
747 if(fDebug > 1)
748 printf(" AliCalorimeterUtils::SetGeometryTransformationMatrices() - Load EMCAL misalignment matrices. \n");
749 if(!strcmp(inputEvent->GetName(),"AliESDEvent")) {
750 for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
751 if(((AliESDEvent*)inputEvent)->GetEMCALMatrix(mod)) {
752 //printf("EMCAL: mod %d, matrix %p\n",mod, ((AliESDEvent*)inputEvent)->GetEMCALMatrix(mod));
753 fEMCALGeo->SetMisalMatrix(((AliESDEvent*)inputEvent)->GetEMCALMatrix(mod),mod) ;
754 fEMCALGeoMatrixSet = kTRUE;//At least one, so good
755 }
756 }// loop over super modules
757 }//ESD as input
758 else {
759 if(fDebug > 1)
760 printf("AliCalorimeterUtils::SetGeometryTransformationMatrices() - Setting of EMCAL transformation matrixes for AODs not implemented yet. \n Import geometry.root file\n");
761 }//AOD as input
762 }//EMCAL geo && no geoManager
763
764 //Get the PHOS transformation geometry matrices from ESD
765 if (!gGeoManager && fPHOSGeo && !fPHOSGeoMatrixSet) {
766 if(fDebug > 1)
767 printf(" AliCalorimeterUtils::SetGeometryTransformationMatrices() - Load PHOS misalignment matrices. \n");
768 if(!strcmp(inputEvent->GetName(),"AliESDEvent")) {
769 for(Int_t mod=0; mod < 5; mod++){
770 if(((AliESDEvent*)inputEvent)->GetPHOSMatrix(mod)) {
771 //printf("PHOS: mod %d, matrix %p\n",mod, ((AliESDEvent*)inputEvent)->GetPHOSMatrix(mod));
772 fPHOSGeo->SetMisalMatrix(((AliESDEvent*)inputEvent)->GetPHOSMatrix(mod),mod) ;
773 fPHOSGeoMatrixSet = kTRUE; //At least one so good
774 }
775 }// loop over modules
776 }//ESD as input
777 else {
778 if(fDebug > 1)
779 printf("AliCalorimeterUtils::SetGeometryTransformationMatrices() - Setting of EMCAL transformation matrixes for AODs not implemented yet. \n Import geometry.root file\n");
780 }//AOD as input
781 }//PHOS geo and geoManager was not set
782
783}
784