CMS 3D CMS Logo
CaloTowersCreationAlgo.cc
Go to the documentation of this file.
1 #include "RecoLocalCalo/CaloTowersCreator/interface/CaloTowersCreationAlgo.h"
2 #include "Geometry/CaloTopology/interface/HcalTopology.h"
3 #include "Geometry/CaloTopology/interface/CaloTowerTopology.h"
4 #include "Geometry/CaloTopology/interface/CaloTowerConstituentsMap.h"
5 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
6 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
7 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
8 #include "FWCore/MessageLogger/interface/MessageLogger.h"
9 #include "Math/Interpolator.h"
10 #include <cmath>
11 
12 //#define EDM_ML_DEBUG
13 
14 CaloTowersCreationAlgo::CaloTowersCreationAlgo()
15  : theEBthreshold(-1000.),
16  theEEthreshold(-1000.),
17 
18  theUseEtEBTresholdFlag(false),
19  theUseEtEETresholdFlag(false),
20  theUseSymEBTresholdFlag(false),
21  theUseSymEETresholdFlag(false),
22 
23  theHcalThreshold(-1000.),
24  theHBthreshold(-1000.),
25  theHBthreshold1(-1000.),
26  theHBthreshold2(-1000.),
27  theHESthreshold(-1000.),
28  theHESthreshold1(-1000.),
29  theHEDthreshold(-1000.),
30  theHEDthreshold1(-1000.),
31  theHOthreshold0(-1000.),
32  theHOthresholdPlus1(-1000.),
33  theHOthresholdMinus1(-1000.),
34  theHOthresholdPlus2(-1000.),
35  theHOthresholdMinus2(-1000.),
36  theHF1threshold(-1000.),
37  theHF2threshold(-1000.),
38  theEBGrid(std::vector<double>(5, 10.)),
39  theEBWeights(std::vector<double>(5, 1.)),
40  theEEGrid(std::vector<double>(5, 10.)),
41  theEEWeights(std::vector<double>(5, 1.)),
42  theHBGrid(std::vector<double>(5, 10.)),
43  theHBWeights(std::vector<double>(5, 1.)),
44  theHESGrid(std::vector<double>(5, 10.)),
45  theHESWeights(std::vector<double>(5, 1.)),
46  theHEDGrid(std::vector<double>(5, 10.)),
47  theHEDWeights(std::vector<double>(5, 1.)),
48  theHOGrid(std::vector<double>(5, 10.)),
49  theHOWeights(std::vector<double>(5, 1.)),
50  theHF1Grid(std::vector<double>(5, 10.)),
51  theHF1Weights(std::vector<double>(5, 1.)),
52  theHF2Grid(std::vector<double>(5, 10.)),
53  theHF2Weights(std::vector<double>(5, 1.)),
54  theEBweight(1.),
55  theEEweight(1.),
56  theHBweight(1.),
57  theHESweight(1.),
58  theHEDweight(1.),
59  theHOweight(1.),
60  theHF1weight(1.),
61  theHF2weight(1.),
62  theEcutTower(-1000.),
63  theEBSumThreshold(-1000.),
64  theEESumThreshold(-1000.),
65  theEBEScale(50.),
66  theEEEScale(50.),
67  theHBEScale(50.),
68  theHESEScale(50.),
69  theHEDEScale(50.),
70  theHOEScale(50.),
71  theHF1EScale(50.),
72  theHF2EScale(50.),
73  theHcalTopology(nullptr),
74  theGeometry(nullptr),
75  theTowerConstituentsMap(nullptr),
76  theHcalAcceptSeverityLevel(0),
77  theRecoveredHcalHitsAreUsed(false),
78  theRecoveredEcalHitsAreUsed(false),
79  useRejectedHitsOnly(false),
80  theHcalAcceptSeverityLevelForRejectedHit(0),
81  useRejectedRecoveredHcalHits(0),
82  useRejectedRecoveredEcalHits(0),
83  missingHcalRescaleFactorForEcal(0.),
84  theHOIsUsed(true),
85  // (for momentum reconstruction algorithm)
86  theMomConstrMethod(0),
87  theMomHBDepth(0.),
88  theMomHEDepth(0.),
89  theMomEBDepth(0.),
90  theMomEEDepth(0.),
91  theHcalPhase(0) {}
92 
93 CaloTowersCreationAlgo::CaloTowersCreationAlgo(double EBthreshold,
94  double EEthreshold,
95 
96  bool useEtEBTreshold,
97  bool useEtEETreshold,
98  bool useSymEBTreshold,
99  bool useSymEETreshold,
100 
101  double HcalThreshold,
102  double HBthreshold,
103  double HBthreshold1,
104  double HBthreshold2,
105  double HESthreshold,
106  double HESthreshold1,
107  double HEDthreshold,
108  double HEDthreshold1,
109  double HOthreshold0,
110  double HOthresholdPlus1,
111  double HOthresholdMinus1,
112  double HOthresholdPlus2,
113  double HOthresholdMinus2,
114  double HF1threshold,
115  double HF2threshold,
116  double EBweight,
117  double EEweight,
118  double HBweight,
119  double HESweight,
120  double HEDweight,
121  double HOweight,
122  double HF1weight,
123  double HF2weight,
124  double EcutTower,
125  double EBSumThreshold,
126  double EESumThreshold,
127  bool useHO,
128  // (momentum reconstruction algorithm)
129  int momConstrMethod,
130  double momHBDepth,
131  double momHEDepth,
132  double momEBDepth,
133  double momEEDepth,
134  int hcalPhase)
135 
136  : theEBthreshold(EBthreshold),
137  theEEthreshold(EEthreshold),
138 
139  theUseEtEBTresholdFlag(useEtEBTreshold),
140  theUseEtEETresholdFlag(useEtEETreshold),
141  theUseSymEBTresholdFlag(useSymEBTreshold),
142  theUseSymEETresholdFlag(useSymEETreshold),
143 
144  theHcalThreshold(HcalThreshold),
145  theHBthreshold(HBthreshold),
146  theHBthreshold1(HBthreshold1),
147  theHBthreshold2(HBthreshold2),
148  theHESthreshold(HESthreshold),
149  theHESthreshold1(HESthreshold1),
150  theHEDthreshold(HEDthreshold),
151  theHEDthreshold1(HEDthreshold1),
152  theHOthreshold0(HOthreshold0),
153  theHOthresholdPlus1(HOthresholdPlus1),
154  theHOthresholdMinus1(HOthresholdMinus1),
155  theHOthresholdPlus2(HOthresholdPlus2),
156  theHOthresholdMinus2(HOthresholdMinus2),
157  theHF1threshold(HF1threshold),
158  theHF2threshold(HF2threshold),
159  theEBGrid(std::vector<double>(5, 10.)),
160  theEBWeights(std::vector<double>(5, 1.)),
161  theEEGrid(std::vector<double>(5, 10.)),
162  theEEWeights(std::vector<double>(5, 1.)),
163  theHBGrid(std::vector<double>(5, 10.)),
164  theHBWeights(std::vector<double>(5, 1.)),
165  theHESGrid(std::vector<double>(5, 10.)),
166  theHESWeights(std::vector<double>(5, 1.)),
167  theHEDGrid(std::vector<double>(5, 10.)),
168  theHEDWeights(std::vector<double>(5, 1.)),
169  theHOGrid(std::vector<double>(5, 10.)),
170  theHOWeights(std::vector<double>(5, 1.)),
171  theHF1Grid(std::vector<double>(5, 10.)),
172  theHF1Weights(std::vector<double>(5, 1.)),
173  theHF2Grid(std::vector<double>(5, 10.)),
174  theHF2Weights(std::vector<double>(5, 1.)),
175  theEBweight(EBweight),
176  theEEweight(EEweight),
177  theHBweight(HBweight),
178  theHESweight(HESweight),
179  theHEDweight(HEDweight),
180  theHOweight(HOweight),
181  theHF1weight(HF1weight),
182  theHF2weight(HF2weight),
183  theEcutTower(EcutTower),
184  theEBSumThreshold(EBSumThreshold),
185  theEESumThreshold(EESumThreshold),
186  theEBEScale(50.),
187  theEEEScale(50.),
188  theHBEScale(50.),
189  theHESEScale(50.),
190  theHEDEScale(50.),
191  theHOEScale(50.),
192  theHF1EScale(50.),
193  theHF2EScale(50.),
194  theHcalTopology(nullptr),
195  theGeometry(nullptr),
196  theTowerConstituentsMap(nullptr),
197  theHcalAcceptSeverityLevel(0),
198  theRecoveredHcalHitsAreUsed(false),
199  theRecoveredEcalHitsAreUsed(false),
200  useRejectedHitsOnly(false),
201  theHcalAcceptSeverityLevelForRejectedHit(0),
202  useRejectedRecoveredHcalHits(0),
203  useRejectedRecoveredEcalHits(0),
204  missingHcalRescaleFactorForEcal(0.),
205  theHOIsUsed(useHO),
206  // (momentum reconstruction algorithm)
207  theMomConstrMethod(momConstrMethod),
208  theMomHBDepth(momHBDepth),
209  theMomHEDepth(momHEDepth),
210  theMomEBDepth(momEBDepth),
211  theMomEEDepth(momEEDepth),
212  theHcalPhase(hcalPhase) {}
213 
214 CaloTowersCreationAlgo::CaloTowersCreationAlgo(double EBthreshold,
215  double EEthreshold,
216  bool useEtEBTreshold,
217  bool useEtEETreshold,
218  bool useSymEBTreshold,
219  bool useSymEETreshold,
220 
221  double HcalThreshold,
222  double HBthreshold,
223  double HBthreshold1,
224  double HBthreshold2,
225  double HESthreshold,
226  double HESthreshold1,
227  double HEDthreshold,
228  double HEDthreshold1,
229  double HOthreshold0,
230  double HOthresholdPlus1,
231  double HOthresholdMinus1,
232  double HOthresholdPlus2,
233  double HOthresholdMinus2,
234  double HF1threshold,
235  double HF2threshold,
236  const std::vector<double>& EBGrid,
237  const std::vector<double>& EBWeights,
238  const std::vector<double>& EEGrid,
239  const std::vector<double>& EEWeights,
240  const std::vector<double>& HBGrid,
241  const std::vector<double>& HBWeights,
242  const std::vector<double>& HESGrid,
243  const std::vector<double>& HESWeights,
244  const std::vector<double>& HEDGrid,
245  const std::vector<double>& HEDWeights,
246  const std::vector<double>& HOGrid,
247  const std::vector<double>& HOWeights,
248  const std::vector<double>& HF1Grid,
249  const std::vector<double>& HF1Weights,
250  const std::vector<double>& HF2Grid,
251  const std::vector<double>& HF2Weights,
252  double EBweight,
253  double EEweight,
254  double HBweight,
255  double HESweight,
256  double HEDweight,
257  double HOweight,
258  double HF1weight,
259  double HF2weight,
260  double EcutTower,
261  double EBSumThreshold,
262  double EESumThreshold,
263  bool useHO,
264  // (for the momentum construction algorithm)
265  int momConstrMethod,
266  double momHBDepth,
267  double momHEDepth,
268  double momEBDepth,
269  double momEEDepth,
270  int hcalPhase)
271 
272  : theEBthreshold(EBthreshold),
273  theEEthreshold(EEthreshold),
274 
275  theUseEtEBTresholdFlag(useEtEBTreshold),
276  theUseEtEETresholdFlag(useEtEETreshold),
277  theUseSymEBTresholdFlag(useSymEBTreshold),
278  theUseSymEETresholdFlag(useSymEETreshold),
279 
280  theHcalThreshold(HcalThreshold),
281  theHBthreshold(HBthreshold),
282  theHBthreshold1(HBthreshold1),
283  theHBthreshold2(HBthreshold2),
284  theHESthreshold(HESthreshold),
285  theHESthreshold1(HESthreshold1),
286  theHEDthreshold(HEDthreshold),
287  theHEDthreshold1(HEDthreshold1),
288  theHOthreshold0(HOthreshold0),
289  theHOthresholdPlus1(HOthresholdPlus1),
290  theHOthresholdMinus1(HOthresholdMinus1),
291  theHOthresholdPlus2(HOthresholdPlus2),
292  theHOthresholdMinus2(HOthresholdMinus2),
293  theHF1threshold(HF1threshold),
294  theHF2threshold(HF2threshold),
295  theEBGrid(EBGrid),
296  theEBWeights(EBWeights),
297  theEEGrid(EEGrid),
298  theEEWeights(EEWeights),
299  theHBGrid(HBGrid),
300  theHBWeights(HBWeights),
301  theHESGrid(HESGrid),
302  theHESWeights(HESWeights),
303  theHEDGrid(HEDGrid),
304  theHEDWeights(HEDWeights),
305  theHOGrid(HOGrid),
306  theHOWeights(HOWeights),
307  theHF1Grid(HF1Grid),
308  theHF1Weights(HF1Weights),
309  theHF2Grid(HF2Grid),
310  theHF2Weights(HF2Weights),
311  theEBweight(EBweight),
312  theEEweight(EEweight),
313  theHBweight(HBweight),
314  theHESweight(HESweight),
315  theHEDweight(HEDweight),
316  theHOweight(HOweight),
317  theHF1weight(HF1weight),
318  theHF2weight(HF2weight),
319  theEcutTower(EcutTower),
320  theEBSumThreshold(EBSumThreshold),
321  theEESumThreshold(EESumThreshold),
322  theEBEScale(50.),
323  theEEEScale(50.),
324  theHBEScale(50.),
325  theHESEScale(50.),
326  theHEDEScale(50.),
327  theHOEScale(50.),
328  theHF1EScale(50.),
329  theHF2EScale(50.),
330  theHcalTopology(nullptr),
331  theGeometry(nullptr),
332  theTowerConstituentsMap(nullptr),
333  theHcalAcceptSeverityLevel(0),
334  theRecoveredHcalHitsAreUsed(false),
335  theRecoveredEcalHitsAreUsed(false),
336  useRejectedHitsOnly(false),
337  theHcalAcceptSeverityLevelForRejectedHit(0),
338  useRejectedRecoveredHcalHits(0),
339  useRejectedRecoveredEcalHits(0),
340  missingHcalRescaleFactorForEcal(0.),
341  theHOIsUsed(useHO),
342  // (momentum reconstruction algorithm)
343  theMomConstrMethod(momConstrMethod),
344  theMomHBDepth(momHBDepth),
345  theMomHEDepth(momHEDepth),
346  theMomEBDepth(momEBDepth),
347  theMomEEDepth(momEEDepth),
348  theHcalPhase(hcalPhase) {
349  // static int N = 0;
350  // std::cout << "VI Algo " << ++N << std::endl;
351  // nalgo=N;
352 }
353 
354 void CaloTowersCreationAlgo::setGeometry(const CaloTowerTopology* cttopo,
355  const CaloTowerConstituentsMap* ctmap,
356  const HcalTopology* htopo,
357  const CaloGeometry* geo) {
358  theTowerTopology = cttopo;
359  theTowerConstituentsMap = ctmap;
360  theHcalTopology = htopo;
361  theGeometry = geo;
362  theTowerGeometry = geo->getSubdetectorGeometry(DetId::Calo, CaloTowerDetId::SubdetId);
363 
364  //initialize ecal bad channel map
365  ecalBadChs.resize(theTowerTopology->sizeForDenseIndexing(), 0);
366 }
367 
368 void CaloTowersCreationAlgo::begin() {
369  theTowerMap.clear();
370  theTowerMapSize = 0;
371  //hcalDropChMap.clear();
372 }
373 
374 void CaloTowersCreationAlgo::process(const HBHERecHitCollection& hbhe) {
375  for (HBHERecHitCollection::const_iterator hbheItr = hbhe.begin(); hbheItr != hbhe.end(); ++hbheItr)
376  assignHitHcal(&(*hbheItr));
377 }
378 
379 void CaloTowersCreationAlgo::process(const HORecHitCollection& ho) {
380  for (HORecHitCollection::const_iterator hoItr = ho.begin(); hoItr != ho.end(); ++hoItr)
381  assignHitHcal(&(*hoItr));
382 }
383 
384 void CaloTowersCreationAlgo::process(const HFRecHitCollection& hf) {
385  for (HFRecHitCollection::const_iterator hfItr = hf.begin(); hfItr != hf.end(); ++hfItr)
386  assignHitHcal(&(*hfItr));
387 }
388 
389 void CaloTowersCreationAlgo::process(const EcalRecHitCollection& ec) {
390  for (EcalRecHitCollection::const_iterator ecItr = ec.begin(); ecItr != ec.end(); ++ecItr)
391  assignHitEcal(&(*ecItr));
392 }
393 
394 // this method should not be used any more as the towers in the changed format
395 // can not be properly rescaled with the "rescale" method.
396 // "rescale was replaced by "rescaleTowers"
397 //
398 void CaloTowersCreationAlgo::process(const CaloTowerCollection& ctc) {
399  for (CaloTowerCollection::const_iterator ctcItr = ctc.begin(); ctcItr != ctc.end(); ++ctcItr) {
400  rescale(&(*ctcItr));
401  }
402 }
403 
404 void CaloTowersCreationAlgo::finish(CaloTowerCollection& result) {
405  // now copy this map into the final collection
406  result.reserve(theTowerMapSize);
407  // auto k=0U;
408  // if (!theEbHandle.isValid()) std::cout << "VI ebHandle not valid" << std::endl;
409  // if (!theEeHandle.isValid()) std::cout << "VI eeHandle not valid" << std::endl;
410 
411  for (auto const& mt : theTowerMap) {
412  // Convert only if there is at least one constituent in the metatower.
413  // The check of constituents size in the coverted tower is still needed!
414  if (!mt.empty()) {
415  convert(mt.id, mt, result);
416  } // ++k;}
417  }
418 
419  // assert(k==theTowerMapSize);
420  // std::cout << "VI TowerMap " << theTowerMapSize << " " << k << std::endl;
421 
422  theTowerMap.clear(); // save the memory
423  theTowerMapSize = 0;
424 }
425 
426 void CaloTowersCreationAlgo::rescaleTowers(const CaloTowerCollection& ctc, CaloTowerCollection& ctcResult) {
427  for (CaloTowerCollection::const_iterator ctcItr = ctc.begin(); ctcItr != ctc.end(); ++ctcItr) {
428  CaloTowerDetId twrId = ctcItr->id();
429  double newE_em = ctcItr->emEnergy();
430  double newE_had = ctcItr->hadEnergy();
431  double newE_outer = ctcItr->outerEnergy();
432 
433  double threshold = 0.0; // not used: we do not change thresholds
434  double weight = 1.0;
435 
436  // HF
437  if (ctcItr->ietaAbs() >= theTowerTopology->firstHFRing()) {
438  double E_short = 0.5 * newE_had; // from the definitions for HF
439  double E_long = newE_em + 0.5 * newE_had; //
440  // scale
441  E_long *= theHF1weight;
442  E_short *= theHF2weight;
443  // convert
444  newE_em = E_long - E_short;
445  newE_had = 2.0 * E_short;
446  }
447 
448  else { // barrel/endcap
449 
450  // find if its in EB, or EE; determine from first ecal constituent found
451  for (unsigned int iConst = 0; iConst < ctcItr->constituentsSize(); ++iConst) {
452  DetId constId = ctcItr->constituent(iConst);
453  if (constId.det() != DetId::Ecal)
454  continue;
455  getThresholdAndWeight(constId, threshold, weight);
456  newE_em *= weight;
457  break;
458  }
459  // HO
460  for (unsigned int iConst = 0; iConst < ctcItr->constituentsSize(); ++iConst) {
461  DetId constId = ctcItr->constituent(iConst);
462  if (constId.det() != DetId::Hcal)
463  continue;
464  if (HcalDetId(constId).subdet() != HcalOuter)
465  continue;
466  getThresholdAndWeight(constId, threshold, weight);
467  newE_outer *= weight;
468  break;
469  }
470  // HB/HE
471  for (unsigned int iConst = 0; iConst < ctcItr->constituentsSize(); ++iConst) {
472  DetId constId = ctcItr->constituent(iConst);
473  if (constId.det() != DetId::Hcal)
474  continue;
475  if (HcalDetId(constId).subdet() == HcalOuter)
476  continue;
477  getThresholdAndWeight(constId, threshold, weight);
478  newE_had *= weight;
479  if (ctcItr->ietaAbs() > theTowerTopology->firstHERing())
480  newE_outer *= weight;
481  break;
482  }
483 
484  } // barrel/endcap region
485 
486  // now make the new tower
487 
488  double newE_hadTot =
489  (theHOIsUsed && twrId.ietaAbs() <= theTowerTopology->lastHORing()) ? newE_had + newE_outer : newE_had;
490 
491  GlobalPoint emPoint = ctcItr->emPosition();
492  GlobalPoint hadPoint = ctcItr->emPosition();
493 
494  double f_em = 1.0 / cosh(emPoint.eta());
495  double f_had = 1.0 / cosh(hadPoint.eta());
496 
497  CaloTower::PolarLorentzVector towerP4;
498 
499  if (ctcItr->ietaAbs() < theTowerTopology->firstHFRing()) {
500  if (newE_em > 0)
501  towerP4 += CaloTower::PolarLorentzVector(newE_em * f_em, emPoint.eta(), emPoint.phi(), 0);
502  if (newE_hadTot > 0)
503  towerP4 += CaloTower::PolarLorentzVector(newE_hadTot * f_had, hadPoint.eta(), hadPoint.phi(), 0);
504  } else {
505  double newE_tot = newE_em + newE_had;
506  // for HF we use common point for ecal, hcal shower positions regardless of the method
507  if (newE_tot > 0)
508  towerP4 += CaloTower::PolarLorentzVector(newE_tot * f_had, hadPoint.eta(), hadPoint.phi(), 0);
509  }
510 
511  CaloTower rescaledTower(twrId, newE_em, newE_had, newE_outer, -1, -1, towerP4, emPoint, hadPoint);
512  // copy the timings, have to convert back to int, 1 unit = 0.01 ns
513  rescaledTower.setEcalTime(int(ctcItr->ecalTime() * 100.0 + 0.5));
514  rescaledTower.setHcalTime(int(ctcItr->hcalTime() * 100.0 + 0.5));
515  //add topology info
516  rescaledTower.setHcalSubdet(theTowerTopology->lastHBRing(),
517  theTowerTopology->lastHERing(),
518  theTowerTopology->lastHFRing(),
519  theTowerTopology->lastHORing());
520 
521  std::vector<DetId> contains;
522  for (unsigned int iConst = 0; iConst < ctcItr->constituentsSize(); ++iConst) {
523  contains.push_back(ctcItr->constituent(iConst));
524  }
525  rescaledTower.addConstituents(contains);
526 
527  rescaledTower.setCaloTowerStatus(ctcItr->towerStatusWord());
528 
529  ctcResult.push_back(rescaledTower);
530 
531  } // end of loop over towers
532 }
533 
534 void CaloTowersCreationAlgo::assignHitHcal(const CaloRecHit* recHit) {
535  DetId detId = recHit->detid();
536  DetId detIdF(detId);
537  if (detId.det() == DetId::Hcal && theHcalTopology->getMergePositionFlag()) {
538  detIdF = theHcalTopology->idFront(HcalDetId(detId));
539 #ifdef EDM_ML_DEBUG
540  std::cout << "AssignHitHcal DetId " << HcalDetId(detId) << " Front " << HcalDetId(detIdF) << std::endl;
541 #endif
542  }
543 
544  unsigned int chStatusForCT = hcalChanStatusForCaloTower(recHit);
545 
546  // this is for skipping channls: mostly needed for the creation of
547  // bad towers from hits i the bad channel collections.
548  if (chStatusForCT == CaloTowersCreationAlgo::IgnoredChan)
549  return;
550 
551  double threshold, weight;
552  getThresholdAndWeight(detId, threshold, weight);
553 
554  double energy = recHit->energy(); // original RecHit energy is used to apply thresholds
555  double e = energy * weight; // energies scaled by user weight: used in energy assignments
556 
557  // SPECIAL handling of tower 28 merged depths --> half into tower 28 and half into tower 29
558  bool merge(false);
559  if (detIdF.det() == DetId::Hcal && HcalDetId(detIdF).subdet() == HcalEndcap &&
560  (theHcalPhase == 0 || theHcalPhase == 1) &&
561  //HcalDetId(detId).depth()==3 &&
562  HcalDetId(detIdF).ietaAbs() == theHcalTopology->lastHERing() - 1) {
563  merge = theHcalTopology->mergedDepth29(HcalDetId(detIdF));
564 #ifdef EDM_ML_DEBUG
565  std::cout << "Merge " << HcalDetId(detIdF) << ":" << merge << std::endl;
566 #endif
567  }
568  if (merge) {
570 
571  // bad channels are counted regardless of energy threshold
572 
573  if (chStatusForCT == CaloTowersCreationAlgo::BadChan) {
574  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
575  if (towerDetId.null())
576  return;
577  MetaTower& tower28 = find(towerDetId);
578  CaloTowerDetId towerDetId29(towerDetId.ieta() + towerDetId.zside(), towerDetId.iphi());
579  MetaTower& tower29 = find(towerDetId29);
580  tower28.numBadHcalCells += 1;
581  tower29.numBadHcalCells += 1;
582  }
583 
584  else if (0.5 * energy >= threshold) { // not bad channel: use energy if above threshold
585 
586  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
587  if (towerDetId.null())
588  return;
589  MetaTower& tower28 = find(towerDetId);
590  CaloTowerDetId towerDetId29(towerDetId.ieta() + towerDetId.zside(), towerDetId.iphi());
591  MetaTower& tower29 = find(towerDetId29);
592 
593  if (chStatusForCT == CaloTowersCreationAlgo::RecoveredChan) {
594  tower28.numRecHcalCells += 1;
595  tower29.numRecHcalCells += 1;
596  } else if (chStatusForCT == CaloTowersCreationAlgo::ProblematicChan) {
597  tower28.numProbHcalCells += 1;
598  tower29.numProbHcalCells += 1;
599  }
600 
601  // NOTE DIVIDE BY 2!!!
602  double e28 = 0.5 * e;
603  double e29 = 0.5 * e;
604 
605  tower28.E_had += e28;
606  tower28.E += e28;
607  std::pair<DetId, float> mc(detId, e28);
608  tower28.metaConstituents.push_back(mc);
609 
610  tower29.E_had += e29;
611  tower29.E += e29;
612  tower29.metaConstituents.push_back(mc);
613 
614  // time info: do not use in averaging if timing error is found: need
615  // full set of status info to implement: use only "good" channels for now
616 
617  if (chStatusForCT == CaloTowersCreationAlgo::GoodChan) {
618  tower28.hadSumTimeTimesE += (e28 * recHit->time());
619  tower28.hadSumEForTime += e28;
620  tower29.hadSumTimeTimesE += (e29 * recHit->time());
621  tower29.hadSumEForTime += e29;
622  }
623 
624  // store the energy in layer 3 also in E_outer
625  tower28.E_outer += e28;
626  tower29.E_outer += e29;
627  } // not a "bad" hit
628  } // end of special case
629 
630  else {
631  HcalDetId hcalDetId(detId);
632 
634 
635  if (hcalDetId.subdet() == HcalOuter) {
636  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
637  if (towerDetId.null())
638  return;
639  MetaTower& tower = find(towerDetId);
640 
641  if (chStatusForCT == CaloTowersCreationAlgo::BadChan) {
642  if (theHOIsUsed)
643  tower.numBadHcalCells += 1;
644  }
645 
646  else if (energy >= threshold) {
647  tower.E_outer += e; // store HO energy even if HO is not used
648  // add energy of the tower and/or flag if theHOIsUsed
649  if (theHOIsUsed) {
650  tower.E += e;
651 
652  if (chStatusForCT == CaloTowersCreationAlgo::RecoveredChan) {
653  tower.numRecHcalCells += 1;
654  } else if (chStatusForCT == CaloTowersCreationAlgo::ProblematicChan) {
655  tower.numProbHcalCells += 1;
656  }
657  } // HO is used
658 
659  // add HO to constituents even if it is not used: JetMET wants to keep these towers
660  std::pair<DetId, float> mc(detId, e);
661  tower.metaConstituents.push_back(mc);
662 
663  } // not a bad channel, energy above threshold
664 
665  } // HO hit
666 
667  // HF calculates EM fraction differently
668  else if (hcalDetId.subdet() == HcalForward) {
669  if (chStatusForCT == CaloTowersCreationAlgo::BadChan) {
670  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
671  if (towerDetId.null())
672  return;
673  MetaTower& tower = find(towerDetId);
674  tower.numBadHcalCells += 1;
675  }
676 
677  else if (energy >= threshold) {
678  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
679  if (towerDetId.null())
680  return;
681  MetaTower& tower = find(towerDetId);
682 
683  if (hcalDetId.depth() == 1) {
684  // long fiber, so E_EM = E(Long) - E(Short)
685  tower.E_em += e;
686  } else {
687  // short fiber, EHAD = 2 * E(Short)
688  tower.E_em -= e;
689  tower.E_had += 2. * e;
690  }
691  tower.E += e;
692  if (chStatusForCT == CaloTowersCreationAlgo::RecoveredChan) {
693  tower.numRecHcalCells += 1;
694  } else if (chStatusForCT == CaloTowersCreationAlgo::ProblematicChan) {
695  tower.numProbHcalCells += 1;
696  }
697 
698  // put the timing in HCAL -> have to check timing errors when available
699  // for now use only good channels
700  if (chStatusForCT == CaloTowersCreationAlgo::GoodChan) {
701  tower.hadSumTimeTimesE += (e * recHit->time());
702  tower.hadSumEForTime += e;
703  }
704 
705  std::pair<DetId, float> mc(detId, e);
706  tower.metaConstituents.push_back(mc);
707 
708  } // not a bad HF channel, energy above threshold
709 
710  } // HF hit
711 
712  else {
713  // HCAL situation normal in HB/HE
714  if (chStatusForCT == CaloTowersCreationAlgo::BadChan) {
715  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
716  if (towerDetId.null())
717  return;
718  MetaTower& tower = find(towerDetId);
719  tower.numBadHcalCells += 1;
720  } else if (energy >= threshold) {
721  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
722  if (towerDetId.null())
723  return;
724  MetaTower& tower = find(towerDetId);
725  tower.E_had += e;
726  tower.E += e;
727  if (chStatusForCT == CaloTowersCreationAlgo::RecoveredChan) {
728  tower.numRecHcalCells += 1;
729  } else if (chStatusForCT == CaloTowersCreationAlgo::ProblematicChan) {
730  tower.numProbHcalCells += 1;
731  }
732 
733  // Timing information: need specific accessors
734  // for now use only good channels
735  if (chStatusForCT == CaloTowersCreationAlgo::GoodChan) {
736  tower.hadSumTimeTimesE += (e * recHit->time());
737  tower.hadSumEForTime += e;
738  }
739  // store energy in highest depth for towers 18-27 (for electron,photon ID in endcap)
740  // also, store energy in HE part of tower 16 (for JetMET cleanup)
741  HcalDetId hcalDetId(detId);
742  if (hcalDetId.subdet() == HcalEndcap) {
743  if (theHcalPhase == 0) {
744  if ((hcalDetId.depth() == 2 && hcalDetId.ietaAbs() >= 18 && hcalDetId.ietaAbs() < 27) ||
745  (hcalDetId.depth() == 3 && hcalDetId.ietaAbs() == 27) ||
746  (hcalDetId.depth() == 3 && hcalDetId.ietaAbs() == 16)) {
747  tower.E_outer += e;
748  }
749  }
750  //combine depths in phase0-like way
751  else if (theHcalPhase == 1) {
752  if ((hcalDetId.depth() >= 3 && hcalDetId.ietaAbs() >= 18 && hcalDetId.ietaAbs() < 26) ||
753  (hcalDetId.depth() >= 4 && (hcalDetId.ietaAbs() == 26 || hcalDetId.ietaAbs() == 27)) ||
754  (hcalDetId.depth() == 3 && hcalDetId.ietaAbs() == 17) ||
755  (hcalDetId.depth() == 4 && hcalDetId.ietaAbs() == 16)) {
756  tower.E_outer += e;
757  }
758  }
759  }
760 
761  std::pair<DetId, float> mc(detId, e);
762  tower.metaConstituents.push_back(mc);
763 
764  } // not a "bad" channel, energy above threshold
765 
766  } // channel in HBHE (excluding twrs 28,29)
767 
768  } // recHit normal case (not in HE towers 28,29)
769 
770 } // end of assignHitHcal method
771 
772 void CaloTowersCreationAlgo::assignHitEcal(const EcalRecHit* recHit) {
773  DetId detId = recHit->detid();
774 
775  unsigned int chStatusForCT;
776  bool ecalIsBad = false;
777  std::tie(chStatusForCT, ecalIsBad) = ecalChanStatusForCaloTower(recHit);
778 
779  // this is for skipping channls: mostly needed for the creation of
780  // bad towers from hits i the bad channel collections.
781  if (chStatusForCT == CaloTowersCreationAlgo::IgnoredChan)
782  return;
783 
784  double threshold, weight;
785  getThresholdAndWeight(detId, threshold, weight);
786 
787  double energy = recHit->energy(); // original RecHit energy is used to apply thresholds
788  double e = energy * weight; // energies scaled by user weight: used in energy assignments
789 
791 
792  // For ECAL we count all bad channels after the metatower is complete
793 
794  // Include options for symmetric thresholds and cut on Et
795  // for ECAL RecHits
796 
797  bool passEmThreshold = false;
798 
799  if (detId.subdetId() == EcalBarrel) {
800  if (theUseEtEBTresholdFlag)
801  energy /= cosh((theGeometry->getGeometry(detId)->getPosition()).eta());
802  if (theUseSymEBTresholdFlag)
803  passEmThreshold = (fabs(energy) >= threshold);
804  else
805  passEmThreshold = (energy >= threshold);
806 
807  } else if (detId.subdetId() == EcalEndcap) {
808  if (theUseEtEETresholdFlag)
809  energy /= cosh((theGeometry->getGeometry(detId)->getPosition()).eta());
810  if (theUseSymEETresholdFlag)
811  passEmThreshold = (fabs(energy) >= threshold);
812  else
813  passEmThreshold = (energy >= threshold);
814  }
815 
816  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(detId);
817  if (towerDetId.null())
818  return;
819  MetaTower& tower = find(towerDetId);
820 
821  // count bad cells and avoid double counting with those from DB (Recovered are counted bad)
822 
823  // somehow misses some
824  // if ( (chStatusForCT == CaloTowersCreationAlgo::BadChan) & (!ecalIsBad) ) ++tower.numBadEcalCells;
825 
826  // a bit slower...
827  if (chStatusForCT == CaloTowersCreationAlgo::BadChan) {
828  auto thisEcalSevLvl = theEcalSevLvlAlgo->severityLevel(detId);
829  // check if the Ecal severity is ok to keep
830  auto sevit = std::find(theEcalSeveritiesToBeExcluded.begin(), theEcalSeveritiesToBeExcluded.end(), thisEcalSevLvl);
831  if (sevit == theEcalSeveritiesToBeExcluded.end())
832  ++tower.numBadEcalCells; // notinDB
833  }
834 
835  // if (chStatusForCT != CaloTowersCreationAlgo::BadChan && energy >= threshold) {
836  if (chStatusForCT != CaloTowersCreationAlgo::BadChan && passEmThreshold) {
837  tower.E_em += e;
838  tower.E += e;
839 
840  if (chStatusForCT == CaloTowersCreationAlgo::RecoveredChan) {
841  tower.numRecEcalCells += 1;
842  } else if (chStatusForCT == CaloTowersCreationAlgo::ProblematicChan) {
843  tower.numProbEcalCells += 1;
844  }
845 
846  // change when full status info is available
847  // for now use only good channels
848 
849  // add e>0 check (new options allow e<0)
850  if (chStatusForCT == CaloTowersCreationAlgo::GoodChan && e > 0) {
851  tower.emSumTimeTimesE += (e * recHit->time());
852  tower.emSumEForTime += e; // see above
853  }
854 
855  std::pair<DetId, float> mc(detId, e);
856  tower.metaConstituents.push_back(mc);
857  }
858 } // end of assignHitEcal method
859 
860 // This method is not flexible enough for the new CaloTower format.
861 // For now make a quick compatibility "fix" : WILL NOT WORK CORRECTLY with anything
862 // except the default simple p4 assignment!!!
863 // Must be rewritten for full functionality.
864 void CaloTowersCreationAlgo::rescale(const CaloTower* ct) {
865  double threshold, weight;
866  CaloTowerDetId towerDetId = theTowerConstituentsMap->towerOf(ct->id());
867  if (towerDetId.null())
868  return;
869  MetaTower& tower = find(towerDetId);
870 
871  tower.E_em = 0.;
872  tower.E_had = 0.;
873  tower.E_outer = 0.;
874  for (unsigned int i = 0; i < ct->constituentsSize(); i++) {
875  DetId detId = ct->constituent(i);
876  getThresholdAndWeight(detId, threshold, weight);
877  DetId::Detector det = detId.det();
878  if (det == DetId::Ecal) {
879  tower.E_em = ct->emEnergy() * weight;
880  } else {
881  HcalDetId hcalDetId(detId);
882  if (hcalDetId.subdet() == HcalForward) {
883  if (hcalDetId.depth() == 1)
884  tower.E_em = ct->emEnergy() * weight;
885  if (hcalDetId.depth() == 2)
886  tower.E_had = ct->hadEnergy() * weight;
887  } else if (hcalDetId.subdet() == HcalOuter) {
888  tower.E_outer = ct->outerEnergy() * weight;
889  } else {
890  tower.E_had = ct->hadEnergy() * weight;
891  }
892  }
893  tower.E = tower.E_had + tower.E_em + tower.E_outer;
894 
895  // this is to be compliant with the new MetaTower setup
896  // used only for the default simple vector assignment
897  std::pair<DetId, float> mc(detId, 0);
898  tower.metaConstituents.push_back(mc);
899  }
900 
901  // preserve time inforamtion
902  tower.emSumTimeTimesE = ct->ecalTime();
903  tower.hadSumTimeTimesE = ct->hcalTime();
904  tower.emSumEForTime = 1.0;
905  tower.hadSumEForTime = 1.0;
906 }
907 
908 CaloTowersCreationAlgo::MetaTower& CaloTowersCreationAlgo::find(const CaloTowerDetId& detId) {
909  if (theTowerMap.empty()) {
910  theTowerMap.resize(theTowerTopology->sizeForDenseIndexing());
911  }
912 
913  auto& mt = theTowerMap[theTowerTopology->denseIndex(detId)];
914 
915  if (mt.empty()) {
916  mt.id = detId;
917  mt.metaConstituents.reserve(detId.ietaAbs() < theTowerTopology->firstHFRing() ? 12 : 2);
918  ++theTowerMapSize;
919  }
920 
921  return mt;
922 }
923 
924 void CaloTowersCreationAlgo::convert(const CaloTowerDetId& id, const MetaTower& mt, CaloTowerCollection& collection) {
925  assert(id.rawId() != 0);
926 
927  double ecalThres = (id.ietaAbs() <= 17) ? (theEBSumThreshold) : (theEESumThreshold);
928  double E = mt.E;
929  double E_em = mt.E_em;
930  double E_had = mt.E_had;
931  double E_outer = mt.E_outer;
932 
933  // Note: E_outer is used to save HO energy OR energy in the outermost depths in endcap region
934  // In the methods with separate treatment of EM and HAD components:
935  // - HO is not used to determine direction, however HO energy is added to get "total had energy"
936  // => Check if the tower is within HO coverage before adding E_outer to the "total had" energy
937  // else the energy will be double counted
938  // When summing up the energy of the tower these checks are performed in the loops over RecHits
939 
940  std::vector<std::pair<DetId, float> > metaContains = mt.metaConstituents;
941  if (id.ietaAbs() < theTowerTopology->firstHFRing() && E_em < ecalThres) { // ignore EM threshold in HF
942  E -= E_em;
943  E_em = 0;
944  std::vector<std::pair<DetId, float> > metaContains_noecal;
945 
946  for (std::vector<std::pair<DetId, float> >::iterator i = metaContains.begin(); i != metaContains.end(); ++i)
947  if (i->first.det() != DetId::Ecal)
948  metaContains_noecal.push_back(*i);
949  metaContains.swap(metaContains_noecal);
950  }
951  if (id.ietaAbs() < theTowerTopology->firstHFRing() && E_had < theHcalThreshold) {
952  E -= E_had;
953 
954  if (theHOIsUsed && id.ietaAbs() <= theTowerTopology->lastHORing())
955  E -= E_outer; // not subtracted before, think it should be done
956 
957  E_had = 0;
958  E_outer = 0;
959  std::vector<std::pair<DetId, float> > metaContains_nohcal;
960 
961  for (std::vector<std::pair<DetId, float> >::iterator i = metaContains.begin(); i != metaContains.end(); ++i)
962  if (i->first.det() != DetId::Hcal)
963  metaContains_nohcal.push_back(*i);
964  metaContains.swap(metaContains_nohcal);
965  }
966 
967  if (metaContains.empty())
968  return;
969 
970  if (missingHcalRescaleFactorForEcal > 0 && E_had == 0 && E_em > 0) {
971  auto match = hcalDropChMap.find(id);
972  if (match != hcalDropChMap.end() && match->second.second) {
973  E_had = missingHcalRescaleFactorForEcal * E_em;
974  E += E_had;
975  }
976  }
977 
978  double E_had_tot = (theHOIsUsed && id.ietaAbs() <= theTowerTopology->lastHORing()) ? E_had + E_outer : E_had;
979 
980  // create CaloTower using the selected algorithm
981 
982  GlobalPoint emPoint, hadPoint;
983 
984  // this is actually a 4D vector
985  Basic3DVectorF towerP4;
986  bool massless = true;
987  // float mass1=0;
988  float mass2 = 0;
989 
990  // conditional assignment of depths for barrel/endcap
991  // Some additional tuning may be required in the transitional region
992  // 14<|iEta|<19
993  double momEmDepth = 0.;
994  double momHadDepth = 0.;
995  if (id.ietaAbs() <= 17) {
996  momHadDepth = theMomHBDepth;
997  momEmDepth = theMomEBDepth;
998  } else {
999  momHadDepth = theMomHEDepth;
1000  momEmDepth = theMomEEDepth;
1001  }
1002 
1003  switch (theMomConstrMethod) {
1004  // FIXME : move to simple cartesian algebra
1005  case 0: { // Simple 4-momentum assignment
1006  GlobalPoint p = theTowerGeometry->getGeometry(id)->getPosition();
1007  towerP4 = p.basicVector().unit();
1008  towerP4[3] = 1.f; // energy
1009  towerP4 *= E;
1010 
1011  // double pf=1.0/cosh(p.eta());
1012  // if (E>0) towerP4 = CaloTower::PolarLorentzVector(E*pf, p.eta(), p.phi(), 0);
1013 
1014  emPoint = p;
1015  hadPoint = p;
1016  } // end case 0
1017  break;
1018 
1019  case 1: { // separate 4-vectors for ECAL, HCAL, add to get the 4-vector of the tower (=>tower has mass!)
1020  if (id.ietaAbs() < theTowerTopology->firstHFRing()) {
1021  Basic3DVectorF emP4;
1022  if (E_em > 0) {
1023  emPoint = emShwrPos(metaContains, momEmDepth, E_em);
1024  emP4 = emPoint.basicVector().unit();
1025  emP4[3] = 1.f; // energy
1026  towerP4 = emP4 * E_em;
1027 
1028  // double emPf = 1.0/cosh(emPoint.eta());
1029  // towerP4 += CaloTower::PolarLorentzVector(E_em*emPf, emPoint.eta(), emPoint.phi(), 0);
1030  }
1031  if ((E_had + E_outer) > 0) {
1032  massless = (E_em <= 0);
1033  hadPoint = hadShwrPos(id, momHadDepth);
1034  auto lP4 = hadPoint.basicVector().unit();
1035  lP4[3] = 1.f; // energy
1036  if (!massless) {
1037  auto diff = lP4 - emP4;
1038  mass2 = std::sqrt(E_em * E_had_tot * diff.mag2());
1039  }
1040  lP4 *= E_had_tot;
1041  towerP4 += lP4;
1042  /*
1043  if (!massless) {
1044  auto p = towerP4;
1045  double m2 = double(p[3]*p[3]) - double(p[0]*p[0])+double(p[1]*p[1])+double(p[2]*p[2]); mass1 = m2>0 ? std::sqrt(m2) : 0;
1046  }
1047  */
1048  // double hadPf = 1.0/cosh(hadPoint.eta());
1049  // if (E_had_tot>0) {
1050  // towerP4 += CaloTower::PolarLorentzVector(E_had_tot*hadPf, hadPoint.eta(), hadPoint.phi(), 0);
1051  // }
1052  }
1053  } else { // forward detector: use the CaloTower position
1054  GlobalPoint p = theTowerGeometry->getGeometry(id)->getPosition();
1055  towerP4 = p.basicVector().unit();
1056  towerP4[3] = 1.f; // energy
1057  towerP4 *= E;
1058  // double pf=1.0/cosh(p.eta());
1059  // if (E>0) towerP4 = CaloTower::PolarLorentzVector(E*pf, p.eta(), p.phi(), 0); // simple momentum assignment, same position
1060  emPoint = p;
1061  hadPoint = p;
1062  }
1063  } // end case 1
1064  break;
1065 
1066  case 2: { // use ECAL position for the tower (when E_cal>0), else default CaloTower position (massless tower)
1067  if (id.ietaAbs() < theTowerTopology->firstHFRing()) {
1068  if (E_em > 0)
1069  emPoint = emShwrLogWeightPos(metaContains, momEmDepth, E_em);
1070  else
1071  emPoint = theTowerGeometry->getGeometry(id)->getPosition();
1072  towerP4 = emPoint.basicVector().unit();
1073  towerP4[3] = 1.f; // energy
1074  towerP4 *= E;
1075 
1076  // double sumPf = 1.0/cosh(emPoint.eta());
1078 
1079  hadPoint = emPoint;
1080  } else { // forward detector: use the CaloTower position
1081  GlobalPoint p = theTowerGeometry->getGeometry(id)->getPosition();
1082  towerP4 = p.basicVector().unit();
1083  towerP4[3] = 1.f; // energy
1084  towerP4 *= E;
1085 
1086  // double pf=1.0/cosh(p.eta());
1087  // if (E>0) towerP4 = CaloTower::PolarLorentzVector(E*pf, p.eta(), p.phi(), 0); // simple momentum assignment, same position
1088  emPoint = p;
1089  hadPoint = p;
1090  }
1091  } // end case 2
1092  break;
1093 
1094  } // end of decision on p4 reconstruction method
1095 
1096  // insert in collection (remove and return if below threshold)
1097  if
1098  UNLIKELY((towerP4[3] == 0) & (E_outer > 0)) {
1099  float val = theHOIsUsed ? 0 : 1E-9; // to keep backwards compatibility for theHOIsUsed == true
1100  collection.emplace_back(id,
1101  E_em,
1102  E_had,
1103  E_outer,
1104  -1,
1105  -1,
1106  CaloTower::PolarLorentzVector(val, hadPoint.eta(), hadPoint.phi(), 0),
1107  emPoint,
1108  hadPoint);
1109  }
1110  else {
1111  collection.emplace_back(
1112  id, E_em, E_had, E_outer, -1, -1, GlobalVector(towerP4), towerP4[3], mass2, emPoint, hadPoint);
1113  }
1114  auto& caloTower = collection.back();
1115 
1116  // if (!massless) std::cout << "massive " << id <<' ' << mass1 <<' ' << mass2 <<' ' << caloTower.mass() << std::endl;
1117  // std::cout << "CaloTowerVI " <<theMomConstrMethod <<' ' << id <<' '<< E_em <<' '<< E_had <<' '<< E_outer <<' '<< GlobalVector(towerP4) <<' '<< towerP4[3] <<' '<< emPoint <<' '<< hadPoint << std::endl;
1118  //if (towerP4[3]==0) std::cout << "CaloTowerVIzero " << theEcutTower << ' ' << collection.back().eta() <<' '<< collection.back().phi() << std::endl;
1119 
1120  if (caloTower.energy() < theEcutTower) {
1121  collection.pop_back();
1122  return;
1123  }
1124 
1125  // set the timings
1126  float ecalTime = (mt.emSumEForTime > 0) ? mt.emSumTimeTimesE / mt.emSumEForTime : -9999;
1127  float hcalTime = (mt.hadSumEForTime > 0) ? mt.hadSumTimeTimesE / mt.hadSumEForTime : -9999;
1128  caloTower.setEcalTime(compactTime(ecalTime));
1129  caloTower.setHcalTime(compactTime(hcalTime));
1130  //add topology info
1131  caloTower.setHcalSubdet(theTowerTopology->lastHBRing(),
1132  theTowerTopology->lastHERing(),
1133  theTowerTopology->lastHFRing(),
1134  theTowerTopology->lastHORing());
1135 
1136  // set the CaloTower status word =====================================
1137  // Channels must be counter exclusively in the defined cathegories
1138  // "Bad" channels (not used in energy assignment) can be flagged during
1139  // CaloTower creation only if specified in the configuration file
1140 
1141  unsigned int numBadHcalChan = mt.numBadHcalCells;
1142  // unsigned int numBadEcalChan = mt.numBadEcalCells;
1143  unsigned int numBadEcalChan = 0; //
1144 
1145  unsigned int numRecHcalChan = mt.numRecHcalCells;
1146  unsigned int numRecEcalChan = mt.numRecEcalCells;
1147  unsigned int numProbHcalChan = mt.numProbHcalCells;
1148  unsigned int numProbEcalChan = mt.numProbEcalCells;
1149 
1150  // now add dead/off/... channels not used in RecHit reconstruction for HCAL
1151  HcalDropChMap::iterator dropChItr = hcalDropChMap.find(id);
1152  if (dropChItr != hcalDropChMap.end())
1153  numBadHcalChan += dropChItr->second.first;
1154 
1155  // for ECAL the number of all bad channels is obtained here -----------------------
1156 
1157  /*
1158  // old hyper slow algorithm
1159  // get all possible constituents of the tower
1160  std::vector<DetId> allConstituents = theTowerConstituentsMap->constituentsOf(id);
1161 
1162  for (std::vector<DetId>::iterator ac_it=allConstituents.begin();
1163  ac_it!=allConstituents.end(); ++ac_it) {
1164 
1165  if (ac_it->det()!=DetId::Ecal) continue;
1166 
1167  int thisEcalSevLvl = -999;
1168 
1169  if (ac_it->subdetId() == EcalBarrel && theEbHandle.isValid()) {
1170  thisEcalSevLvl = theEcalSevLvlAlgo->severityLevel( *ac_it, *theEbHandle);//, *theEcalChStatus);
1171  }
1172  else if (ac_it->subdetId() == EcalEndcap && theEeHandle.isValid()) {
1173  thisEcalSevLvl = theEcalSevLvlAlgo->severityLevel( *ac_it, *theEeHandle);//, *theEcalChStatus);
1174  }
1175 
1176  // check if the Ecal severity is ok to keep
1177  std::vector<int>::const_iterator sevit = std::find(theEcalSeveritiesToBeExcluded.begin(),
1178  theEcalSeveritiesToBeExcluded.end(),
1179  thisEcalSevLvl);
1180  if (sevit!=theEcalSeveritiesToBeExcluded.end()) {
1181  ++numBadEcalChan;
1182  }
1183 
1184  }
1185 
1186  // compare with fast version
1187 
1188  // hcal:
1189  int inEcals[2] = {0,0};
1190  for (std::vector<std::pair<DetId,float> >::iterator i=metaContains.begin(); i!=metaContains.end(); ++i) {
1191  DetId detId = i->first;
1192  if(detId.det() == DetId::Ecal){
1193  if( detId.subdetId()==EcalBarrel ) inEcals[0] =1;
1194  else if( detId.subdetId()==EcalEndcap ) inEcals[1] =1;
1195  }
1196  }
1197 
1198  auto numBadEcalChanNew = ecalBadChs[theTowerTopology->denseIndex(id)]+mt.numBadEcalCells; // - mt.numRecEcalCells
1199  if (int(numBadEcalChanNew)!=int(numBadEcalChan)) {
1200  std::cout << "VI wrong " << ((inEcals[1]==1) ? "EE" : "" ) << id << " " << numBadEcalChanNew << " " << numBadEcalChan
1201  << " " << mt.numBadEcalCells << " " << mt.numRecEcalCells << std::endl;
1202  }
1203  */
1204 
1205  numBadEcalChan = ecalBadChs[theTowerTopology->denseIndex(id)] + mt.numBadEcalCells; // - mt.numRecEcalCells
1206 
1207  //--------------------------------------------------------------------------------------
1208 
1209  caloTower.setCaloTowerStatus(
1210  numBadHcalChan, numBadEcalChan, numRecHcalChan, numRecEcalChan, numProbHcalChan, numProbEcalChan);
1211 
1212  double maxCellE = -999.0; // for storing the hottest cell E in the calotower
1213 
1214  std::vector<DetId> contains;
1215  contains.reserve(metaContains.size());
1216  for (std::vector<std::pair<DetId, float> >::iterator i = metaContains.begin(); i != metaContains.end(); ++i) {
1217  contains.push_back(i->first);
1218 
1219  if (maxCellE < i->second) {
1220  // need an extra check because of the funny towers that are empty except for the presence of an HO
1221  // hit in the constituents (JetMET wanted them saved)
1222  // This constituent is only used for storing the tower, but should not be concidered as a hot cell canditate for
1223  // configurations with useHO = false
1224 
1225  if (i->first.det() == DetId::Ecal) { // ECAL
1226  maxCellE = i->second;
1227  } else { // HCAL
1228  if (HcalDetId(i->first).subdet() != HcalOuter)
1229  maxCellE = i->second;
1230  else if (theHOIsUsed)
1231  maxCellE = i->second;
1232  }
1233 
1234  } // found higher E cell
1235 
1236  } // loop over matacontains
1237 
1238  caloTower.setConstituents(std::move(contains));
1239  caloTower.setHottestCellE(maxCellE);
1240 
1241  // std::cout << "CaloTowerVI " << nalgo << ' ' << caloTower.id() << ((inEcals[1]==1) ? "EE " : " " ) << caloTower.pt() << ' ' << caloTower.et() << ' ' << caloTower.mass() << ' '
1242  // << caloTower.constituentsSize() <<' '<< caloTower.towerStatusWord() << std::endl;
1243 }
1244 
1245 void CaloTowersCreationAlgo::getThresholdAndWeight(const DetId& detId, double& threshold, double& weight) const {
1246  DetId::Detector det = detId.det();
1247  weight = 0; // in case the hit is not identified
1248 
1249  if (det == DetId::Ecal) {
1250  // may or may not be EB. We'll find out.
1251 
1252  EcalSubdetector subdet = (EcalSubdetector)(detId.subdetId());
1253  if (subdet == EcalBarrel) {
1254  threshold = theEBthreshold;
1255  weight = theEBweight;
1256  if (weight <= 0.) {
1257  ROOT::Math::Interpolator my(theEBGrid, theEBWeights, ROOT::Math::Interpolation::kAKIMA);
1258  weight = my.Eval(theEBEScale);
1259  }
1260  } else if (subdet == EcalEndcap) {
1261  threshold = theEEthreshold;
1262  weight = theEEweight;
1263  if (weight <= 0.) {
1264  ROOT::Math::Interpolator my(theEEGrid, theEEWeights, ROOT::Math::Interpolation::kAKIMA);
1265  weight = my.Eval(theEEEScale);
1266  }
1267  }
1268  } else if (det == DetId::Hcal) {
1269  HcalDetId hcalDetId(detId);
1270  HcalSubdetector subdet = hcalDetId.subdet();
1271  int depth = hcalDetId.depth();
1272 
1273  if (subdet == HcalBarrel) {
1274  threshold = (depth == 1) ? theHBthreshold1 : (depth == 2) ? theHBthreshold2 : theHBthreshold;
1275  weight = theHBweight;
1276  if (weight <= 0.) {
1277  ROOT::Math::Interpolator my(theHBGrid, theHBWeights, ROOT::Math::Interpolation::kAKIMA);
1278  weight = my.Eval(theHBEScale);
1279  }
1280  }
1281 
1282  else if (subdet == HcalEndcap) {
1283  // check if it's single or double tower
1284  if (hcalDetId.ietaAbs() < theHcalTopology->firstHEDoublePhiRing()) {
1285  threshold = (depth == 1) ? theHESthreshold1 : theHESthreshold;
1286  weight = theHESweight;
1287  if (weight <= 0.) {
1288  ROOT::Math::Interpolator my(theHESGrid, theHESWeights, ROOT::Math::Interpolation::kAKIMA);
1289  weight = my.Eval(theHESEScale);
1290  }
1291  } else {
1292  threshold = (depth == 1) ? theHEDthreshold1 : theHEDthreshold;
1293  weight = theHEDweight;
1294  if (weight <= 0.) {
1295  ROOT::Math::Interpolator my(theHEDGrid, theHEDWeights, ROOT::Math::Interpolation::kAKIMA);
1296  weight = my.Eval(theHEDEScale);
1297  }
1298  }
1299  }
1300 
1301  else if (subdet == HcalOuter) {
1302  //check if it's ring 0 or +1 or +2 or -1 or -2
1303  if (hcalDetId.ietaAbs() <= 4)
1304  threshold = theHOthreshold0;
1305  else if (hcalDetId.ieta() < 0) {
1306  // set threshold for ring -1 or -2
1307  threshold = (hcalDetId.ietaAbs() <= 10) ? theHOthresholdMinus1 : theHOthresholdMinus2;
1308  } else {
1309  // set threshold for ring +1 or +2
1310  threshold = (hcalDetId.ietaAbs() <= 10) ? theHOthresholdPlus1 : theHOthresholdPlus2;
1311  }
1312  weight = theHOweight;
1313  if (weight <= 0.) {
1314  ROOT::Math::Interpolator my(theHOGrid, theHOWeights, ROOT::Math::Interpolation::kAKIMA);
1315  weight = my.Eval(theHOEScale);
1316  }
1317  }
1318 
1319  else if (subdet == HcalForward) {
1320  if (hcalDetId.depth() == 1) {
1321  threshold = theHF1threshold;
1322  weight = theHF1weight;
1323  if (weight <= 0.) {
1324  ROOT::Math::Interpolator my(theHF1Grid, theHF1Weights, ROOT::Math::Interpolation::kAKIMA);
1325  weight = my.Eval(theHF1EScale);
1326  }
1327  } else {
1328  threshold = theHF2threshold;
1329  weight = theHF2weight;
1330  if (weight <= 0.) {
1331  ROOT::Math::Interpolator my(theHF2Grid, theHF2Weights, ROOT::Math::Interpolation::kAKIMA);
1332  weight = my.Eval(theHF2EScale);
1333  }
1334  }
1335  }
1336  } else {
1337  edm::LogError("CaloTowersCreationAlgo") << "Bad cell: " << det << std::endl;
1338  }
1339 }
1340 
1341 void CaloTowersCreationAlgo::setEBEScale(double scale) {
1342  if (scale > 0.00001)
1343  *&theEBEScale = scale;
1344  else
1345  *&theEBEScale = 50.;
1346 }
1347 
1348 void CaloTowersCreationAlgo::setEEEScale(double scale) {
1349  if (scale > 0.00001)
1350  *&theEEEScale = scale;
1351  else
1352  *&theEEEScale = 50.;
1353 }
1354 
1355 void CaloTowersCreationAlgo::setHBEScale(double scale) {
1356  if (scale > 0.00001)
1357  *&theHBEScale = scale;
1358  else
1359  *&theHBEScale = 50.;
1360 }
1361 
1362 void CaloTowersCreationAlgo::setHESEScale(double scale) {
1363  if (scale > 0.00001)
1364  *&theHESEScale = scale;
1365  else
1366  *&theHESEScale = 50.;
1367 }
1368 
1369 void CaloTowersCreationAlgo::setHEDEScale(double scale) {
1370  if (scale > 0.00001)
1371  *&theHEDEScale = scale;
1372  else
1373  *&theHEDEScale = 50.;
1374 }
1375 
1376 void CaloTowersCreationAlgo::setHOEScale(double scale) {
1377  if (scale > 0.00001)
1378  *&theHOEScale = scale;
1379  else
1380  *&theHOEScale = 50.;
1381 }
1382 
1383 void CaloTowersCreationAlgo::setHF1EScale(double scale) {
1384  if (scale > 0.00001)
1385  *&theHF1EScale = scale;
1386  else
1387  *&theHF1EScale = 50.;
1388 }
1389 
1390 void CaloTowersCreationAlgo::setHF2EScale(double scale) {
1391  if (scale > 0.00001)
1392  *&theHF2EScale = scale;
1393  else
1394  *&theHF2EScale = 50.;
1395 }
1396 
1397 GlobalPoint CaloTowersCreationAlgo::emCrystalShwrPos(DetId detId, float fracDepth) {
1398  auto cellGeometry = theGeometry->getGeometry(detId);
1399  GlobalPoint point = cellGeometry->getPosition(); // face of the cell
1400 
1401  if (fracDepth <= 0)
1402  return point;
1403  if (fracDepth > 1)
1404  fracDepth = 1;
1405 
1406  const GlobalPoint& backPoint = cellGeometry->getBackPoint();
1407  point += fracDepth * (backPoint - point);
1408 
1409  return point;
1410 }
1411 
1412 GlobalPoint CaloTowersCreationAlgo::hadSegmentShwrPos(DetId detId, float fracDepth) {
1413  // same code as above
1414  return emCrystalShwrPos(detId, fracDepth);
1415 }
1416 
1417 GlobalPoint CaloTowersCreationAlgo::hadShwrPos(const std::vector<std::pair<DetId, float> >& metaContains,
1418  float fracDepth,
1419  double hadE) {
1420  // this is based on available RecHits, can lead to different actual depths if
1421  // hits in multi-depth towers are not all there
1422 #ifdef EDM_ML_DEBUG
1423  std::cout << "hadShwrPos called with " << metaContains.size() << " elements and energy " << hadE << ":" << fracDepth
1424  << std::endl;
1425 #endif
1426  if (hadE <= 0)
1427  return GlobalPoint(0, 0, 0);
1428 
1429  double hadX = 0.0;
1430  double hadY = 0.0;
1431  double hadZ = 0.0;
1432 
1433  int nConst = 0;
1434 
1435  std::vector<std::pair<DetId, float> >::const_iterator mc_it = metaContains.begin();
1436  for (; mc_it != metaContains.end(); ++mc_it) {
1437  if (mc_it->first.det() != DetId::Hcal)
1438  continue;
1439  // do not use HO for deirection calculations for now
1440  if (HcalDetId(mc_it->first).subdet() == HcalOuter)
1441  continue;
1442  ++nConst;
1443 
1444  GlobalPoint p = hadSegmentShwrPos(mc_it->first, fracDepth);
1445 
1446  // longitudinal segmentation: do not weight by energy,
1447  // get the geometrical position
1448  hadX += p.x();
1449  hadY += p.y();
1450  hadZ += p.z();
1451  }
1452 
1453  return GlobalPoint(hadX / nConst, hadY / nConst, hadZ / nConst);
1454 }
1455 
1456 GlobalPoint CaloTowersCreationAlgo::hadShwrPos(CaloTowerDetId towerId, float fracDepth) {
1457  // set depth using geometry of cells that are associated with the
1458  // tower (regardless if they have non-zero energies)
1459 
1460  // if (hadE <= 0) return GlobalPoint(0, 0, 0);
1461 #ifdef EDM_ML_DEBUG
1462  std::cout << "hadShwrPos " << towerId << " frac " << fracDepth << std::endl;
1463 #endif
1464  if (fracDepth < 0)
1465  fracDepth = 0;
1466  else if (fracDepth > 1)
1467  fracDepth = 1;
1468 
1469  GlobalPoint point(0, 0, 0);
1470 
1471  int iEta = towerId.ieta();
1472  int iPhi = towerId.iphi();
1473 
1474  HcalDetId frontCellId, backCellId;
1475 
1476  if (towerId.ietaAbs() >= theTowerTopology->firstHFRing()) {
1477  // forward, take the geometry for long fibers
1478  frontCellId = HcalDetId(HcalForward, towerId.zside() * theTowerTopology->convertCTtoHcal(abs(iEta)), iPhi, 1);
1479  backCellId = HcalDetId(HcalForward, towerId.zside() * theTowerTopology->convertCTtoHcal(abs(iEta)), iPhi, 1);
1480  } else {
1481  //use constituents map
1482  std::vector<DetId> items = theTowerConstituentsMap->constituentsOf(towerId);
1483  int frontDepth = 1000;
1484  int backDepth = -1000;
1485  for (unsigned i = 0; i < items.size(); i++) {
1486  if (items[i].det() != DetId::Hcal)
1487  continue;
1488  HcalDetId hid(items[i]);
1489  if (hid.subdet() == HcalOuter)
1490  continue;
1491  if (!theHcalTopology->validHcal(hid, 2))
1492  continue;
1493 
1494  if (theHcalTopology->idFront(hid).depth() < frontDepth) {
1495  frontCellId = hid;
1496  frontDepth = theHcalTopology->idFront(hid).depth();
1497  }
1498  if (theHcalTopology->idBack(hid).depth() > backDepth) {
1499  backCellId = hid;
1500  backDepth = theHcalTopology->idBack(hid).depth();
1501  }
1502  }
1503 #ifdef EDM_ML_DEBUG
1504  std::cout << "Front " << frontCellId << " Back " << backCellId << " Depths " << frontDepth << ":" << backDepth
1505  << std::endl;
1506 #endif
1507  //fix for tower 28/29 - no tower 29 at highest depths
1508  if (towerId.ietaAbs() == theTowerTopology->lastHERing() && (theHcalPhase == 0 || theHcalPhase == 1)) {
1509  CaloTowerDetId towerId28(towerId.ieta() - towerId.zside(), towerId.iphi());
1510  std::vector<DetId> items28 = theTowerConstituentsMap->constituentsOf(towerId28);
1511 #ifdef EDM_ML_DEBUG
1512  std::cout << towerId28 << " with " << items28.size() << " constituents:";
1513  for (unsigned k = 0; k < items28.size(); ++k)
1514  if (items28[k].det() == DetId::Hcal)
1515  std::cout << " " << HcalDetId(items28[k]);
1516  std::cout << std::endl;
1517 #endif
1518  for (unsigned i = 0; i < items28.size(); i++) {
1519  if (items28[i].det() != DetId::Hcal)
1520  continue;
1521  HcalDetId hid(items28[i]);
1522  if (hid.subdet() == HcalOuter)
1523  continue;
1524 
1525  if (theHcalTopology->idBack(hid).depth() > backDepth) {
1526  backCellId = hid;
1527  backDepth = theHcalTopology->idBack(hid).depth();
1528  }
1529  }
1530  }
1531 #ifdef EDM_ML_DEBUG
1532  std::cout << "Back " << backDepth << " ID " << backCellId << std::endl;
1533 #endif
1534  }
1535  point = hadShwPosFromCells(DetId(frontCellId), DetId(backCellId), fracDepth);
1536 
1537  return point;
1538 }
1539 
1540 GlobalPoint CaloTowersCreationAlgo::hadShwPosFromCells(DetId frontCellId, DetId backCellId, float fracDepth) {
1541  // uses the "front" and "back" cells
1542  // to determine the axis. point set by the predefined depth.
1543 
1544  HcalDetId hid1(frontCellId), hid2(backCellId);
1545  if (theHcalTopology->getMergePositionFlag()) {
1546  hid1 = theHcalTopology->idFront(frontCellId);
1547 #ifdef EDM_ML_DEBUG
1548  std::cout << "Front " << HcalDetId(frontCellId) << " " << hid1 << "\n";
1549 #endif
1550  hid2 = theHcalTopology->idBack(backCellId);
1551 #ifdef EDM_ML_DEBUG
1552  std::cout << "Back " << HcalDetId(backCellId) << " " << hid2 << "\n";
1553 #endif
1554  }
1555 
1556  auto frontCellGeometry = theGeometry->getGeometry(DetId(hid1));
1557  auto backCellGeometry = theGeometry->getGeometry(DetId(hid2));
1558 
1559  GlobalPoint point = frontCellGeometry->getPosition();
1560  const GlobalPoint& backPoint = backCellGeometry->getBackPoint();
1561 
1562  point += fracDepth * (backPoint - point);
1563 
1564  return point;
1565 }
1566 
1567 GlobalPoint CaloTowersCreationAlgo::emShwrPos(const std::vector<std::pair<DetId, float> >& metaContains,
1568  float fracDepth,
1569  double emE) {
1570  if (emE <= 0)
1571  return GlobalPoint(0, 0, 0);
1572 
1573  double emX = 0.0;
1574  double emY = 0.0;
1575  double emZ = 0.0;
1576 
1577  double eSum = 0;
1578 
1579  std::vector<std::pair<DetId, float> >::const_iterator mc_it = metaContains.begin();
1580  for (; mc_it != metaContains.end(); ++mc_it) {
1581  if (mc_it->first.det() != DetId::Ecal)
1582  continue;
1583  GlobalPoint p = emCrystalShwrPos(mc_it->first, fracDepth);
1584  double e = mc_it->second;
1585 
1586  if (e > 0) {
1587  emX += p.x() * e;
1588  emY += p.y() * e;
1589  emZ += p.z() * e;
1590  eSum += e;
1591  }
1592  }
1593 
1594  return GlobalPoint(emX / eSum, emY / eSum, emZ / eSum);
1595 }
1596 
1597 GlobalPoint CaloTowersCreationAlgo::emShwrLogWeightPos(const std::vector<std::pair<DetId, float> >& metaContains,
1598  float fracDepth,
1599  double emE) {
1600  double emX = 0.0;
1601  double emY = 0.0;
1602  double emZ = 0.0;
1603 
1604  double weight = 0;
1605  double sumWeights = 0;
1606  double sumEmE = 0; // add crystals with E/E_EM > 1.5%
1607  double crystalThresh = 0.015 * emE;
1608 
1609  std::vector<std::pair<DetId, float> >::const_iterator mc_it = metaContains.begin();
1610  for (; mc_it != metaContains.end(); ++mc_it) {
1611  if (mc_it->second < 0)
1612  continue;
1613  if (mc_it->first.det() == DetId::Ecal && mc_it->second > crystalThresh)
1614  sumEmE += mc_it->second;
1615  }
1616 
1617  for (mc_it = metaContains.begin(); mc_it != metaContains.end(); ++mc_it) {
1618  if (mc_it->first.det() != DetId::Ecal || mc_it->second < crystalThresh)
1619  continue;
1620 
1621  GlobalPoint p = emCrystalShwrPos(mc_it->first, fracDepth);
1622 
1623  weight = 4.2 + log(mc_it->second / sumEmE);
1624  sumWeights += weight;
1625 
1626  emX += p.x() * weight;
1627  emY += p.y() * weight;
1628  emZ += p.z() * weight;
1629  }
1630 
1631  return GlobalPoint(emX / sumWeights, emY / sumWeights, emZ / sumWeights);
1632 }
1633 
1634 int CaloTowersCreationAlgo::compactTime(float time) {
1635  const float timeUnit = 0.01; // discretization (ns)
1636 
1637  if (time > 300.0)
1638  return 30000;
1639  if (time < -300.0)
1640  return -30000;
1641 
1642  return int(time / timeUnit + 0.5);
1643 }
1644 
1645 //========================================================
1646 //
1647 // Bad/anomolous cell handling
1648 
1649 void CaloTowersCreationAlgo::makeHcalDropChMap() {
1650  // This method fills the map of number of dead channels for the calotower,
1651  // The key of the map is CaloTowerDetId.
1652  // By definition these channels are not going to be in the RecHit collections.
1653  hcalDropChMap.clear();
1654  std::vector<DetId> allChanInStatusCont = theHcalChStatus->getAllChannels();
1655 
1656 #ifdef EDM_ML_DEBUG
1657  std::cout << "DropChMap with " << allChanInStatusCont.size() << " channels" << std::endl;
1658 #endif
1659  for (std::vector<DetId>::iterator it = allChanInStatusCont.begin(); it != allChanInStatusCont.end(); ++it) {
1660  const uint32_t dbStatusFlag = theHcalChStatus->getValues(*it)->getValue();
1661  if (theHcalSevLvlComputer->dropChannel(dbStatusFlag)) {
1662  DetId id = theHcalTopology->mergedDepthDetId(HcalDetId(*it));
1663 
1664  CaloTowerDetId twrId = theTowerConstituentsMap->towerOf(id);
1665 
1666  hcalDropChMap[twrId].first += 1;
1667 
1668  HcalDetId hid(*it);
1669 
1670  // special case for tower 29: if HCAL hit is in depth 3 add to twr 29 as well
1671  if (hid.subdet() == HcalEndcap && (theHcalPhase == 0 || theHcalPhase == 1) &&
1672  hid.ietaAbs() == theHcalTopology->lastHERing() - 1) {
1673  bool merge = theHcalTopology->mergedDepth29(hid);
1674  if (merge) {
1675  CaloTowerDetId twrId29(twrId.ieta() + twrId.zside(), twrId.iphi());
1676  hcalDropChMap[twrId29].first += 1;
1677  }
1678  }
1679  }
1680  }
1681  // now I know how many bad channels, but I also need to know if there's any good ones
1682  if (missingHcalRescaleFactorForEcal > 0) {
1683  for (auto& pair : hcalDropChMap) {
1684  if (pair.second.first == 0)
1685  continue; // unexpected, but just in case
1686  int ngood = 0, nbad = 0;
1687  for (DetId id : theTowerConstituentsMap->constituentsOf(pair.first)) {
1688  if (id.det() != DetId::Hcal)
1689  continue;
1690  HcalDetId hid(id);
1691  if (hid.subdet() != HcalBarrel && hid.subdet() != HcalEndcap)
1692  continue;
1693  const uint32_t dbStatusFlag = theHcalChStatus->getValues(id)->getValue();
1694  if (dbStatusFlag == 0 || !theHcalSevLvlComputer->dropChannel(dbStatusFlag)) {
1695  ngood += 1;
1696  } else {
1697  nbad += 1; // recount, since pair.second.first may include HO
1698  }
1699  }
1700  if (nbad > 0 && nbad >= ngood) {
1701  //uncomment for debug (may be useful to tune the criteria above)
1702  //CaloTowerDetId id(pair.first);
1703  //std::cout << "CaloTower at ieta = " << id.ieta() << ", iphi " << id.iphi() << ": set Hcal as not efficient (ngood =" << ngood << ", nbad = " << nbad << ")" << std::endl;
1704  pair.second.second = true;
1705  }
1706  }
1707  }
1708 }
1709 
1710 void CaloTowersCreationAlgo::makeEcalBadChs() {
1711  // std::cout << "VI making EcalBadChs ";
1712 
1713  // for ECAL the number of all bad channels is obtained here -----------------------
1714 
1715  for (auto ind = 0U; ind < theTowerTopology->sizeForDenseIndexing(); ++ind) {
1716  auto& numBadEcalChan = ecalBadChs[ind];
1717  numBadEcalChan = 0;
1718  auto id = theTowerTopology->detIdFromDenseIndex(ind);
1719 
1720  // this is utterly slow... (can be optmized if really needed)
1721 
1722  // get all possible constituents of the tower
1723  std::vector<DetId> allConstituents = theTowerConstituentsMap->constituentsOf(id);
1724 
1725  for (std::vector<DetId>::iterator ac_it = allConstituents.begin(); ac_it != allConstituents.end(); ++ac_it) {
1726  if (ac_it->det() != DetId::Ecal)
1727  continue;
1728 
1729  auto thisEcalSevLvl = theEcalSevLvlAlgo->severityLevel(*ac_it);
1730 
1731  // check if the Ecal severity is ok to keep
1732  std::vector<int>::const_iterator sevit =
1733  std::find(theEcalSeveritiesToBeExcluded.begin(), theEcalSeveritiesToBeExcluded.end(), thisEcalSevLvl);
1734  if (sevit != theEcalSeveritiesToBeExcluded.end()) {
1735  ++numBadEcalChan;
1736  }
1737  }
1738 
1739  // if (0!=numBadEcalChan) std::cout << id << ":" << numBadEcalChan << ", ";
1740  }
1741 
1742  /*
1743  int tot=0;
1744  for (auto ind=0U; ind<theTowerTopology->sizeForDenseIndexing(); ++ind) {
1745  if (ecalBadChs[ind]!=0) ++tot;
1746  }
1747  std::cout << " | " << tot << std::endl;
1748  */
1749 }
1750 
1752 
1753 unsigned int CaloTowersCreationAlgo::hcalChanStatusForCaloTower(const CaloRecHit* hit) {
1754  HcalDetId hid(hit->detid());
1755  DetId id = theHcalTopology->idFront(hid);
1756 #ifdef EDM_ML_DEBUG
1757  std::cout << "ChanStatusForCaloTower for " << hid << " to " << HcalDetId(id) << std::endl;
1758 #endif
1759  const uint32_t recHitFlag = hit->flags();
1760  const uint32_t dbStatusFlag = theHcalChStatus->getValues(id)->getValue();
1761 
1762  int severityLevel = theHcalSevLvlComputer->getSeverityLevel(id, recHitFlag, dbStatusFlag);
1763  bool isRecovered = theHcalSevLvlComputer->recoveredRecHit(id, recHitFlag);
1764 
1765  // For use with hits rejected in the default reconstruction
1766  if (useRejectedHitsOnly) {
1767  if (!isRecovered) {
1768  if (severityLevel <= int(theHcalAcceptSeverityLevel) ||
1769  severityLevel > int(theHcalAcceptSeverityLevelForRejectedHit))
1770  return CaloTowersCreationAlgo::IgnoredChan;
1771  // this hit was either already accepted or is worse than
1772  } else {
1773  if (theRecoveredHcalHitsAreUsed || !useRejectedRecoveredHcalHits) {
1774  // skip recovered hits either because they were already used or because there was an explicit instruction
1775  return CaloTowersCreationAlgo::IgnoredChan;
1776  } else if (useRejectedRecoveredHcalHits) {
1777  return CaloTowersCreationAlgo::RecoveredChan;
1778  }
1779 
1780  } // recovered channels
1781 
1782  // clasify channels as problematic: no good hits are supposed to be present in the
1783  // extra rechit collections
1784  return CaloTowersCreationAlgo::ProblematicChan;
1785 
1786  } // treatment of rejected hits
1787 
1788  // this is for the regular reconstruction sequence
1789 
1790  if (severityLevel == 0)
1791  return CaloTowersCreationAlgo::GoodChan;
1792 
1793  if (isRecovered) {
1794  return (theRecoveredHcalHitsAreUsed) ? CaloTowersCreationAlgo::RecoveredChan : CaloTowersCreationAlgo::BadChan;
1795  } else {
1796  if (severityLevel > int(theHcalAcceptSeverityLevel)) {
1797  return CaloTowersCreationAlgo::BadChan;
1798  } else {
1799  return CaloTowersCreationAlgo::ProblematicChan;
1800  }
1801  }
1802 }
1803 
1804 std::tuple<unsigned int, bool> CaloTowersCreationAlgo::ecalChanStatusForCaloTower(const EcalRecHit* hit) {
1805  // const DetId id = hit->detid();
1806 
1807  // uint16_t dbStatus = theEcalChStatus->find(id)->getStatusCode();
1808  // uint32_t rhFlags = hit->flags();
1809  // int severityLevel = theEcalSevLvlAlgo->severityLevel(rhFlags, dbStatus);
1810  // The methods above will become private and cannot be usef for flagging ecal spikes.
1811  // Use the recommended interface - we leave the parameters for spilke removal to be specified by ECAL.
1812 
1813  // int severityLevel = 999;
1814 
1815  EcalRecHit const& rh = *reinterpret_cast<EcalRecHit const*>(hit);
1816  int severityLevel = theEcalSevLvlAlgo->severityLevel(rh);
1817 
1818  // if (id.subdetId() == EcalBarrel) severityLevel = theEcalSevLvlAlgo->severityLevel( id, *theEbHandle);//, *theEcalChStatus);
1819  // else if (id.subdetId() == EcalEndcap) severityLevel = theEcalSevLvlAlgo->severityLevel( id, *theEeHandle);//, *theEcalChStatus);
1820 
1821  // there should be no other ECAL types used in this reconstruction
1822 
1823  // The definition of ECAL severity levels uses categories that
1824  // are similar to the defined for CaloTower. (However, the categorization
1825  // for CaloTowers depends on the specified maximum acceptabel severity and therefore cannnot
1826  // be exact correspondence between the two. ECAL has additional categories describing modes of failure.)
1827  // This approach is different from the initial idea and from
1828  // the implementation for HCAL. Still make the logic similar to HCAL so that one has the ability to
1829  // exclude problematic channels as defined by ECAL.
1830  // For definitions of ECAL severity levels see RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h
1831 
1832  bool isBad = (severityLevel == EcalSeverityLevel::kBad);
1833 
1834  bool isRecovered = (severityLevel == EcalSeverityLevel::kRecovered);
1835 
1836  // check if the severity is compatible with our configuration
1837  // This applies to the "default" tower cleaning
1838  std::vector<int>::const_iterator sevit =
1839  std::find(theEcalSeveritiesToBeExcluded.begin(), theEcalSeveritiesToBeExcluded.end(), severityLevel);
1840  bool accepted = (sevit == theEcalSeveritiesToBeExcluded.end());
1841 
1842  // For use with hits that were rejected in the regular reconstruction:
1843  // This is for creating calotowers with lower level of cleaning by merging
1844  // the information from the default towers and a collection of towers created from
1845  // bad rechits
1846 
1847  if (useRejectedHitsOnly) {
1848  if (!isRecovered) {
1849  if (accepted || std::find(theEcalSeveritiesToBeUsedInBadTowers.begin(),
1850  theEcalSeveritiesToBeUsedInBadTowers.end(),
1851  severityLevel) == theEcalSeveritiesToBeUsedInBadTowers.end())
1852  return std::make_tuple(CaloTowersCreationAlgo::IgnoredChan, isBad);
1853  // this hit was either already accepted, or is not eligible for inclusion
1854  } else {
1855  if (theRecoveredEcalHitsAreUsed || !useRejectedRecoveredEcalHits) {
1856  // skip recovered hits either because they were already used or because there was an explicit instruction
1857  return std::make_tuple(CaloTowersCreationAlgo::IgnoredChan, isBad);
1858  ;
1859  } else if (useRejectedRecoveredEcalHits) {
1860  return std::make_tuple(CaloTowersCreationAlgo::RecoveredChan, isBad);
1861  }
1862 
1863  } // recovered channels
1864 
1865  // clasify channels as problematic
1866  return std::make_tuple(CaloTowersCreationAlgo::ProblematicChan, isBad);
1867 
1868  } // treatment of rejected hits
1869 
1870  // for normal reconstruction
1871  if (severityLevel == EcalSeverityLevel::kGood)
1872  return std::make_tuple(CaloTowersCreationAlgo::GoodChan, false);
1873 
1874  if (isRecovered) {
1875  return std::make_tuple(
1876  (theRecoveredEcalHitsAreUsed) ? CaloTowersCreationAlgo::RecoveredChan : CaloTowersCreationAlgo::BadChan, true);
1877  } else {
1878  return std::make_tuple(accepted ? CaloTowersCreationAlgo::ProblematicChan : CaloTowersCreationAlgo::BadChan, isBad);
1879  }
1880 }
CaloTower::emEnergy
double emEnergy() const
Definition: CaloTower.h:130
change_name.diff
diff
Definition: change_name.py:13
CaloTowersCreationAlgo::emShwrLogWeightPos
GlobalPoint emShwrLogWeightPos(const std::vector< std::pair< DetId, float >> &metaContains, float fracDepth, double totEmE)
Definition: CaloTowersCreationAlgo.cc:1597
CaloTowersCreationAlgo::IgnoredChan
Definition: CaloTowersCreationAlgo.h:390
CaloTowersCreationAlgo::theUseEtEBTresholdFlag
bool theUseEtEBTresholdFlag
Definition: CaloTowersCreationAlgo.h:290
CaloTowerDetId::ietaAbs
int ietaAbs() const
get the absolute value of the tower ieta
Definition: CaloTowerDetId.h:28
CaloTower::setHcalTime
void setHcalTime(int t)
Definition: CaloTower.h:84
HLT_FULL_cff.HF2Weights
HF2Weights
Definition: HLT_FULL_cff.py:8557
CaloTower::ecalTime
float ecalTime() const
Definition: CaloTower.h:196
CaloTowersCreationAlgo::theMomConstrMethod
int theMomConstrMethod
Definition: CaloTowersCreationAlgo.h:362
EcalRecHit
Definition: EcalRecHit.h:15
mps_fire.i
i
Definition: mps_fire.py:428
CaloTowersCreationAlgo::setHF1EScale
void setHF1EScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1383
EcalSeverityLevel::kGood
Definition: EcalSeverityLevel.h:19
CaloTowersCreationAlgo::theEcalSeveritiesToBeUsedInBadTowers
std::vector< int > theEcalSeveritiesToBeUsedInBadTowers
Definition: CaloTowersCreationAlgo.h:349
CaloTowersCreationAlgo::emCrystalShwrPos
GlobalPoint emCrystalShwrPos(DetId detId, float fracDepth)
Definition: CaloTowersCreationAlgo.cc:1397
edm::SortedCollection::const_iterator
std::vector< T >::const_iterator const_iterator
Definition: SortedCollection.h:80
HcalSeverityLevelComputer::recoveredRecHit
bool recoveredRecHit(const DetId &myid, const uint32_t &myflag) const
Definition: HcalSeverityLevelComputer.cc:363
HcalTopology::idFront
HcalDetId idFront(const HcalDetId &id) const
Definition: HcalTopology.h:170
MessageLogger.h
HcalTopology::validHcal
bool validHcal(const HcalDetId &id) const
Definition: HcalTopology.cc:230
CaloTowersCreationAlgo::assignHitHcal
void assignHitHcal(const CaloRecHit *recHit)
Definition: CaloTowersCreationAlgo.cc:534
CaloTowersCreationAlgo::theEBSumThreshold
double theEBSumThreshold
Definition: CaloTowersCreationAlgo.h:310
funct::false
false
Definition: Factorize.h:29
L1EGammaCrystalsEmulatorProducer_cfi.scale
scale
Definition: L1EGammaCrystalsEmulatorProducer_cfi.py:10
HLT_FULL_cff.EBWeights
EBWeights
Definition: HLT_FULL_cff.py:8571
CaloTowerTopology::firstHFRing
int firstHFRing() const
Definition: CaloTowerTopology.h:40
CaloTowersCreationAlgo::theMomEEDepth
double theMomEEDepth
Definition: CaloTowersCreationAlgo.h:366
CaloTowersCreationAlgo::theHcalTopology
const HcalTopology * theHcalTopology
Definition: CaloTowersCreationAlgo.h:321
EcalSeverityLevelAlgo::severityLevel
EcalSeverityLevel::SeverityLevel severityLevel(const DetId &id) const
Evaluate status from id use channelStatus from DB.
Definition: EcalSeverityLevelAlgo.cc:85
CaloTowersCreationAlgo::useRejectedRecoveredEcalHits
unsigned int useRejectedRecoveredEcalHits
Definition: CaloTowersCreationAlgo.h:352
CaloTowersCreationAlgo::MetaTower::numProbHcalCells
int numProbHcalCells
Definition: CaloTowersCreationAlgo.h:274
CaloTowerConstituentsMap::constituentsOf
std::vector< DetId > constituentsOf(const CaloTowerDetId &id) const
Get the constituent detids for this tower id ( not yet implemented )
Definition: CaloTowerConstituentsMap.cc:69
CaloTowersCreationAlgo::MetaTower::numBadHcalCells
int numBadHcalCells
Definition: CaloTowersCreationAlgo.h:273
CaloTowersCreationAlgo::theEBGrid
std::vector< double > theEBGrid
Definition: CaloTowersCreationAlgo.h:300
HcalTopology::mergedDepth29
bool mergedDepth29(HcalDetId id) const
Definition: HcalTopology.h:111
CaloTowersCreationAlgo::setGeometry
void setGeometry(const CaloTowerTopology *cttopo, const CaloTowerConstituentsMap *ctmap, const HcalTopology *htopo, const CaloGeometry *geo)
Definition: CaloTowersCreationAlgo.cc:354
CaloTowersCreationAlgo::theHESweight
double theHESweight
Definition: CaloTowersCreationAlgo.h:309
CaloTowersCreationAlgo::finish
void finish(CaloTowerCollection &destCollection)
Definition: CaloTowersCreationAlgo.cc:404
CaloTowersParam_cfi.mc
mc
Definition: CaloTowersParam_cfi.py:8
CaloTowersCreationAlgo::setHOEScale
void setHOEScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1376
CaloTowersCreationAlgo::MetaTower::hadSumEForTime
float hadSumEForTime
Definition: CaloTowersCreationAlgo.h:270
HLT_FULL_cff.HF1Grid
HF1Grid
Definition: HLT_FULL_cff.py:8570
CaloTowerConstituentsMap.h
CaloTowersCreationAlgo::theGeometry
const CaloGeometry * theGeometry
Definition: CaloTowersCreationAlgo.h:322
CaloTowersCreationAlgo::theHEDweight
double theHEDweight
Definition: CaloTowersCreationAlgo.h:309
mps_merge.weight
weight
Definition: mps_merge.py:88
AlCaHLTBitMon_ParallelJobs.p
p
Definition: AlCaHLTBitMon_ParallelJobs.py:153
CaloTowersCreationAlgo::RecoveredChan
Definition: CaloTowersCreationAlgo.h:390
reco::Candidate::PolarLorentzVector
math::PtEtaPhiMLorentzVector PolarLorentzVector
Lorentz vector.
Definition: Candidate.h:38
DetId::det
constexpr Detector det() const
get the detector field from this detid
Definition: DetId.h:46
CaloTowersCreationAlgo::rescale
void rescale(const CaloTower *ct)
Definition: CaloTowersCreationAlgo.cc:864
CaloTowersCreationAlgo::hcalChanStatusForCaloTower
unsigned int hcalChanStatusForCaloTower(const CaloRecHit *hit)
Definition: CaloTowersCreationAlgo.cc:1753
gather_cfg.cout
cout
Definition: gather_cfg.py:144
HLT_FULL_cff.HESGrid
HESGrid
Definition: HLT_FULL_cff.py:8588
HcalTopology
Definition: HcalTopology.h:26
DetId::Hcal
Definition: DetId.h:28
HLT_FULL_cff.HF1Weights
HF1Weights
Definition: HLT_FULL_cff.py:8584
hgcalTowerProducer_cfi.tower
tower
Definition: hgcalTowerProducer_cfi.py:4
CaloTowersCreationAlgo::theHF2Grid
std::vector< double > theHF2Grid
Definition: CaloTowersCreationAlgo.h:307
CaloTowersCreationAlgo::convert
void convert(const CaloTowerDetId &id, const MetaTower &mt, CaloTowerCollection &collection)
Definition: CaloTowersCreationAlgo.cc:924
HLT_FULL_cff.EESumThreshold
EESumThreshold
Definition: HLT_FULL_cff.py:8548
CaloTower::id
CaloTowerDetId id() const
Definition: CaloTower.h:123
CaloGeometry::getSubdetectorGeometry
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:34
CaloTowersCreationAlgo::theHF1EScale
double theHF1EScale
Definition: CaloTowersCreationAlgo.h:318
CaloTowersCreationAlgo::theTowerMap
MetaTowerMap theTowerMap
Definition: CaloTowersCreationAlgo.h:375
if
if(0==first)
Definition: CAHitNtupletGeneratorKernelsImpl.h:48
CaloTowersCreationAlgo::theHESthreshold
double theHESthreshold
Definition: CaloTowersCreationAlgo.h:296
CaloTowersCreationAlgo::theEEthreshold
double theEEthreshold
Definition: CaloTowersCreationAlgo.h:289
CaloTowerTopology
Definition: CaloTowerTopology.h:13
cms::cuda::assert
assert(be >=bs)
protons_cff.time
time
Definition: protons_cff.py:39
edm::SortedCollection
Definition: SortedCollection.h:49
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:222
digitizers_cfi.ecalTime
ecalTime
Definition: digitizers_cfi.py:77
CaloTower::outerEnergy
double outerEnergy() const
Definition: CaloTower.h:132
CaloTowersCreationAlgo::theHBEScale
double theHBEScale
Definition: CaloTowersCreationAlgo.h:314
CaloTowersCreationAlgo::setHF2EScale
void setHF2EScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1390
CaloTowerTopology::lastHORing
int lastHORing() const
Definition: CaloTowerTopology.h:43
HLT_FULL_cff.EEGrid
EEGrid
Definition: HLT_FULL_cff.py:8559
HcalDetId::depth
constexpr int depth() const
get the tower depth
Definition: HcalDetId.h:164
CaloTowersCreationAlgo::missingHcalRescaleFactorForEcal
float missingHcalRescaleFactorForEcal
Definition: CaloTowersCreationAlgo.h:355
DetId::null
constexpr bool null() const
is this a null id ?
Definition: DetId.h:59
HLT_FULL_cff.HBGrid
HBGrid
Definition: HLT_FULL_cff.py:8549
CaloTowersCreationAlgo.h
CaloTowersCreationAlgo::theEEWeights
std::vector< double > theEEWeights
Definition: CaloTowersCreationAlgo.h:301
EcalSubdetector
EcalSubdetector
Definition: EcalSubdetector.h:10
CaloTowerDetId::SubdetId
static const int SubdetId
Definition: CaloTowerDetId.h:34
HcalTopology::mergedDepthDetId
HcalDetId mergedDepthDetId(const HcalDetId &id) const
Definition: HcalTopology.h:166
GlobalVector
Global3DVector GlobalVector
Definition: GlobalVector.h:10
CaloRecHit
Definition: CaloRecHit.h:23
CaloTowerTopology::firstHERing
int firstHERing() const
Definition: CaloTowerTopology.h:38
HcalBarrel
Definition: HcalAssistant.h:33
CaloTowersCreationAlgo::theHcalAcceptSeverityLevelForRejectedHit
unsigned int theHcalAcceptSeverityLevelForRejectedHit
Definition: CaloTowersCreationAlgo.h:348
CaloTowersCreationAlgo::useRejectedHitsOnly
bool useRejectedHitsOnly
Definition: CaloTowersCreationAlgo.h:347
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
HcalTopology::getMergePositionFlag
bool getMergePositionFlag() const
Definition: HcalTopology.h:167
CaloTowersCreationAlgo::theHOGrid
std::vector< double > theHOGrid
Definition: CaloTowersCreationAlgo.h:305
edm::SortedCollection::push_back
void push_back(T const &t)
Definition: SortedCollection.h:188
DetId::Calo
Definition: DetId.h:29
CaloTowersCreationAlgo::theTowerConstituentsMap
const CaloTowerConstituentsMap * theTowerConstituentsMap
Definition: CaloTowersCreationAlgo.h:323
HcalCondObjectContainer::getValues
const Item * getValues(DetId fId, bool throwOnFail=true) const
Definition: HcalCondObjectContainer.h:159
EcalBarrel
Definition: EcalSubdetector.h:10
rpcPointValidation_cfi.recHit
recHit
Definition: rpcPointValidation_cfi.py:7
HLT_FULL_cff.HF2Grid
HF2Grid
Definition: HLT_FULL_cff.py:8610
CaloTowersCreationAlgo::theHOEScale
double theHOEScale
Definition: CaloTowersCreationAlgo.h:317
mps_monitormerge.items
list items
Definition: mps_monitormerge.py:29
CaloTowersCreationAlgo::hadSegmentShwrPos
GlobalPoint hadSegmentShwrPos(DetId detId, float fracDepth)
Definition: CaloTowersCreationAlgo.cc:1412
ecaldqm::towerId
unsigned towerId(DetId const &, EcalElectronicsMapping const *)
Definition: EcalDQMCommonUtils.cc:77
HLT_FULL_cff.HBWeights
HBWeights
Definition: HLT_FULL_cff.py:8564
CaloTowersCreationAlgo::find
MetaTower & find(const CaloTowerDetId &id)
looks for a given tower in the internal cache. If it can't find it, it makes it.
Definition: CaloTowersCreationAlgo.cc:908
CaloTowersCreationAlgo::setHBEScale
void setHBEScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1355
CaloTowersCreationAlgo::theHF2weight
double theHF2weight
Definition: CaloTowersCreationAlgo.h:309
CaloTowersCreationAlgo::theEcutTower
double theEcutTower
Definition: CaloTowersCreationAlgo.h:310
CaloTowersCreationAlgo::theHEDWeights
std::vector< double > theHEDWeights
Definition: CaloTowersCreationAlgo.h:304
CaloTowersCreationAlgo::theHF2threshold
double theHF2threshold
Definition: CaloTowersCreationAlgo.h:299
HcalTopology::firstHEDoublePhiRing
int firstHEDoublePhiRing() const
Definition: HcalTopology.h:101
CaloTowersCreationAlgo::theHOIsUsed
bool theHOIsUsed
only affects energy and ET calculation. HO is still recorded in the tower
Definition: CaloTowersCreationAlgo.h:358
DetId
Definition: DetId.h:17
CaloTowersCreationAlgo::compactTime
int compactTime(float time)
Definition: CaloTowersCreationAlgo.cc:1634
UNLIKELY
#define UNLIKELY(x)
Definition: Likely.h:21
CaloTowersCreationAlgo::begin
void begin()
Definition: CaloTowersCreationAlgo.cc:368
CaloTowersCreationAlgo::theHF1weight
double theHF1weight
Definition: CaloTowersCreationAlgo.h:309
CaloGeometry
Definition: CaloGeometry.h:21
CaloTowersCreationAlgo::theHBGrid
std::vector< double > theHBGrid
Definition: CaloTowersCreationAlgo.h:302
CaloTowersCreationAlgo::MetaTower
Definition: CaloTowersCreationAlgo.h:262
photonIsolationHIProducer_cfi.hf
hf
Definition: photonIsolationHIProducer_cfi.py:9
CaloTowersCreationAlgo::theTowerMapSize
unsigned int theTowerMapSize
Definition: CaloTowersCreationAlgo.h:376
CaloTowersCreationAlgo::hadShwrPos
GlobalPoint hadShwrPos(const std::vector< std::pair< DetId, float >> &metaContains, float fracDepth, double hadE)
CaloTowersCreationAlgo::makeEcalBadChs
void makeEcalBadChs()
Definition: CaloTowersCreationAlgo.cc:1710
edm::SortedCollection::begin
const_iterator begin() const
Definition: SortedCollection.h:262
CaloTowersCreationAlgo::assignHitEcal
void assignHitEcal(const EcalRecHit *recHit)
adds a single hit to the tower
Definition: CaloTowersCreationAlgo.cc:772
PVValHelper::eta
Definition: PVValidationHelpers.h:70
CaloTower::constituentsSize
size_t constituentsSize() const
Definition: CaloTower.h:125
reco::LeafCandidate::PolarLorentzVector
math::PtEtaPhiMLorentzVector PolarLorentzVector
Lorentz vector.
Definition: LeafCandidate.h:25
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
CaloTowerTopology::lastHFRing
int lastHFRing() const
Definition: CaloTowerTopology.h:41
CaloTower::hadEnergy
double hadEnergy() const
Definition: CaloTower.h:131
CaloTowersCreationAlgo::theHOthresholdPlus1
double theHOthresholdPlus1
Definition: CaloTowersCreationAlgo.h:298
CaloTowersCreationAlgo::theHcalChStatus
const HcalChannelQuality * theHcalChStatus
Definition: CaloTowersCreationAlgo.h:328
CaloTowersCreationAlgo::theHEDthreshold
double theHEDthreshold
Definition: CaloTowersCreationAlgo.h:297
CaloTowersCreationAlgo::theHESthreshold1
double theHESthreshold1
Definition: CaloTowersCreationAlgo.h:296
CaloTowerConstituentsMap
Definition: CaloTowerConstituentsMap.h:18
CaloTower::setHcalSubdet
void setHcalSubdet(int lastHB, int lastHE, int lastHF, int lastHO)
Definition: CaloTower.h:85
CaloTower::setCaloTowerStatus
void setCaloTowerStatus(unsigned int numBadHcalChan, unsigned int numBadEcalChan, unsigned int numRecHcalChan, unsigned int numRecEcalChan, unsigned int numProbHcalChan, unsigned int numProbEcalChan)
Definition: CaloTower.cc:221
CaloTowersCreationAlgo::CaloTowersCreationAlgo
CaloTowersCreationAlgo()
Definition: CaloTowersCreationAlgo.cc:14
HcalOuter
Definition: HcalAssistant.h:35
HLT_FULL_cff.HcalThreshold
HcalThreshold
Definition: HLT_FULL_cff.py:8556
HCALHighEnergyHPDFilter_cfi.energy
energy
Definition: HCALHighEnergyHPDFilter_cfi.py:5
CaloTowersCreationAlgo::theHcalThreshold
double theHcalThreshold
Definition: CaloTowersCreationAlgo.h:293
CaloTowersCreationAlgo::theEBthreshold
double theEBthreshold
Definition: CaloTowersCreationAlgo.h:289
CaloTowersCreationAlgo::emShwrPos
GlobalPoint emShwrPos(const std::vector< std::pair< DetId, float >> &metaContains, float fracDepth, double totEmE)
Definition: CaloTowersCreationAlgo.cc:1567
GlobalPoint
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
dqmdumpme.k
k
Definition: dqmdumpme.py:60
CaloTower::setEcalTime
void setEcalTime(int t)
Definition: CaloTower.h:83
Point3DBase< float, GlobalTag >
HLT_FULL_cff.EcutTower
EcutTower
Definition: HLT_FULL_cff.py:8589
CaloTowerTopology.h
HLT_FULL_cff.HOWeights
HOWeights
Definition: HLT_FULL_cff.py:8558
LEDCalibrationChannels.depth
depth
Definition: LEDCalibrationChannels.py:65
CaloTowersCreationAlgo::theMomEBDepth
double theMomEBDepth
Definition: CaloTowersCreationAlgo.h:365
HcalSeverityLevelComputer::dropChannel
bool dropChannel(const uint32_t &mystatus) const
Definition: HcalSeverityLevelComputer.cc:395
DetId::subdetId
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector's numbering enum)
Definition: DetId.h:48
mitigatedMETSequence_cff.U
U
Definition: mitigatedMETSequence_cff.py:36
CaloTowersCreationAlgo::GoodChan
Definition: CaloTowersCreationAlgo.h:390
CaloTowersCreationAlgo::theHESWeights
std::vector< double > theHESWeights
Definition: CaloTowersCreationAlgo.h:303
EcalEndcap
Definition: EcalSubdetector.h:10
HcalSeverityLevelComputer::getSeverityLevel
int getSeverityLevel(const DetId &myid, const uint32_t &myflag, const uint32_t &mystatus) const
Definition: HcalSeverityLevelComputer.cc:304
CaloTowersCreationAlgo::setHEDEScale
void setHEDEScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1369
HcalChannelStatus::getValue
uint32_t getValue() const
Definition: HcalChannelStatus.h:60
funct::true
true
Definition: Factorize.h:173
CaloTowersCreationAlgo::theHOweight
double theHOweight
Definition: CaloTowersCreationAlgo.h:309
HLT_FULL_cff.HEDGrid
HEDGrid
Definition: HLT_FULL_cff.py:8602
CaloTowersCreationAlgo::theHEDGrid
std::vector< double > theHEDGrid
Definition: CaloTowersCreationAlgo.h:304
CaloSubdetectorGeometry.h
CaloTowersCreationAlgo::theMomHBDepth
double theMomHBDepth
Definition: CaloTowersCreationAlgo.h:363
CaloTowersCreationAlgo::setEEEScale
void setEEEScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1348
HcalDetId::ieta
constexpr int ieta() const
get the cell ieta
Definition: HcalDetId.h:155
CaloTowersCreationAlgo::theEEGrid
std::vector< double > theEEGrid
Definition: CaloTowersCreationAlgo.h:301
CaloTowersCreationAlgo::theHcalPhase
int theHcalPhase
Definition: CaloTowersCreationAlgo.h:397
CaloTowersCreationAlgo::theHOthreshold0
double theHOthreshold0
Definition: CaloTowersCreationAlgo.h:298
CaloTowerTopology::detIdFromDenseIndex
CaloTowerDetId detIdFromDenseIndex(uint32_t din) const
Definition: CaloTowerTopology.cc:247
CaloTowersCreationAlgo::MetaTower::hadSumTimeTimesE
float hadSumTimeTimesE
Definition: CaloTowersCreationAlgo.h:269
HLT_FULL_cff.EBSumThreshold
EBSumThreshold
Definition: HLT_FULL_cff.py:8541
CaloGeometry::getGeometry
std::shared_ptr< const CaloCellGeometry > getGeometry(const DetId &id) const
Get the cell geometry of a given detector id.
Definition: CaloGeometry.cc:60
CaloTowersCreationAlgo::theEBWeights
std::vector< double > theEBWeights
Definition: CaloTowersCreationAlgo.h:300
edm::SortedCollection::end
const_iterator end() const
Definition: SortedCollection.h:267
HcalCondObjectContainer::getAllChannels
std::vector< DetId > getAllChannels() const
Definition: HcalCondObjectContainer.h:312
CaloTowersCreationAlgo::theHF1Weights
std::vector< double > theHF1Weights
Definition: CaloTowersCreationAlgo.h:306
CaloTowersCreationAlgo::theHcalAcceptSeverityLevel
unsigned int theHcalAcceptSeverityLevel
Definition: CaloTowersCreationAlgo.h:339
match
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:10
CaloTowersCreationAlgo::process
void process(const HBHERecHitCollection &hbhe)
Definition: CaloTowersCreationAlgo.cc:374
HcalTopology::idBack
HcalDetId idBack(const HcalDetId &id) const
Definition: HcalTopology.h:171
CaloTowersCreationAlgo::makeHcalDropChMap
void makeHcalDropChMap()
Definition: CaloTowersCreationAlgo.cc:1649
PV3DBase::eta
T eta() const
Definition: PV3DBase.h:73
HcalDetId::subdet
constexpr HcalSubdetector subdet() const
get the subdetector
Definition: HcalDetId.h:138
HcalDetId
Definition: HcalDetId.h:12
createfilelist.int
int
Definition: createfilelist.py:10
CaloTowerTopology::lastHERing
int lastHERing() const
Definition: CaloTowerTopology.h:39
CaloTower
Definition: CaloTower.h:26
CaloTowersCreationAlgo::MetaTower::metaConstituents
std::vector< std::pair< DetId, float > > metaConstituents
Definition: CaloTowersCreationAlgo.h:266
universalConfigTemplate.collection
collection
Definition: universalConfigTemplate.py:81
CaloTowersCreationAlgo::theHEDEScale
double theHEDEScale
Definition: CaloTowersCreationAlgo.h:316
CaloTowersCreationAlgo::theHBthreshold
double theHBthreshold
Definition: CaloTowersCreationAlgo.h:295
trackerHitRTTI::vector
Definition: trackerHitRTTI.h:21
HLT_FULL_cff.EBGrid
EBGrid
Definition: HLT_FULL_cff.py:8612
photonIsolationHIProducer_cfi.ho
ho
Definition: photonIsolationHIProducer_cfi.py:10
HLT_FULL_cff.HEDWeights
HEDWeights
Definition: HLT_FULL_cff.py:8568
CaloTowersCreationAlgo::useRejectedRecoveredHcalHits
unsigned int useRejectedRecoveredHcalHits
Definition: CaloTowersCreationAlgo.h:351
DetId::Detector
Detector
Definition: DetId.h:24
CaloSubdetectorGeometry::getGeometry
virtual std::shared_ptr< const CaloCellGeometry > getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
Definition: CaloSubdetectorGeometry.cc:36
HcalTopology::lastHERing
int lastHERing() const
Definition: HcalTopology.h:94
CaloTowerConstituentsMap::towerOf
CaloTowerDetId towerOf(const DetId &id) const
Get the tower id for this det id (or null if not known)
Definition: CaloTowerConstituentsMap.cc:26
HLT_FULL_cff.EEWeights
EEWeights
Definition: HLT_FULL_cff.py:8561
edm::contains
bool contains(EventRange const &lh, EventID const &rh)
Definition: EventRange.cc:37
CaloTowersCreationAlgo::theMomHEDepth
double theMomHEDepth
Definition: CaloTowersCreationAlgo.h:364
edm::LogError
Log< level::Error, false > LogError
Definition: MessageLogger.h:123
DetId::Ecal
Definition: DetId.h:27
CaloTowersCreationAlgo::hadShwPosFromCells
GlobalPoint hadShwPosFromCells(DetId frontCell, DetId backCell, float fracDepth)
Definition: CaloTowersCreationAlgo.cc:1540
CaloTowersCreationAlgo::MetaTower::E_outer
float E_outer
Definition: CaloTowersCreationAlgo.h:268
CaloTowersCreationAlgo::theEESumThreshold
double theEESumThreshold
Definition: CaloTowersCreationAlgo.h:310
PV3DBase::basicVector
const BasicVectorType & basicVector() const
Definition: PV3DBase.h:53
CaloTowersCreationAlgo::rescaleTowers
void rescaleTowers(const CaloTowerCollection &ctInput, CaloTowerCollection &ctResult)
Definition: CaloTowersCreationAlgo.cc:426
photonIsolationHIProducer_cfi.hbhe
hbhe
Definition: photonIsolationHIProducer_cfi.py:8
HcalSubdetector
HcalSubdetector
Definition: HcalAssistant.h:31
CaloTowerDetId::iphi
int iphi() const
get the tower iphi
Definition: CaloTowerDetId.cc:30
HcalForward
Definition: HcalAssistant.h:36
CaloTowerTopology::lastHBRing
int lastHBRing() const
Definition: CaloTowerTopology.h:37
CaloTowersCreationAlgo::theHBthreshold2
double theHBthreshold2
Definition: CaloTowersCreationAlgo.h:295
CaloTowersCreationAlgo::theTowerGeometry
const CaloSubdetectorGeometry * theTowerGeometry
Definition: CaloTowersCreationAlgo.h:324
CaloTower::constituent
DetId constituent(size_t i) const
Definition: CaloTower.h:126
HcalTopology.h
CaloCellGeometry.h
CaloTowerTopology::convertCTtoHcal
int convertCTtoHcal(int ct_ieta) const
Definition: CaloTowerTopology.cc:67
CaloTowersCreationAlgo::theHF2Weights
std::vector< double > theHF2Weights
Definition: CaloTowersCreationAlgo.h:307
heppy_batch.val
val
Definition: heppy_batch.py:351
CaloTowersCreationAlgo::theEEweight
double theEEweight
Definition: CaloTowersCreationAlgo.h:308
HLT_FULL_cff.HESWeights
HESWeights
Definition: HLT_FULL_cff.py:8592
eostools.move
def move(src, dest)
Definition: eostools.py:511
CaloTowersCreationAlgo::theHcalSevLvlComputer
const HcalSeverityLevelComputer * theHcalSevLvlComputer
Definition: CaloTowersCreationAlgo.h:331
std
Definition: JetResolutionObject.h:76
interestingDetIdCollectionProducer_cfi.severityLevel
severityLevel
Definition: interestingDetIdCollectionProducer_cfi.py:10
CaloTowersCreationAlgo::theHOthresholdMinus2
double theHOthresholdMinus2
Definition: CaloTowersCreationAlgo.h:299
CaloTowerDetId::ieta
int ieta() const
get the tower ieta
Definition: CaloTowerDetId.h:30
HcalEndcap
Definition: HcalAssistant.h:34
CaloTowersCreationAlgo::theHOthresholdMinus1
double theHOthresholdMinus1
Definition: CaloTowersCreationAlgo.h:298
CaloTowersCreationAlgo::MetaTower::E
float E
Definition: CaloTowersCreationAlgo.h:268
CaloTowersCreationAlgo::theHBweight
double theHBweight
Definition: CaloTowersCreationAlgo.h:309
CaloTowersCreationAlgo::theEBEScale
double theEBEScale
Definition: CaloTowersCreationAlgo.h:312
CaloTowersCreationAlgo::getThresholdAndWeight
void getThresholdAndWeight(const DetId &detId, double &threshold, double &weight) const
helper method to look up the appropriate threshold & weight
Definition: CaloTowersCreationAlgo.cc:1245
EcalSeverityLevel::kRecovered
Definition: EcalSeverityLevel.h:21
CaloTowersCreationAlgo::theEEEScale
double theEEEScale
Definition: CaloTowersCreationAlgo.h:313
CaloTowersCreationAlgo::theHOthresholdPlus2
double theHOthresholdPlus2
Definition: CaloTowersCreationAlgo.h:299
CaloTowersCreationAlgo::theEcalSeveritiesToBeExcluded
std::vector< int > theEcalSeveritiesToBeExcluded
Definition: CaloTowersCreationAlgo.h:340
CaloGeometry.h
CaloTowersCreationAlgo::hcalDropChMap
HcalDropChMap hcalDropChMap
Definition: CaloTowersCreationAlgo.h:381
CaloTowersCreationAlgo::ecalBadChs
std::vector< unsigned short > ecalBadChs
Definition: CaloTowersCreationAlgo.h:385
CaloTowersCreationAlgo::setEBEScale
void setEBEScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1341
CaloTowersCreationAlgo::ProblematicChan
Definition: CaloTowersCreationAlgo.h:390
CaloTowerTopology::sizeForDenseIndexing
uint32_t sizeForDenseIndexing() const
Definition: CaloTowerTopology.h:56
TtSemiLepEvtBuilder_cfi.mt
mt
Definition: TtSemiLepEvtBuilder_cfi.py:47
CaloTowersCreationAlgo::theHESGrid
std::vector< double > theHESGrid
Definition: CaloTowersCreationAlgo.h:303
CaloTowersCreationAlgo::theHESEScale
double theHESEScale
Definition: CaloTowersCreationAlgo.h:315
CaloTowersCreationAlgo::theEcalSevLvlAlgo
const EcalSeverityLevelAlgo * theEcalSevLvlAlgo
Definition: CaloTowersCreationAlgo.h:334
dqm-mbProfile.log
log
Definition: dqm-mbProfile.py:17
CaloTowersCreationAlgo::MetaTower::numRecHcalCells
int numRecHcalCells
Definition: CaloTowersCreationAlgo.h:273
CaloTowersCreationAlgo::theHF1Grid
std::vector< double > theHF1Grid
Definition: CaloTowersCreationAlgo.h:306
L1TowerCalibrationProducer_cfi.iEta
iEta
Definition: L1TowerCalibrationProducer_cfi.py:60
mps_fire.result
result
Definition: mps_fire.py:311
CaloTowersCreationAlgo::theHBWeights
std::vector< double > theHBWeights
Definition: CaloTowersCreationAlgo.h:302
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
CaloTowersCreationAlgo::setHESEScale
void setHESEScale(double scale)
Definition: CaloTowersCreationAlgo.cc:1362
EcalSeverityLevel::kBad
Definition: EcalSeverityLevel.h:24
CaloTowerDetId::zside
int zside() const
get the z-side of the tower (1/-1)
Definition: CaloTowerDetId.h:26
fftjetcommon_cfi.Interpolator
Interpolator
Definition: fftjetcommon_cfi.py:188
HcalDetId::ietaAbs
constexpr int ietaAbs() const
get the absolute value of the cell ieta
Definition: HcalDetId.h:148
CaloTower::hcalTime
float hcalTime() const
Definition: CaloTower.h:197
HLT_FULL_cff.useHO
useHO
Definition: HLT_FULL_cff.py:11497
CaloTowersCreationAlgo::theUseEtEETresholdFlag
bool theUseEtEETresholdFlag
Definition: CaloTowersCreationAlgo.h:290
CaloTowersCreationAlgo::theTowerTopology
const CaloTowerTopology * theTowerTopology
Definition: CaloTowersCreationAlgo.h:320
CaloTowersCreationAlgo::theUseSymEBTresholdFlag
bool theUseSymEBTresholdFlag
Definition: CaloTowersCreationAlgo.h:291
CaloTowersCreationAlgo::MetaTower::E_had
float E_had
Definition: CaloTowersCreationAlgo.h:268
CaloTowersCreationAlgo::theHF1threshold
double theHF1threshold
Definition: CaloTowersCreationAlgo.h:299
CaloTowersCreationAlgo::theHOWeights
std::vector< double > theHOWeights
Definition: CaloTowersCreationAlgo.h:305
point
*vegas h *****************************************************used in the default bin number in original ***version of VEGAS is ***a higher bin number might help to derive a more precise ***grade subtle point
Definition: invegas.h:5
remoteMonitoring_LED_IterMethod_cfg.threshold
threshold
Definition: remoteMonitoring_LED_IterMethod_cfg.py:430
MatrixUtil.merge
def merge(dictlist, TELL=False)
Definition: MatrixUtil.py:201
CaloTowersCreationAlgo::theUseSymEETresholdFlag
bool theUseSymEETresholdFlag
Definition: CaloTowersCreationAlgo.h:291
CaloTowersCreationAlgo::ecalChanStatusForCaloTower
std::tuple< unsigned int, bool > ecalChanStatusForCaloTower(const EcalRecHit *hit)
Definition: CaloTowersCreationAlgo.cc:1804
CaloTowersCreationAlgo::theHEDthreshold1
double theHEDthreshold1
Definition: CaloTowersCreationAlgo.h:297
CaloTowersCreationAlgo::theHBthreshold1
double theHBthreshold1
Definition: CaloTowersCreationAlgo.h:295
CaloTowersCreationAlgo::theRecoveredHcalHitsAreUsed
bool theRecoveredHcalHitsAreUsed
Definition: CaloTowersCreationAlgo.h:342
Basic3DVector::unit
Basic3DVector unit() const
Definition: extBasic3DVector.h:151
CaloTowersCreationAlgo::theHF2EScale
double theHF2EScale
Definition: CaloTowersCreationAlgo.h:319
CaloTower::addConstituents
void addConstituents(const std::vector< DetId > &ids)
Definition: CaloTower.cc:199
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
Basic3DVector< float >
CaloTowersCreationAlgo::theRecoveredEcalHitsAreUsed
bool theRecoveredEcalHitsAreUsed
Definition: CaloTowersCreationAlgo.h:343
weight
Definition: weight.py:1
merge
Definition: merge.py:1
CaloTowersCreationAlgo::theEBweight
double theEBweight
Definition: CaloTowersCreationAlgo.h:308
hit
Definition: SiStripHitEffFromCalibTree.cc:88
HLT_FULL_cff.HOGrid
HOGrid
Definition: HLT_FULL_cff.py:8611
CaloTowerTopology::denseIndex
uint32_t denseIndex(const DetId &id) const
Definition: CaloTowerTopology.cc:235
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
CaloTowerDetId
Definition: CaloTowerDetId.h:12
CaloTowersCreationAlgo::BadChan
Definition: CaloTowersCreationAlgo.h:390
HLT_FULL_cff.missingHcalRescaleFactorForEcal
missingHcalRescaleFactorForEcal
Definition: HLT_FULL_cff.py:8566
Generated for CMSSW Reference Manual by doxygen 1.8.16