CMS 3D CMS Logo
MatrixHCALDetIds.cc
Go to the documentation of this file.
1 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
2 
3 #include "Calibration/IsolatedParticles/interface/MatrixHCALDetIds.h"
4 #include "Calibration/IsolatedParticles/interface/FindDistCone.h"
5 #include "Calibration/IsolatedParticles/interface/DebugInfo.h"
6 #include "Geometry/HcalTowerAlgo/interface/HcalGeometry.h"
7 
8 #include <algorithm>
9 #include <iostream>
10 
11 //#define EDM_ML_DEBUG
12 
13 namespace spr {
14 
15  std::vector<DetId> matrixHCALIds(
16  std::vector<DetId>& dets, const HcalTopology* topology, int ieta, int iphi, bool includeHO, bool debug) {
17 #ifdef EDM_ML_DEBUG
18  if (debug) {
19  std::cout << "matrixHCALIds::Add " << ieta << " rows and " << iphi << " columns of cells for " << dets.size()
20  << " cells" << std::endl;
21  spr::debugHcalDets(0, dets);
22  }
23 #endif
24  std::vector<DetId> vdetN = spr::newHCALIdNS(dets, 0, topology, true, ieta, iphi, debug);
25  std::vector<DetId> vdetS = spr::newHCALIdNS(dets, 0, topology, false, ieta, iphi, debug);
26  for (unsigned int i1 = 0; i1 < vdetS.size(); i1++) {
27  if (std::count(vdetN.begin(), vdetN.end(), vdetS[i1]) == 0)
28  vdetN.push_back(vdetS[i1]);
29  }
30 
31  vdetS = spr::matrixHCALIdsDepth(vdetN, topology, includeHO, debug);
32 
33 #ifdef EDM_ML_DEBUG
34  if (debug) {
35  std::cout << "matrixHCALIds::Total number of cells found is " << vdetS.size() << std::endl;
36  spr::debugHcalDets(0, vdetS);
37  }
38 #endif
39  return vdetS;
40  }
41 
42  std::vector<DetId> matrixHCALIds(const DetId& det,
43  const CaloGeometry* geo,
44  const HcalTopology* topology,
45  double dR,
46  const GlobalVector& trackMom,
47  bool includeHO,
48  bool debug) {
49  HcalDetId hcdet = HcalDetId(det);
50  GlobalPoint core = (static_cast<const HcalGeometry*>(geo->getSubdetectorGeometry(hcdet)))->getPosition(hcdet);
51  std::vector<DetId> dets, vdetx;
52  dets.push_back(det);
53  int ietaphi = (int)(dR / 15.0) + 1;
54  std::vector<DetId> vdets = spr::matrixHCALIds(dets, topology, ietaphi, ietaphi, includeHO, debug);
55  for (unsigned int i = 0; i < vdets.size(); ++i) {
56  HcalDetId hcdet = HcalDetId(vdets[i]);
57  GlobalPoint rpoint = (static_cast<const HcalGeometry*>(geo->getSubdetectorGeometry(hcdet)))->getPosition(hcdet);
58  if (spr::getDistInPlaneTrackDir(core, trackMom, rpoint) < dR) {
59  vdetx.push_back(vdets[i]);
60  }
61  }
62 
63 #ifdef EDM_ML_DEBUG
64  if (debug) {
65  std::cout << "matrixHCALIds::Final List of cells for dR " << dR << " is with " << vdetx.size()
66  << " from original list of " << vdets.size() << " cells" << std::endl;
67  spr::debugHcalDets(0, vdetx);
68  }
69 #endif
70  return vdetx;
71  }
72 
73  std::vector<DetId> matrixHCALIds(std::vector<DetId>& dets,
74  const HcalTopology* topology,
75  int ietaE,
76  int ietaW,
77  int iphiN,
78  int iphiS,
79  bool includeHO,
80  bool debug) {
81 #ifdef EDM_ML_DEBUG
82  if (debug) {
83  std::cout << "matrixHCALIds::Add " << ietaE << "|" << ietaW << " rows and " << iphiN << "|" << iphiS
84  << " columns of cells for " << dets.size() << " cells" << std::endl;
85  spr::debugHcalDets(0, dets);
86  }
87 #endif
88  std::vector<DetId> vdetN = spr::newHCALIdNS(dets, 0, topology, true, ietaE, ietaW, iphiN, iphiS, debug);
89  std::vector<DetId> vdetS = spr::newHCALIdNS(dets, 0, topology, false, ietaE, ietaW, iphiN, iphiS, debug);
90  for (unsigned int i1 = 0; i1 < vdetS.size(); i1++) {
91  if (std::count(vdetN.begin(), vdetN.end(), vdetS[i1]) == 0)
92  vdetN.push_back(vdetS[i1]);
93  }
94 
95  vdetS = spr::matrixHCALIdsDepth(vdetN, topology, includeHO, debug);
96 
97 #ifdef EDM_ML_DEBUG
98  if (debug) {
99  std::cout << "matrixHCALIds::Total number of cells found is " << vdetS.size() << std::endl;
100  spr::debugHcalDets(0, vdetS);
101  }
102 #endif
103  return vdetS;
104  }
105 
106  std::vector<DetId> newHCALIdNS(std::vector<DetId>& dets,
107  unsigned int last,
108  const HcalTopology* topology,
109  bool shiftNorth,
110  int ieta,
111  int iphi,
112  bool debug) {
113 #ifdef EDM_ML_DEBUG
114  if (debug) {
115  std::cout << "newHCALIdNS::Add " << iphi << " columns of cells along " << shiftNorth << " for "
116  << (dets.size() - last) << " cells" << std::endl;
117  spr::debugHcalDets(last, dets);
118  }
119 #endif
120  std::vector<DetId> vdets;
121  vdets.insert(vdets.end(), dets.begin(), dets.end());
122  std::vector<DetId> vdetE, vdetW;
123  if (last == 0) {
124  vdetE = spr::newHCALIdEW(dets, last, topology, true, ieta, debug);
125  vdetW = spr::newHCALIdEW(dets, last, topology, false, ieta, debug);
126  for (unsigned int i1 = 0; i1 < vdetW.size(); i1++) {
127  if (std::count(vdets.begin(), vdets.end(), vdetW[i1]) == 0)
128  vdets.push_back(vdetW[i1]);
129  }
130  for (unsigned int i1 = 0; i1 < vdetE.size(); i1++) {
131  if (std::count(vdets.begin(), vdets.end(), vdetE[i1]) == 0)
132  vdets.push_back(vdetE[i1]);
133  }
134 #ifdef EDM_ML_DEBUG
135  if (debug) {
136  std::cout << "newHCALIdNS::With Added cells along E/W results a set of " << (vdets.size() - dets.size())
137  << " new cells" << std::endl;
138  spr::debugHcalDets(dets.size(), vdets);
139  }
140 #endif
141  }
142  unsigned int last0 = vdets.size();
143  if (iphi > 0) {
144  std::vector<DetId> vdetnew;
145  for (unsigned int i1 = last; i1 < dets.size(); i1++) {
146  std::vector<DetId> vdet;
147  if (shiftNorth)
148  vdet = topology->north(dets[i1]);
149  else
150  vdet = topology->south(dets[i1]);
151  for (unsigned int i2 = 0; i2 < vdet.size(); i2++) {
152  if (std::count(vdets.begin(), vdets.end(), vdet[i2]) == 0)
153  vdetnew.push_back(vdet[i2]);
154  }
155  }
156  iphi--;
157  vdetE = spr::newHCALIdEW(vdetnew, 0, topology, true, ieta, debug);
158  vdetW = spr::newHCALIdEW(vdetnew, 0, topology, false, ieta, debug);
159  for (unsigned int i2 = 0; i2 < vdetW.size(); i2++) {
160  if (std::count(vdets.begin(), vdets.end(), vdetW[i2]) == 0 &&
161  std::count(vdetnew.begin(), vdetnew.end(), vdetW[i2]) == 0)
162  vdets.push_back(vdetW[i2]);
163  }
164  for (unsigned int i2 = 0; i2 < vdetE.size(); i2++) {
165  if (std::count(vdets.begin(), vdets.end(), vdetE[i2]) == 0 &&
166  std::count(vdetnew.begin(), vdetnew.end(), vdetE[i2]) == 0)
167  vdets.push_back(vdetE[i2]);
168  }
169  last = vdets.size();
170  vdets.insert(vdets.end(), vdetnew.begin(), vdetnew.end());
171 #ifdef EDM_ML_DEBUG
172  if (debug) {
173  std::cout << "newHCALIdNS::Addition results a set of " << (vdets.size() - last0) << " new cells" << std::endl;
174  spr::debugHcalDets(last0, vdets);
175  }
176 #endif
177  last0 = last;
178  }
179 
180  if (iphi > 0) {
181  last = last0;
182  return spr::newHCALIdNS(vdets, last, topology, shiftNorth, ieta, iphi, debug);
183  } else {
184 #ifdef EDM_ML_DEBUG
185  if (debug) {
186  std::cout << "newHCALIdNS::Final list consists of " << vdets.size() << " cells" << std::endl;
187  spr::debugHcalDets(0, vdets);
188  }
189 #endif
190  return vdets;
191  }
192  }
193 
194  std::vector<DetId> newHCALIdNS(std::vector<DetId>& dets,
195  unsigned int last,
196  const HcalTopology* topology,
197  bool shiftNorth,
198  int ietaE,
199  int ietaW,
200  int iphiN,
201  int iphiS,
202  bool debug) {
203 #ifdef EDM_ML_DEBUG
204  if (debug) {
205  std::cout << "newHCALIdNS::Add " << iphiN << "|" << iphiS << " columns of cells along " << shiftNorth << " for "
206  << (dets.size() - last) << " cells" << std::endl;
207  spr::debugHcalDets(last, dets);
208  }
209 #endif
210  std::vector<DetId> vdets;
211  vdets.insert(vdets.end(), dets.begin(), dets.end());
212  std::vector<DetId> vdetE, vdetW;
213  if (last == 0) {
214  vdetE = spr::newHCALIdEW(dets, last, topology, true, ietaE, ietaW, debug);
215  vdetW = spr::newHCALIdEW(dets, last, topology, false, ietaE, ietaW, debug);
216  for (unsigned int i1 = 0; i1 < vdetW.size(); i1++) {
217  if (std::count(vdets.begin(), vdets.end(), vdetW[i1]) == 0)
218  vdets.push_back(vdetW[i1]);
219  }
220  for (unsigned int i1 = 0; i1 < vdetE.size(); i1++) {
221  if (std::count(vdets.begin(), vdets.end(), vdetE[i1]) == 0)
222  vdets.push_back(vdetE[i1]);
223  }
224 #ifdef EDM_ML_DEBUG
225  if (debug) {
226  std::cout << "newHCALIdNS::With Added cells along E/W results a set of " << (vdets.size() - dets.size())
227  << " new cells" << std::endl;
228  spr::debugHcalDets(dets.size(), vdets);
229  }
230 #endif
231  }
232  unsigned int last0 = vdets.size();
233  int iphi = iphiS;
234  if (shiftNorth)
235  iphi = iphiN;
236  if (iphi > 0) {
237  std::vector<DetId> vdetnew;
238  for (unsigned int i1 = last; i1 < dets.size(); i1++) {
239  std::vector<DetId> vdet;
240  if (shiftNorth)
241  vdet = topology->north(dets[i1]);
242  else
243  vdet = topology->south(dets[i1]);
244  for (unsigned int i2 = 0; i2 < vdet.size(); i2++) {
245  if (std::count(vdets.begin(), vdets.end(), vdet[i2]) == 0)
246  vdetnew.push_back(vdet[i2]);
247  }
248  }
249  iphi--;
250  vdetE = spr::newHCALIdEW(vdetnew, 0, topology, true, ietaE, ietaW, debug);
251  vdetW = spr::newHCALIdEW(vdetnew, 0, topology, false, ietaE, ietaW, debug);
252  for (unsigned int i2 = 0; i2 < vdetW.size(); i2++) {
253  if (std::count(vdets.begin(), vdets.end(), vdetW[i2]) == 0 &&
254  std::count(vdetnew.begin(), vdetnew.end(), vdetW[i2]) == 0)
255  vdets.push_back(vdetW[i2]);
256  }
257  for (unsigned int i2 = 0; i2 < vdetE.size(); i2++) {
258  if (std::count(vdets.begin(), vdets.end(), vdetE[i2]) == 0 &&
259  std::count(vdetnew.begin(), vdetnew.end(), vdetE[i2]) == 0)
260  vdets.push_back(vdetE[i2]);
261  }
262  last = vdets.size();
263  vdets.insert(vdets.end(), vdetnew.begin(), vdetnew.end());
264 #ifdef EDM_ML_DEBUG
265  if (debug) {
266  std::cout << "newHCALIdNS::Addition results a set of " << (vdets.size() - last0) << " new cells" << std::endl;
267  spr::debugHcalDets(last0, vdets);
268  }
269 #endif
270  last0 = last;
271  }
272  if (shiftNorth)
273  iphiN = iphi;
274  else
275  iphiS = iphi;
276 
277  if (iphi > 0) {
278  last = last0;
279  return spr::newHCALIdNS(vdets, last, topology, shiftNorth, ietaE, ietaW, iphiN, iphiS, debug);
280  } else {
281 #ifdef EDM_ML_DEBUG
282  if (debug) {
283  std::cout << "newHCALIdNS::Final list consists of " << vdets.size() << " cells" << std::endl;
284  spr::debugHcalDets(0, vdets);
285  }
286 #endif
287  return vdets;
288  }
289  }
290 
291  std::vector<DetId> newHCALIdEW(
292  std::vector<DetId>& dets, unsigned int last, const HcalTopology* topology, bool shiftEast, int ieta, bool debug) {
293 #ifdef EDM_ML_DEBUG
294  if (debug) {
295  std::cout << "newHCALIdEW::Add " << ieta << " rows of cells along " << shiftEast << " for "
296  << (dets.size() - last) << " cells" << std::endl;
297  spr::debugHcalDets(last, dets);
298  }
299 #endif
300  std::vector<DetId> vdets;
301  vdets.insert(vdets.end(), dets.begin(), dets.end());
302  if (ieta > 0) {
303  for (unsigned int i1 = last; i1 < dets.size(); i1++) {
304  std::vector<DetId> vdet;
305  if (shiftEast)
306  vdet = topology->east(dets[i1]);
307  else
308  vdet = topology->west(dets[i1]);
309  for (unsigned int i2 = 0; i2 < vdet.size(); i2++) {
310  if (std::count(vdets.begin(), vdets.end(), vdet[i2]) == 0)
311  vdets.push_back(vdet[i2]);
312  }
313  }
314  ieta--;
315  }
316 
317 #ifdef EDM_ML_DEBUG
318  if (debug) {
319  std::cout << "newHCALIdEW::Addition results a set of " << (vdets.size() - dets.size()) << " new cells"
320  << std::endl;
321  spr::debugHcalDets(dets.size(), vdets);
322  }
323 #endif
324  if (ieta > 0) {
325  last = dets.size();
326  return spr::newHCALIdEW(vdets, last, topology, shiftEast, ieta, debug);
327  } else {
328 #ifdef EDM_ML_DEBUG
329  if (debug) {
330  std::cout << "newHCALIdEW::Final list (EW) consists of " << vdets.size() << " cells" << std::endl;
331  spr::debugHcalDets(0, vdets);
332  }
333 #endif
334  return vdets;
335  }
336  }
337 
338  std::vector<DetId> newHCALIdEW(std::vector<DetId>& dets,
339  unsigned int last,
340  const HcalTopology* topology,
341  bool shiftEast,
342  int ietaE,
343  int ietaW,
344  bool debug) {
345 #ifdef EDM_ML_DEBUG
346  if (debug) {
347  std::cout << "newHCALIdEW::Add " << ietaE << "|" << ietaW << " rows of cells along " << shiftEast << " for "
348  << (dets.size() - last) << " cells" << std::endl;
349  spr::debugHcalDets(last, dets);
350  }
351 #endif
352  int ieta = ietaW;
353  if (shiftEast)
354  ieta = ietaE;
355  std::vector<DetId> vdets;
356  vdets.insert(vdets.end(), dets.begin(), dets.end());
357  if (ieta > 0) {
358  for (unsigned int i1 = last; i1 < dets.size(); i1++) {
359  std::vector<DetId> vdet;
360  if (shiftEast)
361  vdet = topology->east(dets[i1]);
362  else
363  vdet = topology->west(dets[i1]);
364  for (unsigned int i2 = 0; i2 < vdet.size(); i2++) {
365  if (std::count(vdets.begin(), vdets.end(), vdet[i2]) == 0)
366  vdets.push_back(vdet[i2]);
367  }
368  }
369  ieta--;
370  }
371  if (shiftEast)
372  ietaE = ieta;
373  else
374  ietaW = ieta;
375 
376 #ifdef EDM_ML_DEBUG
377  if (debug) {
378  std::cout << "newHCALIdEW::Addition results a set of " << (vdets.size() - dets.size()) << " new cells"
379  << std::endl;
380  spr::debugHcalDets(dets.size(), vdets);
381  }
382 #endif
383  if (ieta > 0) {
384  last = dets.size();
385  return spr::newHCALIdEW(vdets, last, topology, shiftEast, ietaE, ietaW, debug);
386  } else {
387 #ifdef EDM_ML_DEBUG
388  if (debug) {
389  std::cout << "newHCALIdEW::Final list (EW) consists of " << vdets.size() << " cells" << std::endl;
390  spr::debugHcalDets(0, vdets);
391  }
392 #endif
393  return vdets;
394  }
395  }
396 
397  std::vector<DetId> matrixHCALIdsDepth(std::vector<DetId>& dets,
398  const HcalTopology* topology,
399  bool includeHO,
400  bool
401 #ifdef EDM_ML_DEBUG
402  debug
403 #endif
404  ) {
405 
406 #ifdef EDM_ML_DEBUG
407  if (debug) {
408  std::cout << "matrixHCALIdsDepth::Add cells with higher depths with HO"
409  << "Flag set to " << includeHO << " to existing " << dets.size() << " cells" << std::endl;
410  spr::debugHcalDets(0, dets);
411  }
412 #endif
413  std::vector<DetId> vdets(dets);
414  for (unsigned int i1 = 0; i1 < dets.size(); i1++) {
415  HcalDetId vdet = dets[i1];
416  for (int idepth = 0; idepth < 3; idepth++) {
417  std::vector<DetId> vUpDetId = topology->up(vdet);
418  if (!vUpDetId.empty()) {
419  if (includeHO || vUpDetId[0].subdetId() != (int)(HcalOuter)) {
420  int n = std::count(vdets.begin(), vdets.end(), vUpDetId[0]);
421  if (n == 0) {
422 #ifdef EDM_ML_DEBUG
423  if (debug)
424  std::cout << "matrixHCALIdsDepth:: Depth " << idepth << " " << vdet << " " << (HcalDetId)vUpDetId[0]
425  << std::endl;
426 #endif
427  vdets.push_back(vUpDetId[0]);
428  }
429  }
430  vdet = vUpDetId[0];
431  }
432  }
433  }
434 
435 #ifdef EDM_ML_DEBUG
436  if (debug) {
437  std::cout << "matrixHCALIdsDepth::Final list contains " << vdets.size() << " cells" << std::endl;
438  spr::debugHcalDets(0, vdets);
439  }
440 #endif
441  return vdets;
442  }
443 
444 } // namespace spr
spr::newHCALIdNS
std::vector< DetId > newHCALIdNS(std::vector< DetId > &dets, unsigned int last, const HcalTopology *topology, bool shiftNorth, int ieta, int iphi, bool debug=false)
Definition: MatrixHCALDetIds.cc:106
CaloGeometry::getSubdetectorGeometry
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:34
spr::getDistInPlaneTrackDir
double getDistInPlaneTrackDir(const GlobalPoint &caloPoint, const GlobalVector &caloVector, const GlobalPoint &rechitPoint, bool debug=false)
Definition: FindDistCone.cc:12
mps_fire.i
i
Definition: mps_fire.py:341
ecaldqm::topology
CaloTopology const * topology(0)
HcalTopology::up
std::vector< DetId > up(const DetId &id) const override
Definition: HcalTopology.cc:446
core
Definition: __init__.py:1
testProducerWithPsetDescEmpty_cfi.i2
i2
Definition: testProducerWithPsetDescEmpty_cfi.py:46
spr::debugHcalDets
void debugHcalDets(unsigned int, std::vector< DetId > &)
Definition: DebugInfo.cc:39
HcalTopology::south
std::vector< DetId > south(const DetId &id) const override
Definition: HcalTopology.cc:435
DetId
Definition: DetId.h:17
debug
#define debug
Definition: HDRShower.cc:19
HGC3DClusterGenMatchSelector_cfi.dR
dR
Definition: HGC3DClusterGenMatchSelector_cfi.py:7
testProducerWithPsetDescEmpty_cfi.i1
i1
Definition: testProducerWithPsetDescEmpty_cfi.py:45
spr::matrixHCALIdsDepth
std::vector< DetId > matrixHCALIdsDepth(std::vector< DetId > &dets, const HcalTopology *topology, bool includeHO=false, bool debug=false)
Definition: MatrixHCALDetIds.cc:397
HcalTopology::west
std::vector< DetId > west(const DetId &id) const override
Definition: HcalTopology.cc:413
HcalTopology::north
std::vector< DetId > north(const DetId &id) const override
Definition: HcalTopology.cc:424
spr::matrixHCALIds
std::vector< DetId > matrixHCALIds(std::vector< DetId > &dets, const HcalTopology *topology, int ieta, int iphi, bool includeHO=false, bool debug=false)
Definition: MatrixHCALDetIds.cc:15
EDM_ML_DEBUG
#define EDM_ML_DEBUG
Definition: HcalHBHEMuonSimAnalyzer.cc:41
gather_cfg.cout
cout
Definition: gather_cfg.py:144
HcalTopology::east
std::vector< DetId > east(const DetId &id) const override
Definition: HcalTopology.cc:402
spr::newHCALIdEW
std::vector< DetId > newHCALIdEW(std::vector< DetId > &dets, unsigned int last, const HcalTopology *topology, bool shiftEast, int ieta, bool debug=false)
Definition: MatrixHCALDetIds.cc:291
Generated for CMSSW Reference Manual by doxygen 1.8.11