1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
|
#ifndef DATAFORMATS_CALOTOWERS_CALOTOWER_H
#define DATAFORMATS_CALOTOWERS_CALOTOWER_H 1
#include "DataFormats/Candidate/interface/LeafCandidate.h"
#include "DataFormats/DetId/interface/DetId.h"
#include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
#include "DataFormats/Math/interface/Vector3D.h"
#include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
#include "Rtypes.h"
#include <vector>
#include <cmath>
#include "DataFormats/GeometryVector/interface/GlobalPoint.h"
/** \class CaloTower
\author J. Mans - Minnesota
*/
// Original author: J. Mans - Minnesota
//
// Modified: Anton Anastassov (Northwestern)
// Make CaloTower inherit from LeafCandidate,
// add new members and accessors.
class CaloTower : public reco::LeafCandidate {
public:
typedef CaloTowerDetId key_type; // for SortedCollection
// Default constructor
CaloTower();
// Constructors from values
CaloTower(const CaloTowerDetId& id,
double emE,
double hadE,
double outerE,
int ecal_tp,
int hcal_tp,
const PolarLorentzVector& p4,
const GlobalPoint& emPosition,
const GlobalPoint& hadPosition);
CaloTower(const CaloTowerDetId& id,
double emE,
double hadE,
double outerE,
int ecal_tp,
int hcal_tp,
const LorentzVector& p4,
const GlobalPoint& emPosition,
const GlobalPoint& hadPosition);
CaloTower(CaloTowerDetId id,
float emE,
float hadE,
float outerE,
int ecal_tp,
int hcal_tp,
GlobalVector p3,
float iEnergy,
bool massless,
GlobalPoint emPosition,
GlobalPoint hadPosition);
CaloTower(CaloTowerDetId id,
float emE,
float hadE,
float outerE,
int ecal_tp,
int hcal_tp,
GlobalVector p3,
float iEnergy,
float imass,
GlobalPoint emPosition,
GlobalPoint hadPosition);
// setters
void addConstituent(DetId id) { constituents_.push_back(id); }
void addConstituents(const std::vector<DetId>& ids);
void setConstituents(std::vector<DetId>&& ids) { constituents_ = std::move(ids); }
void setEcalTime(int t) { ecalTime_ = t; };
void setHcalTime(int t) { hcalTime_ = t; };
void setHcalSubdet(int lastHB, int lastHE, int lastHF, int lastHO) {
int ct_ieta = ietaAbs();
if (ct_ieta <= lastHB)
subdet_ = HcalBarrel;
else if (ct_ieta <= lastHE)
subdet_ = HcalEndcap;
else if (ct_ieta <= lastHF)
subdet_ = HcalForward;
//account for HO separately
if (ct_ieta <= lastHO)
inHO_ = true;
else
inHO_ = false;
//account for gap/crossover tower separately
if (ct_ieta == lastHB)
inHBHEgap_ = true;
else
inHBHEgap_ = false;
}
// set CaloTower status based on the number of
// bad/recovered/problematic cells in ECAL and HCAL
void setCaloTowerStatus(unsigned int numBadHcalChan,
unsigned int numBadEcalChan,
unsigned int numRecHcalChan,
unsigned int numRecEcalChan,
unsigned int numProbHcalChan,
unsigned int numProbEcalChan);
void setCaloTowerStatus(uint32_t s) { twrStatusWord_ = s; }
// set the hottest cell energy in the tower
void setHottestCellE(double e) { hottestCellE_ = e; }
// getters
CaloTowerDetId id() const { return id_; }
const std::vector<DetId>& constituents() const { return constituents_; }
size_t constituentsSize() const { return constituents_.size(); }
DetId constituent(size_t i) const { return constituents_[i]; }
// energy contributions from different detectors
// energy in HO ("outerEnergy")is not included in "hadEnergy"
double emEnergy() const { return emE_; }
double hadEnergy() const { return hadE_; }
double outerEnergy() const { return (inHO_) ? outerE_ : 0.0; }
// transverse energies wrt to vtx (0,0,0)
double emEt() const { return emE_ * sin(theta()); }
double hadEt() const { return hadE_ * sin(theta()); }
double outerEt() const { return (inHO_) ? outerE_ * sin(theta()) : 0.0; }
// preserve the inherited default accessors where applicable
// (user gets default p4 wrt to vtx (0,0,0) using p4(), etc.
using LeafCandidate::et;
using LeafCandidate::p;
using LeafCandidate::p4;
// recalculated wrt user provided vertex Z position;
math::PtEtaPhiMLorentzVector p4(double vtxZ) const;
double p(double vtxZ) const { return p4(vtxZ).P(); }
double et(double vtxZ) const { return p4(vtxZ).Et(); }
double emEt(double vtxZ) const { return emE_ * sin(p4(vtxZ).theta()); }
double hadEt(double vtxZ) const { return hadE_ * sin(p4(vtxZ).theta()); }
double outerEt(double vtxZ) const { return (inHO_) ? outerE_ * sin(p4(vtxZ).theta()) : 0.0; }
// recalculated wrt vertex provided as 3D point
math::PtEtaPhiMLorentzVector p4(const Point& v) const;
double p(const Point& v) const { return p4(v).P(); }
double et(const Point& v) const { return p4(v).Et(); }
double emEt(const Point& v) const { return emE_ * sin(p4(v).theta()); }
double hadEt(const Point& v) const { return hadE_ * sin(p4(v).theta()); }
double outerEt(const Point& v) const { return (inHO_) ? outerE_ * sin(p4(v).theta()) : 0.0; }
double hottestCellE() const { return hottestCellE_; }
// Access to p4 comming from HO alone: requested by JetMET to add/subtract HO contributions
// to the tower for cases when the tower collection was created without/with HO
math::PtEtaPhiMLorentzVector p4_HO() const;
math::PtEtaPhiMLorentzVector p4_HO(double vtxZ) const;
math::PtEtaPhiMLorentzVector p4_HO(const Point& v) const;
// the reference poins in ECAL and HCAL for direction determination
// algorithm and parameters for selecting these points are set in the CaloTowersCreator
const GlobalPoint& emPosition() const { return emPosition_; }
const GlobalPoint& hadPosition() const { return hadPosition_; }
int emLvl1() const { return emLvl1_; }
int hadLv11() const { return hadLvl1_; }
// energy contained in depths>1 in the HE for 18<|iEta|<29
double hadEnergyHeOuterLayer() const { return (subdet_ == HcalEndcap) ? outerE_ : 0; }
double hadEnergyHeInnerLayer() const { return (subdet_ == HcalEndcap) ? hadE_ - outerE_ : 0; }
// energy in the tower by HCAL subdetector
// This is trivial except for tower 16
// needed by JetMET cleanup in AOD.
double energyInHB() const; // { return (id_.ietaAbs()<16)? outerE_ : 0.0; }
double energyInHE() const;
double energyInHF() const;
double energyInHO() const;
// time (ns) in ECAL/HCAL components of the tower based on weigted sum of the times in the contributing RecHits
float ecalTime() const { return float(ecalTime_) * 0.01; }
float hcalTime() const { return float(hcalTime_) * 0.01; }
// position information on the tower
int ieta() const { return id_.ieta(); }
int ietaAbs() const { return id_.ietaAbs(); }
int iphi() const { return id_.iphi(); }
int zside() const { return id_.zside(); }
int numCrystals() const;
// methods to retrieve status information from the CaloTower:
// number of bad/recovered/problematic cells in the tower
// separately for ECAL and HCAL
unsigned int numBadEcalCells() const { return (twrStatusWord_ & 0x1F); }
unsigned int numRecoveredEcalCells() const { return ((twrStatusWord_ >> 5) & 0x1F); }
unsigned int numProblematicEcalCells() const { return ((twrStatusWord_ >> 10) & 0x1F); }
unsigned int numBadHcalCells() const { return ((twrStatusWord_ >> 15) & 0x7); }
unsigned int numRecoveredHcalCells() const { return ((twrStatusWord_ >> 18) & 0x7); }
unsigned int numProblematicHcalCells() const { return ((twrStatusWord_ >> 21) & 0x7); }
// the status word itself
uint32_t towerStatusWord() const { return twrStatusWord_; }
private:
CaloTowerDetId id_;
uint32_t twrStatusWord_;
// positions of assumed EM and HAD shower positions
GlobalPoint emPosition_;
GlobalPoint hadPosition_;
//hcal subdetector info
HcalSubdetector subdet_{HcalEmpty};
bool inHO_{false}, inHBHEgap_{false};
// time
int ecalTime_;
int hcalTime_;
float emE_, hadE_, outerE_;
// for Jet ID use: hottest cell (ECAl or HCAL)
float hottestCellE_;
int emLvl1_, hadLvl1_;
std::vector<DetId> constituents_;
// vertex correction of EM and HAD momentum components:
// internally used in the transformation of the CaloTower p4
// for 3D vertex
math::PtEtaPhiMLorentzVector hadP4(const Point& v) const;
math::PtEtaPhiMLorentzVector emP4(const Point& v) const;
// taking only z-component
math::PtEtaPhiMLorentzVector hadP4(double vtxZ) const;
math::PtEtaPhiMLorentzVector emP4(double vtxZ) const;
};
std::ostream& operator<<(std::ostream& s, const CaloTower& ct);
inline bool operator==(const CaloTower& t1, const CaloTower& t2) { return t1.id() == t2.id(); }
#endif
|