Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:14:19 

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