Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-09-18 05:06:57 

0001 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0002 #include <cmath>
0003 #include <iostream>
0004 #include <cassert>
0005 #include <algorithm>
0006 #include "FWCore/Utilities/interface/Exception.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "DataFormats/Math/interface/GeantUnits.h"
0009 
0010 using namespace geant_units;
0011 using namespace geant_units::operators;
0012 
0013 static const int IPHI_MAX = 72;
0014 
0015 //#define EDM_ML_DEBUG
0016 
0017 HcalTopology::HcalTopology(const HcalDDDRecConstants* hcons, const bool mergePosition)
0018     : hcons_(hcons),
0019       mergePosition_(mergePosition),
0020       excludeHB_(false),
0021       excludeHE_(false),
0022       excludeHO_(false),
0023       excludeHF_(false),
0024       firstHBRing_(1),
0025       firstHERing_(999),
0026       lastHERing_(0),
0027       firstHFRing_(29),
0028       lastHFRing_(41),
0029       firstHORing_(1),
0030       lastHORing_(15),
0031       firstHEDoublePhiRing_(999),
0032       firstHEQuadPhiRing_(999),
0033       firstHFQuadPhiRing_(40),
0034       firstHETripleDepthRing_(999),
0035       singlePhiBins_(IPHI_MAX),
0036       doublePhiBins_(36),
0037       maxPhiHE_(IPHI_MAX) {
0038   mode_ = (HcalTopologyMode::Mode)(hcons_->getTopoMode());
0039   triggerMode_ = (HcalTopologyMode::TriggerMode)(hcons_->getTriggerMode());
0040   maxDepthHB_ = hcons_->getMaxDepth(0);
0041   maxDepthHE_ = hcons_->getMaxDepth(1);
0042   maxDepthHF_ = hcons_->getMaxDepth(2);
0043   etaBinsHB_ = hcons_->getEtaBins(0);
0044   etaBinsHE_ = hcons_->getEtaBins(1);
0045   nEtaHB_ = (hcons_->getEtaRange(0)).second - (hcons_->getEtaRange(0)).first + 1;
0046   lastHBRing_ = firstHBRing_ + nEtaHB_ - 1;
0047   if (hcons_->getNPhi(1) > maxPhiHE_)
0048     maxPhiHE_ = hcons_->getNPhi(1);
0049   for (auto& i : etaBinsHE_) {
0050     if (firstHERing_ > i.ieta)
0051       firstHERing_ = i.ieta;
0052     if (lastHERing_ < i.ieta)
0053       lastHERing_ = i.ieta;
0054     int unit = static_cast<int>((i.dphi / 5.0_deg) + 0.01);
0055     if (unit == 2 && firstHEDoublePhiRing_ > i.ieta)
0056       firstHEDoublePhiRing_ = i.ieta;
0057     if (unit == 4 && firstHEQuadPhiRing_ > i.ieta)
0058       firstHEQuadPhiRing_ = i.ieta;
0059     if (i.layer.size() > 2 && firstHETripleDepthRing_ > i.ieta)
0060       firstHETripleDepthRing_ = i.ieta;
0061   }
0062   if (firstHERing_ > lastHERing_) {
0063     firstHERing_ = lastHERing_ = firstHEDoublePhiRing_ = firstHEQuadPhiRing_ = firstHETripleDepthRing_ = nEtaHE_ = 0;
0064   } else {
0065     nEtaHE_ = (lastHERing_ - firstHERing_ + 1);
0066   }
0067   if (mode_ == HcalTopologyMode::LHC) {
0068     topoVersion_ = 0;         //DL
0069     HBSize_ = kHBSizePreLS1;  // qie-per-fiber * fiber/rm * rm/rbx * rbx/barrel * barrel/hcal
0070     HESize_ = kHESizePreLS1;  // qie-per-fiber * fiber/rm * rm/rbx * rbx/endcap * endcap/hcal
0071     HOSize_ = kHOSizePreLS1;  // ieta * iphi * 2
0072     HFSize_ = kHFSizePreLS1;  // ieta * iphi * depth * 2
0073     CALIBSize_ = kCALIBSizePreLS1;
0074     numberOfShapes_ = 87;
0075   } else if (mode_ == HcalTopologyMode::SLHC) {  // need to know more eventually
0076     topoVersion_ = 10;
0077     HBSize_ = nEtaHB_ * IPHI_MAX * maxDepthHB_ * 2;
0078     HESize_ = nEtaHE_ * maxPhiHE_ * maxDepthHE_ * 2;
0079     HOSize_ = (lastHORing_ - firstHORing_ + 1) * IPHI_MAX * 2;                // ieta * iphi * 2
0080     HFSize_ = (lastHFRing_ - firstHFRing_ + 1) * IPHI_MAX * maxDepthHF_ * 2;  // ieta * iphi * depth * 2
0081     CALIBSize_ = kOffCalibHFX_;
0082     numberOfShapes_ = (maxPhiHE_ > 72) ? 1200 : 500;
0083   }
0084   maxEta_ = (lastHERing_ > lastHFRing_) ? lastHERing_ : lastHFRing_;
0085   if (triggerMode_ == HcalTopologyMode::TriggerMode_2009) {
0086     HTSize_ = kHTSizePreLS1;
0087   } else {
0088     HTSize_ = kHTSizePhase1;
0089   }
0090 
0091 #ifdef EDM_ML_DEBUG
0092   edm::LogVerbatim("HCalGeom") << "Topo sizes " << HBSize_ << ":" << HESize_ << ":" << HOSize_ << ":" << HFSize_ << ":"
0093                                << HTSize_ << ":" << CALIBSize_ << " for mode " << mode_ << ":" << triggerMode_;
0094 #endif
0095 
0096   //The transition between HE/HF in eta
0097   etaTableHF = hcons_->getEtaTableHF();
0098   etaTable = hcons_->getEtaTable();
0099   dPhiTableHF = hcons_->getPhiTableHF();
0100   dPhiTable = hcons_->getPhiTable();
0101   phioff = hcons_->getPhiOffs();
0102   std::pair<int, int> ietaHF = hcons_->getEtaRange(2);
0103   etaHE2HF_ = firstHFRing_;
0104   etaHF2HE_ = lastHERing_;
0105   if (etaBinsHE_.size() > 1) {
0106     double eta = etaBinsHE_[etaBinsHE_.size() - 1].etaMax;
0107     for (unsigned int i = 1; i < etaTableHF.size(); ++i) {
0108       if (eta < etaTableHF[i]) {
0109         etaHE2HF_ = ietaHF.first + i - 1;
0110         break;
0111       }
0112     }
0113     eta = etaTableHF[0];
0114     for (auto& i : etaBinsHE_) {
0115       if (eta < i.etaMax) {
0116         etaHF2HE_ = i.ieta;
0117         break;
0118       }
0119     }
0120   }
0121   const double fiveDegInRad = 5.0_deg;
0122   for (double k : dPhiTable) {
0123     int units = (int)(k / fiveDegInRad + 0.5);
0124     unitPhi.emplace_back(units);
0125   }
0126   for (double k : dPhiTableHF) {
0127     int units = (int)(k / fiveDegInRad + 0.5);
0128     unitPhiHF.emplace_back(units);
0129   }
0130   int nEta = hcons_->getNEta();
0131   for (int ring = 1; ring <= nEta; ++ring) {
0132     std::vector<int> segmentation = hcons_->getDepth(ring - 1, false);
0133     setDepthSegmentation(ring, segmentation, false);
0134 #ifdef EDM_ML_DEBUG
0135     edm::LogVerbatim("HCalGeom") << "Set segmentation for ring " << ring << " with " << segmentation.size()
0136                                  << " elements:";
0137     for (unsigned int k = 0; k < segmentation.size(); ++k)
0138       edm::LogVerbatim("HCalGeom") << "[" << k << "] " << segmentation[k];
0139 #endif
0140     segmentation = hcons_->getDepth(ring - 1, true);
0141     setDepthSegmentation(ring, segmentation, true);
0142 #ifdef EDM_ML_DEBUG
0143     edm::LogVerbatim("HCalGeom") << "Set Plan-1 segmentation for ring " << ring << " with " << segmentation.size()
0144                                  << " elements:";
0145     for (unsigned int k = 0; k < segmentation.size(); ++k)
0146       edm::LogVerbatim("HCalGeom") << "[" << k << "] " << segmentation[k];
0147 #endif
0148   }
0149 
0150 #ifdef EDM_ML_DEBUG
0151   edm::LogVerbatim("HCalGeom") << "Constants in HcalTopology " << firstHBRing_ << ":" << lastHBRing_ << " "
0152                                << firstHERing_ << ":" << lastHERing_ << ":" << firstHEDoublePhiRing_ << ":"
0153                                << firstHEQuadPhiRing_ << ":" << firstHETripleDepthRing_ << " " << firstHFRing_ << ":"
0154                                << lastHFRing_ << ":" << firstHFQuadPhiRing_ << " " << firstHORing_ << ":" << lastHORing_
0155                                << " " << maxDepthHB_ << ":" << maxDepthHE_ << " " << nEtaHB_ << ":" << nEtaHE_ << " "
0156                                << etaHE2HF_ << ":" << etaHF2HE_ << " " << maxPhiHE_;
0157 #endif
0158 }
0159 
0160 HcalTopology::HcalTopology(HcalTopologyMode::Mode mode,
0161                            int maxDepthHB,
0162                            int maxDepthHE,
0163                            HcalTopologyMode::TriggerMode tmode)
0164     : hcons_(nullptr),
0165       mergePosition_(false),
0166       excludeHB_(false),
0167       excludeHE_(false),
0168       excludeHO_(false),
0169       excludeHF_(false),
0170       mode_(mode),
0171       triggerMode_(tmode),
0172       firstHBRing_(1),
0173       lastHBRing_(16),
0174       firstHERing_(16),
0175       lastHERing_(29),
0176       firstHFRing_(29),
0177       lastHFRing_(41),
0178       firstHORing_(1),
0179       lastHORing_(15),
0180       firstHEDoublePhiRing_((mode == HcalTopologyMode::H2 || mode == HcalTopologyMode::H2HE) ? (22) : (21)),
0181       firstHEQuadPhiRing_(999),
0182       firstHFQuadPhiRing_(40),
0183       firstHETripleDepthRing_((mode == HcalTopologyMode::H2 || mode == HcalTopologyMode::H2HE) ? (24) : (27)),
0184       singlePhiBins_(IPHI_MAX),
0185       doublePhiBins_(36),
0186       maxDepthHB_(maxDepthHB),
0187       maxDepthHE_(maxDepthHE),
0188       maxDepthHF_(2),
0189       etaHE2HF_(30),
0190       etaHF2HE_(29),
0191       maxPhiHE_(IPHI_MAX),
0192       HBSize_(kHBSizePreLS1),
0193       HESize_(kHESizePreLS1),
0194       HOSize_(kHOSizePreLS1),
0195       HFSize_(kHFSizePreLS1),
0196       HTSize_(kHTSizePreLS1),
0197       CALIBSize_(kCALIBSizePreLS1),
0198       numberOfShapes_(
0199           ((mode == HcalTopologyMode::SLHC) || (mode_ == HcalTopologyMode::Run3) || (mode_ == HcalTopologyMode::Run4))
0200               ? 500
0201               : 87) {
0202   if (mode_ == HcalTopologyMode::LHC) {
0203     topoVersion_ = 0;         //DL
0204     HBSize_ = kHBSizePreLS1;  // qie-per-fiber * fiber/rm * rm/rbx * rbx/barrel * barrel/hcal
0205     HESize_ = kHESizePreLS1;  // qie-per-fiber * fiber/rm * rm/rbx * rbx/endcap * endcap/hcal
0206     HOSize_ = kHOSizePreLS1;  // ieta * iphi * 2
0207     HFSize_ = kHFSizePreLS1;  // phi * eta * depth * pm
0208   } else if ((mode_ == HcalTopologyMode::SLHC) || (mode_ == HcalTopologyMode::Run3) ||
0209              (mode_ == HcalTopologyMode::Run4)) {  // need to know more eventually
0210     HBSize_ = maxDepthHB * 16 * IPHI_MAX * 2;
0211     HESize_ = maxDepthHE * (29 - 16 + 1) * maxPhiHE_ * 2;
0212     HOSize_ = 15 * IPHI_MAX * 2;                // ieta * iphi * 2
0213     HFSize_ = IPHI_MAX * 13 * maxDepthHF_ * 2;  // phi * eta * depth * pm
0214     CALIBSize_ = kOffCalibHFX_;
0215     topoVersion_ = 10;
0216   }
0217   nEtaHB_ = (lastHBRing_ - firstHBRing_ + 1);
0218   nEtaHE_ = (lastHERing_ - firstHERing_ + 1);
0219   if (triggerMode_ == HcalTopologyMode::TriggerMode_2009) {
0220     HTSize_ = kHTSizePreLS1;
0221   } else {
0222     HTSize_ = kHTSizePhase1;
0223   }
0224 
0225   edm::LogWarning("HCalGeom") << "This is an incomplete constructor of HcalTopology - be warned that many "
0226                               << "functionalities will not be there - revert from this - get from EventSetup";
0227 }
0228 
0229 bool HcalTopology::valid(const DetId& id) const {
0230   assert(id.det() == DetId::Hcal);
0231   return validHcal(id);
0232 }
0233 
0234 bool HcalTopology::validHcal(const HcalDetId& id) const {
0235   // check the raw rules
0236   return validRaw(id) && !isExcluded(id);
0237 }
0238 
0239 bool HcalTopology::validDetId(HcalSubdetector subdet, int ieta, int iphi, int depth) const {
0240   return validHcal(HcalDetId(subdet, ieta, iphi, depth));
0241 }
0242 
0243 bool HcalTopology::validHT(const HcalTrigTowerDetId& id) const {
0244   if (id.iphi() < 1 || id.iphi() > IPHI_MAX || id.ieta() == 0)
0245     return false;
0246   if (id.depth() != 0)
0247     return false;
0248   if (maxDepthHE_ == 0) {
0249     if (id.ietaAbs() > lastHBRing_ && id.ietaAbs() < firstHFRing_)
0250       return false;
0251   }
0252   // Version 2 TPs should be for HBHE when using 1TS filter scheme
0253   if (id.version() == 0 or id.version() == 2) {
0254     if (id.ietaAbs() > 28) {
0255       if (triggerMode_ >= HcalTopologyMode::TriggerMode_2017)
0256         return false;
0257       if (triggerMode_ == HcalTopologyMode::TriggerMode_2018legacy)
0258         return false;
0259       if ((id.iphi() % 4) != 1)
0260         return false;
0261       if (id.ietaAbs() > 32)
0262         return false;
0263     }
0264   } else if (id.version() == 1) {
0265     if (triggerMode_ == HcalTopologyMode::TriggerMode_2009)
0266       return false;
0267     if (id.ietaAbs() < 30 || id.ietaAbs() > 41)
0268       return false;
0269     if (id.ietaAbs() > 29 && ((id.iphi() % 2) == 0))
0270       return false;
0271     if (id.ietaAbs() > 39 && ((id.iphi() % 4) != 3))
0272       return false;
0273   } else if (id.version() > 2) {
0274     // only versions 0, 1, and 2 are supported
0275     return false;
0276   }
0277 
0278   return true;
0279 }
0280 
0281 bool HcalTopology::validCalib(const HcalCalibDetId& tid) const {
0282   bool ok(false);
0283   if (tid.calibFlavor() == HcalCalibDetId::CalibrationBox) {
0284     HcalSubdetector subdet = tid.hcalSubdet();
0285     int ieta = tid.ieta();
0286     int chan = tid.cboxChannel();
0287     unsigned int iphi = static_cast<unsigned int>(tid.iphi());
0288     if (subdet == HcalBarrel) {
0289       if ((std::find(etaCalibHB_, etaCalibHB_ + nEtaCalibHB_, ieta) != (etaCalibHB_ + nEtaCalibHB_)) &&
0290           (std::find(chanCalibHB_, chanCalibHB_ + nchanCalibHB_, chan) != (chanCalibHB_ + nchanCalibHB_)) &&
0291           (iphi >= minPhi_) && (iphi <= maxPhi_))
0292         ok = true;
0293     } else if (subdet == HcalEndcap) {
0294       if ((std::find(etaCalibHE_, etaCalibHE_ + nEtaCalibHE_, ieta) != (etaCalibHE_ + nEtaCalibHE_)) &&
0295           (std::find(chanCalibHE1_, chanCalibHE1_ + nchanCalibHE1_, chan) != (chanCalibHE1_ + nchanCalibHE1_) ||
0296            (chan == chanCalibHE2_)) &&
0297           (iphi >= minPhi_) && (iphi <= maxPhi_))
0298         ok = true;
0299     } else if (subdet == HcalForward) {
0300       if ((std::find(etaCalibHF_, etaCalibHF_ + nEtaCalibHF_, ieta) != (etaCalibHF_ + nEtaCalibHF_)) &&
0301           (std::find(chanCalibHF1_, chanCalibHF1_ + nchanCalibHF1_, chan) != (chanCalibHF1_ + nchanCalibHF1_) ||
0302            (chan == chanCalibHF2_)) &&
0303           (iphi >= minPhi_) && (iphi <= maxPhi_))
0304         ok = true;
0305     } else if (subdet == HcalOuter) {
0306       if ((std::find(etaCalibHO_, etaCalibHO_ + nEtaCalibHO_, ieta) != (etaCalibHO_ + nEtaCalibHO_)) &&
0307           (std::find(chanCalibHO_, chanCalibHO_ + nchanCalibHO_, chan) != (chanCalibHO_ + nchanCalibHO_) ||
0308            (chan == chanCalibHOs_)) &&
0309           (iphi >= minPhi_) && (iphi <= maxPhi_))
0310         ok = true;
0311     }
0312   } else if (tid.calibFlavor() == HcalCalibDetId::HOCrosstalk) {
0313     int ieta = std::abs(tid.ieta());
0314     unsigned int iphi = static_cast<unsigned int>(tid.iphi());
0315     if ((std::find(etaCalibHOX_, etaCalibHOX_ + nEtaCalibHOX_, ieta) != (etaCalibHOX_ + nEtaCalibHOX_)) &&
0316         (iphi >= minPhi_) && (iphi <= maxPhi_))
0317       ok = true;
0318   } else if (tid.calibFlavor() == HcalCalibDetId::HBX) {
0319     int ieta = std::abs(tid.ieta());
0320     unsigned int iphi = static_cast<unsigned int>(tid.iphi());
0321     if ((ieta == etaCalibHBX_) && (iphi >= minPhi_) && (iphi <= maxPhi_))
0322       ok = true;
0323   } else if (tid.calibFlavor() == HcalCalibDetId::HEX) {
0324     int ieta = std::abs(tid.ieta());
0325     unsigned int iphi = static_cast<unsigned int>(tid.iphi());
0326     if ((std::find(etaCalibHEX_, etaCalibHEX_ + nEtaCalibHEX_, ieta) != (etaCalibHEX_ + nEtaCalibHEX_)) &&
0327         (iphi >= minPhi_) && (iphi <= maxPhi_))
0328       ok = true;
0329   } else if ((tid.calibFlavor() == HcalCalibDetId::uMNqie) || (tid.calibFlavor() == HcalCalibDetId::LASERMON) ||
0330              (tid.calibFlavor() == HcalCalibDetId::CastorRadFacility)) {
0331     ok = true;
0332   }
0333   return ok;
0334 }
0335 
0336 bool HcalTopology::validHcal(const HcalDetId& id, const unsigned int flag) const {
0337   /* original logic show here because condensed form below is rather terse
0338   // check the raw rules
0339   bool ok = validHcal(id);
0340   if (flag == 0) { // This is all what is needed
0341   } else if (flag == 1) { // See if it is in the to be merged list and merged list
0342     if (hcons_->isPlan1MergedId(id))          ok = true;
0343     else if (hcons_->isPlan1ToBeMergedId(id)) ok = false;
0344   } else if (!ok) {
0345     ok = hcons_->isPlan1MergedId(id);
0346   }
0347   return ok;
0348   */
0349   return (flag > 0 and hcons_->isPlan1MergedId(id)) or
0350          ((flag != 1 or !hcons_->isPlan1ToBeMergedId(id)) and validHcal(id));
0351 }
0352 
0353 bool HcalTopology::isExcluded(const HcalDetId& id) const {
0354   bool exed = false;
0355   // first, check the full detector exclusions...  (fast)
0356   switch (id.subdet()) {
0357     case (HcalBarrel):
0358       exed = excludeHB_;
0359       break;
0360     case (HcalEndcap):
0361       exed = excludeHE_;
0362       break;
0363     case (HcalOuter):
0364       exed = excludeHO_;
0365       break;
0366     case (HcalForward):
0367       exed = excludeHF_;
0368       break;
0369     default:
0370       exed = false;
0371   }
0372   // next, check the list (slower)
0373   if (!exed && !exclusionList_.empty()) {
0374     std::vector<HcalDetId>::const_iterator i = std::lower_bound(exclusionList_.begin(), exclusionList_.end(), id);
0375     if (i != exclusionList_.end() && *i == id)
0376       exed = true;
0377   }
0378   return exed;
0379 }
0380 
0381 void HcalTopology::exclude(const HcalDetId& id) {
0382   std::vector<HcalDetId>::iterator i = std::lower_bound(exclusionList_.begin(), exclusionList_.end(), id);
0383   if (i == exclusionList_.end() || *i != id) {
0384     exclusionList_.insert(i, id);
0385   }
0386 }
0387 
0388 void HcalTopology::excludeSubdetector(HcalSubdetector subdet) {
0389   switch (subdet) {
0390     case (HcalBarrel):
0391       excludeHB_ = true;
0392       break;
0393     case (HcalEndcap):
0394       excludeHE_ = true;
0395       break;
0396     case (HcalOuter):
0397       excludeHO_ = true;
0398       break;
0399     case (HcalForward):
0400       excludeHF_ = true;
0401       break;
0402     default:
0403       break;
0404   }
0405 }
0406 
0407 std::vector<DetId> HcalTopology::east(const DetId& id) const {
0408   std::vector<DetId> vNeighborsDetId;
0409   HcalDetId neighbors[2];
0410   for (int i = 0; i < decIEta(HcalDetId(id), neighbors); i++) {
0411     if (neighbors[i].oldFormat())
0412       neighbors[i].changeForm();
0413     vNeighborsDetId.emplace_back(DetId(neighbors[i].rawId()));
0414   }
0415   return vNeighborsDetId;
0416 }
0417 
0418 std::vector<DetId> HcalTopology::west(const DetId& id) const {
0419   std::vector<DetId> vNeighborsDetId;
0420   HcalDetId neighbors[2];
0421   for (int i = 0; i < incIEta(HcalDetId(id), neighbors); i++) {
0422     if (neighbors[i].oldFormat())
0423       neighbors[i].changeForm();
0424     vNeighborsDetId.emplace_back(DetId(neighbors[i].rawId()));
0425   }
0426   return vNeighborsDetId;
0427 }
0428 
0429 std::vector<DetId> HcalTopology::north(const DetId& id) const {
0430   std::vector<DetId> vNeighborsDetId;
0431   HcalDetId neighbor;
0432   if (incIPhi(HcalDetId(id), neighbor)) {
0433     if (neighbor.oldFormat())
0434       neighbor.changeForm();
0435     vNeighborsDetId.emplace_back(DetId(neighbor.rawId()));
0436   }
0437   return vNeighborsDetId;
0438 }
0439 
0440 std::vector<DetId> HcalTopology::south(const DetId& id) const {
0441   std::vector<DetId> vNeighborsDetId;
0442   HcalDetId neighbor;
0443   if (decIPhi(HcalDetId(id), neighbor)) {
0444     if (neighbor.oldFormat())
0445       neighbor.changeForm();
0446     vNeighborsDetId.emplace_back(DetId(neighbor.rawId()));
0447   }
0448   return vNeighborsDetId;
0449 }
0450 
0451 std::vector<DetId> HcalTopology::up(const DetId& id) const {
0452   HcalDetId neighbor = id;
0453   std::vector<DetId> vNeighborsDetId;
0454   if (incrementDepth(neighbor)) {
0455     if (neighbor.oldFormat())
0456       neighbor.changeForm();
0457     vNeighborsDetId.emplace_back(neighbor);
0458   }
0459   return vNeighborsDetId;
0460 }
0461 
0462 std::vector<DetId> HcalTopology::down(const DetId& id) const {
0463   HcalDetId neighbor = id;
0464   std::vector<DetId> vNeighborsDetId;
0465   if (decrementDepth(neighbor)) {
0466     if (neighbor.oldFormat())
0467       neighbor.changeForm();
0468     vNeighborsDetId.emplace_back(neighbor);
0469   }
0470   return vNeighborsDetId;
0471 }
0472 
0473 int HcalTopology::exclude(HcalSubdetector subdet, int ieta1, int ieta2, int iphi1, int iphi2, int depth1, int depth2) {
0474   bool exed = false;
0475   // first, check the full detector exclusions...  (fast)
0476   switch (subdet) {
0477     case (HcalBarrel):
0478       exed = excludeHB_;
0479       break;
0480     case (HcalEndcap):
0481       exed = excludeHE_;
0482       break;
0483     case (HcalOuter):
0484       exed = excludeHO_;
0485       break;
0486     case (HcalForward):
0487       exed = excludeHF_;
0488       break;
0489     default:
0490       exed = false;
0491   }
0492   if (exed)
0493     return 0;  // if the whole detector is excluded...
0494 
0495   int ieta_l = std::min(ieta1, ieta2);
0496   int ieta_h = std::max(ieta1, ieta2);
0497   int iphi_l = std::min(iphi1, iphi2);
0498   int iphi_h = std::max(iphi1, iphi2);
0499   int depth_l = std::min(depth1, depth2);
0500   int depth_h = std::max(depth1, depth2);
0501 
0502   int n = 0;
0503   for (int ieta = ieta_l; ieta <= ieta_h; ieta++)
0504     for (int iphi = iphi_l; iphi <= iphi_h; iphi++)
0505       for (int depth = depth_l; depth <= depth_h; depth++) {
0506         HcalDetId id(subdet, ieta, iphi, depth);
0507         if (validRaw(id)) {  // use 'validRaw' to include check validity in "uncut" detector
0508           exclude(id);
0509           n++;
0510         }
0511       }
0512   return n;
0513 }
0514 
0515 /** Basic rules used to derive this code:
0516       
0517   HB has 72 towers in iphi.  Ieta 1-14 have depth=1, Ieta 15-16 have depth=1 or 2.
0518 
0519   HE ieta=16-20 have 72 towers in iphi
0520      ieta=21-29 have 36 towers in iphi
0521      ieta=16 is depth 3 only
0522      ieta=17 is depth 1 only
0523      ieta=18-26 & 29 have depth 1 and 2
0524      ieta=27-28 has depth 1-3
0525 
0526   HF ieta=29-39 have 36 in iphi
0527      ieta=40-41 have 18 in iphi (71,3,7,11...)
0528      all have two depths
0529 
0530 
0531   HO has 15 towers in ieta and 72 in iphi and depth = 4 (one value)
0532 
0533   At H2:
0534 
0535   HE ieta 17 is two depths
0536   HE ieta 22- have 36 towers in iphi (starts one higher)
0537   HE ieta 24- has three depths
0538 
0539   */
0540 
0541 bool HcalTopology::validDetIdPreLS1(const HcalDetId& id) const {
0542   const HcalSubdetector sd(id.subdet());
0543   const int ie(id.ietaAbs());
0544   const int ip(id.iphi());
0545   const int dp(id.depth());
0546 
0547   return ((ip >= 1) && (ip <= IPHI_MAX) && (dp >= 1) && (ie >= 1) &&
0548           (((sd == HcalBarrel) && (((ie <= 14) && (dp == 1)) || (((ie == 15) || (ie == 16)) && (dp <= 2)))) ||
0549            ((sd == HcalEndcap) &&
0550             (((ie == firstHERing()) && (dp == 3)) || ((ie == 17) && (dp == 1)) ||
0551              ((ie >= 18) && (ie <= 20) && (dp <= 2)) || ((ie >= 21) && (ie <= 26) && (dp <= 2) && (ip % 2 == 1)) ||
0552              ((ie >= 27) && (ie <= 28) && (dp <= 3) && (ip % 2 == 1)) || ((ie == 29) && (dp <= 2) && (ip % 2 == 1)))) ||
0553            ((sd == HcalOuter) && (ie <= 15) && (dp == 4)) ||
0554            ((sd == HcalForward) && (dp <= 2) &&
0555             (((ie >= firstHFRing()) && (ie < firstHFQuadPhiRing()) && (ip % 2 == 1)) ||
0556              ((ie >= firstHFQuadPhiRing()) && (ie <= lastHFRing()) && (ip % 4 == 3))))));
0557 }
0558 
0559 /** Is this a valid cell id? */
0560 bool HcalTopology::validRaw(const HcalDetId& id) const {
0561   bool ok = true;
0562   int ieta = id.ieta();
0563   int aieta = id.ietaAbs();
0564   int depth = id.depth();
0565   int iphi = id.iphi();
0566   int zside = id.zside();
0567   HcalSubdetector subdet = id.subdet();
0568   int maxPhi = (subdet == HcalEndcap) ? maxPhiHE_ : IPHI_MAX;
0569   if ((ieta == 0 || iphi <= 0 || iphi > maxPhi) || aieta > maxEta_)
0570     ok = false;  // outer limits
0571 
0572   if (ok) {
0573     if (subdet == HcalBarrel) {
0574       if ((mode_ == HcalTopologyMode::SLHC) || (mode_ == HcalTopologyMode::H2HE) || (mode_ == HcalTopologyMode::Run3) ||
0575           (mode_ == HcalTopologyMode::Run4)) {
0576         if ((aieta > lastHBRing()) || (depth > hcons_->getMaxDepth(0, aieta, iphi, zside)) ||
0577             (depth < hcons_->getMinDepth(0, aieta, iphi, zside)))
0578           ok = false;
0579       } else {
0580         if (aieta > lastHBRing() || depth > 2 || (aieta <= 14 && depth > 1))
0581           ok = false;
0582       }
0583     } else if (subdet == HcalEndcap) {
0584       if ((mode_ == HcalTopologyMode::SLHC) || (mode_ == HcalTopologyMode::H2HE) || (mode_ == HcalTopologyMode::Run3) ||
0585           (mode_ == HcalTopologyMode::Run4)) {
0586         if ((depth > hcons_->getMaxDepth(1, aieta, iphi, zside)) ||
0587             (depth < hcons_->getMinDepth(1, aieta, iphi, zside)) || (aieta < firstHERing()) || (aieta > lastHERing())) {
0588           ok = false;
0589         } else {
0590           for (const auto& i : etaBinsHE_) {
0591             if (aieta == i.ieta) {
0592               if (aieta >= firstHEDoublePhiRing() && (iphi % 2) == 0)
0593                 ok = false;
0594               if (aieta >= firstHEQuadPhiRing() && (iphi % 4) != 3)
0595                 ok = false;
0596               if (aieta + 1 == hcons_->getNoff(1)) {
0597                 if (depth < 1)
0598                   ok = false;
0599               } else {
0600                 if (depth < i.depthStart)
0601                   ok = false;
0602               }
0603               break;
0604             }
0605           }
0606         }
0607       } else {
0608         if (depth > hcons_->getMaxDepth(1, aieta, iphi, zside) || aieta < firstHERing() || aieta > lastHERing() ||
0609             (aieta == firstHERing() && depth != hcons_->getDepthEta16(2, iphi, zside)) ||
0610             (aieta == 17 && depth != 1 && mode_ != HcalTopologyMode::H2) ||  // special case at H2
0611             (((aieta >= 17 && aieta < firstHETripleDepthRing()) || aieta == lastHERing()) && depth > 2) ||
0612             (aieta >= firstHEDoublePhiRing() && (iphi % 2) == 0))
0613           ok = false;
0614       }
0615     } else if (subdet == HcalOuter) {
0616       if (aieta > lastHORing() || iphi > IPHI_MAX || depth != 4)
0617         ok = false;
0618     } else if (subdet == HcalForward) {
0619       if (aieta < firstHFRing() || aieta > lastHFRing() || ((iphi % 2) == 0) ||
0620           (depth > hcons_->maxHFDepth(ieta, iphi)) || (aieta >= firstHFQuadPhiRing() && ((iphi + 1) % 4) != 0))
0621         ok = false;
0622     } else if (subdet == HcalTriggerTower) {
0623       ok = validHT(HcalTrigTowerDetId(id.rawId()));
0624     } else if (subdet == HcalOther) {
0625       ok = validCalib(HcalCalibDetId(id.rawId()));
0626     } else {
0627       ok = false;
0628     }
0629   }
0630   return ok;
0631 }
0632 
0633 bool HcalTopology::incIPhi(const HcalDetId& id, HcalDetId& neighbor) const {
0634   bool ok = valid(id);
0635   if (ok) {
0636     switch (id.subdet()) {
0637       case (HcalBarrel):
0638       case (HcalOuter):
0639         if (id.iphi() == IPHI_MAX)
0640           neighbor = HcalDetId(id.subdet(), id.ieta(), 1, id.depth());
0641         else
0642           neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() + 1, id.depth());
0643         break;
0644       case (HcalEndcap):
0645         if (id.ietaAbs() >= firstHEQuadPhiRing()) {
0646           if (id.iphi() == IPHI_MAX - 1)
0647             neighbor = HcalDetId(id.subdet(), id.ieta(), 3, id.depth());
0648           else
0649             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() + 4, id.depth());
0650         } else if (id.ietaAbs() >= firstHEDoublePhiRing()) {
0651           if (id.iphi() == IPHI_MAX - 1)
0652             neighbor = HcalDetId(id.subdet(), id.ieta(), 1, id.depth());
0653           else
0654             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() + 2, id.depth());
0655         } else {
0656           if (id.iphi() == maxPhiHE_)
0657             neighbor = HcalDetId(id.subdet(), id.ieta(), 1, id.depth());
0658           else
0659             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() + 1, id.depth());
0660         }
0661         break;
0662       case (HcalForward):
0663         if (id.ietaAbs() >= firstHFQuadPhiRing()) {
0664           if (id.iphi() == IPHI_MAX - 1)
0665             neighbor = HcalDetId(id.subdet(), id.ieta(), 3, id.depth());
0666           else
0667             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() + 4, id.depth());
0668         } else {
0669           if (id.iphi() == IPHI_MAX - 1)
0670             neighbor = HcalDetId(id.subdet(), id.ieta(), 1, id.depth());
0671           else
0672             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() + 2, id.depth());
0673         }
0674         if (!validRaw(neighbor))
0675           ok = false;
0676         break;
0677       default:
0678         ok = false;
0679     }
0680   }
0681   return ok;
0682 }
0683 
0684 /** Get the neighbor (if present) of the given cell with lower iphi */
0685 bool HcalTopology::decIPhi(const HcalDetId& id, HcalDetId& neighbor) const {
0686   bool ok = valid(id);
0687   if (ok) {
0688     switch (id.subdet()) {
0689       case (HcalBarrel):
0690       case (HcalOuter):
0691         if (id.iphi() == 1)
0692           neighbor = HcalDetId(id.subdet(), id.ieta(), IPHI_MAX, id.depth());
0693         else
0694           neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() - 1, id.depth());
0695         break;
0696       case (HcalEndcap):
0697         if (id.ietaAbs() >= firstHEQuadPhiRing()) {
0698           if (id.iphi() == 3)
0699             neighbor = HcalDetId(id.subdet(), id.ieta(), IPHI_MAX - 1, id.depth());
0700           else
0701             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() - 4, id.depth());
0702         } else if (id.ietaAbs() >= firstHEDoublePhiRing()) {
0703           if (id.iphi() == 1)
0704             neighbor = HcalDetId(id.subdet(), id.ieta(), IPHI_MAX - 1, id.depth());
0705           else
0706             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() - 2, id.depth());
0707         } else {
0708           if (id.iphi() == 1)
0709             neighbor = HcalDetId(id.subdet(), id.ieta(), maxPhiHE_, id.depth());
0710           else
0711             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() - 1, id.depth());
0712         }
0713         break;
0714       case (HcalForward):
0715         if (id.ietaAbs() >= firstHFQuadPhiRing()) {
0716           if (id.iphi() == 3)
0717             neighbor = HcalDetId(id.subdet(), id.ieta(), IPHI_MAX - 1, id.depth());
0718           else
0719             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() - 4, id.depth());
0720         } else {
0721           if (id.iphi() == 1)
0722             neighbor = HcalDetId(id.subdet(), id.ieta(), IPHI_MAX - 1, id.depth());
0723           else
0724             neighbor = HcalDetId(id.subdet(), id.ieta(), id.iphi() - 2, id.depth());
0725         }
0726         if (!validRaw(neighbor))
0727           ok = false;
0728         break;
0729       default:
0730         ok = false;
0731     }
0732   }
0733   return ok;
0734 }
0735 
0736 int HcalTopology::incIEta(const HcalDetId& id, HcalDetId neighbors[2]) const {
0737   if (id.zside() == 1)
0738     return incAIEta(id, neighbors);
0739   else
0740     return decAIEta(id, neighbors);
0741 }
0742 
0743 int HcalTopology::decIEta(const HcalDetId& id, HcalDetId neighbors[2]) const {
0744   if (id.zside() == 1)
0745     return decAIEta(id, neighbors);
0746   else
0747     return incAIEta(id, neighbors);
0748 }
0749 
0750 /** Increasing in |ieta|, there is always at most one neighbor */
0751 int HcalTopology::incAIEta(const HcalDetId& id, HcalDetId neighbors[2]) const {
0752   int n = 1;
0753   int aieta = id.ietaAbs();
0754 
0755   if (aieta == firstHEDoublePhiRing() - 1 && (id.iphi() % 2) == 0)
0756     neighbors[0] = HcalDetId(id.subdet(), (aieta + 1) * id.zside(), id.iphi() - 1, id.depth());
0757   else if (aieta == firstHFQuadPhiRing() - 1 && ((id.iphi() + 1) % 4) != 0)
0758     neighbors[0] =
0759         HcalDetId(id.subdet(), (aieta + 1) * id.zside(), ((id.iphi() == 1) ? (71) : (id.iphi() - 2)), id.depth());
0760   else if (aieta == firstHEQuadPhiRing() - 1 && ((id.iphi() + 1) % 4) != 0)
0761     neighbors[0] =
0762         HcalDetId(id.subdet(), (aieta + 1) * id.zside(), ((id.iphi() == 1) ? (71) : (id.iphi() - 2)), id.depth());
0763   else if (aieta == lastHBRing() && id.subdet() == HcalBarrel)
0764     neighbors[0] = HcalDetId(HcalEndcap, (aieta + 1) * id.zside(), id.iphi(), 1);
0765   else if (aieta == lastHERing() && id.subdet() == HcalEndcap)
0766     neighbors[0] = HcalDetId(HcalForward, etaHE2HF_ * id.zside(), id.iphi(), 1);
0767   else
0768     neighbors[0] = HcalDetId(id.subdet(), (aieta + 1) * id.zside(), id.iphi(), id.depth());
0769 
0770   if (!valid(neighbors[0]))
0771     n = 0;
0772   return n;
0773 }
0774 
0775 /** Decreasing in |ieta|, there are two neighbors of 40 and 21*/
0776 int HcalTopology::decAIEta(const HcalDetId& id, HcalDetId neighbors[2]) const {
0777   int n = 1;
0778   int aieta = id.ietaAbs();
0779 
0780   if (aieta == firstHEDoublePhiRing()) {
0781     n = 2;
0782     neighbors[0] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), id.iphi(), id.depth());
0783     neighbors[1] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), id.iphi() + 1, id.depth());
0784   } else if (aieta == firstHFQuadPhiRing()) {
0785     n = 2;
0786     neighbors[0] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), id.iphi(), id.depth());
0787     if (id.iphi() == IPHI_MAX - 1)
0788       neighbors[1] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), 1, id.depth());
0789     else
0790       neighbors[1] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), id.iphi() + 2, id.depth());
0791   } else if (aieta == firstHEQuadPhiRing()) {
0792     n = 2;
0793     neighbors[0] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), id.iphi(), id.depth());
0794     if (id.iphi() == IPHI_MAX - 1)
0795       neighbors[1] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), 1, id.depth());
0796     else
0797       neighbors[1] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), id.iphi() + 2, id.depth());
0798   } else if (aieta == 1) {
0799     neighbors[0] = HcalDetId(id.subdet(), -aieta * id.zside(), id.iphi(), id.depth());
0800   } else if (aieta == firstHERing() && id.subdet() == HcalEndcap) {
0801     neighbors[0] = HcalDetId(HcalBarrel, (aieta - 1) * id.zside(), id.iphi(), 1);
0802   } else if (aieta == firstHFRing() && id.subdet() == HcalForward) {
0803     neighbors[0] = HcalDetId(HcalEndcap, etaHF2HE_ * id.zside(), id.iphi(), 1);
0804   } else
0805     neighbors[0] = HcalDetId(id.subdet(), (aieta - 1) * id.zside(), id.iphi(), id.depth());
0806 
0807   if (!valid(neighbors[0]) && n == 2) {
0808     if (!valid(neighbors[1]))
0809       n = 0;
0810     else {
0811       n = 1;
0812       neighbors[0] = neighbors[1];
0813     }
0814   }
0815   if (n == 2 && !valid(neighbors[1]))
0816     n = 1;
0817   if (n == 1 && !valid(neighbors[0]))
0818     n = 0;
0819 
0820   return n;
0821 }
0822 
0823 void HcalTopology::depthBinInformation(
0824     HcalSubdetector subdet, int etaRing, int iphi, int zside, int& nDepthBins, int& startingBin) const {
0825   if (subdet == HcalBarrel) {
0826     if ((mode_ == HcalTopologyMode::SLHC) || (mode_ == HcalTopologyMode::H2HE) || (mode_ == HcalTopologyMode::Run3) ||
0827         (mode_ == HcalTopologyMode::Run4)) {
0828       startingBin = hcons_->getMinDepth(0, etaRing, iphi, zside);
0829       if (etaRing == lastHBRing()) {
0830         nDepthBins = hcons_->getDepthEta16(1, iphi, zside) - startingBin + 1;
0831       } else {
0832         nDepthBins = hcons_->getMaxDepth(0, etaRing, iphi, zside) - startingBin + 1;
0833       }
0834     } else {
0835       if (etaRing <= 14) {
0836         nDepthBins = 1;
0837         startingBin = 1;
0838       } else {
0839         nDepthBins = 2;
0840         startingBin = 1;
0841       }
0842     }
0843   } else if (subdet == HcalEndcap) {
0844     if ((mode_ == HcalTopologyMode::SLHC) || (mode_ == HcalTopologyMode::H2HE) || (mode_ == HcalTopologyMode::Run3) ||
0845         (mode_ == HcalTopologyMode::Run4)) {
0846       if (etaRing == firstHERing()) {
0847         startingBin = hcons_->getDepthEta16(2, iphi, zside);
0848       } else {
0849         startingBin = hcons_->getMinDepth(1, etaRing, iphi, zside);
0850       }
0851       nDepthBins = hcons_->getMaxDepth(1, etaRing, iphi, zside) - startingBin + 1;
0852     } else {
0853       if (etaRing == firstHERing()) {
0854         nDepthBins = 1;
0855         startingBin = 3;
0856       } else if (etaRing == 17) {
0857         nDepthBins = 1;
0858         startingBin = 1;
0859       } else if (etaRing == lastHERing()) {
0860         nDepthBins = 2;
0861         startingBin = 1;
0862       } else {
0863         nDepthBins = (etaRing >= firstHETripleDepthRing()) ? 3 : 2;
0864         startingBin = 1;
0865       }
0866     }
0867   } else if (subdet == HcalForward) {
0868     nDepthBins = maxDepthHF_;
0869     startingBin = 1;
0870   } else if (subdet == HcalOuter) {
0871     nDepthBins = 1;
0872     startingBin = 4;
0873   } else {
0874     edm::LogWarning("HCalGeom") << "Bad HCAL subdetector " << subdet;
0875   }
0876 }
0877 
0878 bool HcalTopology::incrementDepth(HcalDetId& detId) const {
0879   HcalSubdetector subdet = detId.subdet();
0880   int ieta = detId.ieta();
0881   int etaRing = detId.ietaAbs();
0882   int depth = detId.depth();
0883   int iphi = detId.iphi();
0884   int zside = detId.zside();
0885   int nDepthBins(0), startingBin(0);
0886   depthBinInformation(subdet, etaRing, iphi, zside, nDepthBins, startingBin);
0887 
0888   // see if the new depth bin exists
0889   ++depth;
0890   if (depth >= (startingBin + nDepthBins)) {
0891     // handle on a case-by-case basis
0892     if (subdet == HcalBarrel && etaRing < lastHORing()) {
0893       // HO
0894       subdet = HcalOuter;
0895       depth = 4;
0896     } else if (subdet == HcalBarrel && etaRing == lastHBRing()) {
0897       // overlap
0898       subdet = HcalEndcap;
0899       if ((mode_ == HcalTopologyMode::SLHC) || (mode_ == HcalTopologyMode::H2HE) || (mode_ == HcalTopologyMode::Run3) ||
0900           (mode_ == HcalTopologyMode::Run4))
0901         depth = hcons_->getDepthEta16(2, iphi, zside);
0902     } else if ((subdet == HcalEndcap) && (etaRing == lastHERing() - 1) && (mode_ != HcalTopologyMode::SLHC) &&
0903                (mode_ != HcalTopologyMode::Run3) && (mode_ != HcalTopologyMode::Run4)) {
0904       // guard ring HF29 is behind HE 28
0905       subdet = HcalForward;
0906       (ieta > 0) ? ++ieta : --ieta;
0907       depth = 1;
0908     } else if ((subdet == HcalEndcap) && (etaRing == lastHERing()) && (mode_ != HcalTopologyMode::SLHC) &&
0909                (mode_ != HcalTopologyMode::Run3) && (mode_ != HcalTopologyMode::Run4)) {
0910       // split cells go to bigger granularity.  Ring 29 -> 28
0911       (ieta > 0) ? --ieta : ++ieta;
0912     } else {
0913       // no more chances
0914       detId = HcalDetId();
0915       return false;
0916     }
0917   }
0918   detId = HcalDetId(subdet, ieta, iphi, depth);
0919   return validRaw(detId);
0920 }
0921 
0922 bool HcalTopology::decrementDepth(HcalDetId& detId) const {
0923   HcalSubdetector subdet = detId.subdet();
0924   int ieta = detId.ieta();
0925   int etaRing = detId.ietaAbs();
0926   int depth = detId.depth();
0927   int iphi = detId.iphi();
0928   int zside = detId.zside();
0929   int nDepthBins, startingBin;
0930   depthBinInformation(subdet, etaRing, iphi, zside, nDepthBins, startingBin);
0931 
0932   // see if the new depth bin exists
0933   --depth;
0934   if ((subdet == HcalOuter) || (subdet == HcalEndcap && etaRing == firstHERing())) {
0935     subdet = HcalBarrel;
0936     for (int i = 0; i < nEtaHB_; ++i) {
0937       if (etaRing == etaBinsHB_[i].ieta) {
0938         depth = etaBinsHB_[i].depthStart + etaBinsHB_[i].layer.size() - 1;
0939         break;
0940       }
0941     }
0942   } else if ((subdet == HcalEndcap) && (etaRing == lastHERing()) && (depth == hcons_->getDepthEta29(iphi, zside, 0)) &&
0943              (mode_ != HcalTopologyMode::SLHC) && (mode_ != HcalTopologyMode::Run3) &&
0944              (mode_ != HcalTopologyMode::Run4)) {
0945     (ieta > 0) ? --ieta : ++ieta;
0946   } else if (depth <= 0) {
0947     if (subdet == HcalForward && etaRing == firstHFRing()) {
0948       // overlap
0949       subdet = HcalEndcap;
0950       etaRing = etaHF2HE_;
0951       ieta = (ieta > 0) ? etaRing : -etaRing;
0952       for (const auto& i : etaBinsHE_) {
0953         if (etaRing == i.ieta) {
0954           depth = i.depthStart + i.layer.size() - 1;
0955           break;
0956         }
0957       }
0958     } else {
0959       // no more chances
0960       detId = HcalDetId();
0961       return false;
0962     }
0963   }
0964   detId = HcalDetId(subdet, ieta, detId.iphi(), depth);
0965   return validRaw(detId);
0966 }
0967 
0968 int HcalTopology::nPhiBins(int etaRing) const {
0969   int lastPhiBin = singlePhiBins_;
0970   if (etaRing >= firstHFQuadPhiRing() || etaRing >= firstHEQuadPhiRing())
0971     lastPhiBin = doublePhiBins_ / 2;
0972   else if (etaRing >= firstHEDoublePhiRing())
0973     lastPhiBin = doublePhiBins_;
0974   if (hcons_ && etaRing >= hcons_->getEtaRange(1).first && etaRing <= hcons_->getEtaRange(1).second) {
0975     return nPhiBins(HcalBarrel, etaRing);
0976   }
0977   return lastPhiBin;
0978 }
0979 
0980 int HcalTopology::nPhiBins(HcalSubdetector bc, int etaRing) const {
0981   double phiTableVal;
0982   if (bc == HcalForward) {
0983     phiTableVal = dPhiTableHF[etaRing - firstHFRing_];
0984   } else {
0985     phiTableVal = dPhiTable[etaRing - firstHBRing_];
0986   }
0987   int lastPhiBin = 0;
0988   if (phiTableVal != 0.0)
0989     lastPhiBin = static_cast<int>((2._pi / phiTableVal) + 0.001);
0990   return lastPhiBin;
0991 }
0992 
0993 int HcalTopology::maxDepth() const {
0994   int maxd1 = std::max(maxDepthHB_, maxDepthHE_);
0995   int maxd2 = std::max(maxDepthHF_, minMaxDepth_);
0996   return std::max(maxd1, maxd2);
0997 }
0998 
0999 int HcalTopology::maxDepth(HcalSubdetector bc) const {
1000   if (bc == HcalBarrel)
1001     return maxDepthHB_;
1002   else if (bc == HcalEndcap)
1003     return maxDepthHE_;
1004   else if (bc == HcalForward)
1005     return maxDepthHF_;
1006   else
1007     return 4;
1008 }
1009 
1010 int HcalTopology::etaRing(HcalSubdetector bc, double abseta) const {
1011   int etaring = firstHBRing_;
1012   if (bc == HcalForward) {
1013     etaring = firstHFRing_;
1014     for (unsigned int k = 0; k < etaTableHF.size() - 1; ++k) {
1015       if (abseta < etaTableHF[k + 1]) {
1016         etaring += k;
1017         break;
1018       }
1019     }
1020   } else {
1021     for (unsigned int k = 0; k < etaTable.size() - 1; ++k) {
1022       if (abseta < etaTable[k + 1]) {
1023         etaring += k;
1024         break;
1025       }
1026     }
1027     if (abseta >= etaTable[etaTable.size() - 1])
1028       etaring = lastHERing_;
1029   }
1030   return etaring;
1031 }
1032 
1033 int HcalTopology::phiBin(HcalSubdetector bc, int etaring, double phi) const {
1034   //put phi in correct range (0->2pi)
1035   int index(0);
1036   if (bc == HcalBarrel) {
1037     index = (etaring - firstHBRing_);
1038     phi -= phioff[0];
1039   } else if (bc == HcalEndcap) {
1040     index = (etaring - firstHBRing_);
1041     phi -= phioff[1];
1042   } else if (bc == HcalForward) {
1043     index = (etaring - firstHFRing_);
1044     if (index < static_cast<int>(dPhiTableHF.size())) {
1045       if (index >= 0 && unitPhiHF[index] > 2)
1046         phi -= phioff[4];
1047       else
1048         phi -= phioff[2];
1049     }
1050   }
1051   if (index < 0)
1052     index = 0;
1053   if (phi < 0.0)
1054     phi += 2._pi;
1055   else if (phi > 2._pi)
1056     phi -= 2._pi;
1057   int phibin(1), unit(1);
1058   if (bc == HcalForward) {
1059     if (index < (int)(dPhiTableHF.size())) {
1060       unit = unitPhiHF[index];
1061       phibin = static_cast<int>(phi / dPhiTableHF[index]) + 1;
1062     }
1063   } else {
1064     if (index < (int)(dPhiTable.size())) {
1065       phibin = static_cast<int>(phi / dPhiTable[index]) + 1;
1066       unit = unitPhi[index];
1067     }
1068   }
1069   int iphi(phibin);
1070   if (unit == 2)
1071     iphi = (phibin - 1) * 2 + 1;
1072   else if (unit == 4)
1073     iphi = (phibin - 1) * 4 + 3;
1074   return iphi;
1075 }
1076 
1077 void HcalTopology::getDepthSegmentation(const unsigned ring, std::vector<int>& readoutDepths, const bool one) const {
1078   // if it doesn't exist, return the first entry with a lower index.  So if we only
1079   // have entries for 1 and 17, any input from 1-16 should return the entry for ring 1
1080   SegmentationMap::const_iterator pos;
1081   if (!one) {
1082     pos = depthSegmentation_.upper_bound(ring);
1083     if (pos == depthSegmentation_.begin()) {
1084       throw cms::Exception("HcalTopology") << "No depth segmentation found for ring" << ring;
1085     }
1086   } else {
1087     pos = depthSegmentationOne_.upper_bound(ring);
1088     if (pos == depthSegmentationOne_.begin()) {
1089       throw cms::Exception("HcalTopology") << "No depth segmentation found for ring" << ring;
1090     }
1091   }
1092   --pos;
1093   // pos now refers to the last element with key <= ring.
1094   readoutDepths = pos->second;
1095 }
1096 
1097 void HcalTopology::setDepthSegmentation(const unsigned ring, const std::vector<int>& readoutDepths, const bool one) {
1098   if (one) {
1099     depthSegmentationOne_[ring] = readoutDepths;
1100   } else {
1101     depthSegmentation_[ring] = readoutDepths;
1102   }
1103 }
1104 
1105 std::pair<int, int> HcalTopology::segmentBoundaries(const unsigned ring, const unsigned depth, const bool one) const {
1106   std::vector<int> readoutDepths;
1107   getDepthSegmentation(ring, readoutDepths, one);
1108   int d1 = std::lower_bound(readoutDepths.begin(), readoutDepths.end(), depth) - readoutDepths.begin();
1109   int d2 = std::upper_bound(readoutDepths.begin(), readoutDepths.end(), depth) - readoutDepths.begin();
1110   return std::pair<int, int>(d1, d2);
1111 }
1112 
1113 double HcalTopology::etaMax(HcalSubdetector subdet) const {
1114   double eta(0);
1115   switch (subdet) {
1116     case (HcalBarrel):
1117       if (lastHBRing_ < (int)(etaTable.size()))
1118         eta = etaTable[lastHBRing_];
1119       break;
1120     case (HcalEndcap):
1121       if (lastHERing_ < (int)(etaTable.size()) && nEtaHE_ > 0)
1122         eta = etaTable[lastHERing_];
1123       break;
1124     case (HcalOuter):
1125       if (lastHORing_ < (int)(etaTable.size()))
1126         eta = etaTable[lastHORing_];
1127       break;
1128     case (HcalForward):
1129       if (!etaTableHF.empty())
1130         eta = etaTableHF[etaTableHF.size() - 1];
1131       break;
1132     default:
1133       eta = 0;
1134   }
1135   return eta;
1136 }
1137 std::pair<double, double> HcalTopology::etaRange(HcalSubdetector subdet, int keta) const {
1138   int ieta = (keta > 0) ? keta : -keta;
1139   if (subdet == HcalForward) {
1140     if (ieta >= firstHFRing_) {
1141       unsigned int ii = (unsigned int)(ieta - firstHFRing_);
1142       if (ii + 1 < etaTableHF.size())
1143         return std::pair<double, double>(etaTableHF[ii], etaTableHF[ii + 1]);
1144     }
1145   } else {
1146     int ietal = (mode_ == HcalTopologyMode::LHC && ieta == lastHERing_ - 1) ? (ieta + 1) : ieta;
1147     if ((ietal < (int)(etaTable.size())) && (ieta > 0))
1148       return std::pair<double, double>(etaTable[ieta - 1], etaTable[ietal]);
1149   }
1150   return std::pair<double, double>(0, 0);
1151 }
1152 
1153 unsigned int HcalTopology::detId2denseIdPreLS1(const DetId& id) const {
1154   HcalDetId hid(id);
1155   const HcalSubdetector sd(hid.subdet());
1156   const int ip(hid.iphi());
1157   const int ie(hid.ietaAbs());
1158   const int dp(hid.depth());
1159   const int zn(hid.zside() < 0 ? 1 : 0);
1160   unsigned int retval =
1161       ((sd == HcalBarrel)
1162            ? (ip - 1) * 18 + dp - 1 + ie - (ie < 16 ? 1 : 0) + zn * kHBhalf
1163            : ((sd == HcalEndcap)
1164                   ? 2 * kHBhalf + (ip - 1) * 8 + (ip / 2) * 20 +
1165                         ((ie == 16 || ie == 17)
1166                              ? ie - 16
1167                              : ((ie >= 18 && ie <= 20)
1168                                     ? 2 + 2 * (ie - 18) + dp - 1
1169                                     : ((ie >= 21 && ie <= 26)
1170                                            ? 8 + 2 * (ie - 21) + dp - 1
1171                                            : ((ie >= 27 && ie <= 28) ? 20 + 3 * (ie - 27) + dp - 1
1172                                                                      : 26 + 2 * (ie - 29) + dp - 1)))) +
1173                         zn * kHEhalf
1174                   : ((sd == HcalOuter)
1175                          ? 2 * kHBhalf + 2 * kHEhalf + (ip - 1) * 15 + (ie - 1) + zn * kHOhalf
1176                          : ((sd == HcalForward) ? 2 * kHBhalf + 2 * kHEhalf + 2 * kHOhalf + ((ip - 1) / 4) * 4 +
1177                                                       ((ip - 1) / 2) * 22 + 2 * (ie - 29) + (dp - 1) + zn * kHFhalf
1178                                                 : 0xFFFFFFFFu))));
1179   return retval;
1180 }
1181 
1182 unsigned int HcalTopology::detId2denseIdHB(const DetId& id) const {
1183   HcalDetId hid(id);
1184   const int ip(hid.iphi());
1185   const int ie(hid.ietaAbs());
1186   const int dp(hid.depth());
1187   const int zn(hid.zside() < 0 ? 1 : 0);
1188   unsigned int retval = 0xFFFFFFFFu;
1189   if (topoVersion_ == 0) {
1190     retval = (ip - 1) * 18 + dp - 1 + ie - (ie < 16 ? 1 : 0) + zn * kHBhalf;
1191   } else if (topoVersion_ == 10) {
1192     retval = (dp - 1) + maxDepthHB_ * (ip - 1);
1193     if (hid.ieta() > 0)
1194       retval += maxDepthHB_ * IPHI_MAX * (hid.ieta() - firstHBRing());
1195     else
1196       retval += maxDepthHB_ * IPHI_MAX * (hid.ieta() + lastHBRing() + nEtaHB_);
1197   }
1198   return retval;
1199 }
1200 
1201 unsigned int HcalTopology::detId2denseIdHE(const DetId& id) const {
1202   HcalDetId hid(id);
1203   const int ip(hid.iphi());
1204   const int ie(hid.ietaAbs());
1205   const int dp(hid.depth());
1206   const int zn(hid.zside() < 0 ? 1 : 0);
1207   unsigned int retval = 0xFFFFFFFFu;
1208   if (topoVersion_ == 0) {
1209     retval = (ip - 1) * 8 + (ip / 2) * 20 +
1210              ((ie == 16 || ie == 17)
1211                   ? ie - 16
1212                   : ((ie >= 18 && ie <= 20)
1213                          ? 2 + 2 * (ie - 18) + dp - 1
1214                          : ((ie >= 21 && ie <= 26) ? 8 + 2 * (ie - 21) + dp - 1
1215                                                    : ((ie >= 27 && ie <= 28) ? 20 + 3 * (ie - 27) + dp - 1
1216                                                                              : 26 + 2 * (ie - 29) + dp - 1)))) +
1217              zn * kHEhalf;
1218   } else if (topoVersion_ == 10) {
1219     retval = (dp - 1) + maxDepthHE_ * (ip - 1);
1220     if (hid.ieta() > 0)
1221       retval += maxDepthHE_ * maxPhiHE_ * (hid.ieta() - firstHERing());
1222     else
1223       retval += maxDepthHE_ * maxPhiHE_ * (hid.ieta() + lastHERing() + nEtaHE_);
1224   }
1225   return retval;
1226 }
1227 
1228 unsigned int HcalTopology::detId2denseIdHO(const DetId& id) const {
1229   HcalDetId hid(id);
1230   const int ip(hid.iphi());
1231   const int ie(hid.ietaAbs());
1232   const int zn(hid.zside() < 0 ? 1 : 0);
1233 
1234   unsigned int retval = 0xFFFFFFFFu;
1235   if (topoVersion_ == 0) {
1236     retval = (ip - 1) * 15 + (ie - 1) + zn * kHOhalf;
1237   } else if (topoVersion_ == 10) {
1238     if (hid.ieta() > 0)
1239       retval = (ip - 1) + IPHI_MAX * (hid.ieta() - 1);
1240     else
1241       retval = (ip - 1) + IPHI_MAX * (30 + hid.ieta());
1242   }
1243   return retval;
1244 }
1245 
1246 unsigned int HcalTopology::detId2denseIdHF(const DetId& id) const {
1247   HcalDetId hid(id);
1248   const int ip(hid.iphi());
1249   const int ie(hid.ietaAbs());
1250   const int dp(hid.depth());
1251   const int zn(hid.zside() < 0 ? 1 : 0);
1252 
1253   unsigned int retval = 0xFFFFFFFFu;
1254   if (topoVersion_ == 0) {
1255     retval = ((ip - 1) / 4) * 4 + ((ip - 1) / 2) * 22 + 2 * (ie - 29) + (dp - 1) + zn * kHFhalf;
1256   } else if (topoVersion_ == 10) {
1257     retval = dp - 1 + 2 * (ip - 1);
1258     if (hid.ieta() > 0)
1259       retval += maxDepthHF_ * IPHI_MAX * (hid.ieta() - 29);
1260     else
1261       retval += maxDepthHF_ * IPHI_MAX * ((41 + 13) + hid.ieta());
1262   }
1263   return retval;
1264 }
1265 
1266 unsigned int HcalTopology::detId2denseIdHT(const DetId& id) const {
1267   HcalTrigTowerDetId tid(id);
1268   int zside = tid.zside();
1269   unsigned int ietaAbs = tid.ietaAbs();
1270   unsigned int iphi = tid.iphi();
1271   unsigned int ivers = tid.version();
1272 
1273   unsigned int index;
1274   if (ivers == 0) {
1275     if ((iphi - 1) % 4 == 0)
1276       index = (iphi - 1) * 32 + (ietaAbs - 1) - (12 * ((iphi - 1) / 4));
1277     else
1278       index = (iphi - 1) * 28 + (ietaAbs - 1) + (4 * (((iphi - 1) / 4) + 1));
1279     if (zside == -1)
1280       index += kHThalf;
1281   } else {
1282     index = kHTSizePreLS1;
1283     if (zside == -1)
1284       index += ((kHTSizePhase1 - kHTSizePreLS1) / 2);
1285     index += (36 * (ietaAbs - 30) + ((iphi - 1) / 2));
1286   }
1287 
1288   return index;
1289 }
1290 
1291 unsigned int HcalTopology::detId2denseIdCALIB(const DetId& id) const {
1292   HcalCalibDetId tid(id);
1293   int channel = tid.cboxChannel();
1294   int ieta = tid.ieta();
1295   int iphi = tid.iphi();
1296   int zside = tid.zside();
1297   unsigned int index = 0xFFFFFFFFu;
1298 
1299   if (tid.calibFlavor() == HcalCalibDetId::CalibrationBox) {
1300     HcalSubdetector subDet = tid.hcalSubdet();
1301 
1302     if (subDet == HcalBarrel) {
1303 #ifdef EDM_ML_DEBUG
1304       edm::LogVerbatim("HCalGeom") << "CALIB_HB:  dphi = 4 (18 phi values), 3 channel types (0,1,2), eta = -1 or 1\n   "
1305                                       "        total of 18*3*2=108 channels";
1306 #endif
1307       auto indx = std::find(etaCalibHB_, etaCalibHB_ + nEtaCalibHB_, ieta);
1308       auto kndx = std::find(chanCalibHB_, chanCalibHB_ + nchanCalibHB_, channel);
1309       if (indx != etaCalibHB_ + nEtaCalibHB_ && kndx != chanCalibHB_ + nchanCalibHB_) {
1310         int keta = static_cast<int>(indx - etaCalibHB_);
1311         int kchn = static_cast<int>(kndx - chanCalibHB_);
1312         index = ((iphi + 1) / mPhiCalibHB_ - 1) + kPhiCalibHB_ * kchn + keta * kchanCalibHB_ + kOffCalibHB_;
1313       }
1314     } else if (subDet == HcalEndcap) {
1315 #ifdef EDM_ML_DEBUG
1316       edm::LogVerbatim("HCalGeom") << "CALIB_HE:  dphi = 4 (18 phi values), 7 channel types (0,1,2,3,4,5,6), eta = "
1317                                       "-1/+1\n           total of 18*7*2=252 channels    if (channel > 2) channel -= 1";
1318 #endif
1319       auto indx = std::find(etaCalibHE_, etaCalibHE_ + nEtaCalibHE_, ieta);
1320       if (indx != etaCalibHE_ + nEtaCalibHE_) {
1321         int keta = static_cast<int>(indx - etaCalibHE_);
1322         auto kndx = std::find(chanCalibHE1_, chanCalibHE1_ + nchanCalibHE1_, channel);
1323         if (kndx != chanCalibHE1_ + nchanCalibHE1_) {
1324           int kchn = static_cast<int>(kndx - chanCalibHE1_);
1325           index = ((iphi + 1) / mPhiCalibHE_ - 1) + kPhiCalibHE_ * kchn + keta * kchanCalibHE1_ + kOffCalibHE1_;
1326         } else if (channel == chanCalibHE2_) {
1327           index = ((iphi + 1) / mPhiCalibHE_ - 1) + keta * kchanCalibHE2_ + kOffCalibHE2_;
1328         }
1329       }
1330     } else if (subDet == HcalForward) {
1331 #ifdef EDM_ML_DEBUG
1332       edm::LogVerbatim("HCalGeom")
1333           << "CALIB_HF:  dphi = 18 (4 phi values), 3 channel types (0,1,8), eta = -1 or +1\n           or iphi = 1; "
1334              "channel = 9, eta = -1 or +1      total channels 4*3*2=24 + 2";
1335 #endif
1336       auto indx = std::find(etaCalibHF_, etaCalibHF_ + nEtaCalibHF_, ieta);
1337       if (indx != etaCalibHF_ + nEtaCalibHF_) {
1338         int keta = static_cast<int>(indx - etaCalibHF_);
1339         auto kndx = std::find(chanCalibHF1_, chanCalibHF1_ + nchanCalibHF1_, channel);
1340         if (kndx != chanCalibHF1_ + nchanCalibHF1_) {
1341           int kchn = static_cast<int>(kndx - chanCalibHF1_);
1342           index = ((iphi - 1) / mPhiCalibHF1_) + kPhiCalibHF1_ * kchn + keta * kchanCalibHF1_ + kOffCalibHF1_;
1343         } else if (channel == chanCalibHF2_) {
1344           index = keta * kchanCalibHF2_ + kOffCalibHF2_;
1345         }
1346       }
1347     } else if (subDet == HcalOuter) {
1348 #ifdef EDM_ML_DEBUG
1349       edm::LogVerbatim("HCalGeom") << "CALIB_HO:  ";
1350 #endif
1351       auto indx = std::find(etaCalibHO_, etaCalibHO_ + nEtaCalibHO_, ieta);
1352       if (indx != etaCalibHO_ + nEtaCalibHO_) {
1353         int keta = static_cast<int>(indx - etaCalibHO_);
1354         auto kndx = std::find(chanCalibHO_, chanCalibHO_ + nchanCalibHO_, channel);
1355         //there are 5 special calib crosstalk channels, one in each ring
1356         if (channel == chanCalibHOs_) {
1357           index = keta + kOffCalibHO2_;
1358         } else if (kndx != chanCalibHO_ + nchanCalibHO_) {
1359           //for HOM/HOP dphi = 6 (12 phi values),  2 channel types (0,1), eta = -2,-1 or 1,2
1360           //for HO0/YB0 dphi = 12 (6 phi values),  2 channel types (0,1), eta = 0
1361           int kchn = static_cast<int>(kndx - chanCalibHO_);
1362           int kphi = (ieta == 0) ? ((iphi + 1) / mPhiCalibHO0_ - 1) : ((iphi + 1) / mPhiCalibHO1_ - 1);
1363           if (ieta < 0) {
1364             index = kphi + kPhiCalibHO2_ * kchn + kPhiCalibHO1_ * keta + kOffCalibHO1_;
1365           } else if (ieta > 0) {
1366             index = kphi + kPhiCalibHO2_ * kchn + kPhiCalibHO1_ * (keta - 1) + kPhiCalibHO0_ + kOffCalibHO1_;
1367           } else {
1368             index = kphi + kPhiCalibHO2_ * kchn + kPhiCalibHO1_ * keta + kOffCalibHO1_;
1369           }
1370         }
1371       }
1372     } else {
1373       edm::LogWarning("HGCalGeom") << "HCAL Det Id not valid!";
1374     }
1375   } else if (tid.calibFlavor() == HcalCalibDetId::HOCrosstalk) {
1376 #ifdef EDM_ML_DEBUG
1377     edm::LogVerbatim("HCalGeom") << "HX:  for YB0/HO0 phi is grouped in 6 groups of 6 with dphi=2 but the transitions "
1378                                     "are 1 or 3 in such a way that the %36 operation yeilds unique values for every "
1379                                     "iphi\n     ieta = 0 for HO2M/HO1M ieta=2 for HO1P/HO2P; /ieta = 1 YB0/HO0";
1380 #endif
1381     int kphi = ((std::abs(ieta) == etaCalibHOX_[0]) ? ((iphi % 2 == 0) ? (iphi / 2 - 1) : (iphi - 1) / 2) : (iphi - 1));
1382     if (std::abs(ieta) == etaCalibHOX_[0]) {
1383       index = kphi + (((zside + 1) * nPhiCalibHOX_[0]) / 2) + nPhiCalibHOX_[1] + kOffCalibHOX_;
1384     } else if (std::abs(ieta) == etaCalibHOX_[1]) {
1385       index = kphi + ((zside + 1) * nPhiCalibHOX_[1]) + kOffCalibHOX_;
1386     }
1387   } else if (tid.calibFlavor() == HcalCalibDetId::HBX) {
1388     if (std::abs(ieta) == etaCalibHBX_) {
1389       index = kOffCalibHBX_ + (iphi - 1) + (zside + 1) * kPhiCalibHBX_ / 2;
1390     }
1391   } else if (tid.calibFlavor() == HcalCalibDetId::HEX) {
1392     if (std::abs(ieta) == etaCalibHEX_[0]) {
1393       index = kOffCalibHEX_ + (iphi - 1) / mPhiCalibHEX_ + (zside + 1) * kPhiCalibHEX_ / 2;
1394     } else if (std::abs(ieta) == etaCalibHEX_[1]) {
1395       index = kOffCalibHEX_ + (iphi - 1) / mPhiCalibHEX_ + (zside + 1) * kPhiCalibHEX_ / 2 + 2 * kPhiCalibHEX_;
1396     }
1397   }
1398 #ifdef EDM_ML_DEBUG
1399   edm::LogVerbatim("HCalGeom") << "  " << ieta << "  " << zside << "  " << iphi << "  " << index;
1400 #endif
1401   return index;
1402 }
1403 
1404 HcalCalibDetId HcalTopology::denseId2detIdCALIB(const unsigned int& hid) const {
1405   HcalCalibDetId id;
1406   unsigned int hid0(hid);
1407   if (hid0 < kOffCalibHOX_) {
1408     HcalSubdetector subdet(HcalEmpty);
1409     int ieta(0), iphi(0), ichan(0), ctype(0);
1410     if (hid0 < kOffCalibHE1_) {
1411       int id0 = static_cast<int>(hid0);
1412       subdet = HcalBarrel;
1413       iphi = hid0 % kPhiCalibHB_;
1414       int keta = (hid0 < kchanCalibHB_) ? 0 : 1;
1415       ieta = etaCalibHB_[keta];
1416       ichan = (id0 - iphi - keta * kchanCalibHB_) / kPhiCalibHB_;
1417       iphi = mPhiCalibHB_ * (iphi + 1) - 1;
1418       ctype = chanCalibHB_[ichan];
1419     } else if (hid0 < kOffCalibHF1_) {
1420       hid0 -= kOffCalibHE1_;
1421       int id0 = static_cast<int>(hid0);
1422       subdet = HcalEndcap;
1423       iphi = hid0 % kPhiCalibHE_;
1424       int keta = (hid0 < kchanCalibHE1_) ? 0 : 1;
1425       ieta = etaCalibHE_[keta];
1426       ichan = (id0 - iphi - keta * kchanCalibHE1_) / kPhiCalibHE_;
1427       iphi = mPhiCalibHE_ * (iphi + 1) - 1;
1428       ctype = chanCalibHE1_[ichan];
1429     } else if (hid0 < kOffCalibHO1_) {
1430       hid0 -= kOffCalibHF1_;
1431       int id0 = static_cast<int>(hid0);
1432       subdet = HcalForward;
1433       iphi = hid0 % kPhiCalibHF1_;
1434       int keta = (hid0 < kchanCalibHF1_) ? 0 : 1;
1435       ieta = etaCalibHF_[keta];
1436       ichan = (id0 - iphi - keta * kchanCalibHF1_) / kPhiCalibHF1_;
1437       iphi = mPhiCalibHF1_ * iphi + 1;
1438       ctype = chanCalibHF1_[ichan];
1439     } else if (hid0 < kOffCalibHO2_) {
1440       hid0 -= kOffCalibHO1_;
1441       int id0 = static_cast<int>(hid0);
1442       subdet = HcalOuter;
1443       unsigned int kphi = hid0 % kPhiCalibHO2_;
1444       if (kphi < 2 * kPhiCalibHO1_) {
1445         ieta = (kphi >= kPhiCalibHO1_) ? etaCalibHO_[1] : etaCalibHO_[0];
1446         iphi = kphi % kPhiCalibHO1_;
1447         ichan = (id0 - iphi - (ieta + 2) * kPhiCalibHO1_) / kPhiCalibHO2_;
1448         iphi = (iphi + 1) * mPhiCalibHO1_ - 1;
1449       } else if (kphi < (2 * kPhiCalibHO1_ + kPhiCalibHO0_)) {
1450         ieta = etaCalibHO_[2];
1451         iphi = kphi % kPhiCalibHO0_;
1452         ichan = (id0 - iphi - (ieta + 2) * kPhiCalibHO1_) / kPhiCalibHO2_;
1453         iphi = (iphi + 1) * mPhiCalibHO0_ - 1;
1454       } else {
1455         ieta = (kphi >= 3 * kPhiCalibHO1_ + kPhiCalibHO0_) ? etaCalibHO_[4] : etaCalibHO_[3];
1456         iphi = kphi % kPhiCalibHO1_;
1457         ichan = (id0 - iphi - (ieta + 1) * kPhiCalibHO1_ - kPhiCalibHO0_) / kPhiCalibHO2_;
1458         iphi = (iphi + 1) * mPhiCalibHO1_ - 1;
1459       }
1460       ctype = chanCalibHO_[ichan];
1461     } else if (hid0 < kOffCalibHE2_) {
1462       hid0 -= kOffCalibHO2_;
1463       subdet = HcalOuter;
1464       iphi = phiCalibHO_[hid0];
1465       ctype = static_cast<int>(chanCalibHOs_);
1466       ieta = etaCalibHO_[hid0];
1467     } else if (hid0 < kOffCalibHF2_) {
1468       hid0 -= kOffCalibHE2_;
1469       subdet = HcalEndcap;
1470       iphi = hid0 % kPhiCalibHE_;
1471       int keta = (hid0 < kchanCalibHE2_) ? 0 : 1;
1472       ieta = etaCalibHE_[keta];
1473       iphi = mPhiCalibHE_ * (iphi + 1) - 1;
1474       ctype = chanCalibHE2_;
1475     } else {
1476       hid0 -= kOffCalibHF2_;
1477       subdet = HcalForward;
1478       int keta = (hid0 < kchanCalibHF2_) ? 0 : 1;
1479       ieta = etaCalibHF_[keta];
1480       iphi = phiCalibHF2_;
1481       ctype = chanCalibHF2_;
1482     }
1483     id = HcalCalibDetId(subdet, ieta, iphi, ctype);
1484 #ifdef EDM_ML_DEBUG
1485     edm::LogVerbatim("HCalGeom") << "CalibrationBox: " << hid << " o/p " << ieta << ":" << iphi << ":" << ichan << ":"
1486                                  << ctype << " " << id;
1487 #endif
1488   } else if (hid < kOffCalibHBX_) {
1489     hid0 -= kOffCalibHOX_;
1490     int iphi, ieta;
1491     if (hid0 < nPhiCalibHOX_[1]) {
1492       iphi = static_cast<int>(hid0) + 1;
1493       ieta = -etaCalibHOX_[1];
1494     } else if (hid0 < (nPhiCalibHOX_[1] + nPhiCalibHOX_[0])) {
1495       hid0 -= nPhiCalibHOX_[1];
1496       iphi = ((hid0 + phiCalibHOX1_) % phiCalibHOX2_ < phiCalibHOX3_) ? 2 * hid0 + 1 : 2 * hid0 + 2;
1497       ieta = -etaCalibHOX_[0];
1498     } else if (hid0 < (nPhiCalibHOX_[1] + 2 * nPhiCalibHOX_[0])) {
1499       hid0 -= (nPhiCalibHOX_[1] + nPhiCalibHOX_[0]);
1500       iphi = ((hid0 + phiCalibHOX1_) % phiCalibHOX2_ < phiCalibHOX3_) ? 2 * hid0 + 1 : 2 * hid0 + 2;
1501       ieta = etaCalibHOX_[0];
1502     } else {
1503       hid0 -= (nPhiCalibHOX_[1] + 2 * nPhiCalibHOX_[0]);
1504       iphi = static_cast<int>(hid0) + 1;
1505       ieta = etaCalibHOX_[1];
1506     }
1507     id = HcalCalibDetId(HcalCalibDetId::HOCrosstalk, ieta, iphi);
1508 #ifdef EDM_ML_DEBUG
1509     edm::LogVerbatim("HCalGeom") << "HOCrossTalk: " << hid << ":" << hid0 << " o/p " << ieta << ":" << iphi << " "
1510                                  << id;
1511 #endif
1512   } else if (hid < kOffCalibHEX_) {
1513     hid0 -= kOffCalibHBX_;
1514     int ieta = (hid0 >= kPhiCalibHBX_) ? etaCalibHBX_ : -etaCalibHBX_;
1515     int iphi = (hid0 % kPhiCalibHBX_) + 1;
1516     id = HcalCalibDetId(HcalCalibDetId::HBX, ieta, iphi);
1517 #ifdef EDM_ML_DEBUG
1518     edm::LogVerbatim("HCalGeom") << "HBX: " << hid << ":" << hid0 << " o/p " << ieta << ":" << iphi << " " << id;
1519 #endif
1520   } else if (hid < kOffCalibHFX_) {
1521     hid0 -= kOffCalibHEX_;
1522     int iphi = 2 * (hid0 % kPhiCalibHEX_) + 1;
1523     int ieta = ((hid0 < kPhiCalibHEX_)
1524                     ? -etaCalibHEX_[0]
1525                     : ((hid0 < 2 * kPhiCalibHEX_) ? etaCalibHEX_[0]
1526                                                   : ((hid0 < 3 * kPhiCalibHEX_) ? -etaCalibHEX_[1] : etaCalibHEX_[1])));
1527     id = HcalCalibDetId(HcalCalibDetId::HEX, ieta, iphi);
1528 #ifdef EDM_ML_DEBUG
1529     edm::LogVerbatim("HCalGeom") << "HEX: " << hid << ":" << hid0 << " o/p " << ieta << ":" << iphi << " " << id;
1530 #endif
1531   }
1532   return id;
1533 }
1534 
1535 unsigned int HcalTopology::detId2denseId(const DetId& id) const {
1536   unsigned int retval(0);
1537   if (topoVersion_ == 0) {  // pre-LS1
1538     retval = detId2denseIdPreLS1(id);
1539   } else if (topoVersion_ == 10) {
1540     HcalDetId hid(id);
1541     if (hid.subdet() == HcalBarrel) {
1542       retval = (hid.depth() - 1) + maxDepthHB_ * (hid.iphi() - 1);
1543       if (hid.ieta() > 0)
1544         retval += maxDepthHB_ * IPHI_MAX * (hid.ieta() - firstHBRing());
1545       else
1546         retval += maxDepthHB_ * IPHI_MAX * (hid.ieta() + lastHBRing() + nEtaHB_);
1547     } else if (hid.subdet() == HcalEndcap) {
1548       retval = HBSize_;
1549       retval += (hid.depth() - 1) + maxDepthHE_ * (hid.iphi() - 1);
1550       if (hid.ieta() > 0)
1551         retval += maxDepthHE_ * maxPhiHE_ * (hid.ieta() - firstHERing());
1552       else
1553         retval += maxDepthHE_ * maxPhiHE_ * (hid.ieta() + lastHERing() + nEtaHE_);
1554     } else if (hid.subdet() == HcalOuter) {
1555       retval = HBSize_ + HESize_;
1556       if (hid.ieta() > 0)
1557         retval += (hid.iphi() - 1) + IPHI_MAX * (hid.ieta() - 1);
1558       else
1559         retval += (hid.iphi() - 1) + IPHI_MAX * (30 + hid.ieta());
1560     } else if (hid.subdet() == HcalForward) {
1561       retval = HBSize_ + HESize_ + HOSize_;
1562       retval += (hid.depth() - 1) + maxDepthHF_ * (hid.iphi() - 1);
1563       if (hid.ieta() > 0)
1564         retval += maxDepthHF_ * IPHI_MAX * (hid.ieta() - 29);
1565       else
1566         retval += maxDepthHF_ * IPHI_MAX * ((41 + 13) + hid.ieta());
1567     } else {
1568       retval = 0xFFFFFFFu;
1569     }
1570   }
1571 #ifdef EDM_ML_DEBUG
1572   edm::LogVerbatim("HcalGeom") << "DetId2Dense " << topoVersion_ << " ID " << std::hex << id.rawId() << std::dec
1573                                << " | " << HcalDetId(id) << " : " << std::hex << retval << std::dec;
1574 #endif
1575   return retval;
1576 }
1577 
1578 DetId HcalTopology::denseId2detId(unsigned int denseid) const {
1579   HcalSubdetector sd(HcalBarrel);
1580   int ie(0);
1581   int ip(0);
1582   int dp(0);
1583   int in(denseid);
1584   int iz(1);
1585   if (topoVersion_ == 0) {  //DL// pre-LS1
1586     if (denseid < kSizeForDenseIndexingPreLS1) {
1587       if (in > 2 * (kHBhalf + kHEhalf + kHOhalf) - 1) {  // HF
1588         sd = HcalForward;
1589         in -= 2 * (kHBhalf + kHEhalf + kHOhalf);
1590         iz = (in < kHFhalf ? 1 : -1);
1591         in %= kHFhalf;
1592         ip = 4 * (in / 48);
1593         in %= 48;
1594         ip += 1 + (in > 21 ? 2 : 0);
1595         if (3 == ip % 4)
1596           in -= 22;
1597         ie = 29 + in / 2;
1598         dp = 1 + in % 2;
1599       } else if (in > 2 * (kHBhalf + kHEhalf) - 1) {  // HO
1600         sd = HcalOuter;
1601         in -= 2 * (kHBhalf + kHEhalf);
1602         iz = (in < kHOhalf ? 1 : -1);
1603         in %= kHOhalf;
1604         dp = 4;
1605         ip = 1 + in / 15;
1606         ie = 1 + (in - 15 * (ip - 1));
1607       } else if (in > 2 * kHBhalf - 1) {  // Endcap
1608         sd = HcalEndcap;
1609         in -= 2 * kHBhalf;
1610         iz = (in < kHEhalf ? 1 : -1);
1611         in %= kHEhalf;
1612         ip = 2 * (in / 36);
1613         in %= 36;
1614         ip += 1 + in / 28;
1615         if (0 == ip % 2)
1616           in %= 28;
1617         ie = 15 + (in < 2 ? 1 + in : 2 + (in < 20 ? 1 + (in - 2) / 2 : 9 + (in < 26 ? 1 + (in - 20) / 3 : 3)));
1618         dp = (in < 1
1619                   ? 3
1620                   : (in < 2 ? 1 : (in < 20 ? 1 + (in - 2) % 2 : (in < 26 ? 1 + (in - 20) % 3 : (1 + (in - 26) % 2)))));
1621       } else {  // barrel
1622         iz = (in < kHBhalf ? 1 : -1);
1623         in %= kHBhalf;
1624         ip = in / 18 + 1;
1625         in %= 18;
1626         if (in < 14) {
1627           dp = 1;
1628           ie = in + 1;
1629         } else {
1630           in %= 14;
1631           dp = 1 + in % 2;
1632           ie = 15 + in / 2;
1633         }
1634       }
1635     }
1636   } else if (topoVersion_ == 10) {
1637     if (denseid < ncells()) {
1638       if (denseid >= (HBSize_ + HESize_ + HOSize_)) {
1639         sd = HcalForward;
1640         in -= (HBSize_ + HESize_ + HOSize_);
1641         dp = (in % maxDepthHF_) + 1;
1642         ip = (in - dp + 1) % (maxDepthHF_ * IPHI_MAX);
1643         ip = (ip / maxDepthHF_) + 1;
1644         ie = (in - dp + 1 - maxDepthHF_ * (ip - 1)) / (IPHI_MAX * maxDepthHF_);
1645         if (ie > 12) {
1646           ie = 54 - ie;
1647           iz = -1;
1648         } else {
1649           ie += 29;
1650           iz = 1;
1651         }
1652       } else if (denseid >= (HBSize_ + HESize_)) {
1653         sd = HcalOuter;
1654         in -= (HBSize_ + HESize_);
1655         dp = 4;
1656         ip = (in % IPHI_MAX) + 1;
1657         ie = (in - ip + 1) / IPHI_MAX;
1658         if (ie > 14) {
1659           ie = 30 - ie;
1660           iz = -1;
1661         } else {
1662           ie += 1;
1663           iz = 1;
1664         }
1665       } else if (denseid >= (HBSize_)) {
1666         sd = HcalEndcap;
1667         in -= (HBSize_);
1668         dp = (in % maxDepthHE_) + 1;
1669         ip = (in - dp + 1) % (maxDepthHE_ * maxPhiHE_);
1670         ip = (ip / maxDepthHE_) + 1;
1671         ie = (in - dp + 1 - maxDepthHE_ * (ip - 1)) / (maxPhiHE_ * maxDepthHE_);
1672         if (ie >= nEtaHE_) {
1673           ie = lastHERing() + nEtaHE_ - ie;
1674           iz = -1;
1675         } else {
1676           ie = firstHERing() + ie;
1677           iz = 1;
1678         }
1679       } else {
1680         sd = HcalBarrel;
1681         dp = (in % maxDepthHB_) + 1;
1682         ip = (in - dp + 1) % (maxDepthHB_ * IPHI_MAX);
1683         ip = (ip / maxDepthHB_) + 1;
1684         ie = (in - dp + 1 - maxDepthHB_ * (ip - 1)) / (IPHI_MAX * maxDepthHB_);
1685         if (ie >= nEtaHB_) {
1686           ie = lastHBRing() + nEtaHB_ - ie;
1687           iz = -1;
1688         } else {
1689           ie = firstHBRing() + ie;
1690           iz = 1;
1691         }
1692       }
1693     }
1694   }
1695   HcalDetId hid(sd, iz * int(ie), ip, dp);
1696 #ifdef EDM_ML_DEBUG
1697   edm::LogVerbatim("HcalGeom") << "Dens2Det " << topoVersion_ << " i/p " << std::hex << denseid << " : " << hid.rawId()
1698                                << std::dec << " | " << hid;
1699 #endif
1700   return hid;
1701 }
1702 
1703 unsigned int HcalTopology::ncells() const { return HBSize_ + HESize_ + HOSize_ + HFSize_; }
1704 
1705 int HcalTopology::topoVersion() const { return topoVersion_; }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!