AliPhysics  2b88e80 (2b88e80)
AliEMCALRecoUtils.cxx
Go to the documentation of this file.
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 
16 // ROOT includes
17 #include <TH2F.h>
18 #include <TArrayI.h>
19 #include <TArrayF.h>
20 #include <TObjArray.h>
21 
22 // STEER includes
23 #include "AliVCluster.h"
24 #include "AliVCaloCells.h"
25 #include "AliLog.h"
26 #include "AliPID.h"
27 #include "AliESDEvent.h"
28 #include "AliAODEvent.h"
29 #include "AliESDtrack.h"
30 #include "AliAODTrack.h"
31 #include "AliExternalTrackParam.h"
32 #include "AliESDfriendTrack.h"
33 #include "AliTrackerBase.h"
34 #include "AliMCEvent.h"
35 
36 // EMCAL includes
37 #include "AliEMCALRecoUtils.h"
38 #include "AliEMCALGeometry.h"
39 #include "AliTrackerBase.h"
40 #include "AliEMCALPIDUtils.h"
41 
43 ClassImp(AliEMCALRecoUtils) ;
45 
51 //_____________________________________
52 AliEMCALRecoUtils::AliEMCALRecoUtils():
53  fParticleType(0), fPosAlgo(0),
54  fW0(0), fShowerShapeCellLocationType(0),
55  fNonLinearityFunction(0), fNonLinearThreshold(0),
56  fSmearClusterEnergy(kFALSE), fRandom(),
57  fCellsRecalibrated(kFALSE), fRecalibration(kFALSE), fEMCALRecalibrationFactors(),
58  fConstantTimeShift(0), fTimeRecalibration(kFALSE), fEMCALTimeRecalibrationFactors(), fLowGain(kFALSE),
59  fUseL1PhaseInTimeRecalibration(kFALSE), fEMCALL1PhaseInTimeRecalibration(),
60  fUseRunCorrectionFactors(kFALSE),
61  fRemoveBadChannels(kFALSE), fRecalDistToBadChannels(kFALSE), fEMCALBadChannelMap(),
62  fNCellsFromEMCALBorder(0), fNoEMCALBorderAtEta0(kTRUE),
63  fRejectExoticCluster(kFALSE), fRejectExoticCells(kFALSE),
64  fExoticCellFraction(0), fExoticCellDiffTime(0), fExoticCellMinAmplitude(0),
65  fPIDUtils(), fAODFilterMask(0),
66  fAODHybridTracks(0), fAODTPCOnlyTracks(0),
67  fMatchedTrackIndex(0x0), fMatchedClusterIndex(0x0),
68  fResidualEta(0x0), fResidualPhi(0x0), fCutEtaPhiSum(kFALSE), fCutEtaPhiSeparate(kFALSE),
69  fCutR(0), fCutEta(0), fCutPhi(0),
70  fClusterWindow(0), fMass(0),
71  fStepSurface(0), fStepCluster(0),
72  fITSTrackSA(kFALSE), fEMCalSurfaceDistance(440.),
73  fTrackCutsType(0), fCutMinTrackPt(0), fCutMinNClusterTPC(0),
74  fCutMinNClusterITS(0), fCutMaxChi2PerClusterTPC(0), fCutMaxChi2PerClusterITS(0),
75  fCutRequireTPCRefit(kFALSE), fCutRequireITSRefit(kFALSE), fCutAcceptKinkDaughters(kFALSE),
76  fCutMaxDCAToVertexXY(0), fCutMaxDCAToVertexZ(0), fCutDCAToVertex2D(kFALSE),
77  fCutRequireITSStandAlone(kFALSE), fCutRequireITSpureSA(kFALSE),
78  fNMCGenerToAccept(0), fMCGenerToAcceptForTrack(1)
79 {
80  // Init parameters
81  InitParameters();
82 
83  for(Int_t j = 0; j < 5; j++) fMCGenerToAccept[j] = "";
84 
85  // Track matching arrays init
86  fMatchedTrackIndex = new TArrayI();
87  fMatchedClusterIndex = new TArrayI();
88  fResidualPhi = new TArrayF();
89  fResidualEta = new TArrayF();
90  fPIDUtils = new AliEMCALPIDUtils();
91 }
92 
93 //
94 // Copy constructor.
95 //
96 //______________________________________________________________________
97 AliEMCALRecoUtils::AliEMCALRecoUtils(const AliEMCALRecoUtils & reco)
98 : AliEMCALRecoUtilsBase(reco),
99  fParticleType(reco.fParticleType), fPosAlgo(reco.fPosAlgo),
100  fW0(reco.fW0), fShowerShapeCellLocationType(reco.fShowerShapeCellLocationType),
101  fNonLinearityFunction(reco.fNonLinearityFunction), fNonLinearThreshold(reco.fNonLinearThreshold),
102  fSmearClusterEnergy(reco.fSmearClusterEnergy), fRandom(),
103  fCellsRecalibrated(reco.fCellsRecalibrated),
104  fRecalibration(reco.fRecalibration), fEMCALRecalibrationFactors(NULL),
105  fConstantTimeShift(reco.fConstantTimeShift),
106  fTimeRecalibration(reco.fTimeRecalibration), fEMCALTimeRecalibrationFactors(NULL),
107  fLowGain(reco.fLowGain),
108  fUseL1PhaseInTimeRecalibration(reco.fUseL1PhaseInTimeRecalibration),
109  fEMCALL1PhaseInTimeRecalibration(reco.fEMCALL1PhaseInTimeRecalibration),
110  fUseRunCorrectionFactors(reco.fUseRunCorrectionFactors),
111  fRemoveBadChannels(reco.fRemoveBadChannels), fRecalDistToBadChannels(reco.fRecalDistToBadChannels),
112  fEMCALBadChannelMap(NULL),
113  fNCellsFromEMCALBorder(reco.fNCellsFromEMCALBorder), fNoEMCALBorderAtEta0(reco.fNoEMCALBorderAtEta0),
114  fRejectExoticCluster(reco.fRejectExoticCluster), fRejectExoticCells(reco.fRejectExoticCells),
115  fExoticCellFraction(reco.fExoticCellFraction), fExoticCellDiffTime(reco.fExoticCellDiffTime),
116  fExoticCellMinAmplitude(reco.fExoticCellMinAmplitude),
117  fPIDUtils(reco.fPIDUtils), fAODFilterMask(reco.fAODFilterMask),
118  fAODHybridTracks(reco.fAODHybridTracks), fAODTPCOnlyTracks(reco.fAODTPCOnlyTracks),
119  fMatchedTrackIndex( reco.fMatchedTrackIndex? new TArrayI(*reco.fMatchedTrackIndex):0x0),
120  fMatchedClusterIndex(reco.fMatchedClusterIndex?new TArrayI(*reco.fMatchedClusterIndex):0x0),
121  fResidualEta( reco.fResidualEta? new TArrayF(*reco.fResidualEta):0x0),
122  fResidualPhi( reco.fResidualPhi? new TArrayF(*reco.fResidualPhi):0x0),
123  fCutEtaPhiSum(reco.fCutEtaPhiSum), fCutEtaPhiSeparate(reco.fCutEtaPhiSeparate),
124  fCutR(reco.fCutR), fCutEta(reco.fCutEta), fCutPhi(reco.fCutPhi),
125  fClusterWindow(reco.fClusterWindow),
126  fMass(reco.fMass), fStepSurface(reco.fStepSurface), fStepCluster(reco.fStepCluster),
127  fITSTrackSA(reco.fITSTrackSA), fEMCalSurfaceDistance(440.),
128  fTrackCutsType(reco.fTrackCutsType), fCutMinTrackPt(reco.fCutMinTrackPt),
129  fCutMinNClusterTPC(reco.fCutMinNClusterTPC), fCutMinNClusterITS(reco.fCutMinNClusterITS),
130  fCutMaxChi2PerClusterTPC(reco.fCutMaxChi2PerClusterTPC), fCutMaxChi2PerClusterITS(reco.fCutMaxChi2PerClusterITS),
131  fCutRequireTPCRefit(reco.fCutRequireTPCRefit), fCutRequireITSRefit(reco.fCutRequireITSRefit),
132  fCutAcceptKinkDaughters(reco.fCutAcceptKinkDaughters), fCutMaxDCAToVertexXY(reco.fCutMaxDCAToVertexXY),
133  fCutMaxDCAToVertexZ(reco.fCutMaxDCAToVertexZ), fCutDCAToVertex2D(reco.fCutDCAToVertex2D),
134  fCutRequireITSStandAlone(reco.fCutRequireITSStandAlone), fCutRequireITSpureSA(reco.fCutRequireITSpureSA),
135  fNMCGenerToAccept(reco.fNMCGenerToAccept), fMCGenerToAcceptForTrack(reco.fMCGenerToAcceptForTrack)
136 {
137  for (Int_t i = 0; i < 15 ; i++) { fMisalRotShift[i] = reco.fMisalRotShift[i] ;
138  fMisalTransShift[i] = reco.fMisalTransShift[i] ; }
139  for (Int_t i = 0; i < 10 ; i++) { fNonLinearityParams[i] = reco.fNonLinearityParams[i] ; }
140  for (Int_t i = 0; i < 3 ; i++) { fSmearClusterParam[i] = reco.fSmearClusterParam[i] ; }
141  for (Int_t j = 0; j < 5 ; j++) { fMCGenerToAccept[j] = reco.fMCGenerToAccept[j] ; }
142 
143  if(reco.fEMCALBadChannelMap) {
144  // Copy constructor - not taking ownership over calibration histograms
145  fEMCALBadChannelMap = new TObjArray(reco.fEMCALBadChannelMap->GetEntries());
146  fEMCALBadChannelMap->SetOwner(false);
147  for(int ism = 0; ism < reco.fEMCALBadChannelMap->GetEntries(); ism++) fEMCALBadChannelMap->AddAt(reco.fEMCALBadChannelMap->At(ism), ism);
148  }
149 
150  if(reco.fEMCALRecalibrationFactors) {
151  // Copy constructor - not taking ownership over calibration histograms
152  fEMCALRecalibrationFactors = new TObjArray(reco.fEMCALRecalibrationFactors->GetEntries());
153  fEMCALRecalibrationFactors->SetOwner(false);
154  for(int ism = 0; ism < reco.fEMCALRecalibrationFactors->GetEntries(); ism++) fEMCALRecalibrationFactors->AddAt(reco.fEMCALRecalibrationFactors->At(ism), ism);
155  }
156 
157  if(reco.fEMCALTimeRecalibrationFactors) {
158  // Copy constructor - not taking ownership over calibration histograms
159  fEMCALTimeRecalibrationFactors = new TObjArray(reco.fEMCALTimeRecalibrationFactors->GetEntries());
160  fEMCALTimeRecalibrationFactors->SetOwner(false);
161  for(int ism = 0; ism < reco.fEMCALTimeRecalibrationFactors->GetEntries(); ism++) fEMCALTimeRecalibrationFactors->AddAt(reco.fEMCALTimeRecalibrationFactors->At(ism), ism);
162  }
163 }
164 
168 //______________________________________________________________________
169 AliEMCALRecoUtils & AliEMCALRecoUtils::operator = (const AliEMCALRecoUtils & reco)
170 {
171  if (this == &reco)return *this;
172  ((TNamed *)this)->operator=(reco);
173 
174  for (Int_t i = 0; i < 15 ; i++) { fMisalTransShift[i] = reco.fMisalTransShift[i] ;
175  fMisalRotShift[i] = reco.fMisalRotShift[i] ; }
176  for (Int_t i = 0; i < 10 ; i++) { fNonLinearityParams[i] = reco.fNonLinearityParams[i] ; }
177  for (Int_t i = 0; i < 3 ; i++) { fSmearClusterParam[i] = reco.fSmearClusterParam[i] ; }
178 
179  fParticleType = reco.fParticleType;
180  fPosAlgo = reco.fPosAlgo;
181  fW0 = reco.fW0;
182  fShowerShapeCellLocationType = reco.fShowerShapeCellLocationType;
183 
184  fNonLinearityFunction = reco.fNonLinearityFunction;
185  fNonLinearThreshold = reco.fNonLinearThreshold;
186  fSmearClusterEnergy = reco.fSmearClusterEnergy;
187 
188  fCellsRecalibrated = reco.fCellsRecalibrated;
189  fRecalibration = reco.fRecalibration;
190 
191  fConstantTimeShift = reco.fConstantTimeShift;
192  fTimeRecalibration = reco.fTimeRecalibration;
193  fLowGain = reco.fLowGain;
194 
195  fUseL1PhaseInTimeRecalibration = reco.fUseL1PhaseInTimeRecalibration;
196  fEMCALL1PhaseInTimeRecalibration = reco.fEMCALL1PhaseInTimeRecalibration;
197 
198  fUseRunCorrectionFactors = reco.fUseRunCorrectionFactors;
199 
200  fRemoveBadChannels = reco.fRemoveBadChannels;
201  fRecalDistToBadChannels = reco.fRecalDistToBadChannels;
202 
203  fNCellsFromEMCALBorder = reco.fNCellsFromEMCALBorder;
204  fNoEMCALBorderAtEta0 = reco.fNoEMCALBorderAtEta0;
205 
206  fRejectExoticCluster = reco.fRejectExoticCluster;
207  fRejectExoticCells = reco.fRejectExoticCells;
208  fExoticCellFraction = reco.fExoticCellFraction;
209  fExoticCellDiffTime = reco.fExoticCellDiffTime;
210  fExoticCellMinAmplitude = reco.fExoticCellMinAmplitude;
211 
212  fPIDUtils = reco.fPIDUtils;
213 
214  fAODFilterMask = reco.fAODFilterMask;
215  fAODHybridTracks = reco.fAODHybridTracks;
216  fAODTPCOnlyTracks = reco.fAODTPCOnlyTracks;
217 
218  fCutEtaPhiSum = reco.fCutEtaPhiSum;
219  fCutEtaPhiSeparate = reco.fCutEtaPhiSeparate;
220  fCutR = reco.fCutR;
221  fCutEta = reco.fCutEta;
222  fCutPhi = reco.fCutPhi;
223  fClusterWindow = reco.fClusterWindow;
224  fMass = reco.fMass;
225  fStepSurface = reco.fStepSurface;
226  fStepCluster = reco.fStepCluster;
227  fITSTrackSA = reco.fITSTrackSA;
228  fEMCalSurfaceDistance = reco.fEMCalSurfaceDistance;
229 
230  fTrackCutsType = reco.fTrackCutsType;
231  fCutMinTrackPt = reco.fCutMinTrackPt;
232  fCutMinNClusterTPC = reco.fCutMinNClusterTPC;
233  fCutMinNClusterITS = reco.fCutMinNClusterITS;
234  fCutMaxChi2PerClusterTPC = reco.fCutMaxChi2PerClusterTPC;
235  fCutMaxChi2PerClusterITS = reco.fCutMaxChi2PerClusterITS;
236  fCutRequireTPCRefit = reco.fCutRequireTPCRefit;
237  fCutRequireITSRefit = reco.fCutRequireITSRefit;
238  fCutAcceptKinkDaughters = reco.fCutAcceptKinkDaughters;
239  fCutMaxDCAToVertexXY = reco.fCutMaxDCAToVertexXY;
240  fCutMaxDCAToVertexZ = reco.fCutMaxDCAToVertexZ;
241  fCutDCAToVertex2D = reco.fCutDCAToVertex2D;
242  fCutRequireITSStandAlone = reco.fCutRequireITSStandAlone;
243  fCutRequireITSpureSA = reco.fCutRequireITSpureSA;
244 
245  fNMCGenerToAccept = reco.fNMCGenerToAccept;
246  fMCGenerToAcceptForTrack = reco.fMCGenerToAcceptForTrack;
247  for (Int_t j = 0; j < 5 ; j++)
248  fMCGenerToAccept[j] = reco.fMCGenerToAccept[j];
249 
250  //
251  // Assign or copy construct the different TArrays
252  //
253  if (reco.fResidualEta)
254  {
255  if (fResidualEta)
256  *fResidualEta = *reco.fResidualEta;
257  else
258  fResidualEta = new TArrayF(*reco.fResidualEta);
259  }
260  else
261  {
262  if (fResidualEta) delete fResidualEta;
263  fResidualEta = 0;
264  }
265 
266  if (reco.fResidualPhi)
267  {
268  if (fResidualPhi)
269  *fResidualPhi = *reco.fResidualPhi;
270  else
271  fResidualPhi = new TArrayF(*reco.fResidualPhi);
272  }
273  else
274  {
275  if (fResidualPhi) delete fResidualPhi;
276  fResidualPhi = 0;
277  }
278 
279  if (reco.fMatchedTrackIndex)
280  {
281  if (fMatchedTrackIndex)
282  *fMatchedTrackIndex = *reco.fMatchedTrackIndex;
283  else
284  fMatchedTrackIndex = new TArrayI(*reco.fMatchedTrackIndex);
285  }
286  else
287  {
288  if (fMatchedTrackIndex) delete fMatchedTrackIndex;
289  fMatchedTrackIndex = 0;
290  }
291 
292  if (reco.fMatchedClusterIndex)
293  {
294  if (fMatchedClusterIndex)
295  *fMatchedClusterIndex = *reco.fMatchedClusterIndex;
296  else
297  fMatchedClusterIndex = new TArrayI(*reco.fMatchedClusterIndex);
298  }
299  else
300  {
301  if (fMatchedClusterIndex) delete fMatchedClusterIndex;
302  fMatchedClusterIndex = 0;
303  }
304 
305  if(fEMCALBadChannelMap) delete fEMCALBadChannelMap;
306  if(reco.fEMCALBadChannelMap) {
307  // Copy constructor - not taking ownership over calibration histograms
308  fEMCALBadChannelMap = new TObjArray(reco.fEMCALBadChannelMap->GetEntries());
309  fEMCALBadChannelMap->SetOwner(false);
310  for(int ism = 0; ism < reco.fEMCALBadChannelMap->GetEntries(); ism++) fEMCALBadChannelMap->AddAt(reco.fEMCALBadChannelMap->At(ism), ism);
311  }
312 
313  if(fEMCALRecalibrationFactors) delete fEMCALRecalibrationFactors;
314  if(reco.fEMCALRecalibrationFactors) {
315  // Copy constructor - not taking ownership over calibration histograms
316  fEMCALRecalibrationFactors = new TObjArray(reco.fEMCALRecalibrationFactors->GetEntries());
317  fEMCALRecalibrationFactors->SetOwner(false);
318  for(int ism = 0; ism < reco.fEMCALRecalibrationFactors->GetEntries(); ism++) fEMCALRecalibrationFactors->AddAt(reco.fEMCALRecalibrationFactors->At(ism), ism);
319  }
320 
321  if(fEMCALTimeRecalibrationFactors) delete fEMCALTimeRecalibrationFactors;
322  if(reco.fEMCALTimeRecalibrationFactors) {
323  // Copy constructor - not taking ownership over calibration histograms
324  fEMCALTimeRecalibrationFactors = new TObjArray(reco.fEMCALTimeRecalibrationFactors->GetEntries());
325  fEMCALTimeRecalibrationFactors->SetOwner(false);
326  for(int ism = 0; ism < reco.fEMCALTimeRecalibrationFactors->GetEntries(); ism++) fEMCALTimeRecalibrationFactors->AddAt(reco.fEMCALTimeRecalibrationFactors->At(ism), ism);
327  }
328 
329  if(this->fEMCALL1PhaseInTimeRecalibration) delete fEMCALL1PhaseInTimeRecalibration;
330  if(reco.fEMCALL1PhaseInTimeRecalibration) {
331  // Copy constructor - not taking ownership over calibration histograms
332  fEMCALL1PhaseInTimeRecalibration = new TObjArray(reco.fEMCALL1PhaseInTimeRecalibration->GetEntries());
333  fEMCALL1PhaseInTimeRecalibration->SetOwner(false);
334  for(int ism = 0; ism < reco.fEMCALL1PhaseInTimeRecalibration->GetEntries(); ism++) fEMCALL1PhaseInTimeRecalibration->AddAt(reco.fEMCALL1PhaseInTimeRecalibration->At(ism), ism);
335  }
336 
337  return *this;
338 }
339 
343 //_____________________________________
344 AliEMCALRecoUtils::~AliEMCALRecoUtils()
345 {
346  if (fEMCALRecalibrationFactors)
347  {
348  delete fEMCALRecalibrationFactors;
349  }
350 
351  if (fEMCALTimeRecalibrationFactors)
352  {
353  delete fEMCALTimeRecalibrationFactors;
354  }
355 
356  if(fEMCALL1PhaseInTimeRecalibration)
357  {
358  delete fEMCALL1PhaseInTimeRecalibration;
359  }
360 
361  if (fEMCALBadChannelMap)
362  {
363  delete fEMCALBadChannelMap;
364  }
365 
366  delete fMatchedTrackIndex ;
367  delete fMatchedClusterIndex ;
368  delete fResidualEta ;
369  delete fResidualPhi ;
370  delete fPIDUtils ;
371 
372  InitTrackCuts();
373 }
374 
387 //_______________________________________________________________________________
388 Bool_t AliEMCALRecoUtils::AcceptCalibrateCell(Int_t absID, Int_t bc,
389  Float_t & amp, Double_t & time,
390  AliVCaloCells* cells)
391 {
392  AliEMCALGeometry* geom = AliEMCALGeometry::GetInstance();
393 
394  if(!geom)
395  {
396  AliError("No instance of the geometry is available");
397  return kFALSE;
398  }
399 
400  if ( absID < 0 || absID >= 24*48*geom->GetNumberOfSuperModules() )
401  return kFALSE;
402 
403  Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
404 
405  if (!geom->GetCellIndex(absID,imod,iTower,iIphi,iIeta))
406  {
407  // cell absID does not exist
408  amp=0; time = 1.e9;
409  return kFALSE;
410  }
411 
412  geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
413 
414  // Do not include bad channels found in analysis,
415  if (IsBadChannelsRemovalSwitchedOn() && GetEMCALChannelStatus(imod, ieta, iphi))
416  return kFALSE;
417 
418  //Recalibrate energy
419  amp = cells->GetCellAmplitude(absID);
420  if (!fCellsRecalibrated && IsRecalibrationOn())
421  amp *= GetEMCALChannelRecalibrationFactor(imod,ieta,iphi);
422 
423  // Recalibrate time
424  time = cells->GetCellTime(absID);
425  time-=fConstantTimeShift*1e-9; // only in case of old Run1 simulation
426  Bool_t isLowGain = !(cells->GetCellHighGain(absID));//HG = false -> LG = true
427 
428  RecalibrateCellTime(absID,bc,time,isLowGain);
429 
430  //Recalibrate time with L1 phase
431  RecalibrateCellTimeL1Phase(imod, bc, time);
432 
433  return kTRUE;
434 }
435 
445 //_____________________________________________________________________________
446 Bool_t AliEMCALRecoUtils::CheckCellFiducialRegion(const AliEMCALGeometry* geom,
447  const AliVCluster* cluster,
448  AliVCaloCells* cells)
449 {
450  if (!cluster)
451  {
452  AliInfo("Cluster pointer null!");
453  return kFALSE;
454  }
455 
456  // If the distance to the border is 0 or negative just exit accept all clusters
457  if (cells->GetType()==AliVCaloCells::kEMCALCell && fNCellsFromEMCALBorder <= 0 )
458  return kTRUE;
459 
460  Int_t absIdMax = -1, iSM =-1, ieta = -1, iphi = -1;
461  Bool_t shared = kFALSE;
462  GetMaxEnergyCell(geom, cells, cluster, absIdMax, iSM, ieta, iphi, shared);
463 
464  AliDebug(2,Form("Cluster Max AbsId %d, Cell Energy %2.2f, Cluster Energy %2.2f, Ncells from border %d, EMCAL eta=0 %d\n",
465  absIdMax, cells->GetCellAmplitude(absIdMax), cluster->E(), fNCellsFromEMCALBorder, fNoEMCALBorderAtEta0));
466 
467  if (absIdMax==-1) return kFALSE;
468 
469  // Check if the cell is close to the borders:
470  Bool_t okrow = kFALSE;
471  Bool_t okcol = kFALSE;
472 
473  if (iSM < 0 || iphi < 0 || ieta < 0 )
474  {
475  AliFatal(Form("Negative value for super module: %d, or cell ieta: %d, or cell iphi: %d, check EMCAL geometry name\n",
476  iSM,ieta,iphi));
477  return kFALSE; // trick coverity
478  }
479 
480  // Check rows/phi
481  Int_t iPhiLast = 24;
482  if ( geom->GetSMType(iSM) == AliEMCALGeometry::kEMCAL_Half ) iPhiLast /= 2;
483  else if ( geom->GetSMType(iSM) == AliEMCALGeometry::kEMCAL_3rd ) iPhiLast /= 3;// 1/3 sm case
484  else if ( geom->GetSMType(iSM) == AliEMCALGeometry::kDCAL_Ext ) iPhiLast /= 3;// 1/3 sm case
485 
486  if(iphi >= fNCellsFromEMCALBorder && iphi < iPhiLast - fNCellsFromEMCALBorder) okrow = kTRUE;
487 
488  // Check columns/eta
489  Int_t iEtaLast = 48;
490  if ( !fNoEMCALBorderAtEta0 || geom->GetSMType(iSM) == AliEMCALGeometry::kDCAL_Standard )
491  {
492  // consider inner border
493  if ( geom->IsDCALSM(iSM) ) iEtaLast = iEtaLast*2/3;
494 
495  if ( ieta > fNCellsFromEMCALBorder && ieta < iEtaLast-fNCellsFromEMCALBorder ) okcol = kTRUE;
496  }
497  else
498  {
499  if (iSM%2==0)
500  {
501  if (ieta >= fNCellsFromEMCALBorder) okcol = kTRUE;
502  }
503  else
504  {
505  if (ieta < iEtaLast-fNCellsFromEMCALBorder) okcol = kTRUE;
506  }
507  }//eta 0 not checked
508 
509  AliDebug(2,Form("EMCAL Cluster in %d cells fiducial volume: ieta %d, iphi %d, SM %d: column? %d, row? %d\nq",
510  fNCellsFromEMCALBorder, ieta, iphi, iSM, okcol, okrow));
511 
512  if (okcol && okrow)
513  {
514  return kTRUE;
515  }
516  else
517  {
518  AliDebug(2,Form("Reject cluster in border, max cell : ieta %d, iphi %d, SM %d\n",ieta, iphi, iSM));
519 
520  return kFALSE;
521  }
522 }
523 
534 //_______________________________________________________________________________
535 Bool_t AliEMCALRecoUtils::ClusterContainsBadChannel(const AliEMCALGeometry* geom,
536  const UShort_t* cellList,
537  Int_t nCells)
538 {
539  if (!fRemoveBadChannels) return kFALSE;
540  if (!fEMCALBadChannelMap) return kFALSE;
541 
542  Int_t icol = -1;
543  Int_t irow = -1;
544  Int_t imod = -1;
545  for (Int_t iCell = 0; iCell<nCells; iCell++)
546  {
547  //Get the column and row
548  Int_t iTower = -1, iIphi = -1, iIeta = -1;
549  geom->GetCellIndex(cellList[iCell],imod,iTower,iIphi,iIeta);
550 
551  if (fEMCALBadChannelMap->GetEntries() <= imod) continue;
552 
553  geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol);
554 
555  if (GetEMCALChannelStatus(imod, icol, irow))
556  {
557  AliDebug(2,Form("Cluster with bad channel: SM %d, col %d, row %d\n",imod, icol, irow));
558  return kTRUE;
559  }
560  }// cell cluster loop
561 
562  return kFALSE;
563 }
564 
576 //___________________________________________________________________________
577 Float_t AliEMCALRecoUtils::GetECross(Int_t absID, Double_t tcell,
578  AliVCaloCells* cells, Int_t bc)
579 {
580  AliEMCALGeometry * geom = AliEMCALGeometry::GetInstance();
581 
582  if(!geom)
583  {
584  AliError("No instance of the geometry is available");
585  return -1;
586  }
587 
588  Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
589  geom->GetCellIndex(absID,imod,iTower,iIphi,iIeta);
590  geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
591 
592  // Get close cells index, energy and time, not in corners
593 
594  Int_t absID1 = -1;
595  Int_t absID2 = -1;
596 
597  if ( iphi < AliEMCALGeoParams::fgkEMCALRows-1) absID1 = geom->GetAbsCellIdFromCellIndexes(imod, iphi+1, ieta);
598  if ( iphi > 0 ) absID2 = geom->GetAbsCellIdFromCellIndexes(imod, iphi-1, ieta);
599 
600  // In case of cell in eta = 0 border, depending on SM shift the cross cell index
601 
602  Int_t absID3 = -1;
603  Int_t absID4 = -1;
604 
605  if ( ieta == AliEMCALGeoParams::fgkEMCALCols-1 && !(imod%2) )
606  {
607  absID3 = geom-> GetAbsCellIdFromCellIndexes(imod+1, iphi, 0);
608  absID4 = geom-> GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
609  }
610  else if ( ieta == 0 && imod%2 )
611  {
612  absID3 = geom-> GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
613  absID4 = geom-> GetAbsCellIdFromCellIndexes(imod-1, iphi, AliEMCALGeoParams::fgkEMCALCols-1);
614  }
615  else
616  {
617  if ( ieta < AliEMCALGeoParams::fgkEMCALCols-1 )
618  absID3 = geom-> GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
619  if ( ieta > 0 )
620  absID4 = geom-> GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
621  }
622 
623  //printf("IMOD %d, AbsId %d, a %d, b %d, c %d e %d \n",imod,absID,absID1,absID2,absID3,absID4);
624 
625  Float_t ecell1 = 0, ecell2 = 0, ecell3 = 0, ecell4 = 0;
626  Double_t tcell1 = 0, tcell2 = 0, tcell3 = 0, tcell4 = 0;
627 
628  AcceptCalibrateCell(absID1,bc, ecell1,tcell1,cells);
629  AcceptCalibrateCell(absID2,bc, ecell2,tcell2,cells);
630  AcceptCalibrateCell(absID3,bc, ecell3,tcell3,cells);
631  AcceptCalibrateCell(absID4,bc, ecell4,tcell4,cells);
632 
633  if (TMath::Abs(tcell-tcell1)*1.e9 > fExoticCellDiffTime) ecell1 = 0 ;
634  if (TMath::Abs(tcell-tcell2)*1.e9 > fExoticCellDiffTime) ecell2 = 0 ;
635  if (TMath::Abs(tcell-tcell3)*1.e9 > fExoticCellDiffTime) ecell3 = 0 ;
636  if (TMath::Abs(tcell-tcell4)*1.e9 > fExoticCellDiffTime) ecell4 = 0 ;
637 
638  return ecell1+ecell2+ecell3+ecell4;
639 }
640 
651 //_____________________________________________________________________________________________
652 Bool_t AliEMCALRecoUtils::IsExoticCell(Int_t absID, AliVCaloCells* cells, Int_t bc)
653 {
654  if (!fRejectExoticCells) return kFALSE;
655 
656  Float_t ecell = 0;
657  Double_t tcell = 0;
658  Bool_t accept = AcceptCalibrateCell(absID, bc, ecell ,tcell ,cells);
659 
660  if (!accept) return kTRUE; // reject this cell
661 
662  if (ecell < fExoticCellMinAmplitude) return kFALSE; // do not reject low energy cells
663 
664  Float_t eCross = GetECross(absID,tcell,cells,bc);
665 
666  if (1-eCross/ecell > fExoticCellFraction)
667  {
668  AliDebug(2,Form("AliEMCALRecoUtils::IsExoticCell() - EXOTIC CELL id %d, eCell %f, eCross %f, 1-eCross/eCell %f\n",
669  absID,ecell,eCross,1-eCross/ecell));
670  return kTRUE;
671  }
672 
673  return kFALSE;
674 }
675 
685 //___________________________________________________________________
686 Bool_t AliEMCALRecoUtils::IsExoticCluster(const AliVCluster *cluster,
687  AliVCaloCells *cells,
688  Int_t bc)
689 {
690  if (!cluster)
691  {
692  AliInfo("Cluster pointer null!");
693  return kFALSE;
694  }
695 
696  if (!fRejectExoticCluster) return kFALSE;
697 
698  // Get highest energy tower
699  AliEMCALGeometry* geom = AliEMCALGeometry::GetInstance();
700 
701  if(!geom)
702  {
703  AliError("No instance of the geometry is available");
704  return kFALSE;
705  }
706 
707  Int_t iSupMod = -1, absId = -1, ieta = -1, iphi = -1;
708  Bool_t shared = kFALSE;
709  GetMaxEnergyCell(geom, cells, cluster, absId, iSupMod, ieta, iphi, shared);
710 
711  return IsExoticCell(absId,cells,bc);
712 }
713 
722 //_______________________________________________________________________
723 Float_t AliEMCALRecoUtils::SmearClusterEnergy(const AliVCluster* cluster)
724 {
725  if (!cluster)
726  {
727  AliInfo("Cluster pointer null!");
728  return 0;
729  }
730 
731  Float_t energy = cluster->E() ;
732  Float_t rdmEnergy = energy ;
733 
734  if (fSmearClusterEnergy)
735  {
736  rdmEnergy = fRandom.Gaus(energy,fSmearClusterParam[0] * TMath::Sqrt(energy) +
737  fSmearClusterParam[1] * energy +
738  fSmearClusterParam[2] );
739  AliDebug(2, Form("Energy: original %f, smeared %f\n", energy, rdmEnergy));
740  }
741 
742  return rdmEnergy;
743 }
744 
753 //____________________________________________________________________________
754 Float_t AliEMCALRecoUtils::CorrectClusterEnergyLinearity(AliVCluster* cluster)
755 {
756  if (!cluster)
757  {
758  AliInfo("Cluster pointer null!");
759  return 0;
760  }
761 
762  Float_t energy = cluster->E();
763  if (energy < 0.100)
764  {
765  // Clusters with less than 100 MeV or negative are not possible
766  AliInfo(Form("Too low cluster energy!, E = %f < 0.100 GeV",energy));
767  return 0;
768  }
769  else if(energy > 300.)
770  {
771  AliInfo(Form("Too high cluster energy!, E = %f GeV, do not apply non linearity",energy));
772  return energy;
773  }
774 
775  switch (fNonLinearityFunction)
776  {
777  case kPi0MC:
778  {
779  //Non-Linearity correction (from MC with function ([0]*exp(-[1]/E))+(([2]/([3]*2.*TMath::Pi())*exp(-(E-[4])^2/(2.*[3]^2)))))
780  //fNonLinearityParams[0] = 1.014;
781  //fNonLinearityParams[1] =-0.03329;
782  //fNonLinearityParams[2] =-0.3853;
783  //fNonLinearityParams[3] = 0.5423;
784  //fNonLinearityParams[4] =-0.4335;
785  energy *= (fNonLinearityParams[0]*exp(-fNonLinearityParams[1]/energy))+
786  ((fNonLinearityParams[2]/(fNonLinearityParams[3]*2.*TMath::Pi())*
787  exp(-(energy-fNonLinearityParams[4])*(energy-fNonLinearityParams[4])/(2.*fNonLinearityParams[3]*fNonLinearityParams[3]))));
788  break;
789  }
790 
791  case kPi0MCv2:
792  {
793  //Non-Linearity correction (from MC with function [0]/((x+[1])^[2]))+1;
794  //fNonLinearityParams[0] = 3.11111e-02;
795  //fNonLinearityParams[1] =-5.71666e-02;
796  //fNonLinearityParams[2] = 5.67995e-01;
797 
798  energy *= fNonLinearityParams[0]/TMath::Power(energy+fNonLinearityParams[1],fNonLinearityParams[2])+1;
799  break;
800  }
801 
802  case kPi0MCv3:
803  {
804  //Same as beam test corrected, change parameters
805  //fNonLinearityParams[0] = 9.81039e-01
806  //fNonLinearityParams[1] = 1.13508e-01;
807  //fNonLinearityParams[2] = 1.00173e+00;
808  //fNonLinearityParams[3] = 9.67998e-02;
809  //fNonLinearityParams[4] = 2.19381e+02;
810  //fNonLinearityParams[5] = 6.31604e+01;
811  //fNonLinearityParams[6] = 1;
812  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
813 
814  break;
815  }
816 
817 
818  case kPi0GammaGamma:
819  {
820  //Non-Linearity correction (from Olga Data with function p0+p1*exp(-p2*E))
821  //fNonLinearityParams[0] = 1.04;
822  //fNonLinearityParams[1] = -0.1445;
823  //fNonLinearityParams[2] = 1.046;
824  energy /= (fNonLinearityParams[0]+fNonLinearityParams[1]*exp(-fNonLinearityParams[2]*energy)); //Olga function
825  break;
826  }
827 
828  case kPi0GammaConversion:
829  {
830  //Non-Linearity correction (Nicolas from Dimitri Data with function C*[1-a*exp(-b*E)])
831  //fNonLinearityParams[0] = 0.139393/0.1349766;
832  //fNonLinearityParams[1] = 0.0566186;
833  //fNonLinearityParams[2] = 0.982133;
834  energy /= fNonLinearityParams[0]*(1-fNonLinearityParams[1]*exp(-fNonLinearityParams[2]*energy));
835 
836  break;
837  }
838 
839  case kBeamTest:
840  {
841  //From beam test, Alexei's results, for different ZS thresholds
842  // th=30 MeV; th = 45 MeV; th = 75 MeV
843  //fNonLinearityParams[0] = 1.007; 1.003; 1.002
844  //fNonLinearityParams[1] = 0.894; 0.719; 0.797
845  //fNonLinearityParams[2] = 0.246; 0.334; 0.358
846  //Rescale the param[0] with 1.03
847  energy /= fNonLinearityParams[0]/(1+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]));
848 
849  break;
850  }
851 
852  case kBeamTestCorrected:
853  {
854  // From beam test, corrected for material between beam and EMCAL
855  //fNonLinearityParams[0] = 0.99078
856  //fNonLinearityParams[1] = 0.161499;
857  //fNonLinearityParams[2] = 0.655166;
858  //fNonLinearityParams[3] = 0.134101;
859  //fNonLinearityParams[4] = 163.282;
860  //fNonLinearityParams[5] = 23.6904;
861  //fNonLinearityParams[6] = 0.978;
862  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
863 
864  break;
865  }
866 
867  case kBeamTestCorrectedv2:
868  {
869  // From beam test, corrected for material between beam and EMCAL
870  // Different function to kBeamTestCorrected
871  //fNonLinearityParams[0] = 0.983504;
872  //fNonLinearityParams[1] = 0.210106;
873  //fNonLinearityParams[2] = 0.897274;
874  //fNonLinearityParams[3] = 0.0829064;
875  //fNonLinearityParams[4] = 152.299;
876  //fNonLinearityParams[5] = 31.5028;
877  //fNonLinearityParams[6] = 0.968;
878  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
879 
880  break;
881  }
882 
883  case kBeamTestCorrectedv3:
884  {
885  // Same function as kBeamTestCorrectedv2, different default parametrization.
886  //fNonLinearityParams[0] = 0.976941;
887  //fNonLinearityParams[1] = 0.162310;
888  //fNonLinearityParams[2] = 1.08689;
889  //fNonLinearityParams[3] = 0.0819592;
890  //fNonLinearityParams[4] = 152.338;
891  //fNonLinearityParams[5] = 30.9594;
892  //fNonLinearityParams[6] = 0.9615;
893  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
894 
895  break;
896  }
897 
898  case kSDMv5:
899  {
900  //Based on fit to the MC/data using kNoCorrection on the data - utilizes symmetric decay method and kPi0MCv5(MC) - 28 Oct 2013
901  //fNonLinearityParams[0] = 1.0;
902  //fNonLinearityParams[1] = 6.64778e-02;
903  //fNonLinearityParams[2] = 1.570;
904  //fNonLinearityParams[3] = 9.67998e-02;
905  //fNonLinearityParams[4] = 2.19381e+02;
906  //fNonLinearityParams[5] = 6.31604e+01;
907  //fNonLinearityParams[6] = 1.01286;
908  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5])))) * (0.964 + exp(-3.132-0.435*energy*2.0));
909 
910  break;
911  }
912 
913  case kPi0MCv5:
914  {
915  //Based on comparing MC truth information to the reconstructed energy of clusters.
916  //fNonLinearityParams[0] = 1.0;
917  //fNonLinearityParams[1] = 6.64778e-02;
918  //fNonLinearityParams[2] = 1.570;
919  //fNonLinearityParams[3] = 9.67998e-02;
920  //fNonLinearityParams[4] = 2.19381e+02;
921  //fNonLinearityParams[5] = 6.31604e+01;
922  //fNonLinearityParams[6] = 1.01286;
923  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
924 
925  break;
926  }
927 
928  case kSDMv6:
929  {
930  //Based on fit to the MC/data using kNoCorrection on the data
931  // - utilizes symmetric decay method and kPi0MCv6(MC) - 09 Dec 2014
932  // - parameters constrained by the test beam data as well
933  // described in the note: https://aliceinfo.cern.ch/Notes/node/211 - Sec 3.1.2 (Test Beam Constrained SDM).
934  //fNonLinearityParams[0] = 1.0;
935  //fNonLinearityParams[1] = 0.237767;
936  //fNonLinearityParams[2] = 0.651203;
937  //fNonLinearityParams[3] = 0.183741;
938  //fNonLinearityParams[4] = 155.427;
939  //fNonLinearityParams[5] = 17.0335;
940  //fNonLinearityParams[6] = 0.987054;
941  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
942 
943  break;
944  }
945 
946  case kPi0MCv6:
947  {
948  //Based on comparing MC truth information to the reconstructed energy of clusters.
949  // described in the note: https://aliceinfo.cern.ch/Notes/node/211 - Sec 3.1.2 (Test Beam Constrained SDM).
950  //fNonLinearityParams[0] = 1.0;
951  //fNonLinearityParams[1] = 0.0797873;
952  //fNonLinearityParams[2] = 1.68322;
953  //fNonLinearityParams[3] = 0.0806098;
954  //fNonLinearityParams[4] = 244.586;
955  //fNonLinearityParams[5] = 116.938;
956  //fNonLinearityParams[6] = 1.00437;
957  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
958 
959  break;
960  }
961 
962  case kPCMv1:
963  {
964  //based on symmetric decays of pi0 meson
965  // described in the note: https://aliceinfo.cern.ch/Notes/node/211 - Sec 3.1.2 (Test Beam Constrained SDM).
966  // parameters vary from MC to MC
967  //fNonLinearityParams[0] = 0.984876;
968  //fNonLinearityParams[1] = -9.999609;
969  //fNonLinearityParams[2] = -4.999891;
970  //fNonLinearityParams[3] = 0.;
971  //fNonLinearityParams[4] = 0.;
972  //fNonLinearityParams[5] = 0.;
973  //fNonLinearityParams[6] = 0.;
974  energy /= TMath::Power(fNonLinearityParams[0] + exp(fNonLinearityParams[1] + fNonLinearityParams[2]*energy),2);//result coming from calo-conv method
975 
976  break;
977  }
978 
979  case kPCMplusBTCv1:
980  {
981  //convolution of TestBeamCorrectedv3 with PCM method
982  //Based on comparing MC truth information to the reconstructed energy of clusters.
983  // described in the note: https://aliceinfo.cern.ch/Notes/node/211 - Sec 3.1.2 (Test Beam Constrained SDM).
984  // parameters vary from MC to MC
985  //fNonLinearityParams[0] = 0.976941;
986  //fNonLinearityParams[1] = 0.162310;
987  //fNonLinearityParams[2] = 1.08689;
988  //fNonLinearityParams[3] = 0.0819592;
989  //fNonLinearityParams[4] = 152.338;
990  //fNonLinearityParams[5] = 30.9594;
991  //fNonLinearityParams[6] = 0.9615;
992  //fNonLinearityParams[7] = 0.984876;
993  //fNonLinearityParams[8] = -9.999609;
994  //fNonLinearityParams[9] = -4.999891;
995  energy *= fNonLinearityParams[6]/(fNonLinearityParams[0]*(1./(1.+fNonLinearityParams[1]*exp(-energy/fNonLinearityParams[2]))*1./(1.+fNonLinearityParams[3]*exp((energy-fNonLinearityParams[4])/fNonLinearityParams[5]))));
996  energy /= TMath::Power(fNonLinearityParams[7] + exp(fNonLinearityParams[8] + fNonLinearityParams[9]*energy),2);//result coming from calo-conv method
997 
998  break;
999  }
1000 
1001  case kPCMsysv1:
1002  {
1003  // Systematic variation of kPCMv1
1004  //Based on comparing MC truth information to the reconstructed energy of clusters.
1005  // described in the note: https://aliceinfo.cern.ch/Notes/node/211 - Sec 3.1.2 (Test Beam Constrained SDM).
1006  // parameters vary from MC to MC
1007  //fNonLinearityParams[0] = 0.0;
1008  //fNonLinearityParams[1] = 1.0;
1009  //fNonLinearityParams[2] = 1.0;
1010  //fNonLinearityParams[3] = 0.0;
1011  //fNonLinearityParams[4] = 1.0;
1012  //fNonLinearityParams[5] = 0.0;
1013  //fNonLinearityParams[6] = 0.0;
1014  energy /= TMath::Power( (fNonLinearityParams[0] + fNonLinearityParams[1] * TMath::Power(energy,fNonLinearityParams[2]) ) /
1015  (fNonLinearityParams[3] + fNonLinearityParams[4] * TMath::Power(energy,fNonLinearityParams[5]) ) + fNonLinearityParams[6], 2);//result coming from calo-conv method
1016 
1017  break;
1018  }
1019 
1020 
1021  case kNoCorrection:
1022  AliDebug(2,"No correction on the energy\n");
1023  break;
1024 
1025  }
1026 
1027  return energy;
1028 }
1029 
1034 //__________________________________________________
1035 void AliEMCALRecoUtils::InitNonLinearityParam()
1036 {
1037  if (fNonLinearityFunction == kPi0MC) {
1038  fNonLinearityParams[0] = 1.014;
1039  fNonLinearityParams[1] = -0.03329;
1040  fNonLinearityParams[2] = -0.3853;
1041  fNonLinearityParams[3] = 0.5423;
1042  fNonLinearityParams[4] = -0.4335;
1043  }
1044 
1045  if (fNonLinearityFunction == kPi0MCv2) {
1046  fNonLinearityParams[0] = 3.11111e-02;
1047  fNonLinearityParams[1] =-5.71666e-02;
1048  fNonLinearityParams[2] = 5.67995e-01;
1049  }
1050 
1051  if (fNonLinearityFunction == kPi0MCv3) {
1052  fNonLinearityParams[0] = 9.81039e-01;
1053  fNonLinearityParams[1] = 1.13508e-01;
1054  fNonLinearityParams[2] = 1.00173e+00;
1055  fNonLinearityParams[3] = 9.67998e-02;
1056  fNonLinearityParams[4] = 2.19381e+02;
1057  fNonLinearityParams[5] = 6.31604e+01;
1058  fNonLinearityParams[6] = 1;
1059  }
1060 
1061  if (fNonLinearityFunction == kPi0GammaGamma) {
1062  fNonLinearityParams[0] = 1.04;
1063  fNonLinearityParams[1] = -0.1445;
1064  fNonLinearityParams[2] = 1.046;
1065  }
1066 
1067  if (fNonLinearityFunction == kPi0GammaConversion) {
1068  fNonLinearityParams[0] = 0.139393;
1069  fNonLinearityParams[1] = 0.0566186;
1070  fNonLinearityParams[2] = 0.982133;
1071  }
1072 
1073  if (fNonLinearityFunction == kBeamTest) {
1074  if (fNonLinearThreshold == 30) {
1075  fNonLinearityParams[0] = 1.007;
1076  fNonLinearityParams[1] = 0.894;
1077  fNonLinearityParams[2] = 0.246;
1078  }
1079  if (fNonLinearThreshold == 45) {
1080  fNonLinearityParams[0] = 1.003;
1081  fNonLinearityParams[1] = 0.719;
1082  fNonLinearityParams[2] = 0.334;
1083  }
1084  if (fNonLinearThreshold == 75) {
1085  fNonLinearityParams[0] = 1.002;
1086  fNonLinearityParams[1] = 0.797;
1087  fNonLinearityParams[2] = 0.358;
1088  }
1089  }
1090 
1091  if (fNonLinearityFunction == kBeamTestCorrected) {
1092  fNonLinearityParams[0] = 0.99078;
1093  fNonLinearityParams[1] = 0.161499;
1094  fNonLinearityParams[2] = 0.655166;
1095  fNonLinearityParams[3] = 0.134101;
1096  fNonLinearityParams[4] = 163.282;
1097  fNonLinearityParams[5] = 23.6904;
1098  fNonLinearityParams[6] = 0.978;
1099  }
1100 
1101  if (fNonLinearityFunction == kBeamTestCorrectedv2) {
1102  // Parameters until November 2015, use now kBeamTestCorrectedv3
1103  fNonLinearityParams[0] = 0.983504;
1104  fNonLinearityParams[1] = 0.210106;
1105  fNonLinearityParams[2] = 0.897274;
1106  fNonLinearityParams[3] = 0.0829064;
1107  fNonLinearityParams[4] = 152.299;
1108  fNonLinearityParams[5] = 31.5028;
1109  fNonLinearityParams[6] = 0.968;
1110  }
1111 
1112  if (fNonLinearityFunction == kBeamTestCorrectedv3) {
1113 
1114  // New parametrization of kBeamTestCorrectedv2
1115  // excluding point at 0.5 GeV from Beam Test Data
1116  // https://indico.cern.ch/event/438805/contribution/1/attachments/1145354/1641875/emcalPi027August2015.pdf
1117 
1118  fNonLinearityParams[0] = 0.976941;
1119  fNonLinearityParams[1] = 0.162310;
1120  fNonLinearityParams[2] = 1.08689;
1121  fNonLinearityParams[3] = 0.0819592;
1122  fNonLinearityParams[4] = 152.338;
1123  fNonLinearityParams[5] = 30.9594;
1124  fNonLinearityParams[6] = 0.9615;
1125  }
1126 
1127  if (fNonLinearityFunction == kSDMv5) {
1128  fNonLinearityParams[0] = 1.0;
1129  fNonLinearityParams[1] = 6.64778e-02;
1130  fNonLinearityParams[2] = 1.570;
1131  fNonLinearityParams[3] = 9.67998e-02;
1132  fNonLinearityParams[4] = 2.19381e+02;
1133  fNonLinearityParams[5] = 6.31604e+01;
1134  fNonLinearityParams[6] = 1.01286;
1135  }
1136 
1137  if (fNonLinearityFunction == kPi0MCv5) {
1138  fNonLinearityParams[0] = 1.0;
1139  fNonLinearityParams[1] = 6.64778e-02;
1140  fNonLinearityParams[2] = 1.570;
1141  fNonLinearityParams[3] = 9.67998e-02;
1142  fNonLinearityParams[4] = 2.19381e+02;
1143  fNonLinearityParams[5] = 6.31604e+01;
1144  fNonLinearityParams[6] = 1.01286;
1145  }
1146 
1147  if (fNonLinearityFunction == kSDMv6) {
1148  fNonLinearityParams[0] = 1.0;
1149  fNonLinearityParams[1] = 0.237767;
1150  fNonLinearityParams[2] = 0.651203;
1151  fNonLinearityParams[3] = 0.183741;
1152  fNonLinearityParams[4] = 155.427;
1153  fNonLinearityParams[5] = 17.0335;
1154  fNonLinearityParams[6] = 0.987054;
1155  }
1156 
1157  if (fNonLinearityFunction == kPi0MCv6) {
1158  fNonLinearityParams[0] = 1.0;
1159  fNonLinearityParams[1] = 0.0797873;
1160  fNonLinearityParams[2] = 1.68322;
1161  fNonLinearityParams[3] = 0.0806098;
1162  fNonLinearityParams[4] = 244.586;
1163  fNonLinearityParams[5] = 116.938;
1164  fNonLinearityParams[6] = 1.00437;
1165  }
1166 
1167 if (fNonLinearityFunction == kPCMv1) {
1168  //parameters change from MC production to MC production, they need to set for each period
1169  fNonLinearityParams[0] = 0.984876;
1170  fNonLinearityParams[1] = -9.999609;
1171  fNonLinearityParams[2] = -4.999891;
1172  fNonLinearityParams[3] = 0.;
1173  fNonLinearityParams[4] = 0.;
1174  fNonLinearityParams[5] = 0.;
1175  fNonLinearityParams[6] = 0.;
1176  }
1177 
1178  if (fNonLinearityFunction == kPCMplusBTCv1) {
1179  // test beam corrected values convoluted with symmetric meson decays values
1180  // for test beam:
1181  // https://indico.cern.ch/event/438805/contribution/1/attachments/1145354/1641875/emcalPi027August2015.pdf
1182  // for PCM method:
1183  // https://aliceinfo.cern.ch/Notes/node/211
1184  fNonLinearityParams[0] = 0.976941;
1185  fNonLinearityParams[1] = 0.162310;
1186  fNonLinearityParams[2] = 1.08689;
1187  fNonLinearityParams[3] = 0.0819592;
1188  fNonLinearityParams[4] = 152.338;
1189  fNonLinearityParams[5] = 30.9594;
1190  fNonLinearityParams[6] = 0.9615;
1191  fNonLinearityParams[7] = 0.984876;
1192  fNonLinearityParams[8] = -9.999609;
1193  fNonLinearityParams[9] = -4.999891;
1194  }
1195 
1196  if (fNonLinearityFunction == kPCMsysv1) {
1197  //systematics for kPCMv1
1198  // for PCM method:
1199  // https://aliceinfo.cern.ch/Notes/node/211
1200  fNonLinearityParams[0] = 0.0;
1201  fNonLinearityParams[1] = 1.0;
1202  fNonLinearityParams[2] = 1.0;
1203  fNonLinearityParams[3] = 0.0;
1204  fNonLinearityParams[4] = 1.0;
1205  fNonLinearityParams[5] = 0.0;
1206  fNonLinearityParams[6] = 0.0;
1207  }
1208 }
1209 
1210 
1224 //____________________________________________________________________
1225 void AliEMCALRecoUtils::GetMaxEnergyCell(const AliEMCALGeometry *geom,
1226  AliVCaloCells* cells,
1227  const AliVCluster* clu,
1228  Int_t & absId,
1229  Int_t & iSupMod,
1230  Int_t & ieta,
1231  Int_t & iphi,
1232  Bool_t & shared)
1233 {
1234  Double_t eMax = -1.;
1235  Double_t eCell = -1.;
1236  Float_t fraction = 1.;
1237  Float_t recalFactor = 1.;
1238  Int_t cellAbsId = -1 ;
1239 
1240  Int_t iTower = -1;
1241  Int_t iIphi = -1;
1242  Int_t iIeta = -1;
1243  Int_t iSupMod0= -1;
1244 
1245  if (!clu)
1246  {
1247  AliInfo("Cluster pointer null!");
1248  absId=-1; iSupMod0=-1, ieta = -1; iphi = -1; shared = -1;
1249  return;
1250  }
1251 
1252  for (Int_t iDig=0; iDig< clu->GetNCells(); iDig++)
1253  {
1254  cellAbsId = clu->GetCellAbsId(iDig);
1255  fraction = clu->GetCellAmplitudeFraction(iDig);
1256  //printf("a Cell %d, id, %d, amp %f, fraction %f\n",iDig,cellAbsId,cells->GetCellAmplitude(cellAbsId),fraction);
1257  if (fraction < 1e-4) fraction = 1.; // in case unfolding is off
1258 
1259  geom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
1260  geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
1261 
1262  if (iDig==0)
1263  {
1264  iSupMod0=iSupMod;
1265  } else if (iSupMod0!=iSupMod)
1266  {
1267  shared = kTRUE;
1268  //printf("AliEMCALRecoUtils::GetMaxEnergyCell() - SHARED CLUSTER\n");
1269  }
1270 
1271  if (!fCellsRecalibrated && IsRecalibrationOn())
1272  recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
1273 
1274  eCell = cells->GetCellAmplitude(cellAbsId)*fraction*recalFactor;
1275  //printf("b Cell %d, id, %d, amp %f, fraction %f\n",iDig,cellAbsId,eCell,fraction);
1276  if (eCell > eMax)
1277  {
1278  eMax = eCell;
1279  absId = cellAbsId;
1280  //printf("\t new max: cell %d, e %f, ecell %f\n",maxId, eMax,eCell);
1281  }
1282  }// cell loop
1283 
1284  //Get from the absid the supermodule, tower and eta/phi numbers
1285  geom->GetCellIndex(absId,iSupMod,iTower,iIphi,iIeta);
1286  //Gives SuperModule and Tower numbers
1287  geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
1288  //printf("Max id %d, iSM %d, col %d, row %d\n",absId,iSupMod,ieta,iphi);
1289  //printf("Max end---\n");
1290 }
1291 
1300 //_________________________________________________________
1301 Float_t AliEMCALRecoUtils::GetCellWeight(Float_t eCell, Float_t eCluster) const
1302 {
1303  if (eCell > 0 && eCluster > 0)
1304  {
1305  if ( fW0 > 0 ) return TMath::Max( 0., fW0 + TMath::Log( eCell / eCluster ) ) ;
1306  else return TMath::Log( eCluster / eCell ) ;
1307  }
1308  else
1309  return 0. ;
1310 }
1311 
1315 //______________________________________
1316 void AliEMCALRecoUtils::InitParameters()
1317 {
1318  fParticleType = kPhoton;
1319  fPosAlgo = kUnchanged;
1320  fW0 = 4.5;
1321 
1322  fNonLinearityFunction = kNoCorrection;
1323  fNonLinearThreshold = 30;
1324 
1325  fExoticCellFraction = 0.97;
1326  fExoticCellDiffTime = 1e6;
1327  fExoticCellMinAmplitude = 4.0;
1328 
1329  fAODFilterMask = 128;
1330  fAODHybridTracks = kFALSE;
1331  fAODTPCOnlyTracks = kTRUE;
1332 
1333  fCutEtaPhiSum = kTRUE;
1334  fCutEtaPhiSeparate = kFALSE;
1335 
1336  fCutR = 0.05;
1337  fCutEta = 0.025;
1338  fCutPhi = 0.05;
1339 
1340  fClusterWindow = 100;
1341  fMass = 0.139;
1342 
1343  fStepSurface = 20.;
1344  fStepCluster = 5.;
1345  fTrackCutsType = kLooseCut;
1346 
1347  fCutMinTrackPt = 0;
1348  fCutMinNClusterTPC = -1;
1349  fCutMinNClusterITS = -1;
1350 
1351  fCutMaxChi2PerClusterTPC = 1e10;
1352  fCutMaxChi2PerClusterITS = 1e10;
1353 
1354  fCutRequireTPCRefit = kFALSE;
1355  fCutRequireITSRefit = kFALSE;
1356  fCutAcceptKinkDaughters = kFALSE;
1357 
1358  fCutMaxDCAToVertexXY = 1e10;
1359  fCutMaxDCAToVertexZ = 1e10;
1360  fCutDCAToVertex2D = kFALSE;
1361 
1362  fCutRequireITSStandAlone = kFALSE; //MARCEL
1363  fCutRequireITSpureSA = kFALSE; //Marcel
1364 
1365  // Misalignment matrices
1366  for (Int_t i = 0; i < 15 ; i++)
1367  {
1368  fMisalTransShift[i] = 0.;
1369  fMisalRotShift[i] = 0.;
1370  }
1371 
1372  // Non linearity
1373  for (Int_t i = 0; i < 10 ; i++) fNonLinearityParams[i] = 0.;
1374 
1375  // For kBeamTestCorrectedv2 case, but default is no correction
1376  fNonLinearityParams[0] = 0.983504;
1377  fNonLinearityParams[1] = 0.210106;
1378  fNonLinearityParams[2] = 0.897274;
1379  fNonLinearityParams[3] = 0.0829064;
1380  fNonLinearityParams[4] = 152.299;
1381  fNonLinearityParams[5] = 31.5028;
1382  fNonLinearityParams[6] = 0.968;
1383 
1384  // Cluster energy smearing
1385  fSmearClusterEnergy = kFALSE;
1386  fSmearClusterParam[0] = 0.07; // * sqrt E term
1387  fSmearClusterParam[1] = 0.00; // * E term
1388  fSmearClusterParam[2] = 0.00; // constant
1389 }
1390 
1394 //_____________________________________________________
1395 void AliEMCALRecoUtils::InitEMCALRecalibrationFactors()
1396 {
1397  AliDebug(2,"AliCalorimeterUtils::InitEMCALRecalibrationFactors()");
1398 
1399  // In order to avoid rewriting the same histograms
1400  Bool_t oldStatus = TH1::AddDirectoryStatus();
1401  TH1::AddDirectory(kFALSE);
1402 
1403  fEMCALRecalibrationFactors = new TObjArray(22);
1404  for (int i = 0; i < 22; i++)
1405  fEMCALRecalibrationFactors->Add(new TH2F(Form("EMCALRecalFactors_SM%d",i),
1406  Form("EMCALRecalFactors_SM%d",i), 48, 0, 48, 24, 0, 24));
1407  //Init the histograms with 1
1408  for (Int_t sm = 0; sm < 22; sm++)
1409  {
1410  for (Int_t i = 0; i < 48; i++)
1411  {
1412  for (Int_t j = 0; j < 24; j++)
1413  {
1414  SetEMCALChannelRecalibrationFactor(sm,i,j,1.);
1415  }
1416  }
1417  }
1418 
1419  fEMCALRecalibrationFactors->SetOwner(kTRUE);
1420  fEMCALRecalibrationFactors->Compress();
1421 
1422  // In order to avoid rewriting the same histograms
1423  TH1::AddDirectory(oldStatus);
1424 }
1425 
1429 //_________________________________________________________
1430 void AliEMCALRecoUtils::InitEMCALTimeRecalibrationFactors()
1431 {
1432  AliDebug(2,"AliCalorimeterUtils::InitEMCALRecalibrationFactors()");
1433 
1434  // In order to avoid rewriting the same histograms
1435  Bool_t oldStatus = TH1::AddDirectoryStatus();
1436  TH1::AddDirectory(kFALSE);
1437 
1438  if(fLowGain) fEMCALTimeRecalibrationFactors = new TObjArray(8);
1439  else fEMCALTimeRecalibrationFactors = new TObjArray(4);
1440 
1441  for (int i = 0; i < 4; i++)
1442  fEMCALTimeRecalibrationFactors->Add(new TH1F(Form("hAllTimeAvBC%d",i),
1443  Form("hAllTimeAvBC%d",i),
1444  48*24*22,0.,48*24*22) );
1445  // Init the histograms with 0
1446  for (Int_t iBC = 0; iBC < 4; iBC++)
1447  {
1448  for (Int_t iCh = 0; iCh < 48*24*22; iCh++)
1449  SetEMCALChannelTimeRecalibrationFactor(iBC,iCh,0.,kFALSE);
1450  }
1451 
1452  if(fLowGain) {
1453  for (int iBC = 0; iBC < 4; iBC++) {
1454  fEMCALTimeRecalibrationFactors->Add(new TH1F(Form("hAllTimeAvLGBC%d",iBC),
1455  Form("hAllTimeAvLGBC%d",iBC),
1456  48*24*22,0.,48*24*22) );
1457  for (Int_t iCh = 0; iCh < 48*24*22; iCh++)
1458  SetEMCALChannelTimeRecalibrationFactor(iBC,iCh,0.,kTRUE);
1459  }
1460  }
1461 
1462 
1463 
1464  fEMCALTimeRecalibrationFactors->SetOwner(kTRUE);
1465  fEMCALTimeRecalibrationFactors->Compress();
1466 
1467  // In order to avoid rewriting the same histograms
1468  TH1::AddDirectory(oldStatus);
1469 }
1470 
1474 //____________________________________________________
1475 void AliEMCALRecoUtils::InitEMCALBadChannelStatusMap()
1476 {
1477  AliDebug(2,"AliEMCALRecoUtils::InitEMCALBadChannelStatusMap()");
1478 
1479  // In order to avoid rewriting the same histograms
1480  Bool_t oldStatus = TH1::AddDirectoryStatus();
1481  TH1::AddDirectory(kFALSE);
1482 
1483  fEMCALBadChannelMap = new TObjArray(22);
1484  //TH2F * hTemp = new TH2I("EMCALBadChannelMap","EMCAL SuperModule bad channel map", 48, 0, 48, 24, 0, 24);
1485 
1486  for (int i = 0; i < 22; i++)
1487  fEMCALBadChannelMap->Add(new TH2I(Form("EMCALBadChannelMap_Mod%d",i),Form("EMCALBadChannelMap_Mod%d",i), 48, 0, 48, 24, 0, 24));
1488 
1489  fEMCALBadChannelMap->SetOwner(kTRUE);
1490  fEMCALBadChannelMap->Compress();
1491 
1492  // In order to avoid rewriting the same histograms
1493  TH1::AddDirectory(oldStatus);
1494 }
1495 
1499 //___________________________________________________________
1500 void AliEMCALRecoUtils::InitEMCALL1PhaseInTimeRecalibration()
1501 {
1502  AliDebug(2,"AliEMCALRecoUtils::InitEMCALL1PhaseInTimeRecalibrationFactors()");
1503 
1504  // In order to avoid rewriting the same histograms
1505  Bool_t oldStatus = TH1::AddDirectoryStatus();
1506  TH1::AddDirectory(kFALSE);
1507 
1508  fEMCALL1PhaseInTimeRecalibration = new TObjArray(1);
1509 
1510  fEMCALL1PhaseInTimeRecalibration->Add(new TH1C("h0","EMCALL1phaseForSM", 22, 0, 22));
1511 
1512  for (Int_t i = 0; i < 22; i++) //loop over SMs, default value = 0
1513  SetEMCALL1PhaseInTimeRecalibrationForSM(i,0);
1514 
1515  fEMCALL1PhaseInTimeRecalibration->SetOwner(kTRUE);
1516  fEMCALL1PhaseInTimeRecalibration->Compress();
1517 
1518  // In order to avoid rewriting the same histograms
1519  TH1::AddDirectory(oldStatus);
1520 }
1521 
1531 //____________________________________________________________________________
1532 void AliEMCALRecoUtils::RecalibrateClusterEnergy(const AliEMCALGeometry* geom,
1533  AliVCluster * cluster,
1534  AliVCaloCells * cells,
1535  Int_t bc)
1536 {
1537  if (!cluster)
1538  {
1539  AliInfo("Cluster pointer null!");
1540  return;
1541  }
1542 
1543  // Get the cluster number of cells and list of absId, check what kind of cluster do we have.
1544  UShort_t * index = cluster->GetCellsAbsId() ;
1545  Double_t * fraction = cluster->GetCellsAmplitudeFraction() ;
1546  Int_t ncells = cluster->GetNCells();
1547 
1548  // Initialize some used variables
1549  Float_t energy = 0;
1550  Int_t absId =-1;
1551  Int_t icol =-1, irow =-1, imod=1;
1552  Float_t factor = 1, frac = 0;
1553  Int_t absIdMax = -1;
1554  Float_t emax = 0;
1555 
1556  // Loop on the cells, get the cell amplitude and recalibration factor, multiply and and to the new energy
1557  for (Int_t icell = 0; icell < ncells; icell++)
1558  {
1559  absId = index[icell];
1560  frac = fraction[icell];
1561  if (frac < 1e-5) frac = 1; //in case of EMCAL, this is set as 0 since unfolding is off
1562 
1563  if (!fCellsRecalibrated && IsRecalibrationOn())
1564  {
1565  // Energy
1566  Int_t iTower = -1, iIphi = -1, iIeta = -1;
1567  geom->GetCellIndex(absId,imod,iTower,iIphi,iIeta);
1568  if (fEMCALRecalibrationFactors->GetEntries() <= imod)
1569  continue;
1570  geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol);
1571  factor = GetEMCALChannelRecalibrationFactor(imod,icol,irow);
1572 
1573  AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - recalibrate cell: module %d, col %d, row %d, cell fraction %f,recalibration factor %f, cell energy %f\n",
1574  imod,icol,irow,frac,factor,cells->GetCellAmplitude(absId)));
1575 
1576  }
1577 
1578  energy += cells->GetCellAmplitude(absId)*factor*frac;
1579 
1580  if (emax < cells->GetCellAmplitude(absId)*factor*frac)
1581  {
1582  emax = cells->GetCellAmplitude(absId)*factor*frac;
1583  absIdMax = absId;
1584  }
1585  }
1586 
1587  AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - Energy before %f, after %f \n",cluster->E(),energy));
1588 
1589  cluster->SetE(energy);
1590 
1591  // Recalculate time of cluster
1592  Double_t timeorg = cluster->GetTOF();
1593  Bool_t isLowGain = !(cells->GetCellHighGain(absIdMax));//HG = false -> LG = true
1594 
1595  Double_t time = cells->GetCellTime(absIdMax);
1596  if (!fCellsRecalibrated && IsTimeRecalibrationOn())
1597  RecalibrateCellTime(absIdMax,bc,time,isLowGain);
1598  time-=fConstantTimeShift*1e-9; // only in case of Run1 old simulations
1599 
1600  // Recalibrate time with L1 phase
1601  if (!fCellsRecalibrated && IsL1PhaseInTimeRecalibrationOn())
1602  RecalibrateCellTimeL1Phase(imod, bc, time);
1603 
1604  cluster->SetTOF(time);
1605 
1606  AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - Time before %f, after %f \n",timeorg,cluster->GetTOF()));
1607 }
1608 
1616 //_______________________________________________________________________
1617 void AliEMCALRecoUtils::RecalibrateCells(AliVCaloCells * cells, Int_t bc)
1618 {
1619  if (!IsRecalibrationOn() && !IsTimeRecalibrationOn() && !IsBadChannelsRemovalSwitchedOn())
1620  return;
1621 
1622  if (!cells)
1623  {
1624  AliInfo("Cells pointer null!");
1625  return;
1626  }
1627 
1628  Short_t absId =-1;
1629  Bool_t accept = kFALSE;
1630  Float_t ecell = 0;
1631  Double_t tcell = 0;
1632  Double_t ecellin = 0;
1633  Double_t tcellin = 0;
1634  Int_t mclabel = -1;
1635  Double_t efrac = 0;
1636 
1637  Int_t nEMcell = cells->GetNumberOfCells() ;
1638  for (Int_t iCell = 0; iCell < nEMcell; iCell++)
1639  {
1640  cells->GetCell( iCell, absId, ecellin, tcellin, mclabel, efrac );
1641 
1642  accept = AcceptCalibrateCell(absId, bc, ecell ,tcell ,cells);
1643  if (!accept)
1644  {
1645  ecell = 0;
1646  tcell = -1;
1647  }
1648 
1649  // Set new values
1650  cells->SetCell(iCell,absId,ecell, tcell, mclabel, efrac);
1651  }
1652 
1653  fCellsRecalibrated = kTRUE;
1654 }
1655 
1664 //_______________________________________________________________________________________________________
1665 void AliEMCALRecoUtils::RecalibrateCellTime(Int_t absId, Int_t bc, Double_t & celltime, Bool_t isLGon) const
1666 {
1667  if (!fCellsRecalibrated && IsTimeRecalibrationOn() && bc >= 0) {
1668  if(fLowGain)
1669  celltime -= GetEMCALChannelTimeRecalibrationFactor(bc%4,absId,isLGon)*1.e-9;
1670  else
1671  celltime -= GetEMCALChannelTimeRecalibrationFactor(bc%4,absId,kFALSE)*1.e-9;
1672  }
1673 }
1674 
1682 //_______________________________________________________________________________________________________
1683 void AliEMCALRecoUtils::RecalibrateCellTimeL1Phase(Int_t iSM, Int_t bc, Double_t & celltime) const
1684 {
1685  if (!fCellsRecalibrated && IsL1PhaseInTimeRecalibrationOn() && bc >= 0)
1686  {
1687  bc=bc%4;
1688 
1689  Float_t offsetPerSM=0.;
1690  Int_t l1PhaseShift = GetEMCALL1PhaseInTimeRecalibrationForSM(iSM);
1691  Int_t l1Phase=l1PhaseShift & 3; //bit operation
1692 
1693  if(bc >= l1Phase)
1694  offsetPerSM = (bc - l1Phase)*25;
1695  else
1696  offsetPerSM = (bc - l1Phase + 4)*25;
1697 
1698  Int_t l1shiftOffset=l1PhaseShift>>2; //bit operation
1699  l1shiftOffset*=25;
1700 
1701  celltime -= offsetPerSM*1.e-9;
1702  celltime -= l1shiftOffset*1.e-9;
1703  }
1704 }
1705 
1721 //______________________________________________________________________________
1722 void AliEMCALRecoUtils::RecalculateCellLabelsRemoveAddedGenerator( Int_t absID, AliVCluster* clus, AliMCEvent* mc,
1723  Float_t & amp, TArrayI & labeArr, TArrayF & eDepArr) const
1724 {
1725  TString genName ;
1726  Float_t eDepFrac[4];
1727 
1728  Float_t edepTotFrac = 1;
1729  Bool_t found = kFALSE;
1730  Float_t ampOrg = amp;
1731 
1732  //
1733  // Get the energy deposition fraction from cluster.
1734  //
1735  for(Int_t icluscell = 0; icluscell < clus->GetNCells(); icluscell++ )
1736  {
1737  if(absID == clus->GetCellAbsId(icluscell))
1738  {
1739  clus->GetCellMCEdepFractionArray(icluscell,eDepFrac);
1740 
1741  found = kTRUE;
1742 
1743  break;
1744  }
1745  }
1746 
1747  if ( !found )
1748  {
1749  AliWarning(Form("Cell abs ID %d NOT found in cluster",absID));
1750  return;
1751  }
1752 
1753  //
1754  // Check if there is a particle contribution from a given generator name.
1755  // If it is not one of the selected generators,
1756  // remove the constribution from the cell.
1757  //
1758  if ( mc && fNMCGenerToAccept > 0 )
1759  {
1760  //printf("Accept contribution from generator? \n");
1761  for(UInt_t imc = 0; imc < 4; imc++)
1762  {
1763  if ( eDepFrac[imc] > 0 && clus->GetNLabels() > imc )
1764  {
1765  mc->GetCocktailGenerator(clus->GetLabelAt(imc),genName);
1766 
1767  Bool_t generOK = kFALSE;
1768  for(Int_t ig = 0; ig < fNMCGenerToAccept; ig++)
1769  {
1770  if ( genName.Contains(fMCGenerToAccept[ig]) ) generOK = kTRUE;
1771  }
1772 
1773  if ( !generOK )
1774  {
1775  amp-=ampOrg*eDepFrac[imc];
1776 
1777  edepTotFrac-=eDepFrac[imc];
1778 
1779  eDepFrac[imc] = 0;
1780  }
1781 
1782  } // eDep > 0
1783  } // 4 possible loop
1784 
1785  } // accept at least one kind of generator
1786 
1787  //
1788  // Add MC label and Edep to corresponding array (to be used later in digits)
1789  //
1790  Int_t nLabels = 0;
1791  for(UInt_t imc = 0; imc < 4; imc++)
1792  {
1793  if ( eDepFrac[imc] > 0 && clus->GetNLabels() > imc && edepTotFrac > 0 )
1794  {
1795  nLabels++;
1796 
1797  labeArr.Set(nLabels);
1798  labeArr.AddAt(clus->GetLabelAt(imc), nLabels-1);
1799 
1800  eDepArr.Set(nLabels);
1801  eDepArr.AddAt( (eDepFrac[imc]/edepTotFrac) * amp, nLabels-1);
1802  // use as deposited energy a fraction of the simulated energy (smeared and with noise)
1803  } // edep > 0
1804 
1805  } // mc cell label loop
1806 
1807  //
1808  // If no label found, reject cell
1809  // It can happen to have this case (4 MC labels per cell is not enough for some cases)
1810  // better to remove. To be treated carefully.
1811  //
1812  if ( nLabels == 0 ) amp = 0;
1813 
1814 }
1815 
1823 //______________________________________________________________________________
1824 void AliEMCALRecoUtils::RecalculateClusterPosition(const AliEMCALGeometry *geom,
1825  AliVCaloCells* cells,
1826  AliVCluster* clu)
1827 {
1828  if (!clu)
1829  {
1830  AliInfo("Cluster pointer null!");
1831  return;
1832  }
1833 
1834  if (fPosAlgo==kPosTowerGlobal) RecalculateClusterPositionFromTowerGlobal( geom, cells, clu);
1835  else if (fPosAlgo==kPosTowerIndex) RecalculateClusterPositionFromTowerIndex ( geom, cells, clu);
1836  else AliDebug(2,"Algorithm to recalculate position not selected, do nothing.");
1837 }
1838 
1847 //_____________________________________________________________________________________________
1848 void AliEMCALRecoUtils::RecalculateClusterPositionFromTowerGlobal(const AliEMCALGeometry *geom,
1849  AliVCaloCells* cells,
1850  AliVCluster* clu)
1851 {
1852  Double_t eCell = 0.;
1853  Float_t fraction = 1.;
1854  Float_t recalFactor = 1.;
1855 
1856  Int_t absId = -1;
1857  Int_t iTower = -1, iIphi = -1, iIeta = -1;
1858  Int_t iSupModMax = -1, iSM=-1, iphi = -1, ieta = -1;
1859  Float_t weight = 0., totalWeight=0.;
1860  Float_t newPos[3] = {-1.,-1.,-1.};
1861  Double_t pLocal[3], pGlobal[3];
1862  Bool_t shared = kFALSE;
1863 
1864  Float_t clEnergy = clu->E(); //Energy already recalibrated previously
1865  if (clEnergy <= 0) return;
1866 
1867  GetMaxEnergyCell(geom, cells, clu, absId, iSupModMax, ieta, iphi,shared);
1868  Double_t depth = GetDepth(clEnergy,fParticleType,iSupModMax) ;
1869 
1870  //printf("** Cluster energy %f, ncells %d, depth %f\n",clEnergy,clu->GetNCells(),depth);
1871 
1872  for (Int_t iDig=0; iDig< clu->GetNCells(); iDig++)
1873  {
1874  absId = clu->GetCellAbsId(iDig);
1875  fraction = clu->GetCellAmplitudeFraction(iDig);
1876  if (fraction < 1e-4) fraction = 1.; // in case unfolding is off
1877 
1878  if (!fCellsRecalibrated)
1879  {
1880  geom->GetCellIndex(absId,iSM,iTower,iIphi,iIeta);
1881  geom->GetCellPhiEtaIndexInSModule(iSM,iTower,iIphi, iIeta,iphi,ieta);
1882  if (IsRecalibrationOn()) {
1883  recalFactor = GetEMCALChannelRecalibrationFactor(iSM,ieta,iphi);
1884  }
1885  }
1886 
1887  eCell = cells->GetCellAmplitude(absId)*fraction*recalFactor;
1888 
1889  weight = GetCellWeight(eCell,clEnergy);
1890  totalWeight += weight;
1891 
1892  geom->RelPosCellInSModule(absId,depth,pLocal[0],pLocal[1],pLocal[2]);
1893  //printf("pLocal (%f,%f,%f), SM %d, absId %d\n",pLocal[0],pLocal[1],pLocal[2],iSupModMax,absId);
1894  geom->GetGlobal(pLocal,pGlobal,iSupModMax);
1895  //printf("pLocal (%f,%f,%f)\n",pGlobal[0],pGlobal[1],pGlobal[2]);
1896 
1897  for (int i=0; i<3; i++ ) newPos[i] += (weight*pGlobal[i]);
1898  }// cell loop
1899 
1900  if (totalWeight>0)
1901  {
1902  for (int i=0; i<3; i++ ) newPos[i] /= totalWeight;
1903  }
1904 
1905  //Float_t pos[]={0,0,0};
1906  //clu->GetPosition(pos);
1907  //printf("OldPos : %2.3f,%2.3f,%2.3f\n",pos[0],pos[1],pos[2]);
1908  //printf("NewPos : %2.3f,%2.3f,%2.3f\n",newPos[0],newPos[1],newPos[2]);
1909 
1910  if (iSupModMax > 1) { //sector 1
1911  newPos[0] +=fMisalTransShift[3];//-=3.093;
1912  newPos[1] +=fMisalTransShift[4];//+=6.82;
1913  newPos[2] +=fMisalTransShift[5];//+=1.635;
1914  //printf(" + : %2.3f,%2.3f,%2.3f\n",fMisalTransShift[3],fMisalTransShift[4],fMisalTransShift[5]);
1915  } else { //sector 0
1916  newPos[0] +=fMisalTransShift[0];//+=1.134;
1917  newPos[1] +=fMisalTransShift[1];//+=8.2;
1918  newPos[2] +=fMisalTransShift[2];//+=1.197;
1919  //printf(" + : %2.3f,%2.3f,%2.3f\n",fMisalTransShift[0],fMisalTransShift[1],fMisalTransShift[2]);
1920  }
1921  //printf("NewPos : %2.3f,%2.3f,%2.3f\n",newPos[0],newPos[1],newPos[2]);
1922 
1923  clu->SetPosition(newPos);
1924 }
1925 
1934 //____________________________________________________________________________________________
1935 void AliEMCALRecoUtils::RecalculateClusterPositionFromTowerIndex(const AliEMCALGeometry *geom,
1936  AliVCaloCells* cells,
1937  AliVCluster* clu)
1938 {
1939  Double_t eCell = 1.;
1940  Float_t fraction = 1.;
1941  Float_t recalFactor = 1.;
1942 
1943  Int_t absId = -1;
1944  Int_t iTower = -1;
1945  Int_t iIphi = -1, iIeta = -1;
1946  Int_t iSupMod = -1, iSupModMax = -1;
1947  Int_t iphi = -1, ieta =-1;
1948  Bool_t shared = kFALSE;
1949 
1950  Float_t clEnergy = clu->E(); //Energy already recalibrated previously.
1951 
1952  if (clEnergy <= 0)
1953  return;
1954  GetMaxEnergyCell(geom, cells, clu, absId, iSupModMax, ieta, iphi,shared);
1955  Float_t depth = GetDepth(clEnergy,fParticleType,iSupMod) ;
1956 
1957  Float_t weight = 0., weightedCol = 0., weightedRow = 0., totalWeight=0.;
1958  Bool_t areInSameSM = kTRUE; //exclude clusters with cells in different SMs for now
1959  Int_t startingSM = -1;
1960 
1961  for (Int_t iDig=0; iDig< clu->GetNCells(); iDig++)
1962  {
1963  absId = clu->GetCellAbsId(iDig);
1964  fraction = clu->GetCellAmplitudeFraction(iDig);
1965  if (fraction < 1e-4) fraction = 1.; // in case unfolding is off
1966 
1967  if (iDig==0) startingSM = iSupMod;
1968  else if (iSupMod != startingSM) areInSameSM = kFALSE;
1969 
1970  eCell = cells->GetCellAmplitude(absId);
1971 
1972  geom->GetCellIndex(absId,iSupMod,iTower,iIphi,iIeta);
1973  geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
1974 
1975  if (!fCellsRecalibrated)
1976  {
1977  if (IsRecalibrationOn()) {
1978  recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
1979  }
1980  }
1981 
1982  eCell = cells->GetCellAmplitude(absId)*fraction*recalFactor;
1983 
1984  weight = GetCellWeight(eCell,clEnergy);
1985  if (weight < 0) weight = 0;
1986  totalWeight += weight;
1987  weightedCol += ieta*weight;
1988  weightedRow += iphi*weight;
1989 
1990  //printf("Max cell? cell %d, amplitude org %f, fraction %f, recalibration %f, amplitude new %f \n",cellAbsId, cells->GetCellAmplitude(cellAbsId), fraction, recalFactor, eCell) ;
1991  }// cell loop
1992 
1993  Float_t xyzNew[]={-1.,-1.,-1.};
1994  if (areInSameSM == kTRUE)
1995  {
1996  //printf("In Same SM\n");
1997  weightedCol = weightedCol/totalWeight;
1998  weightedRow = weightedRow/totalWeight;
1999  geom->RecalculateTowerPosition(weightedRow, weightedCol, iSupModMax, depth, fMisalTransShift, fMisalRotShift, xyzNew);
2000  }
2001  else
2002  {
2003  //printf("In Different SM\n");
2004  geom->RecalculateTowerPosition(iphi, ieta, iSupModMax, depth, fMisalTransShift, fMisalRotShift, xyzNew);
2005  }
2006 
2007  clu->SetPosition(xyzNew);
2008 }
2009 
2017 //___________________________________________________________________________________________
2018 void AliEMCALRecoUtils::RecalculateClusterDistanceToBadChannel(const AliEMCALGeometry * geom,
2019  AliVCaloCells* cells,
2020  AliVCluster * cluster)
2021 {
2022  if (!fRecalDistToBadChannels) return;
2023 
2024  if (!cluster)
2025  {
2026  AliInfo("Cluster pointer null!");
2027  return;
2028  }
2029 
2030  // Get channels map of the supermodule where the cluster is.
2031  Int_t absIdMax = -1, iSupMod =-1, icolM = -1, irowM = -1;
2032  Bool_t shared = kFALSE;
2033  GetMaxEnergyCell(geom, cells, cluster, absIdMax, iSupMod, icolM, irowM, shared);
2034  TH2D* hMap = (TH2D*)fEMCALBadChannelMap->At(iSupMod);
2035 
2036  Int_t dRrow, dRcol;
2037  Float_t minDist = 10000.;
2038  Float_t dist = 0.;
2039 
2040  // Loop on tower status map
2041  for (Int_t irow = 0; irow < AliEMCALGeoParams::fgkEMCALRows; irow++)
2042  {
2043  for (Int_t icol = 0; icol < AliEMCALGeoParams::fgkEMCALCols; icol++)
2044  {
2045  // Check if tower is bad.
2046  if (hMap->GetBinContent(icol,irow)==0) continue;
2047  //printf("AliEMCALRecoUtils::RecalculateDistanceToBadChannels() - \n \t Bad channel in SM %d, col %d, row %d, \n \t Cluster max in col %d, row %d\n",
2048  // iSupMod,icol, irow, icolM,irowM);
2049 
2050  dRrow=TMath::Abs(irowM-irow);
2051  dRcol=TMath::Abs(icolM-icol);
2052  dist=TMath::Sqrt(dRrow*dRrow+dRcol*dRcol);
2053  if (dist < minDist)
2054  {
2055  //printf("MIN DISTANCE TO BAD %2.2f\n",dist);
2056  minDist = dist;
2057  }
2058  }
2059  }
2060 
2061  // In case the cluster is shared by 2 SuperModules, need to check the map of the second Super Module
2062  if (shared)
2063  {
2064  TH2D* hMap2 = 0;
2065  Int_t iSupMod2 = -1;
2066 
2067  // The only possible combinations are (0,1), (2,3) ... (8,9)
2068  if (iSupMod%2) iSupMod2 = iSupMod-1;
2069  else iSupMod2 = iSupMod+1;
2070  hMap2 = (TH2D*)fEMCALBadChannelMap->At(iSupMod2);
2071 
2072  // Loop on tower status map of second super module
2073  for (Int_t irow = 0; irow < AliEMCALGeoParams::fgkEMCALRows; irow++)
2074  {
2075  for (Int_t icol = 0; icol < AliEMCALGeoParams::fgkEMCALCols; icol++)
2076  {
2077  // Check if tower is bad.
2078  if (hMap2->GetBinContent(icol,irow)==0) continue;
2079 
2080  //printf("AliEMCALRecoUtils::RecalculateDistanceToBadChannels(shared) - \n \t Bad channel in SM %d, col %d, row %d \n \t Cluster max in SM %d, col %d, row %d\n",
2081  // iSupMod2,icol, irow,iSupMod,icolM,irowM);
2082  dRrow=TMath::Abs(irow-irowM);
2083 
2084  if (iSupMod%2)
2085  dRcol=TMath::Abs(icol-(AliEMCALGeoParams::fgkEMCALCols+icolM));
2086  else
2087  dRcol=TMath::Abs(AliEMCALGeoParams::fgkEMCALCols+icol-icolM);
2088 
2089  dist=TMath::Sqrt(dRrow*dRrow+dRcol*dRcol);
2090  if (dist < minDist) minDist = dist;
2091  }
2092  }
2093  }// shared cluster in 2 SuperModules
2094 
2095  AliDebug(2,Form("Max cluster cell (SM,col,row)=(%d %d %d) - Distance to Bad Channel %2.2f",iSupMod, icolM, irowM, minDist));
2096  cluster->SetDistanceToBadChannel(minDist);
2097 }
2098 
2105 //__________________________________________________________________
2106 void AliEMCALRecoUtils::RecalculateClusterPID(AliVCluster * cluster)
2107 {
2108  if (!cluster)
2109  {
2110  AliInfo("Cluster pointer null!");
2111  return;
2112  }
2113 
2114  if (cluster->GetM02() != 0)
2115  fPIDUtils->ComputePID(cluster->E(),cluster->GetM02());
2116 
2117  Float_t pidlist[AliPID::kSPECIESCN+1];
2118  for (Int_t i = 0; i < AliPID::kSPECIESCN+1; i++) pidlist[i] = fPIDUtils->GetPIDFinal(i);
2119 
2120  cluster->SetPID(pidlist);
2121 }
2122 
2143 //___________________________________________________________________________________________________________________
2144 void AliEMCALRecoUtils::RecalculateClusterShowerShapeParametersWithCellCuts(const AliEMCALGeometry * geom,
2145  AliVCaloCells* cells, AliVCluster * cluster,
2146  Float_t cellEcut, Float_t cellTimeCut, Int_t bc,
2147  Float_t & enAfterCuts, Float_t & l0, Float_t & l1,
2148  Float_t & disp, Float_t & dEta, Float_t & dPhi,
2149  Float_t & sEta, Float_t & sPhi, Float_t & sEtaPhi)
2150 {
2151  if (!cluster)
2152  {
2153  AliInfo("Cluster pointer null!");
2154  return;
2155  }
2156 
2157  Double_t eCell = 0.;
2158  Double_t tCell = 0.;
2159  Float_t fraction = 1.;
2160  Float_t recalFactor = 1.;
2161  Bool_t isLowGain = kFALSE;
2162 
2163  Int_t iSupMod = -1;
2164  Int_t iTower = -1;
2165  Int_t iIphi = -1;
2166  Int_t iIeta = -1;
2167  Int_t iphi = -1;
2168  Int_t ieta = -1;
2169  Double_t etai = -1.;
2170  Double_t phii = -1.;
2171 
2172  Int_t nstat = 0 ;
2173  Float_t wtot = 0.;
2174  Double_t w = 0.;
2175  Double_t etaMean = 0.;
2176  Double_t phiMean = 0.;
2177 
2178  Double_t pGlobal[3];
2179 
2180  // Loop on cells, calculate the cluster energy, in case a cut on cell energy is added,
2181  // or the non linearity correction was applied
2182  // and to check if the cluster is between 2 SM in eta
2183  Int_t iSM0 = -1;
2184  Bool_t shared = kFALSE;
2185  Float_t energy = 0;
2186 
2187  enAfterCuts = 0;
2188  l0 = 0; l1 = 0;
2189  disp = 0; dEta = 0; dPhi = 0;
2190  sEta = 0; sPhi = 0; sEtaPhi = 0;
2191 
2192  for (Int_t iDigit=0; iDigit < cluster->GetNCells(); iDigit++)
2193  {
2194  // Get from the absid the supermodule, tower and eta/phi numbers
2195 
2196  Int_t absId = cluster->GetCellAbsId(iDigit);
2197 
2198  geom->GetCellIndex(absId,iSupMod,iTower,iIphi,iIeta);
2199  geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi,iIeta, iphi,ieta);
2200 
2201  // Check if there are cells of different SM
2202  if (iDigit == 0 ) iSM0 = iSupMod;
2203  else if (iSupMod!= iSM0) shared = kTRUE;
2204 
2205  // Get the cell energy, if recalibration is on, apply factors
2206  fraction = cluster->GetCellAmplitudeFraction(iDigit);
2207  if (fraction < 1e-4) fraction = 1.; // in case unfolding is off
2208 
2209  if (IsRecalibrationOn())
2210  recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
2211 
2212  eCell = cells->GetCellAmplitude(absId)*fraction*recalFactor;
2213  tCell = cells->GetCellTime (absId);
2214  isLowGain = !(cells->GetCellHighGain(absId));//HG = false -> LG = true
2215 
2216  RecalibrateCellTime(absId, bc, tCell,isLowGain);
2217  tCell*=1e9;
2218  tCell-=fConstantTimeShift;
2219 
2220  if(eCell > 0.05) // at least a minimum 50 MeV cut
2221  energy += eCell;
2222 
2223  if(eCell > cellEcut && TMath::Abs(tCell) < cellTimeCut)
2224  enAfterCuts += eCell;
2225  } // cell loop
2226 
2227 
2228  // Loop on cells to calculate weights and shower shape terms parameters
2229  for (Int_t iDigit=0; iDigit < cluster->GetNCells(); iDigit++)
2230  {
2231  // Get from the absid the supermodule, tower and eta/phi numbers
2232  Int_t absId = cluster->GetCellAbsId(iDigit);
2233 
2234  geom->GetCellIndex(absId,iSupMod,iTower,iIphi,iIeta);
2235  geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi,iIeta, iphi,ieta);
2236 
2237  //Get the cell energy, if recalibration is on, apply factors
2238  fraction = cluster->GetCellAmplitudeFraction(iDigit);
2239  if (fraction < 1e-4) fraction = 1.; // in case unfolding is off
2240 
2241  if (!fCellsRecalibrated && IsRecalibrationOn())
2242  recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
2243 
2244  eCell = cells->GetCellAmplitude(absId)*fraction*recalFactor;
2245  tCell = cells->GetCellTime (absId);
2246  isLowGain = !(cells->GetCellHighGain(absId));//HG = false -> LG = true
2247 
2248  RecalibrateCellTime(absId, bc, tCell,isLowGain);
2249  tCell*=1e9;
2250  tCell-=fConstantTimeShift;
2251 
2252  // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2
2253  // C Side impair SM, nSupMod%2=1; A side pair SM, nSupMod%2=0
2254  if (shared && iSupMod%2) ieta+=AliEMCALGeoParams::fgkEMCALCols;
2255 
2256  if ( energy > 0 && eCell > cellEcut && TMath::Abs(tCell) < cellTimeCut )
2257  {
2258  w = GetCellWeight(eCell, energy);
2259 
2260  // Cell index
2261  if ( fShowerShapeCellLocationType == 0 )
2262  {
2263  etai=(Double_t)ieta;
2264  phii=(Double_t)iphi;
2265  }
2266  // Cell angle location
2267  else if( fShowerShapeCellLocationType == 1 )
2268  {
2269  geom->EtaPhiFromIndex(absId, etai, phii);
2270  etai *= TMath::RadToDeg(); // change units to degrees instead of radians
2271  phii *= TMath::RadToDeg(); // change units to degrees instead of radians
2272  }
2273  else
2274  {
2275  geom->GetGlobal(absId,pGlobal);
2276 
2277  // Cell x-z location
2278  if( fShowerShapeCellLocationType == 2 )
2279  {
2280  etai = pGlobal[2];
2281  phii = pGlobal[0];
2282  }
2283  // Cell r-z location
2284  else
2285  {
2286  etai = pGlobal[2];
2287  phii = TMath::Sqrt(pGlobal[0]*pGlobal[0]+pGlobal[1]*pGlobal[1]);
2288  }
2289  }
2290 
2291  if (w > 0.0)
2292  {
2293  wtot += w ;
2294  nstat++;
2295 
2296  // Shower shape
2297  sEta += w * etai * etai ;
2298  etaMean += w * etai ;
2299  sPhi += w * phii * phii ;
2300  phiMean += w * phii ;
2301  sEtaPhi += w * etai * phii ;
2302  }
2303  }
2304  else if(eCell > 0.05)
2305  AliDebug(2,Form("Wrong energy in cell %f and/or cluster %f\n", eCell, cluster->E()));
2306  } // cell loop
2307 
2308  //printf("sEta %f sPhi %f etaMean %f phiMean %f sEtaPhi %f wtot %f\n",sEta,sPhi,etaMean,phiMean,sEtaPhi, wtot);
2309 
2310  // Normalize to the weight
2311  if (wtot > 0)
2312  {
2313  etaMean /= wtot ;
2314  phiMean /= wtot ;
2315  }
2316  else
2317  AliDebug(2,Form("Wrong weight %f\n", wtot));
2318 
2319  // Loop on cells to calculate dispersion
2320  for (Int_t iDigit=0; iDigit < cluster->GetNCells(); iDigit++)
2321  {
2322  // Get from the absid the supermodule, tower and eta/phi numbers
2323  Int_t absId = cluster->GetCellAbsId(iDigit);
2324 
2325  geom->GetCellIndex(absId,iSupMod,iTower,iIphi,iIeta);
2326  geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi,iIeta, iphi,ieta);
2327 
2328  //Get the cell energy, if recalibration is on, apply factors
2329  fraction = cluster->GetCellAmplitudeFraction(iDigit);
2330  if (fraction < 1e-4) fraction = 1.; // in case unfolding is off
2331  if (IsRecalibrationOn())
2332  recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
2333 
2334  eCell = cells->GetCellAmplitude(absId)*fraction*recalFactor;
2335  tCell = cells->GetCellTime (absId);
2336  isLowGain = !(cells->GetCellHighGain(absId));//HG = false -> LG = true
2337 
2338  RecalibrateCellTime(absId, bc, tCell,isLowGain);
2339  tCell*=1e9;
2340  tCell-=fConstantTimeShift;
2341 
2342  // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2
2343  // C Side impair SM, nSupMod%2=1; A side pair SM, nSupMod%2=0
2344  if (shared && iSupMod%2) ieta+=AliEMCALGeoParams::fgkEMCALCols;
2345 
2346  if ( energy > 0 && eCell > cellEcut && TMath::Abs(tCell) < cellTimeCut )
2347  {
2348  w = GetCellWeight(eCell,cluster->E());
2349 
2350  // Cell index
2351  if ( fShowerShapeCellLocationType == 0 )
2352  {
2353  etai=(Double_t)ieta;
2354  phii=(Double_t)iphi;
2355  }
2356  // Cell angle location
2357  else if( fShowerShapeCellLocationType == 1 )
2358  {
2359  geom->EtaPhiFromIndex(absId, etai, phii);
2360  etai *= TMath::RadToDeg(); // change units to degrees instead of radians
2361  phii *= TMath::RadToDeg(); // change units to degrees instead of radians
2362  }
2363  else
2364  {
2365  geom->GetGlobal(absId,pGlobal);
2366 
2367  // Cell x-z location
2368  if( fShowerShapeCellLocationType == 2 )
2369  {
2370  etai = pGlobal[2];
2371  phii = pGlobal[0];
2372  }
2373  // Cell r-z location
2374  else
2375  {
2376  etai = pGlobal[2];
2377  phii = TMath::Sqrt(pGlobal[0]*pGlobal[0]+pGlobal[1]*pGlobal[1]);
2378  }
2379  }
2380 
2381  if (w > 0.0)
2382  {
2383  disp += w *((etai-etaMean)*(etai-etaMean)+(phii-phiMean)*(phii-phiMean));
2384  dEta += w * (etai-etaMean)*(etai-etaMean) ;
2385  dPhi += w * (phii-phiMean)*(phii-phiMean) ;
2386  }
2387  }
2388  else
2389  AliDebug(2,Form("Wrong energy in cell %f and/or cluster %f\n", eCell, cluster->E()));
2390  }// cell loop
2391 
2392  // Normalize to the weigth and set shower shape parameters
2393  if (wtot > 0 && nstat > 1)
2394  {
2395  disp /= wtot ;
2396  dEta /= wtot ;
2397  dPhi /= wtot ;
2398  sEta /= wtot ;
2399  sPhi /= wtot ;
2400  sEtaPhi /= wtot ;
2401 
2402  sEta -= etaMean * etaMean ;
2403  sPhi -= phiMean * phiMean ;
2404  sEtaPhi -= etaMean * phiMean ;
2405 
2406  l0 = (0.5 * (sEta + sPhi) + TMath::Sqrt( 0.25 * (sEta - sPhi) * (sEta - sPhi) + sEtaPhi * sEtaPhi ));
2407  l1 = (0.5 * (sEta + sPhi) - TMath::Sqrt( 0.25 * (sEta - sPhi) * (sEta - sPhi) + sEtaPhi * sEtaPhi ));
2408 
2409  //printf("sEta %f sPhi %f etaMean %f phiMean %f sEtaPhi %f wtot %f l0 %f l1 %f\n",sEta,sPhi,etaMean,phiMean,sEtaPhi, wtot,l0,l1);
2410 
2411  }
2412  else
2413  {
2414  l0 = 0. ;
2415  l1 = 0. ;
2416  dEta = 0. ; dPhi = 0. ; disp = 0. ;
2417  sEta = 0. ; sPhi = 0. ; sEtaPhi = 0. ;
2418  }
2419 }
2420 
2439 //___________________________________________________________________________________________________________________
2440 void AliEMCALRecoUtils::RecalculateClusterShowerShapeParameters(const AliEMCALGeometry * geom,
2441  AliVCaloCells* cells, AliVCluster * cluster,
2442  Float_t & l0, Float_t & l1,
2443  Float_t & disp, Float_t & dEta, Float_t & dPhi,
2444  Float_t & sEta, Float_t & sPhi, Float_t & sEtaPhi)
2445 {
2446  Float_t newEnergy = 0;
2447  Float_t cellEmin = 0.05; // 50 MeV
2448  Float_t cellTimeWindow = 1000; // open cut
2449  Int_t bc = 0;
2450  AliEMCALRecoUtils::RecalculateClusterShowerShapeParametersWithCellCuts(geom, cells, cluster,
2451  cellEmin, cellTimeWindow, bc,
2452  newEnergy, l0, l1, disp,
2453  dEta, dPhi, sEta, sPhi, sEtaPhi);
2454 }
2455 
2464 //____________________________________________________________________________________________
2465 void AliEMCALRecoUtils::RecalculateClusterShowerShapeParameters(const AliEMCALGeometry * geom,
2466  AliVCaloCells* cells,
2467  AliVCluster * cluster)
2468 {
2469  Float_t l0 = 0., l1 = 0.;
2470  Float_t disp = 0., dEta = 0., dPhi = 0.;
2471  Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
2472 
2473  AliEMCALRecoUtils::RecalculateClusterShowerShapeParameters(geom,cells,cluster,l0,l1,disp,
2474  dEta, dPhi, sEta, sPhi, sEtaPhi);
2475 
2476  cluster->SetM02(l0);
2477  cluster->SetM20(l1);
2478  if (disp > 0. ) cluster->SetDispersion(TMath::Sqrt(disp)) ;
2479 
2480 }
2481 
2494 //____________________________________________________________________________________________
2495 void AliEMCALRecoUtils::RecalculateClusterShowerShapeParametersWithCellCuts(const AliEMCALGeometry * geom,
2496  AliVCaloCells* cells, AliVCluster * cluster,
2497  Float_t cellEcut, Float_t cellTimeCut, Int_t bc,
2498  Float_t & enAfterCuts)
2499 {
2500  Float_t l0 = 0., l1 = 0.;
2501  Float_t disp = 0., dEta = 0., dPhi = 0.;
2502  Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
2503 
2504  AliEMCALRecoUtils::RecalculateClusterShowerShapeParametersWithCellCuts(geom, cells, cluster,
2505  cellEcut, cellTimeCut, bc,
2506  enAfterCuts, l0, l1, disp,
2507  dEta, dPhi, sEta, sPhi, sEtaPhi);
2508 
2509  cluster->SetM02(l0);
2510  cluster->SetM20(l1);
2511  if (disp > 0. ) cluster->SetDispersion(TMath::Sqrt(disp)) ;
2512 
2513 }
2514 
2528 //____________________________________________________________________________
2529 void AliEMCALRecoUtils::FindMatches(AliVEvent *event,
2530  TObjArray * clusterArr,
2531  const AliEMCALGeometry *geom,
2532  AliMCEvent * mc)
2533 {
2534  fMatchedTrackIndex ->Reset();
2535  fMatchedClusterIndex->Reset();
2536  fResidualPhi->Reset();
2537  fResidualEta->Reset();
2538 
2539  fMatchedTrackIndex ->Set(1000);
2540  fMatchedClusterIndex->Set(1000);
2541  fResidualPhi->Set(1000);
2542  fResidualEta->Set(1000);
2543 
2544  AliESDEvent* esdevent = dynamic_cast<AliESDEvent*> (event);
2545  AliAODEvent* aodevent = dynamic_cast<AliAODEvent*> (event);
2546 
2547  // Init the magnetic field if not already on
2548  if (!TGeoGlobalMagField::Instance()->GetField())
2549  {
2550  if (!event->InitMagneticField())
2551  {
2552  AliInfo("Mag Field not initialized, null esd/aod evetn pointers");
2553  }
2554  } // Init mag field
2555 
2556  if (esdevent)
2557  {
2558  UInt_t mask1 = esdevent->GetESDRun()->GetDetectorsInDAQ();
2559  UInt_t mask2 = esdevent->GetESDRun()->GetDetectorsInReco();
2560  Bool_t desc1 = (mask1 >> 3) & 0x1;
2561  Bool_t desc2 = (mask2 >> 3) & 0x1;
2562  if (desc1==0 || desc2==0) {
2563  // AliError(Form("TPC not in DAQ/RECO: %u (%u)/%u (%u)",
2564  // mask1, esdevent->GetESDRun()->GetDetectorsInReco(),
2565  // mask2, esdevent->GetESDRun()->GetDetectorsInDAQ()));
2566  fITSTrackSA=kTRUE;
2567  }
2568  }
2569 
2570  TObjArray *clusterArray = 0x0;
2571  if (!clusterArr)
2572  {
2573  clusterArray = new TObjArray(event->GetNumberOfCaloClusters());
2574  for (Int_t icl=0; icl<event->GetNumberOfCaloClusters(); icl++)
2575  {
2576  AliVCluster *cluster = (AliVCluster*) event->GetCaloCluster(icl);
2577  if (!IsGoodCluster(cluster,geom,(AliVCaloCells*)event->GetEMCALCells()))
2578  continue;
2579  clusterArray->AddAt(cluster,icl);
2580  }
2581  }
2582 
2583  Int_t matched=0;
2584  Double_t cv[21];
2585  TString genName;
2586  for (Int_t i=0; i<21;i++) cv[i]=0;
2587  for (Int_t itr=0; itr<event->GetNumberOfTracks(); itr++)
2588  {
2589  AliExternalTrackParam *trackParam = 0;
2590  Int_t mcLabel = -1;
2591  // If the input event is ESD, the starting point for extrapolation is TPCOut, if available, or TPCInner
2592  AliESDtrack *esdTrack = 0;
2593  AliAODTrack *aodTrack = 0;
2594  if (esdevent)
2595  {
2596  esdTrack = esdevent->GetTrack(itr);
2597  if (!esdTrack) continue;
2598  if (!IsAccepted(esdTrack)) continue;
2599  if (esdTrack->Pt()<fCutMinTrackPt) continue;
2600 
2601  if ( TMath::Abs(esdTrack->Eta()) > 0.9 ) continue;
2602 
2603  // Save some time and memory in case of no DCal present
2604  if( geom->GetNumberOfSuperModules() < 13 ) // Run1 10 (12, 2 not active but present)
2605  {
2606  Double_t phi = esdTrack->Phi()*TMath::RadToDeg();
2607  if ( phi <= 10 || phi >= 250 ) continue;
2608  }
2609 
2610  if (!fITSTrackSA)
2611  trackParam = const_cast<AliExternalTrackParam*>(esdTrack->GetInnerParam()); // if TPC Available
2612  else
2613  trackParam = new AliExternalTrackParam(*esdTrack); // If ITS Track Standing alone
2614 
2615  mcLabel = TMath::Abs(esdTrack->GetLabel());
2616  }
2617 
2618  // If the input event is AOD, the starting point for extrapolation is at vertex
2619  // AOD tracks are selected according to its filterbit.
2620  else if (aodevent)
2621  {
2622  aodTrack = dynamic_cast<AliAODTrack*>(aodevent->GetTrack(itr));
2623 
2624  if (!aodTrack) AliFatal("Not a standard AOD");
2625 
2626  if (!aodTrack) continue;
2627 
2628  if (fAODTPCOnlyTracks)
2629  { // Match with TPC only tracks, default from May 2013, before filter bit 32
2630  //printf("Match with TPC only tracks, accept? %d, test bit 128 <%d> \n", aodTrack->IsTPCOnly(), aodTrack->TestFilterMask(128));
2631  if (!aodTrack->IsTPCConstrained()) continue ;
2632  }
2633  else if (fAODHybridTracks)
2634  { // Match with hybrid tracks
2635  //printf("Match with Hybrid tracks, accept? %d \n", aodTrack->IsHybridGlobalConstrainedGlobal());
2636  if (!aodTrack->IsHybridGlobalConstrainedGlobal()) continue ;
2637  } else
2638  { // Match with tracks on a mask
2639  //printf("Match with tracks having filter bit mask %d, accept? %d \n",fAODFilterMask,aodTrack->TestFilterMask(fAODFilterMask));
2640  if (!aodTrack->TestFilterMask(fAODFilterMask) ) continue; //Select AOD tracks
2641  }
2642 
2643  if (aodTrack->Pt() < fCutMinTrackPt) continue;
2644 
2645  if ( TMath::Abs(aodTrack->Eta()) > 0.9 ) continue;
2646 
2647  // Save some time and memory in case of no DCal present
2648  if( geom->GetNumberOfSuperModules() < 13 ) // Run1 10 (12, 2 not active but present)
2649  {
2650  Double_t phi = aodTrack->Phi()*TMath::RadToDeg();
2651  if ( phi <= 10 || phi >= 250 ) continue;
2652  }
2653 
2654  Double_t pos[3],mom[3];
2655  aodTrack->GetXYZ(pos);
2656  aodTrack->GetPxPyPz(mom);
2657  AliDebug(5,Form("aod track: i=%d | pos=(%5.4f,%5.4f,%5.4f) | mom=(%5.4f,%5.4f,%5.4f) | charge=%d\n",
2658  itr,pos[0],pos[1],pos[2],mom[0],mom[1],mom[2],aodTrack->Charge()));
2659 
2660  trackParam= new AliExternalTrackParam(pos,mom,cv,aodTrack->Charge());
2661 
2662  mcLabel = TMath::Abs(aodTrack->GetLabel());
2663  }
2664 
2665  //Return if the input data is not "AOD" or "ESD"
2666  else
2667  {
2668  AliWarning("Wrong input data type! Should be \"AOD\" or \"ESD\" ");
2669  if (clusterArray)
2670  {
2671  clusterArray->Clear();
2672  delete clusterArray;
2673  }
2674  return;
2675  }
2676 
2677  if (!trackParam) continue;
2678 
2679  //
2680  // Check if track comes from a particular MC generator, do not include it
2681  // if it is not a selected one
2682  //
2683  if( mc && fMCGenerToAcceptForTrack && fNMCGenerToAccept > 0 )
2684  {
2685  mc->GetCocktailGenerator(mcLabel,genName);
2686 
2687  Bool_t generOK = kFALSE;
2688  for(Int_t ig = 0; ig < fNMCGenerToAccept; ig++)
2689  {
2690  if ( genName.Contains(fMCGenerToAccept[ig]) ) generOK = kTRUE;
2691  }
2692 
2693  if ( !generOK ) continue;
2694  }
2695 
2696  // Extrapolate the track to EMCal surface, see AliEMCALRecoUtilsBase
2697  AliExternalTrackParam emcalParam(*trackParam);
2698  Float_t eta, phi, pt;
2699  if (!ExtrapolateTrackToEMCalSurface(&emcalParam, fEMCalSurfaceDistance, fMass, fStepSurface, eta, phi, pt))
2700  {
2701  if (aodevent && trackParam) delete trackParam;
2702  if (fITSTrackSA && trackParam) delete trackParam;
2703  continue;
2704  }
2705 
2706  if ( TMath::Abs(eta) > 0.75 )
2707  {
2708  if ( trackParam && (aodevent || fITSTrackSA) ) delete trackParam;
2709  continue;
2710  }
2711 
2712  // Save some time and memory in case of no DCal present
2713  if ( geom->GetNumberOfSuperModules() < 13 && // Run1 10 (12, 2 not active but present)
2714  ( phi < 70*TMath::DegToRad() || phi > 190*TMath::DegToRad()))
2715  {
2716  if ( trackParam && (aodevent || fITSTrackSA) ) delete trackParam;
2717  continue;
2718  }
2719 
2720  //Find matched clusters
2721  Int_t index = -1;
2722  Float_t dEta = -999, dPhi = -999;
2723  if (!clusterArr)
2724  index = FindMatchedClusterInClusterArr(&emcalParam, &emcalParam, clusterArray, dEta, dPhi);
2725  else
2726  index = FindMatchedClusterInClusterArr(&emcalParam, &emcalParam, clusterArr , dEta, dPhi);
2727 
2728 
2729  if (index>-1)
2730  {
2731  fMatchedTrackIndex ->AddAt(itr,matched);
2732  fMatchedClusterIndex ->AddAt(index,matched);
2733  fResidualEta ->AddAt(dEta,matched);
2734  fResidualPhi ->AddAt(dPhi,matched);
2735  matched++;
2736  }
2737 
2738  if (aodevent && trackParam) delete trackParam;
2739  if (fITSTrackSA && trackParam) delete trackParam;
2740  }//track loop
2741 
2742  if (clusterArray)
2743  {
2744  clusterArray->Clear();
2745  delete clusterArray;
2746  }
2747 
2748  AliDebug(2,Form("Number of matched pairs = %d !\n",matched));
2749 
2750  fMatchedTrackIndex ->Set(matched);
2751  fMatchedClusterIndex ->Set(matched);
2752  fResidualPhi ->Set(matched);
2753  fResidualEta ->Set(matched);
2754 }
2755 
2767 //________________________________________________________________________________
2768 Int_t AliEMCALRecoUtils::FindMatchedClusterInEvent(const AliESDtrack *track,
2769  const AliVEvent *event,
2770  const AliEMCALGeometry *geom,
2771  Float_t &dEta, Float_t &dPhi)
2772 {
2773  Int_t index = -1;
2774 
2775  if ( TMath::Abs(track->Eta()) > 0.9 ) return index;
2776 
2777  // Save some time and memory in case of no DCal present
2778  if( geom->GetNumberOfSuperModules() < 13 ) // Run1 10 (12, 2 not active but present)
2779  {
2780  Double_t phiV = track->Phi()*TMath::RadToDeg();
2781  if ( phiV <= 10 || phiV >= 250 ) return index;
2782  }
2783 
2784  AliExternalTrackParam *trackParam = 0;
2785  if (!fITSTrackSA)
2786  trackParam = const_cast<AliExternalTrackParam*>(track->GetInnerParam()); // If TPC
2787  else
2788  trackParam = new AliExternalTrackParam(*track);
2789 
2790  if (!trackParam) return index;
2791  AliExternalTrackParam emcalParam(*trackParam);
2792 
2793  Float_t eta, phi, pt;
2794  if (!AliEMCALRecoUtilsBase::ExtrapolateTrackToEMCalSurface(&emcalParam, fEMCalSurfaceDistance, fMass, fStepSurface, eta, phi, pt))
2795  {
2796  if (fITSTrackSA) delete trackParam;
2797  return index;
2798  }
2799 
2800  if ( TMath::Abs(eta) > 0.75 )
2801  {
2802  if (fITSTrackSA) delete trackParam;
2803  return index;
2804  }
2805 
2806  // Save some time and memory in case of no DCal present
2807  if ( geom->GetNumberOfSuperModules() < 13 && // Run1 10 (12, 2 not active but present)
2808  ( phi < 70*TMath::DegToRad() || phi > 190*TMath::DegToRad()))
2809  {
2810  if (fITSTrackSA) delete trackParam;
2811  return index;
2812  }
2813 
2814  TObjArray *clusterArr = new TObjArray(event->GetNumberOfCaloClusters());
2815 
2816  for (Int_t icl=0; icl<event->GetNumberOfCaloClusters(); icl++)
2817  {
2818  AliVCluster *cluster = (AliVCluster*) event->GetCaloCluster(icl);
2819  if (!IsGoodCluster(cluster,geom,(AliVCaloCells*)event->GetEMCALCells())) continue;
2820  clusterArr->AddAt(cluster,icl);
2821  }
2822 
2823  index = FindMatchedClusterInClusterArr(&emcalParam, &emcalParam, clusterArr, dEta, dPhi);
2824  clusterArr->Clear();
2825  delete clusterArr;
2826  if (fITSTrackSA) delete trackParam;
2827 
2828  return index;
2829 }
2830 
2841 //_______________________________________________________________________________________________
2842 Int_t AliEMCALRecoUtils::FindMatchedClusterInClusterArr(const AliExternalTrackParam *emcalParam,
2843  AliExternalTrackParam *trkParam,
2844  const TObjArray * clusterArr,
2845  Float_t &dEta, Float_t &dPhi)
2846 {
2847  dEta=-999, dPhi=-999;
2848  Float_t dRMax = fCutR, dEtaMax=fCutEta, dPhiMax=fCutPhi;
2849  Int_t index = -1;
2850  Float_t tmpEta=-999, tmpPhi=-999;
2851 
2852  Double_t exPos[3] = {0.,0.,0.};
2853  if (!emcalParam->GetXYZ(exPos)) return index;
2854 
2855  Float_t clsPos[3] = {0.,0.,0.};
2856  for (Int_t icl=0; icl<clusterArr->GetEntriesFast(); icl++)
2857  {
2858  AliVCluster *cluster = dynamic_cast<AliVCluster*> (clusterArr->At(icl)) ;
2859 
2860  if (!cluster || !cluster->IsEMCAL()) continue;
2861 
2862  cluster->GetPosition(clsPos);
2863 
2864  Double_t dR = TMath::Sqrt(TMath::Power(exPos[0]-clsPos[0],2)+TMath::Power(exPos[1]-clsPos[1],2)+TMath::Power(exPos[2]-clsPos[2],2));
2865  if (dR > fClusterWindow) continue;
2866 
2867  AliExternalTrackParam trkPamTmp (*trkParam);//Retrieve the starting point every time before the extrapolation
2868 
2869  if (!AliEMCALRecoUtilsBase::ExtrapolateTrackToCluster(&trkPamTmp, cluster, fMass, fStepCluster, tmpEta, tmpPhi)) continue;
2870 
2871  if (fCutEtaPhiSum)
2872  {
2873  Float_t tmpR=TMath::Sqrt(tmpEta*tmpEta + tmpPhi*tmpPhi);
2874  if (tmpR<dRMax)
2875  {
2876  dRMax=tmpR;
2877  dEtaMax=tmpEta;
2878  dPhiMax=tmpPhi;
2879  index=icl;
2880  }
2881  }
2882  else if (fCutEtaPhiSeparate)
2883  {
2884  if (TMath::Abs(tmpEta)<TMath::Abs(dEtaMax) && TMath::Abs(tmpPhi)<TMath::Abs(dPhiMax))
2885  {
2886  dEtaMax = tmpEta;
2887  dPhiMax = tmpPhi;
2888  index=icl;
2889  }
2890  }
2891  else
2892  {
2893  AliError("Please specify your cut criteria");
2894  AliError("To cut on sqrt(dEta^2+dPhi^2), use: SwitchOnCutEtaPhiSum()");
2895  AliError("To cut on dEta and dPhi separately, use: SwitchOnCutEtaPhiSeparate()");
2896  return index;
2897  }
2898  }
2899 
2900  dEta=dEtaMax;
2901  dPhi=dPhiMax;
2902 
2903  return index;
2904 }
2905 
2918 //---------------------------------------------------------------------------------
2919 Bool_t AliEMCALRecoUtils::ExtrapolateTrackToCluster(AliExternalTrackParam *trkParam,
2920  const AliVCluster *cluster,
2921  Float_t &tmpEta,
2922  Float_t &tmpPhi)
2923 {
2924  return AliEMCALRecoUtilsBase::ExtrapolateTrackToCluster(trkParam, cluster, fMass, fStepCluster, tmpEta, tmpPhi);
2925 }
2926 
2936 //_______________________________________________________________________
2937 void AliEMCALRecoUtils::GetMatchedResiduals(Int_t clsIndex,
2938  Float_t &dEta, Float_t &dPhi)
2939 {
2940  if (FindMatchedPosForCluster(clsIndex) >= 999)
2941  {
2942  AliDebug(2,"No matched tracks found!\n");
2943  dEta=999.;
2944  dPhi=999.;
2945  return;
2946  }
2947 
2948  dEta = fResidualEta->At(FindMatchedPosForCluster(clsIndex));
2949  dPhi = fResidualPhi->At(FindMatchedPosForCluster(clsIndex));
2950 }
2951 
2961 //______________________________________________________________________________________________
2962 void AliEMCALRecoUtils::GetMatchedClusterResiduals(Int_t trkIndex, Float_t &dEta, Float_t &dPhi)
2963 {
2964  if (FindMatchedPosForTrack(trkIndex) >= 999)
2965  {
2966  AliDebug(2,"No matched cluster found!\n");
2967  dEta=999.;
2968  dPhi=999.;
2969  return;
2970  }
2971 
2972  dEta = fResidualEta->At(FindMatchedPosForTrack(trkIndex));
2973  dPhi = fResidualPhi->At(FindMatchedPosForTrack(trkIndex));
2974 }
2975 
2982 //__________________________________________________________
2983 Int_t AliEMCALRecoUtils::GetMatchedTrackIndex(Int_t clsIndex)
2984 {
2985  if (IsClusterMatched(clsIndex))
2986  return fMatchedTrackIndex->At(FindMatchedPosForCluster(clsIndex));
2987  else
2988  return -1;
2989 }
2990 
2997 //__________________________________________________________
2998 Int_t AliEMCALRecoUtils::GetMatchedClusterIndex(Int_t trkIndex)
2999 {
3000  if (IsTrackMatched(trkIndex))
3001  return fMatchedClusterIndex->At(FindMatchedPosForTrack(trkIndex));
3002  else
3003  return -1;
3004 }
3005 
3013 //______________________________________________________________
3014 Bool_t AliEMCALRecoUtils::IsClusterMatched(Int_t clsIndex) const
3015 {
3016  if (FindMatchedPosForCluster(clsIndex) < 999)
3017  return kTRUE;
3018  else
3019  return kFALSE;
3020 }
3021 
3029 //____________________________________________________________
3030 Bool_t AliEMCALRecoUtils::IsTrackMatched(Int_t trkIndex) const
3031 {
3032  if (FindMatchedPosForTrack(trkIndex) < 999)
3033  return kTRUE;
3034  else
3035  return kFALSE;
3036 }
3037 
3045 //______________________________________________________________________
3046 UInt_t AliEMCALRecoUtils::FindMatchedPosForCluster(Int_t clsIndex) const
3047 {
3048  Float_t tmpR = fCutR;
3049  UInt_t pos = 999;
3050 
3051  for (Int_t i=0; i<fMatchedClusterIndex->GetSize(); i++)
3052  {
3053  if (fMatchedClusterIndex->At(i)==clsIndex)
3054  {
3055  Float_t r = TMath::Sqrt(fResidualEta->At(i)*fResidualEta->At(i) + fResidualPhi->At(i)*fResidualPhi->At(i));
3056  if (r<tmpR)
3057  {
3058  pos=i;
3059  tmpR=r;
3060  AliDebug(3,Form("Matched cluster index: index: %d, dEta: %2.4f, dPhi: %2.4f.\n",
3061  fMatchedClusterIndex->At(i),fResidualEta->At(i),fResidualPhi->At(i)));
3062  }
3063  }
3064  }
3065  return pos;
3066 }
3067 
3075 //____________________________________________________________________
3076 UInt_t AliEMCALRecoUtils::FindMatchedPosForTrack(Int_t trkIndex) const
3077 {
3078  Float_t tmpR = fCutR;
3079  UInt_t pos = 999;
3080 
3081  for (Int_t i=0; i<fMatchedTrackIndex->GetSize(); i++)
3082  {
3083  if (fMatchedTrackIndex->At(i)==trkIndex)
3084  {
3085  Float_t r = TMath::Sqrt(fResidualEta->At(i)*fResidualEta->At(i) + fResidualPhi->At(i)*fResidualPhi->At(i));
3086  if (r<tmpR)
3087  {
3088  pos=i;
3089  tmpR=r;
3090  AliDebug(3,Form("Matched track index: index: %d, dEta: %2.4f, dPhi: %2.4f.\n",
3091  fMatchedTrackIndex->At(i),fResidualEta->At(i),fResidualPhi->At(i)));
3092  }
3093  }
3094  }
3095  return pos;
3096 }
3097 
3112 //__________________________________________________________________________
3113 Bool_t AliEMCALRecoUtils::IsGoodCluster(AliVCluster *cluster,
3114  const AliEMCALGeometry *geom,
3115  AliVCaloCells* cells, Int_t bc)
3116 {
3117  Bool_t isGood=kTRUE;
3118 
3119  if (!cluster || !cluster->IsEMCAL()) return kFALSE;
3120  if (ClusterContainsBadChannel(geom,cluster->GetCellsAbsId(),cluster->GetNCells())) return kFALSE;
3121  if (!CheckCellFiducialRegion(geom,cluster,cells)) return kFALSE;
3122  if (IsExoticCluster(cluster, cells,bc)) return kFALSE;
3123 
3124  return isGood;
3125 }
3126 
3137 //__________________________________________________________
3138 Bool_t AliEMCALRecoUtils::IsAccepted(AliESDtrack *esdTrack)
3139 {
3140  UInt_t status = esdTrack->GetStatus();
3141 
3142  Int_t nClustersITS = esdTrack->GetITSclusters(0);
3143  Int_t nClustersTPC = esdTrack->GetTPCclusters(0);
3144 
3145  Float_t chi2PerClusterITS = -1;
3146  Float_t chi2PerClusterTPC = -1;
3147  if (nClustersITS!=0)
3148  chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
3149  if (nClustersTPC!=0)
3150  chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
3151 
3152  //
3153  // DCA cuts
3154  // Only to be used for primary
3155  //
3156  if ( fTrackCutsType == kGlobalCut )
3157  {
3158  Float_t maxDCAToVertexXYPtDep = 0.0182 + 0.0350/TMath::Power(esdTrack->Pt(),1.01);
3159  // This expression comes from AliESDtrackCuts::GetStandardITSTPCTrackCuts2010()
3160 
3161  //AliDebug(3,Form("Track pT = %f, DCAtoVertexXY = %f",esdTrack->Pt(),MaxDCAToVertexXYPtDep));
3162 
3163  SetMaxDCAToVertexXY(maxDCAToVertexXYPtDep); //Set pT dependent DCA cut to vertex in x-y plane
3164  }
3165  else if( fTrackCutsType == kGlobalCut2011 )
3166  {
3167  Float_t maxDCAToVertexXYPtDep = 0.0105 + 0.0350/TMath::Power(esdTrack->Pt(),1.1);
3168  // This expression comes from AliESDtrackCuts::GetStandardITSTPCTrackCuts2011()
3169 
3170  //AliDebug(3,Form("Track pT = %f, DCAtoVertexXY = %f",esdTrack->Pt(),MaxDCAToVertexXYPtDep));
3171 
3172  SetMaxDCAToVertexXY(maxDCAToVertexXYPtDep); //Set pT dependent DCA cut to vertex in x-y plane
3173  }
3174 
3175  Float_t b[2];
3176  Float_t bCov[3];
3177  esdTrack->GetImpactParameters(b,bCov);
3178  if (bCov[0]<=0 || bCov[2]<=0)
3179  {
3180  AliDebug(1, "Estimated b resolution lower or equal zero!");
3181  bCov[0]=0; bCov[2]=0;
3182  }
3183 
3184  Float_t dcaToVertexXY = b[0];
3185  Float_t dcaToVertexZ = b[1];
3186  Float_t dcaToVertex = -1;
3187 
3188  if (fCutDCAToVertex2D)
3189  dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY / fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY +
3190  dcaToVertexZ*dcaToVertexZ / fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ );
3191  else
3192  dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
3193 
3194  // cut the track?
3195 
3196  Bool_t cuts[kNCuts];
3197  for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
3198 
3199  // track quality cuts
3200  if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
3201  cuts[0]=kTRUE;
3202  if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
3203  cuts[1]=kTRUE;
3204  if (nClustersTPC<fCutMinNClusterTPC)
3205  cuts[2]=kTRUE;
3206  if (nClustersITS<fCutMinNClusterITS)
3207  cuts[3]=kTRUE;
3208  if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
3209  cuts[4]=kTRUE;
3210  if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
3211  cuts[5]=kTRUE;
3212  if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
3213  cuts[6]=kTRUE;
3214  if (fCutDCAToVertex2D && dcaToVertex > 1)
3215  cuts[7] = kTRUE;
3216  if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
3217  cuts[8] = kTRUE;
3218  if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
3219  cuts[9] = kTRUE;
3220 
3221  if (fTrackCutsType == kGlobalCut || fTrackCutsType == kGlobalCut2011)
3222  {
3223  //Require at least one SPD point + anything else in ITS
3224  if ( (esdTrack->HasPointOnITSLayer(0) || esdTrack->HasPointOnITSLayer(1)) == kFALSE)
3225  cuts[10] = kTRUE;
3226  }
3227 
3228  // ITS
3229  if (fCutRequireITSStandAlone || fCutRequireITSpureSA)
3230  {
3231  if ((status & AliESDtrack::kITSin) == 0 || (status & AliESDtrack::kTPCin))
3232  {
3233  // TPC tracks
3234  cuts[11] = kTRUE;
3235  }
3236  else
3237  {
3238  // ITS standalone tracks
3239  if (fCutRequireITSStandAlone && !fCutRequireITSpureSA)
3240  {
3241  if (status & AliESDtrack::kITSpureSA) cuts[11] = kTRUE;
3242  }
3243  else if (fCutRequireITSpureSA)
3244  {
3245  if (!(status & AliESDtrack::kITSpureSA)) cuts[11] = kTRUE;
3246  }
3247  }
3248  }
3249 
3250  Bool_t cut=kFALSE;
3251  for (Int_t i=0; i<kNCuts; i++)
3252  if (cuts[i]) { cut = kTRUE ; }
3253 
3254  // cut the track
3255  if (cut)
3256  return kFALSE;
3257  else
3258  return kTRUE;
3259 }
3260 
3266 //_____________________________________
3267 void AliEMCALRecoUtils::InitTrackCuts()
3268 {
3269  switch (fTrackCutsType)
3270  {
3271  case kTPCOnlyCut:
3272  {
3273  AliInfo(Form("Track cuts for matching: AliESDtrackCuts::GetStandardTPCOnlyTrackCuts()"));
3274  //TPC
3275  SetMinNClustersTPC(70);
3276  SetMaxChi2PerClusterTPC(4);
3277  SetAcceptKinkDaughters(kFALSE);
3278  SetRequireTPCRefit(kFALSE);
3279 
3280  //ITS
3281  SetRequireITSRefit(kFALSE);
3282  SetMaxDCAToVertexZ(3.2);
3283  SetMaxDCAToVertexXY(2.4);
3284  SetDCAToVertex2D(kTRUE);
3285 
3286  break;
3287  }
3288 
3289  case kGlobalCut:
3290  {
3291  AliInfo(Form("Track cuts for matching: AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE)"));
3292  //TPC
3293  SetMinNClustersTPC(70);
3294  SetMaxChi2PerClusterTPC(4);
3295  SetAcceptKinkDaughters(kFALSE);
3296  SetRequireTPCRefit(kTRUE);
3297 
3298  //ITS
3299  SetRequireITSRefit(kTRUE);
3300  SetMaxDCAToVertexZ(2);
3301  SetMaxDCAToVertexXY();
3302  SetDCAToVertex2D(kFALSE);
3303 
3304  break;
3305  }
3306 
3307  case kLooseCut:
3308  {
3309  AliInfo(Form("Track cuts for matching: Loose cut w/o DCA cut"));
3310  SetMinNClustersTPC(50);
3311  SetAcceptKinkDaughters(kTRUE);
3312 
3313  break;
3314  }
3315 
3316  case kITSStandAlone:
3317  {
3318  AliInfo(Form("Track cuts for matching: ITS Stand Alone tracks cut w/o DCA cut"));
3319  SetRequireITSRefit(kTRUE);
3320  SetRequireITSStandAlone(kTRUE);
3321  SetITSTrackSA(kTRUE);
3322  break;
3323  }
3324 
3325  case kGlobalCut2011:
3326  {
3327  AliInfo(Form("Track cuts for matching: AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kTRUE)"));
3328  //TPC
3329  SetMinNClustersTPC(50);
3330  SetMaxChi2PerClusterTPC(4);
3331  SetAcceptKinkDaughters(kFALSE);
3332  SetRequireTPCRefit(kTRUE);
3333 
3334  //ITS
3335  SetRequireITSRefit(kTRUE);
3336  SetMaxDCAToVertexZ(2);
3337  SetMaxDCAToVertexXY();
3338  SetDCAToVertex2D(kFALSE);
3339 
3340  break;
3341  }
3342 
3343  case kLooseCutWithITSrefit:
3344  {
3345  AliInfo(Form("Track cuts for matching: Loose cut w/o DCA cut plus ITSrefit"));
3346  SetMinNClustersTPC(50);
3347  SetAcceptKinkDaughters(kTRUE);
3348  SetRequireITSRefit(kTRUE);
3349 
3350  break;
3351  }
3352  }
3353 }
3354 
3360 //________________________________________________________________________
3361 void AliEMCALRecoUtils::SetClusterMatchedToTrack(const AliVEvent *event)
3362 {
3363  Int_t nTracks = event->GetNumberOfTracks();
3364  for (Int_t iTrack = 0; iTrack < nTracks; ++iTrack)
3365  {
3366  AliVTrack* track = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
3367  if (!track)
3368  {
3369  AliWarning(Form("Could not receive track %d", iTrack));
3370  continue;
3371  }
3372 
3373  Int_t matchClusIndex = GetMatchedClusterIndex(iTrack);
3374  track->SetEMCALcluster(matchClusIndex); //sets -1 if track not matched within residual
3375 
3376  /*the following can be done better if AliVTrack::SetStatus will be there. Patch pending with Andreas/Peter*/
3377  AliESDtrack* esdtrack = dynamic_cast<AliESDtrack*>(track);
3378  if (esdtrack)
3379  {
3380  if (matchClusIndex != -1)
3381  esdtrack->SetStatus(AliESDtrack::kEMCALmatch);
3382  else
3383  esdtrack->ResetStatus(AliESDtrack::kEMCALmatch);
3384  }
3385  else
3386  {
3387  AliAODTrack* aodtrack = dynamic_cast<AliAODTrack*>(track);
3388  if (matchClusIndex != -1)
3389  aodtrack->SetStatus(AliESDtrack::kEMCALmatch);
3390  else
3391  aodtrack->ResetStatus(AliESDtrack::kEMCALmatch);
3392  }
3393  }
3394  AliDebug(2,"Track matched to closest cluster");
3395 }
3396 
3403 //_________________________________________________________________________
3404 void AliEMCALRecoUtils::SetTracksMatchedToCluster(const AliVEvent *event)
3405 {
3406  for (Int_t iClus=0; iClus < event->GetNumberOfCaloClusters(); ++iClus)
3407  {
3408  AliVCluster *cluster = event->GetCaloCluster(iClus);
3409  if (!cluster->IsEMCAL())
3410  continue;
3411 
3412  //
3413  // Remove old matches in cluster
3414  //
3415  if(cluster->GetNTracksMatched() > 0)
3416  {
3417  if(!strcmp("AliESDCaloCluster",Form("%s",cluster->ClassName())))
3418  {
3419  TArrayI arrayTrackMatched(0);
3420  ((AliESDCaloCluster*)cluster)->AddTracksMatched(arrayTrackMatched);
3421  }
3422  else
3423  {
3424  for(Int_t iTrack = 0; iTrack < cluster->GetNTracksMatched(); iTrack++)
3425  {
3426  ((AliAODCaloCluster*)cluster)->RemoveTrackMatched((TObject*)((AliAODCaloCluster*)cluster)->GetTrackMatched(iTrack));
3427  }
3428  }
3429  }
3430 
3431  //
3432  // Find new matches and put them in the cluster
3433  //
3434  Int_t nTracks = event->GetNumberOfTracks();
3435  TArrayI arrayTrackMatched(nTracks);
3436 
3437  // Get the closest track matched to the cluster
3438  Int_t nMatched = 0;
3439  Int_t matchTrackIndex = GetMatchedTrackIndex(iClus);
3440  if (matchTrackIndex != -1)
3441  {
3442  arrayTrackMatched[nMatched] = matchTrackIndex;
3443  nMatched++;
3444  }
3445 
3446  // Get all other tracks matched to the cluster
3447  for (Int_t iTrk=0; iTrk<nTracks; ++iTrk)
3448  {
3449  AliVTrack* track = dynamic_cast<AliVTrack*>(event->GetTrack(iTrk));
3450 
3451  if( !track ) continue;
3452 
3453  if ( iTrk == matchTrackIndex ) continue;
3454 
3455  if ( track->GetEMCALcluster() == iClus )
3456  {
3457  arrayTrackMatched[nMatched] = iTrk;
3458  ++nMatched;
3459  }
3460  }
3461 
3462  AliDebug(2,Form("cluster E %f, N matches %d, first match %d\n",cluster->E(),nMatched,arrayTrackMatched[0]));
3463 
3464  arrayTrackMatched.Set(nMatched);
3465  AliESDCaloCluster *esdcluster = dynamic_cast<AliESDCaloCluster*>(cluster);
3466  if (esdcluster)
3467  esdcluster->AddTracksMatched(arrayTrackMatched);
3468  else if ( nMatched > 0 )
3469  {
3470  AliAODCaloCluster *aodcluster = dynamic_cast<AliAODCaloCluster*>(cluster);
3471  if (aodcluster)
3472  {
3473  aodcluster->AddTrackMatched(event->GetTrack(arrayTrackMatched.At(0)));
3474  //AliAODTrack *aodtrack=dynamic_cast<AliAODTrack*>(event->GetTrack(arrayTrackMatched.At(0)));
3475  //printf("Is the closest matching track with ID %d a 128? %d what's its full filter map? %u\n",aodtrack->GetID(), aodtrack->TestFilterBit(128),aodtrack->GetFilterMap());
3476  //printf("With specs: pt %.4f, eta %.4f, phi %.4f\n",aodtrack->Pt(),aodtrack->Eta(), aodtrack->Phi());
3477  }
3478  }
3479 
3480  Float_t eta= -999, phi = -999;
3481  if (matchTrackIndex != -1)
3482  GetMatchedResiduals(iClus, eta, phi);
3483 
3484  cluster->SetTrackDistance(phi, eta);
3485  }
3486 
3487  AliDebug(2,"Cluster matched to tracks");
3488 }
3489 
3490 void AliEMCALRecoUtils::SetEMCALChannelRecalibrationFactors(const TObjArray *map) {
3491  if(fEMCALRecalibrationFactors) fEMCALRecalibrationFactors->Clear();
3492  else {
3493  fEMCALRecalibrationFactors = new TObjArray(map->GetEntries());
3494  fEMCALRecalibrationFactors->SetOwner(true);
3495  }
3496  if(!fEMCALRecalibrationFactors->IsOwner()){
3497  // Must claim ownership since the new objects are owend by this instance
3498  fEMCALRecalibrationFactors->SetOwner(kTRUE);
3499  }
3500  for(int i = 0; i < map->GetEntries(); i++){
3501  TH2F *hist = dynamic_cast<TH2F *>(map->At(i));
3502  if(!hist) continue;
3503  this->SetEMCALChannelRecalibrationFactors(i, hist);
3504  }
3505 }
3506 
3507 void AliEMCALRecoUtils::SetEMCALChannelRecalibrationFactors(Int_t iSM , const TH2F* h) {
3508  if(!fEMCALRecalibrationFactors){
3509  fEMCALRecalibrationFactors = new TObjArray(iSM);
3510  fEMCALRecalibrationFactors->SetOwner(true);
3511  }
3512  if(fEMCALRecalibrationFactors->GetEntries() <= iSM) fEMCALRecalibrationFactors->Expand(iSM+1);
3513  if(fEMCALRecalibrationFactors->At(iSM)) fEMCALRecalibrationFactors->RemoveAt(iSM);
3514  TH2F *clone = new TH2F(*h);
3515  clone->SetDirectory(NULL);
3516  fEMCALRecalibrationFactors->AddAt(clone,iSM);
3517 }
3518 
3519 void AliEMCALRecoUtils::SetEMCALChannelStatusMap(const TObjArray *map) {
3520  if(fEMCALBadChannelMap) fEMCALBadChannelMap->Clear();
3521  else {
3522  fEMCALBadChannelMap = new TObjArray(map->GetEntries());
3523  fEMCALBadChannelMap->SetOwner(true);
3524  }
3525  if(!fEMCALBadChannelMap->IsOwner()) {
3526  // Must claim ownership since the new objects are owend by this instance
3527  fEMCALBadChannelMap->SetOwner(true);
3528  }
3529  for(int i = 0; i < map->GetEntries(); i++){
3530  TH2I *hist = dynamic_cast<TH2I *>(map->At(i));
3531  if(!hist) continue;
3532  this->SetEMCALChannelStatusMap(i, hist);
3533  }
3534 }
3535 
3536 void AliEMCALRecoUtils::SetEMCALChannelStatusMap(Int_t iSM , const TH2I* h) {
3537  if(!fEMCALBadChannelMap){
3538  fEMCALBadChannelMap = new TObjArray(iSM);
3539  fEMCALBadChannelMap->SetOwner(true);
3540  }
3541  if(fEMCALBadChannelMap->GetEntries() <= iSM) fEMCALBadChannelMap->Expand(iSM+1);
3542  if(fEMCALBadChannelMap->At(iSM)) fEMCALBadChannelMap->RemoveAt(iSM);
3543  TH2I *clone = new TH2I(*h);
3544  clone->SetDirectory(NULL);
3545  fEMCALBadChannelMap->AddAt(clone,iSM);
3546 }
3547 
3548 void AliEMCALRecoUtils::SetEMCALChannelTimeRecalibrationFactors(const TObjArray *map) {
3549  if(fEMCALTimeRecalibrationFactors) fEMCALTimeRecalibrationFactors->Clear();
3550  else {
3551  fEMCALTimeRecalibrationFactors = new TObjArray(map->GetEntries());
3552  fEMCALTimeRecalibrationFactors->SetOwner(true);
3553  }
3554  if(!fEMCALTimeRecalibrationFactors->IsOwner()) {
3555  // Must claim ownership since the new objects are owend by this instance
3556  fEMCALTimeRecalibrationFactors->SetOwner(kTRUE);
3557  }
3558  for(int i = 0; i < map->GetEntries(); i++){
3559  TH1F *hist = dynamic_cast<TH1F *>(map->At(i));
3560  if(!hist) continue;
3561  this->SetEMCALChannelTimeRecalibrationFactors(i, hist);
3562  }
3563 }
3564 
3565 void AliEMCALRecoUtils::SetEMCALChannelTimeRecalibrationFactors(Int_t bc, const TH1F* h){
3566  if(!fEMCALTimeRecalibrationFactors){
3567  fEMCALTimeRecalibrationFactors = new TObjArray(bc);
3568  fEMCALTimeRecalibrationFactors->SetOwner(true);
3569  }
3570  if(fEMCALTimeRecalibrationFactors->GetEntries() <= bc) fEMCALTimeRecalibrationFactors->Expand(bc+1);
3571  if(fEMCALTimeRecalibrationFactors->At(bc)) fEMCALTimeRecalibrationFactors->RemoveAt(bc);
3572  TH1F *clone = new TH1F(*h);
3573  clone->SetDirectory(NULL);
3574  fEMCALTimeRecalibrationFactors->AddAt(clone,bc);
3575 }
3576 
3577 void AliEMCALRecoUtils::SetEMCALL1PhaseInTimeRecalibrationForAllSM(const TObjArray *map) {
3578  if(fEMCALL1PhaseInTimeRecalibration) fEMCALL1PhaseInTimeRecalibration->Clear();
3579  else {
3580  fEMCALL1PhaseInTimeRecalibration = new TObjArray(map->GetEntries());
3581  fEMCALL1PhaseInTimeRecalibration->SetOwner(true);
3582  }
3583  if(!fEMCALL1PhaseInTimeRecalibration->IsOwner()) {
3584  // Must claim ownership since the new objects are owend by this instance
3585  fEMCALL1PhaseInTimeRecalibration->SetOwner(true);
3586  }
3587  for(int i = 0; i < map->GetEntries(); i++) {
3588  TH1C *hist = dynamic_cast<TH1C *>(map->At(i));
3589  if(!hist) continue;
3590  SetEMCALL1PhaseInTimeRecalibrationForAllSM(hist);
3591  }
3592 }
3593 
3594 void AliEMCALRecoUtils::SetEMCALL1PhaseInTimeRecalibrationForAllSM(const TH1C* h) {
3595  if(!fEMCALL1PhaseInTimeRecalibration){
3596  fEMCALL1PhaseInTimeRecalibration = new TObjArray(1);
3597  fEMCALL1PhaseInTimeRecalibration->SetOwner(true);
3598  }
3599  if(fEMCALL1PhaseInTimeRecalibration->GetEntries()<1) fEMCALL1PhaseInTimeRecalibration->Expand(1);
3600  TH1C *clone = new TH1C(*h);
3601  clone->SetDirectory(NULL);
3602  fEMCALL1PhaseInTimeRecalibration->AddAt(clone,0);
3603 }
3604 
3608 //___________________________________________________
3609 void AliEMCALRecoUtils::Print(const Option_t *) const
3610 {
3611  printf("-------------------------------------------------------------------------------------------------------------------------------------- \n");
3612  printf("AliEMCALRecoUtils Settings: \n");
3613  printf("\tMisalignment shifts\n");
3614  for (Int_t i=0; i<5; i++) printf("\t\t sector %d, traslation (x,y,z)=(%f,%f,%f), rotation (x,y,z)=(%f,%f,%f)\n",i,
3615  fMisalTransShift[i*3],fMisalTransShift[i*3+1],fMisalTransShift[i*3+2],
3616  fMisalRotShift[i*3], fMisalRotShift[i*3+1], fMisalRotShift[i*3+2] );
3617  printf("\tNon linearity function %d, parameters:\n", fNonLinearityFunction);
3618  if (fNonLinearityFunction != 3) // print only if not kNoCorrection
3619  for (Int_t i=0; i<10; i++) printf("param[%d]=%f\n",i, fNonLinearityParams[i]);
3620 
3621  printf("\tPosition Recalculation option %d, Particle Type %d, fW0 %2.2f, Recalibrate Data %d \n",fPosAlgo,fParticleType,fW0, fRecalibration);
3622 
3623  printf("\tMatching criteria: ");
3624  if (fCutEtaPhiSum) {
3625  printf("\t\tsqrt(dEta^2+dPhi^2)<%4.3f\n",fCutR);
3626  } else if (fCutEtaPhiSeparate) {
3627  printf("\t\tdEta<%4.3f, dPhi<%4.3f\n",fCutEta,fCutPhi);
3628  } else {
3629  printf("\t\tError\n");
3630  printf("\t\tplease specify your cut criteria\n");
3631  printf("\t\tTo cut on sqrt(dEta^2+dPhi^2), use: SwitchOnCutEtaPhiSum()\n");
3632  printf("\t\tTo cut on dEta and dPhi separately, use: SwitchOnCutEtaPhiSeparate()\n");
3633  }
3634 
3635  printf("\tMass hypothesis = %2.3f [GeV/c^2], extrapolation step to surface = %2.2f[cm], step to cluster = %2.2f[cm]\n",fMass,fStepSurface, fStepCluster);
3636  printf("\tCluster selection window: dR < %2.0f\n",fClusterWindow);
3637 
3638  printf("\tTrack cuts: \n");
3639  printf("\t\tMinimum track pT: %1.2f\n",fCutMinTrackPt);
3640  printf("\t\tAOD track selection: tpc only %d, or hybrid %d, or mask: %d\n",fAODTPCOnlyTracks,fAODHybridTracks, fAODFilterMask);
3641  printf("\t\tTPCRefit = %d, ITSRefit = %d\n",fCutRequireTPCRefit,fCutRequireITSRefit);
3642  printf("\t\tAcceptKinks = %d\n",fCutAcceptKinkDaughters);
3643  printf("\t\tMinNCulsterTPC = %d, MinNClusterITS = %d\n",fCutMinNClusterTPC,fCutMinNClusterITS);
3644  printf("\t\tMaxChi2TPC = %2.2f, MaxChi2ITS = %2.2f\n",fCutMaxChi2PerClusterTPC,fCutMaxChi2PerClusterITS);
3645  printf("\t\tDCSToVertex2D = %d, MaxDCAToVertexXY = %2.2f, MaxDCAToVertexZ = %2.2f\n",fCutDCAToVertex2D,fCutMaxDCAToVertexXY,fCutMaxDCAToVertexZ);
3646  printf("-------------------------------------------------------------------------------------------------------------------------------------- \n");
3647 }
AliEMCALRecoUtils::RecalculateCellLabelsRemoveAddedGenerator
void RecalculateCellLabelsRemoveAddedGenerator(Int_t absID, AliVCluster *clus, AliMCEvent *mc, Float_t &amp, TArrayI &labeArr, TArrayF &eDepArr) const
Definition: AliEMCALRecoUtils.cxx:1722
AliEMCALRecoUtils::IsL1PhaseInTimeRecalibrationOn
Bool_t IsL1PhaseInTimeRecalibrationOn() const
Definition: AliEMCALRecoUtils.h:216
AliEMCALRecoUtils::fMass
Double_t fMass
Mass hypothesis of the track.
Definition: AliEMCALRecoUtils.h:513
AliEMCALRecoUtils::fMCGenerToAcceptForTrack
Bool_t fMCGenerToAcceptForTrack
Activate the removal of tracks entering the track matching that come from a particular generator...
Definition: AliEMCALRecoUtils.h:538
Double_t
double Double_t
Definition: External.C:58
AliEMCALRecoUtils::CheckCellFiducialRegion
Bool_t CheckCellFiducialRegion(const AliEMCALGeometry *geom, const AliVCluster *cluster, AliVCaloCells *cells)
Definition: AliEMCALRecoUtils.cxx:446
AliEMCALRecoUtils::fAODTPCOnlyTracks
Bool_t fAODTPCOnlyTracks
Match with TPC only tracks.
Definition: AliEMCALRecoUtils.h:501
AliEMCALRecoUtils.h
AliEMCALRecoUtils::fCutRequireTPCRefit
Bool_t fCutRequireTPCRefit
Require TPC refit.
Definition: AliEMCALRecoUtils.h:526
TH2F
Definition: External.C:236
AliEMCALRecoUtils::SetRequireITSStandAlone
void SetRequireITSStandAlone(Bool_t b=kFALSE)
Definition: AliEMCALRecoUtils.h:373
AliEMCALRecoUtils::IsClusterMatched
Bool_t IsClusterMatched(Int_t clsIndex) const
Definition: AliEMCALRecoUtils.cxx:3014
AliEMCALRecoUtils::GetEMCALChannelTimeRecalibrationFactor
Float_t GetEMCALChannelTimeRecalibrationFactor(Int_t bc, Int_t absID, Bool_t isLGon=kFALSE) const
Definition: AliEMCALRecoUtils.h:198
AliEMCALRecoUtils::fMatchedClusterIndex
TArrayI * fMatchedClusterIndex
Array that stores indexes of matched clusters.
Definition: AliEMCALRecoUtils.h:504
AliEMCALRecoUtils::fEMCALTimeRecalibrationFactors
TObjArray * fEMCALTimeRecalibrationFactors
Array of histograms with map of time recalibration factors, EMCAL.
Definition: AliEMCALRecoUtils.h:469
AliEMCALRecoUtils::SetMaxDCAToVertexXY
void SetMaxDCAToVertexXY(Float_t dist=1e10)
Definition: AliEMCALRecoUtils.h:370
AliEMCALRecoUtils::fRandom
TRandom3 fRandom
Random generator for cluster energy smearing.
Definition: AliEMCALRecoUtils.h:459
AliEMCALRecoUtils::fTimeRecalibration
Bool_t fTimeRecalibration
Switch on or off the time recalibration.
Definition: AliEMCALRecoUtils.h:468
AliEMCALRecoUtils::fCutEtaPhiSum
Bool_t fCutEtaPhiSum
Place cut on sqrt(dEta^2+dPhi^2)
Definition: AliEMCALRecoUtils.h:507
AliEMCALRecoUtils::FindMatchedClusterInClusterArr
Int_t FindMatchedClusterInClusterArr(const AliExternalTrackParam *emcalParam, AliExternalTrackParam *trkParam, const TObjArray *clusterArr, Float_t &dEta, Float_t &dPhi)
Definition: AliEMCALRecoUtils.cxx:2842
AliEMCALRecoUtils::InitEMCALL1PhaseInTimeRecalibration
void InitEMCALL1PhaseInTimeRecalibration()
Definition: AliEMCALRecoUtils.cxx:1500
AliEMCALRecoUtils::fCutMinNClusterTPC
Int_t fCutMinNClusterTPC
Min number of tpc clusters.
Definition: AliEMCALRecoUtils.h:522
AliEMCALRecoUtils::fCutMaxDCAToVertexZ
Float_t fCutMaxDCAToVertexZ
Track-to-vertex cut in max absolute distance in z-plane.
Definition: AliEMCALRecoUtils.h:530
AliEMCALRecoUtils::SetRequireITSRefit
void SetRequireITSRefit(Bool_t b=kFALSE)
Definition: AliEMCALRecoUtils.h:368
energy
energy
Definition: HFPtSpectrum.C:44
AliEMCALRecoUtils::SetEMCALL1PhaseInTimeRecalibrationForSM
void SetEMCALL1PhaseInTimeRecalibrationForSM(Int_t iSM, Int_t c=0)
Definition: AliEMCALRecoUtils.h:228
AliEMCALRecoUtils::SetDCAToVertex2D
void SetDCAToVertex2D(Bool_t b=kFALSE)
Definition: AliEMCALRecoUtils.h:372
AliEMCALRecoUtils::SetEMCALChannelRecalibrationFactor
void SetEMCALChannelRecalibrationFactor(Int_t iSM, Int_t iCol, Int_t iRow, Double_t c=1)
Definition: AliEMCALRecoUtils.h:177
AliEMCALRecoUtils::fNonLinearThreshold
Int_t fNonLinearThreshold
Non linearity threshold value for kBeamTest non linearity function.
Definition: AliEMCALRecoUtils.h:454
AliEMCALRecoUtils::fITSTrackSA
Bool_t fITSTrackSA
If track matching is to be done with ITS tracks standing alone.
Definition: AliEMCALRecoUtils.h:516
AliEMCALRecoUtils::SetITSTrackSA
void SetITSTrackSA(Bool_t isITS)
Definition: AliEMCALRecoUtils.h:344
AliEMCALRecoUtils::IsExoticCluster
Bool_t IsExoticCluster(const AliVCluster *cluster, AliVCaloCells *cells, Int_t bc=0)
Definition: AliEMCALRecoUtils.cxx:686
AliEMCALRecoUtils::fMisalRotShift
Float_t fMisalRotShift[15]
Cluster position rotation shift parameters.
Definition: AliEMCALRecoUtils.h:441
AliEMCALRecoUtils::SetMinNClustersTPC
void SetMinNClustersTPC(Int_t min=-1)
Definition: AliEMCALRecoUtils.h:363
AliEMCALRecoUtils::fRecalDistToBadChannels
Bool_t fRecalDistToBadChannels
Calculate distance from highest energy tower of cluster to closes bad channel.
Definition: AliEMCALRecoUtils.h:481
AliEMCALRecoUtils::fRecalibration
Bool_t fRecalibration
Switch on or off the recalibration.
Definition: AliEMCALRecoUtils.h:463
AliEMCALRecoUtils::fMatchedTrackIndex
TArrayI * fMatchedTrackIndex
Array that stores indexes of matched tracks.
Definition: AliEMCALRecoUtils.h:503
AliEMCALRecoUtils::SetMaxDCAToVertexZ
void SetMaxDCAToVertexZ(Float_t dist=1e10)
Definition: AliEMCALRecoUtils.h:371
AliEMCALRecoUtils::IsBadChannelsRemovalSwitchedOn
Bool_t IsBadChannelsRemovalSwitchedOn() const
Definition: AliEMCALRecoUtils.h:252
AliEMCALRecoUtils::IsExoticCell
Bool_t IsExoticCell(Int_t absId, AliVCaloCells *cells, Int_t bc=-1)
Definition: AliEMCALRecoUtils.cxx:652
AliEMCALRecoUtils::RecalculateClusterPositionFromTowerGlobal
void RecalculateClusterPositionFromTowerGlobal(const AliEMCALGeometry *geom, AliVCaloCells *cells, AliVCluster *clu)
Definition: AliEMCALRecoUtils.cxx:1848
AliEMCALRecoUtils::fExoticCellDiffTime
Float_t fExoticCellDiffTime
If time of candidate to exotic and close cell is too different (in ns), it must be noisy...
Definition: AliEMCALRecoUtils.h:492
AliEMCALRecoUtils::fResidualPhi
TArrayF * fResidualPhi
Array that stores the residual phi.
Definition: AliEMCALRecoUtils.h:506
AliEMCALRecoUtils::RecalibrateCells
void RecalibrateCells(AliVCaloCells *cells, Int_t bc)
Definition: AliEMCALRecoUtils.cxx:1617
AliEMCALRecoUtils::fCutMinTrackPt
Double_t fCutMinTrackPt
Cut on track pT.
Definition: AliEMCALRecoUtils.h:521
AliEMCALRecoUtils::IsRecalibrationOn
Bool_t IsRecalibrationOn() const
Definition: AliEMCALRecoUtils.h:164
AliEMCALRecoUtils
Some utilities for cluster and cell treatment.
Definition: AliEMCALRecoUtils.h:59
AliEMCALRecoUtils::fCutMaxChi2PerClusterITS
Float_t fCutMaxChi2PerClusterITS
Max its fit chi2 per its cluster.
Definition: AliEMCALRecoUtils.h:525
AliEMCALRecoUtils::RecalculateClusterShowerShapeParametersWithCellCuts
void RecalculateClusterShowerShapeParametersWithCellCuts(const AliEMCALGeometry *geom, AliVCaloCells *cells, AliVCluster *cluster, Float_t cellEcut, Float_t cellTimeCut, Int_t bc, Float_t &enAfterCuts)
Definition: AliEMCALRecoUtils.cxx:2495
AliEMCALRecoUtils::fRemoveBadChannels
Bool_t fRemoveBadChannels
Check the channel status provided and remove clusters with bad channels.
Definition: AliEMCALRecoUtils.h:480
AliEMCALRecoUtils::IsAccepted
Bool_t IsAccepted(AliESDtrack *track)
Definition: AliEMCALRecoUtils.cxx:3138
AliEMCALRecoUtils::fMisalTransShift
Float_t fMisalTransShift[15]
Cluster position translation shift parameters.
Definition: AliEMCALRecoUtils.h:440
AliEMCALRecoUtils::SetAcceptKinkDaughters
void SetAcceptKinkDaughters(Bool_t b=kTRUE)
Definition: AliEMCALRecoUtils.h:369
AliEMCALRecoUtils::RecalculateClusterDistanceToBadChannel
void RecalculateClusterDistanceToBadChannel(const AliEMCALGeometry *geom, AliVCaloCells *cells, AliVCluster *cluster)
Definition: AliEMCALRecoUtils.cxx:2018
AliEMCALRecoUtils::fClusterWindow
Double_t fClusterWindow
Select clusters in the window to be matched.
Definition: AliEMCALRecoUtils.h:512
AliEMCALRecoUtils::fCutRequireITSRefit
Bool_t fCutRequireITSRefit
Require ITS refit.
Definition: AliEMCALRecoUtils.h:527
AliEMCALRecoUtils::fNonLinearityParams
Float_t fNonLinearityParams[10]
Parameters for the non linearity function.
Definition: AliEMCALRecoUtils.h:453
AliEMCALRecoUtils::GetCellWeight
Float_t GetCellWeight(Float_t eCell, Float_t eCluster) const
Definition: AliEMCALRecoUtils.cxx:1301
AliEMCALRecoUtils::fUseL1PhaseInTimeRecalibration
Bool_t fUseL1PhaseInTimeRecalibration
Switch on or off the L1 phase in time recalibration.
Definition: AliEMCALRecoUtils.h:473
AliEMCALRecoUtils::GetEMCALChannelRecalibrationFactor
Float_t GetEMCALChannelRecalibrationFactor(Int_t iSM, Int_t iCol, Int_t iRow) const
Definition: AliEMCALRecoUtils.h:173
AliEMCALRecoUtils::SetClusterMatchedToTrack
void SetClusterMatchedToTrack(const AliVEvent *event)
Definition: AliEMCALRecoUtils.cxx:3361
AliEMCALRecoUtils::fCutDCAToVertex2D
Bool_t fCutDCAToVertex2D
If true a 2D DCA cut is made.
Definition: AliEMCALRecoUtils.h:531
AliEMCALRecoUtils::GetMatchedTrackIndex
Int_t GetMatchedTrackIndex(Int_t clsIndex)
Definition: AliEMCALRecoUtils.cxx:2983
AliEMCALRecoUtils::RecalculateClusterPositionFromTowerIndex
void RecalculateClusterPositionFromTowerIndex(const AliEMCALGeometry *geom, AliVCaloCells *cells, AliVCluster *clu)
Definition: AliEMCALRecoUtils.cxx:1935
Int_t
int Int_t
Definition: External.C:63
AliEMCALRecoUtils::Print
void Print(const Option_t *) const
Definition: AliEMCALRecoUtils.cxx:3609
AliEMCALRecoUtils::FindMatchedPosForCluster
UInt_t FindMatchedPosForCluster(Int_t clsIndex) const
Definition: AliEMCALRecoUtils.cxx:3046
AliEMCALRecoUtils::AcceptCalibrateCell
Bool_t AcceptCalibrateCell(Int_t absId, Int_t bc, Float_t &amp, Double_t &time, AliVCaloCells *cells)
Definition: AliEMCALRecoUtils.cxx:388
UInt_t
unsigned int UInt_t
Definition: External.C:33
AliEMCALRecoUtils::RecalculateClusterShowerShapeParameters
void RecalculateClusterShowerShapeParameters(const AliEMCALGeometry *geom, AliVCaloCells *cells, AliVCluster *cluster)
Definition: AliEMCALRecoUtils.cxx:2465
Float_t
float Float_t
Definition: External.C:68
AliEMCALRecoUtils::fCutR
Float_t fCutR
sqrt(dEta^2+dPhi^2) cut on matching
Definition: AliEMCALRecoUtils.h:509
AliEMCALRecoUtils::fNCellsFromEMCALBorder
Int_t fNCellsFromEMCALBorder
Number of cells from EMCAL border the cell with maximum amplitude has to be.
Definition: AliEMCALRecoUtils.h:485
AliEMCALRecoUtils::fCellsRecalibrated
Bool_t fCellsRecalibrated
Internal bool to check if cells (time/energy) where recalibrated and not recalibrate them when recalc...
Definition: AliEMCALRecoUtils.h:462
AliEMCALRecoUtils::fRejectExoticCluster
Bool_t fRejectExoticCluster
Switch on or off exotic cluster rejection.
Definition: AliEMCALRecoUtils.h:489
AliEMCALRecoUtils::IsTimeRecalibrationOn
Bool_t IsTimeRecalibrationOn() const
Definition: AliEMCALRecoUtils.h:191
AliEMCALRecoUtils::fParticleType
Int_t fParticleType
Particle type for depth calculation, see enum ParticleType.
Definition: AliEMCALRecoUtils.h:442
AliEMCALRecoUtils::GetEMCALChannelStatus
Int_t GetEMCALChannelStatus(Int_t iSM, Int_t iCol, Int_t iRow) const
Definition: AliEMCALRecoUtils.h:262
TH2D
Definition: External.C:228
AliEMCALRecoUtils::SetMaxChi2PerClusterTPC
void SetMaxChi2PerClusterTPC(Float_t max=1e10)
Definition: AliEMCALRecoUtils.h:365
AliEMCALRecoUtils::fConstantTimeShift
Float_t fConstantTimeShift
Apply a 600 ns (+15.8) time shift in case of simulation, shift in ns.
Definition: AliEMCALRecoUtils.h:467
AliEMCALRecoUtils::GetMatchedClusterIndex
Int_t GetMatchedClusterIndex(Int_t trkIndex)
Definition: AliEMCALRecoUtils.cxx:2998
AliEMCALRecoUtils::fNonLinearityFunction
Int_t fNonLinearityFunction
Non linearity function choice, see enum NonlinearityFunctions.
Definition: AliEMCALRecoUtils.h:452
AliEMCALRecoUtils::fCutRequireITSStandAlone
Bool_t fCutRequireITSStandAlone
Require ITSStandAlone.
Definition: AliEMCALRecoUtils.h:532
AliEMCALRecoUtils::fEMCALL1PhaseInTimeRecalibration
TObjArray * fEMCALL1PhaseInTimeRecalibration
Histogram with map of L1 phase per SM, EMCAL.
Definition: AliEMCALRecoUtils.h:474
AliEMCALRecoUtils::fEMCALRecalibrationFactors
TObjArray * fEMCALRecalibrationFactors
Array of histograms with map of recalibration factors, EMCAL.
Definition: AliEMCALRecoUtils.h:464
AliEMCALRecoUtils::fUseRunCorrectionFactors
Bool_t fUseRunCorrectionFactors
Use Run Dependent Correction.
Definition: AliEMCALRecoUtils.h:477
AliEMCALRecoUtils::SetTracksMatchedToCluster
void SetTracksMatchedToCluster(const AliVEvent *event)
Definition: AliEMCALRecoUtils.cxx:3404
AliEMCALRecoUtils::fCutMaxChi2PerClusterTPC
Float_t fCutMaxChi2PerClusterTPC
Max tpc fit chi2 per tpc cluster.
Definition: AliEMCALRecoUtils.h:524
AliEMCALRecoUtils::fExoticCellMinAmplitude
Float_t fExoticCellMinAmplitude
Check for exotic only if amplitud is larger than this value.
Definition: AliEMCALRecoUtils.h:493
AliEMCALRecoUtils::RecalibrateClusterEnergy
void RecalibrateClusterEnergy(const AliEMCALGeometry *geom, AliVCluster *cluster, AliVCaloCells *cells, Int_t bc=-1)
Definition: AliEMCALRecoUtils.cxx:1532
AliEMCALRecoUtils::fW0
Float_t fW0
Energy weight used in cluster position and shower shape calculations.
Definition: AliEMCALRecoUtils.h:445
AliEMCALRecoUtils::fAODHybridTracks
Bool_t fAODHybridTracks
Match with hybrid.
Definition: AliEMCALRecoUtils.h:500
AliEMCALRecoUtils::GetMatchedResiduals
void GetMatchedResiduals(Int_t clsIndex, Float_t &dEta, Float_t &dPhi)
Definition: AliEMCALRecoUtils.cxx:2937
AliEMCALRecoUtils::operator=
AliEMCALRecoUtils & operator=(const AliEMCALRecoUtils &)
Definition: AliEMCALRecoUtils.cxx:169
AliEMCALRecoUtils::SetRequireTPCRefit
void SetRequireTPCRefit(Bool_t b=kFALSE)
Definition: AliEMCALRecoUtils.h:367
AliEMCALRecoUtils::fResidualEta
TArrayF * fResidualEta
Array that stores the residual eta.
Definition: AliEMCALRecoUtils.h:505
Short_t
short Short_t
Definition: External.C:23
AliEMCALRecoUtils::fMCGenerToAccept
TString fMCGenerToAccept[5]
List with name of generators that should not be included.
Definition: AliEMCALRecoUtils.h:537
AliEMCALRecoUtils::fSmearClusterEnergy
Bool_t fSmearClusterEnergy
Smear cluster energy, to be done only for simulated data to match real data.
Definition: AliEMCALRecoUtils.h:457
AliEMCALRecoUtils::IsGoodCluster
Bool_t IsGoodCluster(AliVCluster *cluster, const AliEMCALGeometry *geom, AliVCaloCells *cells, Int_t bc=-1)
Definition: AliEMCALRecoUtils.cxx:3113
AliEMCALRecoUtils::GetMaxEnergyCell
void GetMaxEnergyCell(const AliEMCALGeometry *geom, AliVCaloCells *cells, const AliVCluster *clu, Int_t &absId, Int_t &iSupMod, Int_t &ieta, Int_t &iphi, Bool_t &shared)
Definition: AliEMCALRecoUtils.cxx:1225
AliEMCALRecoUtils::InitEMCALTimeRecalibrationFactors
void InitEMCALTimeRecalibrationFactors()
Definition: AliEMCALRecoUtils.cxx:1430
AliEMCALRecoUtils::fStepSurface
Double_t fStepSurface
Length of step to extrapolate tracks to EMCal surface.
Definition: AliEMCALRecoUtils.h:514
AliEMCALRecoUtils::fAODFilterMask
UInt_t fAODFilterMask
Filter mask to select AOD tracks. Refer to $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C.
Definition: AliEMCALRecoUtils.h:499
AliEMCALRecoUtils::ExtrapolateTrackToCluster
static Bool_t ExtrapolateTrackToCluster(AliExternalTrackParam *trkParam, const AliVCluster *cluster, Double_t mass, Double_t step, Float_t &tmpEta, Float_t &tmpPhi)
Definition: AliEMCALRecoUtils.h:306
AliEMCALRecoUtils::GetMatchedClusterResiduals
void GetMatchedClusterResiduals(Int_t trkIndex, Float_t &dEta, Float_t &dPhi)
Definition: AliEMCALRecoUtils.cxx:2962
AliEMCALRecoUtils::fSmearClusterParam
Float_t fSmearClusterParam[3]
Energy smearing parameters.
Definition: AliEMCALRecoUtils.h:458
AliEMCALRecoUtils::fCutAcceptKinkDaughters
Bool_t fCutAcceptKinkDaughters
Accepting kink daughters?
Definition: AliEMCALRecoUtils.h:528
AliEMCALRecoUtils::RecalculateClusterPosition
void RecalculateClusterPosition(const AliEMCALGeometry *geom, AliVCaloCells *cells, AliVCluster *clu)
Definition: AliEMCALRecoUtils.cxx:1824
AliEMCALRecoUtils::fExoticCellFraction
Float_t fExoticCellFraction
Good cell if fraction < 1-ecross/ecell.
Definition: AliEMCALRecoUtils.h:491
AliEMCALRecoUtils::RecalculateClusterPID
void RecalculateClusterPID(AliVCluster *cluster)
Definition: AliEMCALRecoUtils.cxx:2106
AliEMCALRecoUtils::fPosAlgo
Int_t fPosAlgo
Position recalculation algorithm, see enum PositionAlgorithms.
Definition: AliEMCALRecoUtils.h:443
AliEMCALRecoUtils::SetEMCALL1PhaseInTimeRecalibrationForAllSM
void SetEMCALL1PhaseInTimeRecalibrationForAllSM(const TObjArray *map)
Definition: AliEMCALRecoUtils.cxx:3577
AliEMCALRecoUtils::fPIDUtils
AliEMCALPIDUtils * fPIDUtils
Recalculate PID parameters.
Definition: AliEMCALRecoUtils.h:496
AliEMCALRecoUtils::fRejectExoticCells
Bool_t fRejectExoticCells
Remove exotic cells.
Definition: AliEMCALRecoUtils.h:490
AliEMCALRecoUtils::FindMatchedClusterInEvent
Int_t FindMatchedClusterInEvent(const AliESDtrack *track, const AliVEvent *event, const AliEMCALGeometry *geom, Float_t &dEta, Float_t &dPhi)
Definition: AliEMCALRecoUtils.cxx:2768
AliEMCALRecoUtils::IsTrackMatched
Bool_t IsTrackMatched(Int_t trkIndex) const
Definition: AliEMCALRecoUtils.cxx:3030
AliEMCALRecoUtils::RecalibrateCellTimeL1Phase
void RecalibrateCellTimeL1Phase(Int_t iSM, Int_t bc, Double_t &time) const
Definition: AliEMCALRecoUtils.cxx:1683
AliEMCALRecoUtils::fCutMaxDCAToVertexXY
Float_t fCutMaxDCAToVertexXY
Track-to-vertex cut in max absolute distance in xy-plane.
Definition: AliEMCALRecoUtils.h:529
AliEMCALRecoUtils::FindMatches
void FindMatches(AliVEvent *event, TObjArray *clusterArr=0x0, const AliEMCALGeometry *geom=0x0, AliMCEvent *mc=0x0)
Definition: AliEMCALRecoUtils.cxx:2529
UShort_t
unsigned short UShort_t
Definition: External.C:28
AliEMCALRecoUtils::fEMCALBadChannelMap
TObjArray * fEMCALBadChannelMap
Array of histograms with map of bad channels, EMCAL.
Definition: AliEMCALRecoUtils.h:482
AliEMCALRecoUtils::SetEMCALChannelTimeRecalibrationFactors
void SetEMCALChannelTimeRecalibrationFactors(const TObjArray *map)
Definition: AliEMCALRecoUtils.cxx:3548
AliEMCALRecoUtils::fStepCluster
Double_t fStepCluster
Length of step to extrapolate tracks to clusters.
Definition: AliEMCALRecoUtils.h:515
AliEMCALRecoUtils::fNoEMCALBorderAtEta0
Bool_t fNoEMCALBorderAtEta0
Do fiducial cut in EMCAL region eta = 0?
Definition: AliEMCALRecoUtils.h:486
Option_t
const char Option_t
Definition: External.C:48
AliEMCALRecoUtils::fNMCGenerToAccept
Int_t fNMCGenerToAccept
Number of MC generators that should not be included in analysis.
Definition: AliEMCALRecoUtils.h:536
AliEMCALRecoUtils::SetEMCALChannelTimeRecalibrationFactor
void SetEMCALChannelTimeRecalibrationFactor(Int_t bc, Int_t absID, Double_t c=0, Bool_t isLGon=kFALSE)
Definition: AliEMCALRecoUtils.h:202
AliEMCALRecoUtils::fCutRequireITSpureSA
Bool_t fCutRequireITSpureSA
ITS pure standalone tracks.
Definition: AliEMCALRecoUtils.h:533
AliEMCALRecoUtils::GetEMCALL1PhaseInTimeRecalibrationForSM
Int_t GetEMCALL1PhaseInTimeRecalibrationForSM(Int_t iSM) const
Definition: AliEMCALRecoUtils.h:224
Bool_t
bool Bool_t
Definition: External.C:53
AliEMCALRecoUtils::fCutMinNClusterITS
Int_t fCutMinNClusterITS
Min number of its clusters.
Definition: AliEMCALRecoUtils.h:523
AliEMCALRecoUtils::CorrectClusterEnergyLinearity
Float_t CorrectClusterEnergyLinearity(AliVCluster *clu)
Definition: AliEMCALRecoUtils.cxx:754
AliEMCALRecoUtils::FindMatchedPosForTrack
UInt_t FindMatchedPosForTrack(Int_t trkIndex) const
Definition: AliEMCALRecoUtils.cxx:3076
AliEMCALRecoUtils::ClusterContainsBadChannel
Bool_t ClusterContainsBadChannel(const AliEMCALGeometry *geom, const UShort_t *cellList, Int_t nCells)
Definition: AliEMCALRecoUtils.cxx:535
AliEMCALRecoUtils::fEMCalSurfaceDistance
Double_t fEMCalSurfaceDistance
EMCal surface distance (= 430 by default, the last 10 cm are propagated on a cluster-track pair basis...
Definition: AliEMCALRecoUtils.h:517
AliEMCALRecoUtils::SetEMCALChannelRecalibrationFactors
void SetEMCALChannelRecalibrationFactors(const TObjArray *map)
Definition: AliEMCALRecoUtils.cxx:3490
AliEMCALRecoUtils::GetECross
Float_t GetECross(Int_t absID, Double_t tcell, AliVCaloCells *cells, Int_t bc)
Definition: AliEMCALRecoUtils.cxx:577
AliEMCALRecoUtils::fLowGain
Bool_t fLowGain
Switch on or off calibration with low gain channels.
Definition: AliEMCALRecoUtils.h:470
AliEMCALRecoUtils::SetEMCALChannelStatusMap
void SetEMCALChannelStatusMap(const TObjArray *map)
Definition: AliEMCALRecoUtils.cxx:3519
AliEMCALRecoUtils::fCutPhi
Float_t fCutPhi
dPhi cut on matching
Definition: AliEMCALRecoUtils.h:511
AliEMCALRecoUtils::fTrackCutsType
Int_t fTrackCutsType
ESD track cuts type for matching, see enum TrackCutsType.
Definition: AliEMCALRecoUtils.h:520
AliEMCALRecoUtils::RecalibrateCellTime
void RecalibrateCellTime(Int_t absId, Int_t bc, Double_t &time, Bool_t isLGon=kFALSE) const
Definition: AliEMCALRecoUtils.cxx:1665
AliEMCALRecoUtils::fCutEtaPhiSeparate
Bool_t fCutEtaPhiSeparate
Cut on dEta and dPhi separately.
Definition: AliEMCALRecoUtils.h:508
AliEMCALRecoUtils::fCutEta
Float_t fCutEta
dEta cut on matching
Definition: AliEMCALRecoUtils.h:510
AliEMCALRecoUtils::SmearClusterEnergy
Float_t SmearClusterEnergy(const AliVCluster *clu)
Definition: AliEMCALRecoUtils.cxx:723