]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EMCAL/AliEMCALRecoUtils.h
use eta-phi cuts instead of R-z cuts for track matching, add track momentum cut ...
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALRecoUtils.h
CommitLineData
d9b3567c 1#ifndef ALIEMCALRECOUTILS_H
2#define ALIEMCALRECOUTILS_H
3
4/* $Id: AliEMCALRecoUtils.h 33808 2009-07-15 09:48:08Z gconesab $ */
5
6///////////////////////////////////////////////////////////////////////////////
7//
8// Class AliEMCALRecoUtils
9// Some utilities to recalculate the cluster position or energy linearity
10//
11//
12// Author: Gustavo Conesa (LPSC- Grenoble)
b540d03f 13// Track matching part: Rongrong Ma (Yale)
d9b3567c 14///////////////////////////////////////////////////////////////////////////////
15
16//Root includes
17#include "TNamed.h"
094786cc 18#include "TMath.h"
19#include "TObjArray.h"
bd8c7aef 20#include "TArrayI.h"
21#include "TArrayF.h"
17688f67 22#include "TH2F.h"
d9b3567c 23
24//AliRoot includes
25class AliVCluster;
26class AliVCaloCells;
bd8c7aef 27class AliVEvent;
d9b3567c 28#include "AliLog.h"
b540d03f 29
30// EMCAL includes
094786cc 31class AliEMCALGeometry;
83bfd77a 32class AliEMCALPIDUtils;
bd8c7aef 33class AliESDtrack;
bb6f5f0b 34class AliExternalTrackParam;
d9b3567c 35
36class AliEMCALRecoUtils : public TNamed {
37
38public:
39
40 AliEMCALRecoUtils();
41 AliEMCALRecoUtils(const AliEMCALRecoUtils&);
42 AliEMCALRecoUtils& operator=(const AliEMCALRecoUtils&);
b540d03f 43 virtual ~AliEMCALRecoUtils() ;
44 void Print(const Option_t*) const;
45
46 //enums
4b58ac4f 47 enum NonlinearityFunctions{kPi0MC=0,kPi0GammaGamma=1,kPi0GammaConversion=2,kNoCorrection=3,kBeamTest=4,kBeamTestCorrected=5};
fd6df01c 48 enum PositionAlgorithms{kUnchanged=-1,kPosTowerIndex=0, kPosTowerGlobal=1};
094786cc 49 enum ParticleType{kPhoton=0, kElectron=1,kHadron =2, kUnknown=-1};
b540d03f 50 enum { kNCuts = 11 }; //track matching
51
52 //-----------------------------------------------------
d9b3567c 53 //Position recalculation
b540d03f 54 //-----------------------------------------------------
55
094786cc 56 void RecalculateClusterPosition(AliEMCALGeometry *geom, AliVCaloCells* cells, AliVCluster* clu);
57 void RecalculateClusterPositionFromTowerIndex (AliEMCALGeometry *geom, AliVCaloCells* cells, AliVCluster* clu);
58 void RecalculateClusterPositionFromTowerGlobal(AliEMCALGeometry *geom, AliVCaloCells* cells, AliVCluster* clu);
59
60 Float_t GetCellWeight(const Float_t eCell, const Float_t eCluster) const { return TMath::Max( 0., fW0 + TMath::Log( eCell / eCluster ));}
61
62 Float_t GetDepth(const Float_t eCluster, const Int_t iParticle, const Int_t iSM) const ;
63
64 void GetMaxEnergyCell(AliEMCALGeometry *geom, AliVCaloCells* cells, AliVCluster* clu,
cb231979 65 Int_t & absId, Int_t& iSupMod, Int_t& ieta, Int_t& iphi, Bool_t &shared);
d9b3567c 66
2a71e873 67 Float_t GetMisalTransShift(const Int_t i) const {
68 if(i < 15 ){return fMisalTransShift[i]; }
d9b3567c 69 else { AliInfo(Form("Index %d larger than 15, do nothing\n",i)); return 0.;}
70 }
094786cc 71 Float_t* GetMisalTransShiftArray() {return fMisalTransShift; }
d9b3567c 72
2a71e873 73 void SetMisalTransShift(const Int_t i, const Float_t shift) {
74 if(i < 15 ){fMisalTransShift[i] = shift; }
d9b3567c 75 else { AliInfo(Form("Index %d larger than 15, do nothing\n",i));}
76 }
2a71e873 77 void SetMisalTransShiftArray(Float_t * misal)
78 { for(Int_t i = 0; i < 15; i++)fMisalTransShift[i] = misal[i]; }
d9b3567c 79
2a71e873 80 Float_t GetMisalRotShift(const Int_t i) const {
81 if(i < 15 ){return fMisalRotShift[i]; }
82 else { AliInfo(Form("Index %d larger than 15, do nothing\n",i)); return 0.;}
83 }
094786cc 84 Float_t* GetMisalRotShiftArray() {return fMisalRotShift; }
2a71e873 85
86 void SetMisalRotShift(const Int_t i, const Float_t shift) {
87 if(i < 15 ){fMisalRotShift[i] = shift; }
88 else { AliInfo(Form("Index %d larger than 15, do nothing\n",i));}
89 }
90 void SetMisalRotShiftArray(Float_t * misal)
91 { for(Int_t i = 0; i < 15; i++)fMisalRotShift[i] = misal[i]; }
92
96957075 93 Int_t GetParticleType() const { return fParticleType ;}
94 void SetParticleType(Int_t particle) { fParticleType = particle ;}
2a71e873 95
96957075 96 Int_t GetPositionAlgorithm() const { return fPosAlgo ;}
97 void SetPositionAlgorithm(Int_t alg) { fPosAlgo = alg ;}
2a71e873 98
96957075 99 Float_t GetW0() const { return fW0 ;}
100 void SetW0(Float_t w0) { fW0 = w0 ;}
094786cc 101
b540d03f 102 //-----------------------------------------------------
d9b3567c 103 //Non Linearity
b540d03f 104 //-----------------------------------------------------
105
d9b3567c 106 Float_t CorrectClusterEnergyLinearity(AliVCluster* clu);
107
108 Float_t GetNonLinearityParam(const Int_t i) const {
dff9e2e3 109 if(i < 7 ){return fNonLinearityParams[i]; }
110 else { AliInfo(Form("Index %d larger than 7, do nothing\n",i)); return 0.;}
d9b3567c 111 }
112 void SetNonLinearityParam(const Int_t i, const Float_t param) {
dff9e2e3 113 if(i < 7 ){fNonLinearityParams[i] = param; }
114 else { AliInfo(Form("Index %d larger than 7, do nothing\n",i));}
d9b3567c 115 }
7e0ecb89 116 void InitNonLinearityParam();
117
118 Int_t GetNonLinearityFunction() const { return fNonLinearityFunction ; }
119 void SetNonLinearityFunction(Int_t fun) { fNonLinearityFunction = fun ; InitNonLinearityParam() ; }
120
121 void SetNonLinearityThreshold(Int_t threshold) {fNonLinearThreshold = threshold ;} //only for Alexie's non linearity correction
122 Int_t GetNonLinearityThreshold() const {return fNonLinearThreshold ;}
123
b540d03f 124
125 //-----------------------------------------------------
094786cc 126 //Recalibration
b540d03f 127 //-----------------------------------------------------
128
094786cc 129 void RecalibrateClusterEnergy(AliEMCALGeometry* geom, AliVCluster* cluster, AliVCaloCells * cells);
130
b540d03f 131 Bool_t IsRecalibrationOn() const { return fRecalibration ; }
132 void SwitchOnRecalibration() { fRecalibration = kTRUE ; if(!fEMCALRecalibrationFactors)InitEMCALRecalibrationFactors();}
133 void SwitchOffRecalibration() { fRecalibration = kFALSE ; }
134 void InitEMCALRecalibrationFactors() ;
96957075 135
136 //Recalibrate channels with time dependent corrections
b540d03f 137 void SwitchOnTimeDepCorrection() { fUseTimeCorrectionFactors = kTRUE ; SwitchOnRecalibration();}
138 void SwitchOffTimeDepCorrection() { fUseTimeCorrectionFactors = kFALSE;}
96957075 139 void SetTimeDependentCorrections(Int_t runnumber);
140
094786cc 141 Float_t GetEMCALChannelRecalibrationFactor(Int_t iSM , Int_t iCol, Int_t iRow) const {
142 if(fEMCALRecalibrationFactors) return (Float_t) ((TH2F*)fEMCALRecalibrationFactors->At(iSM))->GetBinContent(iCol,iRow);
143 else return 1;}
144
145 void SetEMCALChannelRecalibrationFactor(Int_t iSM , Int_t iCol, Int_t iRow, Double_t c = 1) {
146 if(!fEMCALRecalibrationFactors) InitEMCALRecalibrationFactors();
147 ((TH2F*)fEMCALRecalibrationFactors->At(iSM))->SetBinContent(iCol,iRow,c);}
148
96957075 149 TH2F * GetEMCALChannelRecalibrationFactors(Int_t iSM) const { return (TH2F*)fEMCALRecalibrationFactors->At(iSM) ;}
150 void SetEMCALChannelRecalibrationFactors(TObjArray *map) { fEMCALRecalibrationFactors = map ;}
151 void SetEMCALChannelRecalibrationFactors(Int_t iSM , TH2F* h) { fEMCALRecalibrationFactors->AddAt(h,iSM) ;}
094786cc 152
b540d03f 153 //-----------------------------------------------------
fd6df01c 154 //Modules fiducial region, remove clusters in borders
b540d03f 155 //-----------------------------------------------------
156
fd6df01c 157 Bool_t CheckCellFiducialRegion(AliEMCALGeometry* geom, AliVCluster* cluster, AliVCaloCells* cells) ;
96957075 158 void SetNumberOfCellsFromEMCALBorder(Int_t n) { fNCellsFromEMCALBorder = n ;}
159 Int_t GetNumberOfCellsFromEMCALBorder() const { return fNCellsFromEMCALBorder ;}
fd6df01c 160
96957075 161 void SwitchOnNoFiducialBorderInEMCALEta0() { fNoEMCALBorderAtEta0 = kTRUE ;}
162 void SwitchOffNoFiducialBorderInEMCALEta0() { fNoEMCALBorderAtEta0 = kFALSE ;}
163 Bool_t IsEMCALNoBorderAtEta0() { return fNoEMCALBorderAtEta0 ;}
fd6df01c 164
b540d03f 165 //-----------------------------------------------------
fd6df01c 166 // Bad channels
b540d03f 167 //-----------------------------------------------------
168
169 Bool_t IsBadChannelsRemovalSwitchedOn() const { return fRemoveBadChannels ;}
170 void SwitchOnBadChannelsRemoval () { fRemoveBadChannels = kTRUE ; if(!fEMCALBadChannelMap)InitEMCALBadChannelStatusMap();}
171 void SwitchOffBadChannelsRemoval() { fRemoveBadChannels = kFALSE ;}
fd6df01c 172
b540d03f 173 Bool_t IsDistanceToBadChannelRecalculated() const { return fRecalDistToBadChannels ;}
174 void SwitchOnDistToBadChannelRecalculation() { fRecalDistToBadChannels = kTRUE ; if(!fEMCALBadChannelMap)InitEMCALBadChannelStatusMap();}
175 void SwitchOffDistToBadChannelRecalculation() { fRecalDistToBadChannels = kFALSE ;}
78467229 176
fd6df01c 177 void InitEMCALBadChannelStatusMap() ;
178
179 Int_t GetEMCALChannelStatus(Int_t iSM , Int_t iCol, Int_t iRow) const {
180 if(fEMCALBadChannelMap) return (Int_t) ((TH2I*)fEMCALBadChannelMap->At(iSM))->GetBinContent(iCol,iRow);
181 else return 0;}//Channel is ok by default
182
183 void SetEMCALChannelStatus(Int_t iSM , Int_t iCol, Int_t iRow, Double_t c = 1) {
184 if(!fEMCALBadChannelMap)InitEMCALBadChannelStatusMap() ;
185 ((TH2I*)fEMCALBadChannelMap->At(iSM))->SetBinContent(iCol,iRow,c);}
186
187 TH2I * GetEMCALChannelStatusMap(Int_t iSM) const {return (TH2I*)fEMCALBadChannelMap->At(iSM);}
188 void SetEMCALChannelStatusMap(TObjArray *map) {fEMCALBadChannelMap = map;}
6fe0e6d0 189 void SetEMCALChannelStatusMap(Int_t iSM , TH2I* h) {fEMCALBadChannelMap->AddAt(h,iSM);}
190
fd6df01c 191 Bool_t ClusterContainsBadChannel(AliEMCALGeometry* geom, UShort_t* cellList, Int_t nCells);
192
b540d03f 193 //-----------------------------------------------------
194 // Recalculate other cluster parameters
195 //-----------------------------------------------------
196
cb231979 197 void RecalculateClusterDistanceToBadChannel(AliEMCALGeometry * geom, AliVCaloCells* cells, AliVCluster * cluster);
83bfd77a 198 void RecalculateClusterPID(AliVCluster * cluster);
cb231979 199
83bfd77a 200 AliEMCALPIDUtils * GetPIDUtils() { return fPIDUtils;}
201
202 void RecalculateClusterShowerShapeParameters(AliEMCALGeometry * geom, AliVCaloCells* cells, AliVCluster * cluster);
203
b540d03f 204 //----------------------------------------------------
205 // Track matching
206 //----------------------------------------------------
bd8c7aef 207
9741c6a0 208 void FindMatches(AliVEvent *event, TObjArray * clusterArr=0x0, AliEMCALGeometry *geom=0x0);
209 Int_t FindMatchedCluster(AliESDtrack *track, AliVEvent *event, AliEMCALGeometry *geom);
fa4287a2 210 Bool_t ExtrapolateTrackToCluster(AliExternalTrackParam *trkParam, AliVCluster *cluster, Float_t &tmpEta, Float_t &tmpPhi);
211 void GetMatchedResiduals(Int_t clsIndex, Float_t &dEta, Float_t &dPhi);
212 void GetMatchedClusterResiduals(Int_t trkIndex, Float_t &dEta, Float_t &dPhi);
bb6f5f0b 213 Int_t GetMatchedTrackIndex(Int_t clsIndex);
214 Int_t GetMatchedClusterIndex(Int_t trkIndex);
215 Bool_t IsClusterMatched(Int_t clsIndex);
216 Bool_t IsTrackMatched(Int_t trkIndex);
217 UInt_t FindMatchedPosForCluster(Int_t clsIndex) const;
218 UInt_t FindMatchedPosForTrack(Int_t trkIndex) const;
219
fa4287a2 220 void SwitchOnCutEtaPhiSum() {fCutEtaPhiSum=kTRUE;fCutEtaPhiSeparate=kFALSE;}
221 void SwitchOnCutEtaPhiSeparate() {fCutEtaPhiSeparate=kTRUE;fCutEtaPhiSum=kFALSE;}
222
223 Float_t GetCutR() const { return fCutR ;}
224 Float_t GetCutEta() const { return fCutEta ;}
225 Float_t GetCutPhi() const { return fCutPhi ;}
226 void SetCutR(Float_t cutR) { fCutR=cutR ;}
227 void SetCutEta(Float_t cutEta) { fCutEta=cutEta ;}
228 void SetCutPhi(Float_t cutPhi) { fCutPhi=cutPhi ;}
229 void SetCutZ(Float_t cutZ) {printf("Obsolete fucntion of cutZ=%1.1f\n",cutZ);} //Obsolete
bb6f5f0b 230
231 Double_t GetMass() const { return fMass ;}
232 Double_t GetStep() const { return fStep ;}
233 void SetMass(Double_t mass){ fMass=mass ;}
234 void SetStep(Double_t step){ fStep=step ;}
235
9741c6a0 236 //Cluster cut
237 Bool_t IsGoodCluster(AliVCluster *cluster, AliEMCALGeometry *geom, AliVCaloCells* cells);
238 Bool_t IsExoticCluster(AliVCluster *cluster);
239
240 void SwitchOnRejectExoticCluster() { fRejectExoticCluster=kTRUE ;}
241 void SwitchOffRejectExoticCluster() { fRejectExoticCluster=kFALSE ;}
242 Bool_t IsRejectExoticCluster() { return fRejectExoticCluster ;}
243
bd8c7aef 244
245 //Track Cuts
b540d03f 246 Bool_t IsAccepted(AliESDtrack *track);
247 void InitTrackCuts();
bd8c7aef 248
249 // track quality cut setters
fa4287a2 250 void SetMinTrackPt(Double_t pt=0) { fCutMinTrackPt = pt ;}
b540d03f 251 void SetMinNClustersTPC(Int_t min=-1) { fCutMinNClusterTPC = min ;}
252 void SetMinNClustersITS(Int_t min=-1) { fCutMinNClusterITS = min ;}
253 void SetMaxChi2PerClusterTPC(Float_t max=1e10) { fCutMaxChi2PerClusterTPC = max ;}
254 void SetMaxChi2PerClusterITS(Float_t max=1e10) { fCutMaxChi2PerClusterITS = max ;}
255 void SetRequireTPCRefit(Bool_t b=kFALSE) { fCutRequireTPCRefit = b ;}
256 void SetRequireITSRefit(Bool_t b=kFALSE) { fCutRequireITSRefit = b ;}
257 void SetAcceptKinkDaughters(Bool_t b=kTRUE) { fCutAcceptKinkDaughters = b ;}
258 void SetMaxDCAToVertexXY(Float_t dist=1e10) { fCutMaxDCAToVertexXY = dist ;}
259 void SetMaxDCAToVertexZ(Float_t dist=1e10) { fCutMaxDCAToVertexZ = dist ;}
260 void SetDCAToVertex2D(Bool_t b=kFALSE) { fCutDCAToVertex2D = b ;}
bd8c7aef 261
fa4287a2 262 // getters
263 Double_t GetMinTrackPt() const { return fCutMinTrackPt ;}
b540d03f 264 Int_t GetMinNClusterTPC() const { return fCutMinNClusterTPC ;}
265 Int_t GetMinNClustersITS() const { return fCutMinNClusterITS ;}
266 Float_t GetMaxChi2PerClusterTPC() const { return fCutMaxChi2PerClusterTPC ;}
267 Float_t GetMaxChi2PerClusterITS() const { return fCutMaxChi2PerClusterITS ;}
268 Bool_t GetRequireTPCRefit() const { return fCutRequireTPCRefit ;}
269 Bool_t GetRequireITSRefit() const { return fCutRequireITSRefit ;}
270 Bool_t GetAcceptKinkDaughters() const { return fCutAcceptKinkDaughters ;}
271 Float_t GetMaxDCAToVertexXY() const { return fCutMaxDCAToVertexXY ;}
272 Float_t GetMaxDCAToVertexZ() const { return fCutMaxDCAToVertexZ ;}
273 Bool_t GetDCAToVertex2D() const { return fCutDCAToVertex2D ;}
bd8c7aef 274
fd6df01c 275
d9b3567c 276private:
277
b540d03f 278 //Position recalculation
96957075 279 Float_t fMisalTransShift[15]; // Shift parameters
280 Float_t fMisalRotShift[15]; // Shift parameters
281 Int_t fNonLinearityFunction; // Non linearity function choice
dff9e2e3 282 Float_t fNonLinearityParams[7]; // Parameters for the non linearity function
96957075 283 Int_t fParticleType; // Particle type for depth calculation
284 Int_t fPosAlgo; // Position recalculation algorithm
285 Float_t fW0; // Weight0
7e0ecb89 286 Int_t fNonLinearThreshold; // Non linearity threshold value for kBeamTesh non linearity function
fd6df01c 287
b540d03f 288 // Recalibration
fd6df01c 289 Bool_t fRecalibration; // Switch on or off the recalibration
290 TObjArray* fEMCALRecalibrationFactors; // Array of histograms with map of recalibration factors, EMCAL
b540d03f 291
292 // Bad Channels
fd6df01c 293 Bool_t fRemoveBadChannels; // Check the channel status provided and remove clusters with bad channels
78467229 294 Bool_t fRecalDistToBadChannels; // Calculate distance from highest energy tower of cluster to closes bad channel
fd6df01c 295 TObjArray* fEMCALBadChannelMap; // Array of histograms with map of bad channels, EMCAL
b540d03f 296
297 // Border cells
fd6df01c 298 Int_t fNCellsFromEMCALBorder; // Number of cells from EMCAL border the cell with maximum amplitude has to be.
299 Bool_t fNoEMCALBorderAtEta0; // Do fiducial cut in EMCAL region eta = 0?
b540d03f 300
bb6f5f0b 301 //Track matching
302 UInt_t fAODFilterMask; // Filter mask to select AOD tracks. Refer to $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C
b540d03f 303 TArrayI * fMatchedTrackIndex; // Array that stores indexes of matched tracks
96957075 304 TArrayI * fMatchedClusterIndex; // Array that stores indexes of matched clusters
fa4287a2 305 TArrayF * fResidualEta; // Array that stores the residual eta
306 TArrayF * fResidualPhi; // Array that stores the residual phi
307 Bool_t fCutEtaPhiSum; // Place cut on sqrt(dEta^2+dPhi^2)
308 Bool_t fCutEtaPhiSeparate; // Cut on dEta and dPhi separately
309 Float_t fCutR; // sqrt(dEta^2+dPhi^2) cut on matching
310 Float_t fCutEta; // dEta cut on matching
311 Float_t fCutPhi; // dPhi cut on matching
bb6f5f0b 312 Double_t fMass; // Mass hypothesis of the track
313 Double_t fStep; // Length of each step used in extrapolation in the unit of cm.
9741c6a0 314
315 // Cluster cuts
316 Bool_t fRejectExoticCluster; // Switch on or off exotic cluster rejection
317
318 // Track cuts
fa4287a2 319 Double_t fCutMinTrackPt; // Cut on track pT
96957075 320 Int_t fCutMinNClusterTPC; // Min number of tpc clusters
321 Int_t fCutMinNClusterITS; // Min number of its clusters
322 Float_t fCutMaxChi2PerClusterTPC; // Max tpc fit chi2 per tpc cluster
323 Float_t fCutMaxChi2PerClusterITS; // Max its fit chi2 per its cluster
324 Bool_t fCutRequireTPCRefit; // Require TPC refit
325 Bool_t fCutRequireITSRefit; // Require ITS refit
326 Bool_t fCutAcceptKinkDaughters; // Accepting kink daughters?
327 Float_t fCutMaxDCAToVertexXY; // Track-to-vertex cut in max absolute distance in xy-plane
328 Float_t fCutMaxDCAToVertexZ; // Track-to-vertex cut in max absolute distance in z-plane
329 Bool_t fCutDCAToVertex2D; // If true a 2D DCA cut is made. Tracks are accepted if sqrt((DCAXY / fCutMaxDCAToVertexXY)^2 + (DCAZ / fCutMaxDCAToVertexZ)^2) < 1 AND sqrt((DCAXY / fCutMinDCAToVertexXY)^2 + (DCAZ / fCutMinDCAToVertexZ)^2) > 1
bd8c7aef 330
b540d03f 331 //PID
96957075 332 AliEMCALPIDUtils * fPIDUtils; // Recalculate PID parameters
333
334 //Time Correction
335 Bool_t fUseTimeCorrectionFactors; // Use Time Dependent Correction
336 Bool_t fTimeCorrectionFactorsSet; // Time Correction set at leat once
83bfd77a 337
fa4287a2 338 ClassDef(AliEMCALRecoUtils, 12)
d9b3567c 339
340};
341
342#endif // ALIEMCALRECOUTILS_H
343
344