Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-07-04 01:51:01

0001 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0002 
0003 #include "DataFormats/Math/interface/GeantUnits.h"
0004 #include "DataFormats/ForwardDetId/interface/HFNoseDetId.h"
0005 #include "DataFormats/ForwardDetId/interface/HGCScintillatorDetId.h"
0006 #include "DataFormats/ForwardDetId/interface/HGCSiliconDetId.h"
0007 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Utilities/interface/Exception.h"
0010 #include "Geometry/HGCalCommonData/interface/HGCalGeomParameters.h"
0011 #include "Geometry/HGCalCommonData/interface/HGCalGeometryMode.h"
0012 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0013 #include "Geometry/HGCalCommonData/interface/HGCalWaferMask.h"
0014 #include "Geometry/HGCalCommonData/interface/HGCalWaferType.h"
0015 
0016 #include <algorithm>
0017 #include <bitset>
0018 #include <iterator>
0019 #include <functional>
0020 #include <numeric>
0021 
0022 //#define EDM_ML_DEBUG
0023 using namespace geant_units::operators;
0024 
0025 HGCalDDDConstants::HGCalDDDConstants(const HGCalParameters* hp, const std::string& name)
0026     : hgpar_(hp), sqrt3_(std::sqrt(3.0)), mode_(hgpar_->mode_), fullAndPart_(waferHexagon8File()) {
0027 #ifdef EDM_ML_DEBUG
0028   edm::LogVerbatim("HGCalGeom") << "Mode " << mode_ << " FullAndPart " << fullAndPart_;
0029 #endif
0030   if (waferHexagon6() || waferHexagon8()) {
0031     rmax_ = (HGCalParameters::k_ScaleFromDDD * (hgpar_->waferR_) * std::cos(30._deg));
0032     rmaxT_ = rmax_ + 0.5 * hgpar_->sensorSeparation_;
0033     hexside_ = 2.0 * rmax_ * tan30deg_;
0034     hexsideT_ = 2.0 * rmaxT_ * tan30deg_;
0035     hgcell_ = std::make_unique<HGCalCell>(2.0 * rmaxT_, hgpar_->nCellsFine_, hgpar_->nCellsCoarse_);
0036     hgcellUV_ = std::make_unique<HGCalCellUV>(
0037         2.0 * rmax_, hgpar_->sensorSeparation_, hgpar_->nCellsFine_, hgpar_->nCellsCoarse_);
0038 #ifdef EDM_ML_DEBUG
0039     edm::LogVerbatim("HGCalGeom") << "rmax_ " << rmax_ << ":" << rmaxT_ << ":" << hexside_ << ":" << hexsideT_
0040                                   << " CellSize " << 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[0] << ":"
0041                                   << 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[1];
0042 #endif
0043   }
0044   if ((mode_ == HGCalGeometryMode::Hexagon8Cassette) || (mode_ == HGCalGeometryMode::TrapezoidCassette)) {
0045     hgcassette_.setParameter(hgpar_->cassettes_, hgpar_->cassetteShift_);
0046 #ifdef EDM_ML_DEBUG
0047     edm::LogVerbatim("HGCalGeom") << "Setup HGCalCassette for " << hgpar_->cassettes_ << " cassettes";
0048 #endif
0049   }
0050   // init maps and constants
0051   modHalf_ = 0;
0052   maxWafersPerLayer_ = 0;
0053   for (int simreco = 0; simreco < 2; ++simreco) {
0054     tot_layers_[simreco] = layersInit((bool)simreco);
0055     max_modules_layer_[simreco].resize(tot_layers_[simreco] + 1);
0056     for (unsigned int layer = 1; layer <= tot_layers_[simreco]; ++layer) {
0057       max_modules_layer_[simreco][layer] = modulesInit(layer, (bool)simreco);
0058       if (simreco == 1) {
0059         modHalf_ += max_modules_layer_[simreco][layer];
0060         maxWafersPerLayer_ = std::max(maxWafersPerLayer_, max_modules_layer_[simreco][layer]);
0061 #ifdef EDM_ML_DEBUG
0062         edm::LogVerbatim("HGCalGeom") << "Layer " << layer << " with " << max_modules_layer_[simreco][layer] << ":"
0063                                       << modHalf_ << " modules in RECO";
0064       } else {
0065         edm::LogVerbatim("HGCalGeom") << "Layer " << layer << " with " << max_modules_layer_[simreco][layer]
0066                                       << " modules in SIM";
0067 #endif
0068       }
0069     }
0070 #ifdef EDM_ML_DEBUG
0071     edm::LogVerbatim("HGCalGeom") << "SimReco " << simreco << " with " << tot_layers_[simreco] << " Layers";
0072 #endif
0073   }
0074   tot_wafers_ = wafers();
0075 
0076 #ifdef EDM_ML_DEBUG
0077   edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants initialized for " << name << " with " << layers(false) << ":"
0078                                 << layers(true) << " layers, " << wafers() << ":" << 2 * modHalf_
0079                                 << " wafers with maximum " << maxWafersPerLayer_ << " per layer and "
0080                                 << "maximum of " << maxCells(false) << ":" << maxCells(true) << " cells";
0081 #endif
0082   if (waferHexagon6() || waferHexagon8()) {
0083     int wminT(9999999), wmaxT(-9999999), kount1(0), kount2(0);
0084     for (unsigned int i = 0; i < getTrFormN(); ++i) {
0085       int lay0 = getTrForm(i).lay;
0086       int wmin(9999999), wmax(-9999999), kount(0);
0087       for (int wafer = 0; wafer < sectors(); ++wafer) {
0088         bool waferIn = waferInLayer(wafer, lay0, true);
0089         if (waferHexagon8()) {
0090           int kndx = HGCalWaferIndex::waferIndex(lay0,
0091                                                  HGCalWaferIndex::waferU(hgpar_->waferCopy_[wafer]),
0092                                                  HGCalWaferIndex::waferV(hgpar_->waferCopy_[wafer]));
0093           waferIn_[kndx] = waferIn;
0094         }
0095         if (waferIn) {
0096           int waferU = ((waferHexagon6()) ? wafer : HGCalWaferIndex::waferU(hgpar_->waferCopy_[wafer]));
0097           if (waferU < wmin)
0098             wmin = waferU;
0099           if (waferU > wmax)
0100             wmax = waferU;
0101           ++kount;
0102         }
0103       }
0104       if (wminT > wmin)
0105         wminT = wmin;
0106       if (wmaxT < wmax)
0107         wmaxT = wmax;
0108       if (kount1 < kount)
0109         kount1 = kount;
0110       kount2 += kount;
0111 #ifdef EDM_ML_DEBUG
0112       int lay1 = getIndex(lay0, true).first;
0113       edm::LogVerbatim("HGCalGeom") << "Index " << i << " Layer " << lay0 << ":" << lay1 << " Wafer " << wmin << ":"
0114                                     << wmax << ":" << kount;
0115 #endif
0116       HGCWaferParam a1{{wmin, wmax, kount}};
0117       waferLayer_[lay0] = a1;
0118     }
0119     waferMax_ = std::array<int, 4>{{wminT, wmaxT, kount1, kount2}};
0120 #ifdef EDM_ML_DEBUG
0121     edm::LogVerbatim("HGCalGeom") << "Overall wafer statistics: " << wminT << ":" << wmaxT << ":" << kount1 << ":"
0122                                   << kount2;
0123 #endif
0124   }
0125 }
0126 
0127 std::pair<int, int> HGCalDDDConstants::assignCell(float x, float y, int lay, int subSec, bool reco) const {
0128   const auto& index = getIndex(lay, reco);
0129   if (index.first < 0)
0130     return std::make_pair(-1, -1);
0131   if (waferHexagon6()) {
0132     float xx = (reco) ? x : HGCalParameters::k_ScaleFromDDD * x;
0133     float yy = (reco) ? y : HGCalParameters::k_ScaleFromDDD * y;
0134 
0135     // First the wafer
0136     int wafer = cellHex(xx, yy, rmax_, hgpar_->waferPosX_, hgpar_->waferPosY_);
0137     if (wafer < 0 || wafer >= static_cast<int>(hgpar_->waferTypeT_.size())) {
0138       edm::LogWarning("HGCalGeom") << "Wafer no. out of bound for " << wafer << ":" << (hgpar_->waferTypeT_).size()
0139                                    << ":" << (hgpar_->waferPosX_).size() << ":" << (hgpar_->waferPosY_).size()
0140                                    << " ***** ERROR *****";
0141       return std::make_pair(-1, -1);
0142     } else {
0143       // Now the cell
0144       xx -= hgpar_->waferPosX_[wafer];
0145       yy -= hgpar_->waferPosY_[wafer];
0146       if (hgpar_->waferTypeT_[wafer] == 1)
0147         return std::make_pair(wafer,
0148                               cellHex(xx,
0149                                       yy,
0150                                       0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[0],
0151                                       hgpar_->cellFineX_,
0152                                       hgpar_->cellFineY_));
0153       else
0154         return std::make_pair(wafer,
0155                               cellHex(xx,
0156                                       yy,
0157                                       0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[1],
0158                                       hgpar_->cellCoarseX_,
0159                                       hgpar_->cellCoarseY_));
0160     }
0161   } else {
0162     return std::make_pair(-1, -1);
0163   }
0164 }
0165 
0166 std::array<int, 5> HGCalDDDConstants::assignCellHex(
0167     float x, float y, int lay, bool reco, bool extend, bool debug) const {
0168   int waferU(0), waferV(0), waferType(-1), cellU(0), cellV(0);
0169   if (waferHexagon8()) {
0170     double xx = (reco) ? HGCalParameters::k_ScaleToDDD * x : x;
0171     double yy = (reco) ? HGCalParameters::k_ScaleToDDD * y : y;
0172     double wt(1.0);
0173 #ifdef EDM_ML_DEBUG
0174     edm::LogVerbatim("HGCalGeom") << "assignCellHex x " << x << ":" << xx << " y " << y << ":" << yy << " Lay " << lay;
0175 #endif
0176     waferFromPosition(xx, yy, lay, waferU, waferV, cellU, cellV, waferType, wt, extend, debug);
0177   }
0178   return std::array<int, 5>{{waferU, waferV, waferType, cellU, cellV}};
0179 }
0180 
0181 std::array<int, 3> HGCalDDDConstants::assignCellTrap(float x, float y, float z, int layer, bool reco) const {
0182   int irad(-1), iphi(-1), type(-1);
0183   const auto& indx = getIndex(layer, reco);
0184   if (indx.first < 0)
0185     return std::array<int, 3>{{irad, iphi, type}};
0186   double xx = (z > 0) ? x : -x;
0187   double phi = (((y == 0.0) && (x == 0.0)) ? 0. : std::atan2(y, xx));
0188   if (phi < 0)
0189     phi += (2.0 * M_PI);
0190   if (indx.second != 0)
0191     iphi = 1 + static_cast<int>(phi / indx.second);
0192   if (mode_ == HGCalGeometryMode::TrapezoidCassette) {
0193     int cassette = HGCalTileIndex::tileCassette(iphi, hgpar_->phiOffset_, hgpar_->nphiCassette_, hgpar_->cassettes_);
0194     auto cshift = hgcassette_.getShift(layer, 1, cassette);
0195 #ifdef EDM_ML_DEBUG
0196     edm::LogVerbatim("HGCalGeom") << "Cassette " << cassette << " Shift " << cshift.first << ":" << cshift.second;
0197 #endif
0198     if (reco) {
0199       x -= cshift.first;
0200       y -= cshift.second;
0201     } else {
0202       x = HGCalParameters::k_ScaleFromDDD * x - cshift.first;
0203       y = HGCalParameters::k_ScaleFromDDD * y - cshift.second;
0204     }
0205   }
0206   type = hgpar_->scintType(layer);
0207   double r = std::sqrt(x * x + y * y);
0208   auto ir = std::lower_bound(hgpar_->radiusLayer_[type].begin(), hgpar_->radiusLayer_[type].end(), r);
0209   irad = static_cast<int>(ir - hgpar_->radiusLayer_[type].begin());
0210   irad = std::clamp(irad, hgpar_->iradMinBH_[indx.first], hgpar_->iradMaxBH_[indx.first]);
0211 #ifdef EDM_ML_DEBUG
0212   edm::LogVerbatim("HGCalGeom") << "assignCellTrap Input " << x << ":" << y << ":" << z << ":" << layer << ":" << reco
0213                                 << " x|r " << xx << ":" << r << " phi " << phi << " o/p " << irad << ":" << iphi << ":"
0214                                 << type;
0215 #endif
0216   return std::array<int, 3>{{irad, iphi, type}};
0217 }
0218 
0219 std::pair<double, double> HGCalDDDConstants::cellEtaPhiTrap(int type, int irad) const {
0220   double dr(0), df(0);
0221   if (tileTrapezoid()) {
0222     double r = 0.5 * ((hgpar_->radiusLayer_[type][irad - 1] + hgpar_->radiusLayer_[type][irad]));
0223     dr = (hgpar_->radiusLayer_[type][irad] - hgpar_->radiusLayer_[type][irad - 1]);
0224     df = r * hgpar_->cellSize_[type];
0225   }
0226   return std::make_pair(dr, df);
0227 }
0228 
0229 bool HGCalDDDConstants::cellInLayer(int waferU, int waferV, int cellU, int cellV, int lay, bool reco) const {
0230   const auto& indx = getIndex(lay, true);
0231   if (indx.first >= 0) {
0232     if (mode_ == HGCalGeometryMode::Hexagon8Cassette) {
0233       int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0234       auto ktr = hgpar_->waferInfoMap_.find(indx);
0235       int part = (ktr != hgpar_->waferInfoMap_.end()) ? (ktr->second).part : HGCalTypes::WaferFull;
0236       return HGCalWaferMask::goodCell(cellU, cellV, part);
0237     } else if (mode_ == HGCalGeometryMode::Hexagon8Module) {
0238       int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0239       auto ktr = hgpar_->waferInfoMap_.find(indx);
0240       int thck(HGCalTypes::WaferFineThin), part(HGCalTypes::WaferFull), rotn(HGCalTypes::WaferOrient0);
0241       if (ktr != hgpar_->waferInfoMap_.end()) {
0242         thck = (ktr->second).type;
0243         part = (ktr->second).part;
0244         rotn = (ktr->second).orient;
0245       }
0246       int ncell = (thck == HGCalTypes::WaferFineThin) ? hgpar_->nCellsFine_ : hgpar_->nCellsCoarse_;
0247       return HGCalWaferMask::goodCell(cellU, cellV, ncell, part, rotn);
0248     } else if (waferHexagon8() || waferHexagon6()) {
0249       const auto& xy = ((waferHexagon8()) ? locateCell(lay, waferU, waferV, cellU, cellV, reco, true, false, false)
0250                                           : locateCell(cellU, lay, waferU, reco));
0251       double rpos = sqrt(xy.first * xy.first + xy.second * xy.second);
0252       return ((rpos >= hgpar_->rMinLayHex_[indx.first]) && (rpos <= hgpar_->rMaxLayHex_[indx.first]));
0253     } else {
0254       return true;
0255     }
0256   } else {
0257     return false;
0258   }
0259 }
0260 
0261 double HGCalDDDConstants::cellThickness(int layer, int waferU, int waferV) const {
0262   double thick(-1);
0263   int type = waferType(layer, waferU, waferV, false);
0264   if (type >= 0) {
0265     if (waferHexagon8()) {
0266       thick = 10000.0 * hgpar_->cellThickness_[type];  // cm to micron
0267     } else if (waferHexagon6()) {
0268       thick = 100.0 * (type + 1);  // type = 1,2,3 for 100,200,300 micron
0269     }
0270   }
0271   return thick;
0272 }
0273 
0274 double HGCalDDDConstants::cellSizeHex(int type) const {
0275   int indx = ((waferHexagon8()) ? ((type >= 1) ? 1 : 0) : ((type == 1) ? 1 : 0));
0276   double cell = (tileTrapezoid() ? 0.5 * hgpar_->cellSize_[indx]
0277                                  : 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[indx]);
0278   return cell;
0279 }
0280 
0281 int32_t HGCalDDDConstants::cellType(int type, int cellU, int cellV, int iz, int fwdBack, int orient) const {
0282   int placement = (orient < 0) ? HGCalCell::cellPlacementOld : HGCalCell::cellPlacementIndex(iz, fwdBack, orient);
0283   int ncell = (type == 0) ? hgpar_->nCellsFine_ : hgpar_->nCellsCoarse_;
0284   auto cellType = HGCalCell::cellType(cellU, cellV, ncell, placement);
0285   return cellType.first;
0286 }
0287 
0288 double HGCalDDDConstants::distFromEdgeHex(double x, double y, double z) const {
0289   // Assming the point is within a hexagonal plane of the wafer, calculate
0290   // the shortest distance from the edge
0291   if (z < 0)
0292     x = -x;
0293   double dist(0);
0294   // Input x, y in Geant4 unit and transformed to CMSSW standard
0295   double xx = HGCalParameters::k_ScaleFromDDD * x;
0296   double yy = HGCalParameters::k_ScaleFromDDD * y;
0297   if (waferHexagon8()) {
0298     int ll = layerIndex(getLayer(z, false), false);
0299     xx -= hgpar_->xLayerHex_[ll];
0300     yy -= hgpar_->yLayerHex_[ll];
0301   }
0302   int sizew = static_cast<int>(hgpar_->waferPosX_.size());
0303   int wafer = sizew;
0304   // Transform to the local coordinate frame of the wafer first
0305   for (int k = 0; k < sizew; ++k) {
0306     double dx = std::abs(xx - hgpar_->waferPosX_[k]);
0307     double dy = std::abs(yy - hgpar_->waferPosY_[k]);
0308     if ((dx <= rmax_) && (dy <= hexside_) && ((dy <= 0.5 * hexside_) || (dx * tan30deg_ <= (hexside_ - dy)))) {
0309       wafer = k;
0310       xx -= hgpar_->waferPosX_[k];
0311       yy -= hgpar_->waferPosY_[k];
0312       break;
0313     }
0314   }
0315   // Look at only one quarter (both x,y are positive)
0316   if (wafer < sizew) {
0317     if (std::abs(yy) < 0.5 * hexside_) {
0318       dist = rmax_ - std::abs(xx);
0319     } else {
0320       dist = 0.5 * ((rmax_ - std::abs(xx)) - sqrt3_ * (std::abs(yy) - 0.5 * hexside_));
0321     }
0322   } else {
0323     dist = 0;
0324   }
0325   dist *= HGCalParameters::k_ScaleToDDD;
0326 #ifdef EDM_ML_DEBUG
0327   edm::LogVerbatim("HGCalGeom") << "DistFromEdgeHex: Local " << xx << ":" << yy << " wafer " << wafer << " flag "
0328                                 << (wafer < sizew) << " Distance " << rmax_ << ":" << (rmax_ - std::abs(xx)) << ":"
0329                                 << (std::abs(yy) - 0.5 * hexside_) << ":" << 0.5 * hexside_ << ":" << dist;
0330 #endif
0331   return dist;
0332 }
0333 
0334 double HGCalDDDConstants::distFromEdgeTrap(double x, double y, double z) const {
0335   // Assming the point is within the eta-phi plane of the scintillator tile,
0336   // calculate the shortest distance from the edge
0337   int lay = getLayer(z, false);
0338   double xx = (z < 0) ? -x : x;
0339   int indx = layerIndex(lay, false);
0340   double r = HGCalParameters::k_ScaleFromDDD * std::sqrt(x * x + y * y);
0341   double phi = (r == 0. ? 0. : std::atan2(y, xx));
0342   if (phi < 0)
0343     phi += (2.0 * M_PI);
0344   int type = hgpar_->scintType(lay);
0345   double cell = hgpar_->scintCellSize(lay);
0346   // Compare with the center of the tile find distances along R and also phi
0347   // Take the smaller value
0348   auto ir = std::lower_bound(hgpar_->radiusLayer_[type].begin(), hgpar_->radiusLayer_[type].end(), r);
0349   int irad = static_cast<int>(ir - hgpar_->radiusLayer_[type].begin());
0350   irad = std::clamp(irad, hgpar_->iradMinBH_[indx], hgpar_->iradMaxBH_[indx]);
0351   int iphi = 1 + static_cast<int>(phi / cell);
0352   double dphi = std::max(0.0, (0.5 * cell - std::abs(phi - (iphi - 0.5) * cell)));
0353   double dist = std::min((r - hgpar_->radiusLayer_[type][irad - 1]), (hgpar_->radiusLayer_[type][irad] - r));
0354 #ifdef EDM_ML_DEBUG
0355   edm::LogVerbatim("HGCalGeom") << "DistFromEdgeTrap: Global " << x << ":" << y << ":" << z << " Layer " << lay
0356                                 << " Index " << indx << ":" << type << " xx " << xx << " R " << r << ":" << irad << ":"
0357                                 << hgpar_->radiusLayer_[type][irad - 1] << ":" << hgpar_->radiusLayer_[type][irad]
0358                                 << " Phi " << phi << ":" << iphi << ":" << (iphi - 0.5) * cell << " cell " << cell
0359                                 << " Dphi " << dphi << " Dist " << dist << ":" << r * dphi;
0360 #endif
0361   return HGCalParameters::k_ScaleToDDD * std::min(r * dphi, dist);
0362 }
0363 
0364 int HGCalDDDConstants::getLayer(double z, bool reco) const {
0365   // Get the layer # from the gloabl z coordinate
0366   unsigned int k = 0;
0367   double zz = (reco ? std::abs(z) : HGCalParameters::k_ScaleFromDDD * std::abs(z));
0368   const auto& zLayerHex = hgpar_->zLayerHex_;
0369   auto itr = std::find_if(zLayerHex.begin() + 1, zLayerHex.end(), [&k, &zz, &zLayerHex](double zLayer) {
0370     ++k;
0371     return zz < 0.5 * (zLayerHex[k - 1] + zLayerHex[k]);
0372   });
0373   int lay = (itr == zLayerHex.end()) ? static_cast<int>(zLayerHex.size()) : k;
0374   if (waferHexagon6() && reco) {
0375     int indx = layerIndex(lay, false);
0376     if (indx >= 0)
0377       lay = hgpar_->layerGroupO_[indx];
0378   } else {
0379     lay += (hgpar_->firstLayer_ - 1);
0380   }
0381   return lay;
0382 }
0383 
0384 HGCalParameters::hgtrap HGCalDDDConstants::getModule(unsigned int indx, bool hexType, bool reco) const {
0385   HGCalParameters::hgtrap mytr;
0386   if (hexType) {
0387     if (indx >= hgpar_->waferTypeL_.size())
0388       edm::LogWarning("HGCalGeom") << "Wafer no. out bound for index " << indx << ":" << (hgpar_->waferTypeL_).size()
0389                                    << ":" << (hgpar_->waferPosX_).size() << ":" << (hgpar_->waferPosY_).size()
0390                                    << " ***** ERROR *****";
0391     unsigned int type =
0392         ((indx < hgpar_->waferTypeL_.size()) ? hgpar_->waferTypeL_[indx] - 1 : HGCSiliconDetId::HGCalCoarseThick);
0393     mytr = hgpar_->getModule(type, reco);
0394   } else {
0395     mytr = hgpar_->getModule(indx, reco);
0396   }
0397   return mytr;
0398 }
0399 
0400 std::vector<HGCalParameters::hgtrap> HGCalDDDConstants::getModules() const {
0401   std::vector<HGCalParameters::hgtrap> mytrs;
0402   for (unsigned int k = 0; k < hgpar_->moduleLayR_.size(); ++k)
0403     mytrs.emplace_back(hgpar_->getModule(k, true));
0404   return mytrs;
0405 }
0406 
0407 int HGCalDDDConstants::getPhiBins(int lay) const { return (tileTrapezoid() ? hgpar_->scintCells(lay) : 0); }
0408 
0409 std::pair<int, int> HGCalDDDConstants::getREtaRange(int lay) const {
0410   int irmin(0), irmax(0);
0411   if (tileTrapezoid()) {
0412     int indx = layerIndex(lay, false);
0413     if ((indx >= 0) && (indx < static_cast<int>(hgpar_->iradMinBH_.size()))) {
0414       irmin = hgpar_->iradMinBH_[indx];
0415       irmax = hgpar_->iradMaxBH_[indx];
0416     }
0417   }
0418   return std::make_pair(irmin, irmax);
0419 }
0420 
0421 std::vector<HGCalParameters::hgtrform> HGCalDDDConstants::getTrForms() const {
0422   std::vector<HGCalParameters::hgtrform> mytrs;
0423   for (unsigned int k = 0; k < hgpar_->trformIndex_.size(); ++k)
0424     mytrs.emplace_back(hgpar_->getTrForm(k));
0425   return mytrs;
0426 }
0427 
0428 int HGCalDDDConstants::getTypeTrap(int layer) const {
0429   // Get the module type for scinitllator
0430   if (tileTrapezoid()) {
0431     return hgpar_->scintType(layer);
0432   } else {
0433     return -1;
0434   }
0435 }
0436 
0437 int HGCalDDDConstants::getTypeHex(int layer, int waferU, int waferV) const {
0438   // Get the module type for a silicon wafer
0439   if (waferHexagon8()) {
0440     auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
0441     return ((itr == hgpar_->typesInLayers_.end() ? 2 : hgpar_->waferTypeL_[itr->second]));
0442   } else {
0443     return -1;
0444   }
0445 }
0446 
0447 std::pair<double, double> HGCalDDDConstants::getXY(int layer, double x, double y, bool forwd) const {
0448   int ll = layer - hgpar_->firstLayer_;
0449   double x0(x), y0(y);
0450   if ((!hgpar_->layerType_.empty()) && (ll < static_cast<int>(hgpar_->layerRotV_.size()))) {
0451     if (forwd) {
0452       x0 = x * hgpar_->layerRotV_[ll].first + y * hgpar_->layerRotV_[ll].second;
0453       y0 = y * hgpar_->layerRotV_[ll].first - x * hgpar_->layerRotV_[ll].second;
0454     } else {
0455       x0 = x * hgpar_->layerRotV_[ll].first - y * hgpar_->layerRotV_[ll].second;
0456       y0 = y * hgpar_->layerRotV_[ll].first + x * hgpar_->layerRotV_[ll].second;
0457     }
0458   }
0459 #ifdef EDM_ML_DEBUG
0460   double x1(x0), y1(y0);
0461   if (ll < static_cast<int>(hgpar_->layerRotV_.size())) {
0462     if (forwd) {
0463       x1 = x0 * hgpar_->layerRotV_[ll].first - y0 * hgpar_->layerRotV_[ll].second;
0464       y1 = y0 * hgpar_->layerRotV_[ll].first + x0 * hgpar_->layerRotV_[ll].second;
0465     } else {
0466       x1 = x0 * hgpar_->layerRotV_[ll].first + y0 * hgpar_->layerRotV_[ll].second;
0467       y1 = y0 * hgpar_->layerRotV_[ll].first - x0 * hgpar_->layerRotV_[ll].second;
0468     }
0469   }
0470   edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants: Layer " << layer << ":" << ll << " mode " << forwd << " x " << x
0471                                 << ":" << x0 << ":" << x1 << " y " << y << ":" << y0 << ":" << y1;
0472 #endif
0473   return std::make_pair(x0, y0);
0474 }
0475 
0476 bool HGCalDDDConstants::isHalfCell(int waferType, int cell) const {
0477   if (waferType < 1 || cell < 0)
0478     return false;
0479   return waferType == 2 ? hgpar_->cellCoarseHalf_[cell] : hgpar_->cellFineHalf_[cell];
0480 }
0481 
0482 bool HGCalDDDConstants::isValidHex(int lay, int mod, int cell, bool reco) const {
0483   // Check validity for a layer|wafer|cell of pre-TDR version
0484   bool result(false), resultMod(false);
0485   int cellmax(0);
0486   if (waferHexagon6()) {
0487     int32_t copyNumber = hgpar_->waferCopy_[mod];
0488     result = ((lay > 0 && lay <= static_cast<int>(layers(reco))));
0489     if (result) {
0490       const int32_t lay_idx = reco ? (lay - 1) * 3 + 1 : lay;
0491       const auto& the_modules = hgpar_->copiesInLayers_[lay_idx];
0492       auto moditr = the_modules.find(copyNumber);
0493       result = resultMod = (moditr != the_modules.end());
0494 #ifdef EDM_ML_DEBUG
0495       if (!result)
0496         edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants: Layer " << lay << ":" << lay_idx << " Copy " << copyNumber
0497                                       << ":" << mod << " Flag " << result;
0498 #endif
0499       if (result) {
0500         if (moditr->second >= 0) {
0501           if (mod >= static_cast<int>(hgpar_->waferTypeT_.size()))
0502             edm::LogWarning("HGCalGeom") << "Module no. out of bound for " << mod << " to be compared with "
0503                                          << (hgpar_->waferTypeT_).size() << " ***** ERROR *****";
0504           cellmax = ((hgpar_->waferTypeT_[mod] - 1 == HGCSiliconDetId::HGCalFine)
0505                          ? static_cast<int>(hgpar_->cellFineX_.size())
0506                          : static_cast<int>(hgpar_->cellCoarseX_.size()));
0507           result = (cell >= 0 && cell <= cellmax);
0508         } else {
0509           result = isValidCell(lay_idx, mod, cell);
0510         }
0511       }
0512     }
0513   }
0514 
0515 #ifdef EDM_ML_DEBUG
0516   if (!result)
0517     edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants: Layer " << lay << ":"
0518                                   << (lay > 0 && (lay <= static_cast<int>(layers(reco)))) << " Module " << mod << ":"
0519                                   << resultMod << " Cell " << cell << ":" << cellmax << ":"
0520                                   << (cell >= 0 && cell <= cellmax) << ":" << maxCells(reco);
0521 #endif
0522   return result;
0523 }
0524 
0525 bool HGCalDDDConstants::isValidHex8(int layer, int modU, int modV, bool fullAndPart) const {
0526   // Check validity for a layer|wafer|cell of post-TDR version
0527   int indx = HGCalWaferIndex::waferIndex(layer, modU, modV);
0528   auto itr = hgpar_->typesInLayers_.find(indx);
0529 #ifdef EDM_ML_DEBUG
0530   edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:WaferType " << layer << ":" << modU << ":" << modV
0531                                 << ":" << indx << " Test " << (itr != hgpar_->typesInLayers_.end());
0532 #endif
0533   if (itr == hgpar_->typesInLayers_.end())
0534     return false;
0535 
0536   if (fullAndPart_) {
0537     auto ktr = hgpar_->waferInfoMap_.find(indx);
0538 #ifdef EDM_ML_DEBUG
0539     edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:WaferInfoMap " << layer << ":" << modU << ":"
0540                                   << modV << ":" << indx << " Test " << (ktr != hgpar_->waferInfoMap_.end());
0541 #endif
0542     if (ktr == hgpar_->waferInfoMap_.end())
0543       return false;
0544   } else {
0545     auto jtr = waferIn_.find(indx);
0546 #ifdef EDM_ML_DEBUG
0547     edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:WaferIn " << jtr->first << ":" << jtr->second;
0548 #endif
0549     if (!(jtr->second))
0550       return false;
0551   }
0552 
0553   if (fullAndPart || fullAndPart_) {
0554     auto ktr = hgpar_->waferTypes_.find(indx);
0555     if (ktr != hgpar_->waferTypes_.end()) {
0556       if (hgpar_->waferMaskMode_ > 0) {
0557         if (ktr->second.first == HGCalTypes::WaferOut)
0558           return false;
0559       } else {
0560         if (ktr->second.first < HGCalTypes::WaferCornerMin)
0561           return false;
0562       }
0563     }
0564   }
0565   return true;
0566 }
0567 
0568 bool HGCalDDDConstants::isValidHex8(int layer, int modU, int modV, int cellU, int cellV, bool fullAndPart) const {
0569   // First check validity for a layer|wafer| of post TDR version
0570   if (!isValidHex8(layer, modU, modV, fullAndPart))
0571     return false;
0572   int indx = HGCalWaferIndex::waferIndex(layer, modU, modV);
0573   auto itr = hgpar_->typesInLayers_.find(indx);
0574   int type = hgpar_->waferTypeL_[itr->second];
0575   int N = ((hgpar_->waferTypeL_[itr->second] == 0) ? hgpar_->nCellsFine_ : hgpar_->nCellsCoarse_);
0576 #ifdef EDM_ML_DEBUG
0577   edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:Cell " << cellU << ":" << cellV << ":" << N
0578                                 << " Tests " << (cellU >= 0) << ":" << (cellU < 2 * N) << ":" << (cellV >= 0) << ":"
0579                                 << (cellV < 2 * N) << ":" << ((cellV - cellU) < N) << ":" << ((cellU - cellV) <= N);
0580 #endif
0581   if ((cellU < 0) || (cellU >= 2 * N) || (cellV < 0) || (cellV >= 2 * N))
0582     return false;
0583   if (((cellV - cellU) >= N) || ((cellU - cellV) > N))
0584     return false;
0585 
0586   return isValidCell8(layer, modU, modV, cellU, cellV, type);
0587 }
0588 
0589 bool HGCalDDDConstants::isValidTrap(int layer, int irad, int iphi) const {
0590   // Check validity for a layer|eta|phi of scintillator
0591   const auto& indx = getIndex(layer, true);
0592   if (indx.first < 0)
0593     return false;
0594   return ((irad >= hgpar_->iradMinBH_[indx.first]) && (irad <= (hgpar_->iradMaxBH_[indx.first] + 1)) && (iphi > 0) &&
0595           (iphi <= hgpar_->scintCells(layer)));
0596 }
0597 
0598 int HGCalDDDConstants::lastLayer(bool reco) const {
0599   return (hgpar_->firstLayer_ + tot_layers_[static_cast<int>(reco)] - 1);
0600 }
0601 
0602 unsigned int HGCalDDDConstants::layers(bool reco) const { return tot_layers_[static_cast<int>(reco)]; }
0603 
0604 int HGCalDDDConstants::layerIndex(int lay, bool reco) const {
0605   int ll = lay - hgpar_->firstLayer_;
0606   if (ll < 0 || ll >= static_cast<int>(hgpar_->layerIndex_.size()))
0607     return -1;
0608   if (waferHexagon6()) {
0609     if (reco && ll >= static_cast<int>(hgpar_->depthIndex_.size()))
0610       return -1;
0611     return (reco ? hgpar_->depthLayerF_[ll] : hgpar_->layerIndex_[ll]);
0612   } else {
0613     return (hgpar_->layerIndex_[ll]);
0614   }
0615 }
0616 
0617 unsigned int HGCalDDDConstants::layersInit(bool reco) const {
0618   return (reco ? hgpar_->depthIndex_.size() : hgpar_->layerIndex_.size());
0619 }
0620 
0621 std::pair<float, float> HGCalDDDConstants::localToGlobal8(
0622     int lay, int waferU, int waferV, double localX, double localY, bool reco, bool debug) const {
0623   double x(localX), y(localY);
0624   bool rotx =
0625       ((!hgpar_->layerType_.empty()) && (hgpar_->layerType_[lay - hgpar_->firstLayer_] == HGCalTypes::WaferCenterR));
0626   if (debug)
0627     edm::LogVerbatim("HGCalGeom") << "LocalToGlobal8 " << lay << ":" << (lay - hgpar_->firstLayer_) << ":" << rotx
0628                                   << " Local (" << x << ":" << y << ") Reco " << reco;
0629   if (!reco) {
0630     x *= HGCalParameters::k_ScaleToDDD;
0631     y *= HGCalParameters::k_ScaleToDDD;
0632   }
0633   const auto& xy = waferPositionNoRot(lay, waferU, waferV, reco, debug);
0634   x += xy.first;
0635   y += xy.second;
0636   int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0637   auto ktr = hgpar_->waferInfoMap_.find(indx);
0638   if ((mode_ == HGCalGeometryMode::Hexagon8Cassette) && (ktr != hgpar_->waferInfoMap_.end())) {
0639     auto cshift = hgcassette_.getShift(lay, 1, (ktr->second).cassette);
0640     if (debug)
0641       edm::LogVerbatim("HGCalGeom") << "Cassette " << (ktr->second).cassette << " Shift " << cshift.first << ":"
0642                                     << cshift.second;
0643     if (!reco) {
0644       x += (HGCalParameters::k_ScaleToDDD)*cshift.first;
0645       y += (HGCalParameters::k_ScaleToDDD)*cshift.second;
0646     } else {
0647       x += cshift.first;
0648       y += cshift.second;
0649     }
0650   }
0651   if (debug)
0652     edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by adding " << xy.first << ":" << xy.second;
0653   return (rotx ? getXY(lay, x, y, false) : std::make_pair(x, y));
0654 }
0655 
0656 std::pair<float, float> HGCalDDDConstants::locateCell(int cell, int lay, int type, bool reco) const {
0657   // type refers to wafer # for hexagon cell
0658   float x(999999.), y(999999.);
0659   const auto& index = getIndex(lay, reco);
0660   int i = index.first;
0661   if (i < 0)
0662     return std::make_pair(x, y);
0663   if (waferHexagon6()) {
0664     x = hgpar_->waferPosX_[type];
0665     y = hgpar_->waferPosY_[type];
0666 #ifdef EDM_ML_DEBUG
0667     float x0(x), y0(y);
0668 #endif
0669     if (hgpar_->waferTypeT_[type] - 1 == HGCSiliconDetId::HGCalFine) {
0670       x += hgpar_->cellFineX_[cell];
0671       y += hgpar_->cellFineY_[cell];
0672     } else {
0673       x += hgpar_->cellCoarseX_[cell];
0674       y += hgpar_->cellCoarseY_[cell];
0675     }
0676 #ifdef EDM_ML_DEBUG
0677     edm::LogVerbatim("HGCalGeom") << "LocateCell (Wafer) " << x0 << ":" << y0 << " Final " << x << ":" << y;
0678 #endif
0679     if (!reco) {
0680       x *= HGCalParameters::k_ScaleToDDD;
0681       y *= HGCalParameters::k_ScaleToDDD;
0682     }
0683   }
0684   return std::make_pair(x, y);
0685 }
0686 
0687 std::pair<float, float> HGCalDDDConstants::locateCell(
0688     int lay, int waferU, int waferV, int cellU, int cellV, bool reco, bool all, bool norot, bool debug) const {
0689   double x(0), y(0);
0690   int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0691   auto itr = hgpar_->typesInLayers_.find(indx);
0692   int type = ((itr == hgpar_->typesInLayers_.end()) ? 2 : hgpar_->waferTypeL_[itr->second]);
0693   int layertype = layerType(lay);
0694   bool rotx = (norot) ? false : (layertype == HGCalTypes::WaferCenterR);
0695   if (debug) {
0696     edm::LogVerbatim("HGCalGeom") << "LocateCell " << lay << ":" << (lay - hgpar_->firstLayer_) << ":" << layertype
0697                                   << ":" << rotx << ":" << waferU << ":" << waferV << ":" << indx << ":"
0698                                   << (itr == hgpar_->typesInLayers_.end()) << ":" << type << " Flags " << reco << ":"
0699                                   << all;
0700   }
0701   auto ktr = hgpar_->waferInfoMap_.end();
0702   int place(HGCalCell::cellPlacementOld);
0703   if (waferHexagon8File()) {
0704     if (mode_ == HGCalGeometryMode::Hexagon8Cassette) {
0705       ktr = hgpar_->waferInfoMap_.find(indx);
0706       if (ktr != hgpar_->waferInfoMap_.end())
0707         place = HGCalCell::cellPlacementIndex(1, HGCalTypes::layerFrontBack(layertype), (ktr->second).orient);
0708     }
0709     auto xy = hgcell_->cellUV2XY2(cellU, cellV, place, type);
0710     x = xy.first;
0711     y = xy.second;
0712     if (debug)
0713       edm::LogVerbatim("HGCalGeom") << "Type " << type << " Place " << place << " Cell " << cellU << ":" << cellV
0714                                     << " Position " << x << ":" << y;
0715   } else {
0716     int kndx = cellV * 100 + cellU;
0717     if (type == 0) {
0718       auto jtr = hgpar_->cellFineIndex_.find(kndx);
0719       if (jtr != hgpar_->cellFineIndex_.end()) {
0720         x = hgpar_->cellFineX_[jtr->second];
0721         y = hgpar_->cellFineY_[jtr->second];
0722       }
0723       if (debug)
0724         edm::LogVerbatim("HGCalGeom") << "Fine " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
0725                                       << (jtr != hgpar_->cellFineIndex_.end());
0726     } else {
0727       auto jtr = hgpar_->cellCoarseIndex_.find(kndx);
0728       if (jtr != hgpar_->cellCoarseIndex_.end()) {
0729         x = hgpar_->cellCoarseX_[jtr->second];
0730         y = hgpar_->cellCoarseY_[jtr->second];
0731       }
0732       if (debug)
0733         edm::LogVerbatim("HGCalGeom") << "Coarse " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y
0734                                       << ":" << (jtr != hgpar_->cellCoarseIndex_.end());
0735     }
0736   }
0737   if (!reco) {
0738     x *= HGCalParameters::k_ScaleToDDD;
0739     y *= HGCalParameters::k_ScaleToDDD;
0740   }
0741   if (all) {
0742     const auto& xy = waferPositionNoRot(lay, waferU, waferV, reco, debug);
0743     x += xy.first;
0744     y += xy.second;
0745     if ((mode_ == HGCalGeometryMode::Hexagon8Cassette) && (ktr != hgpar_->waferInfoMap_.end())) {
0746       auto cshift = hgcassette_.getShift(lay, 1, (ktr->second).cassette);
0747       if (debug)
0748         edm::LogVerbatim("HGCalGeom") << "Cassette " << (ktr->second).cassette << " Shift " << cshift.first << ":"
0749                                       << cshift.second;
0750       if (!reco) {
0751         x += (HGCalParameters::k_ScaleToDDD)*cshift.first;
0752         y += (HGCalParameters::k_ScaleToDDD)*cshift.second;
0753       } else {
0754         x += cshift.first;
0755         y += cshift.second;
0756       }
0757     }
0758     if (debug)
0759       edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by adding " << xy.first << ":" << xy.second;
0760   }
0761   return (rotx ? getXY(lay, x, y, false) : std::make_pair(x, y));
0762 }
0763 
0764 std::pair<float, float> HGCalDDDConstants::locateCell(const HGCSiliconDetId& id, bool debug) const {
0765   return locateCell(id.layer(), id.waferU(), id.waferV(), id.cellU(), id.cellV(), true, true, false, debug);
0766 }
0767 
0768 std::pair<float, float> HGCalDDDConstants::locateCell(const HGCScintillatorDetId& id, bool debug) const {
0769   return locateCellTrap(id.layer(), id.iradius(), id.iphi(), true, debug);
0770 }
0771 
0772 std::pair<float, float> HGCalDDDConstants::locateCellHex(int cell, int wafer, bool reco) const {
0773   float x(0), y(0);
0774   if (hgpar_->waferTypeT_[wafer] - 1 == HGCSiliconDetId::HGCalFine) {
0775     x = hgpar_->cellFineX_[cell];
0776     y = hgpar_->cellFineY_[cell];
0777   } else {
0778     x = hgpar_->cellCoarseX_[cell];
0779     y = hgpar_->cellCoarseY_[cell];
0780   }
0781   if (!reco) {
0782     x *= HGCalParameters::k_ScaleToDDD;
0783     y *= HGCalParameters::k_ScaleToDDD;
0784   }
0785   return std::make_pair(x, y);
0786 }
0787 
0788 std::pair<float, float> HGCalDDDConstants::locateCellTrap(int lay, int irad, int iphi, bool reco, bool debug) const {
0789   float x(0), y(0);
0790   const auto& indx = getIndex(lay, reco);
0791   if (indx.first >= 0) {
0792     int ir = std::abs(irad);
0793     int type = hgpar_->scintType(lay);
0794     double phi = (iphi - 0.5) * indx.second;
0795     double z = hgpar_->zLayerHex_[indx.first];
0796     double r = 0.5 * (hgpar_->radiusLayer_[type][ir - 1] + hgpar_->radiusLayer_[type][ir]);
0797     std::pair<double, double> range = rangeR(z, true);
0798     if (debug)
0799       edm::LogVerbatim("HGCalGeom") << "locateCellTrap:: Input " << lay << ":" << irad << ":" << iphi << ":" << reco
0800                                     << " IR " << ir << ":" << hgpar_->iradMinBH_[indx.first] << ":"
0801                                     << hgpar_->iradMaxBH_[indx.first] << " Type " << type << " Z " << indx.first << ":"
0802                                     << z << " phi " << phi << " R " << r << ":" << range.first << ":" << range.second;
0803     if ((mode_ != HGCalGeometryMode::TrapezoidFile) && (mode_ != HGCalGeometryMode::TrapezoidModule) &&
0804         (mode_ != HGCalGeometryMode::TrapezoidCassette))
0805       r = std::max(range.first, std::min(r, range.second));
0806     x = r * std::cos(phi);
0807     y = r * std::sin(phi);
0808     if (mode_ == HGCalGeometryMode::TrapezoidCassette) {
0809       int cassette = HGCalTileIndex::tileCassette(iphi, hgpar_->phiOffset_, hgpar_->nphiCassette_, hgpar_->cassettes_);
0810       auto cshift = hgcassette_.getShift(lay, 1, cassette);
0811       if (debug)
0812         edm::LogVerbatim("HGCalGeom") << "Cassette " << cassette << " Shift " << cshift.first << ":" << cshift.second;
0813       x += cshift.first;
0814       y += cshift.second;
0815     }
0816     if (irad < 0)
0817       x = -x;
0818   }
0819   if (!reco) {
0820     x *= HGCalParameters::k_ScaleToDDD;
0821     y *= HGCalParameters::k_ScaleToDDD;
0822   }
0823   return std::make_pair(x, y);
0824 }
0825 
0826 bool HGCalDDDConstants::maskCell(const DetId& detId, int corners) const {
0827   bool mask(false);
0828   if (corners > 2 && corners <= static_cast<int>(HGCalParameters::k_CornerSize)) {
0829     if (waferHexagon8()) {
0830       int N(0), layer(0), waferU(0), waferV(0), u(0), v(0);
0831       if (detId.det() == DetId::Forward) {
0832         HFNoseDetId id(detId);
0833         N = getUVMax(id.type());
0834         layer = id.layer();
0835         waferU = id.waferU();
0836         waferV = id.waferV();
0837         u = id.cellU();
0838         v = id.cellV();
0839       } else {
0840         HGCSiliconDetId id(detId);
0841         N = getUVMax(id.type());
0842         layer = id.layer();
0843         waferU = id.waferU();
0844         waferV = id.waferV();
0845         u = id.cellU();
0846         v = id.cellV();
0847       }
0848       int wl = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
0849       auto itr = hgpar_->waferTypes_.find(wl);
0850       auto ktr = hgpar_->waferInfoMap_.find(wl);
0851 #ifdef EDM_ML_DEBUG
0852       edm::LogVerbatim("HGCalGeom") << "MaskCell: Layer " << layer << " Wafer " << waferU << ":" << waferV << " Index "
0853                                     << wl << ":" << (itr != hgpar_->waferTypes_.end()) << ":"
0854                                     << (ktr != hgpar_->waferInfoMap_.end());
0855 #endif
0856       if (mode_ == HGCalGeometryMode::Hexagon8Cassette) {
0857         int part = (ktr != hgpar_->waferInfoMap_.end()) ? (ktr->second).part : HGCalTypes::WaferFull;
0858         mask = !(HGCalWaferMask::goodCell(u, v, part));
0859       } else if (itr != hgpar_->waferTypes_.end()) {
0860         if ((itr->second).second <= HGCalTypes::k_OffsetRotation)
0861           mask = HGCalWaferMask::maskCell(u, v, N, (itr->second).first, (itr->second).second, corners);
0862         else
0863           mask = !(HGCalWaferMask::goodCell(
0864               u, v, N, (itr->second).first, ((itr->second).second - HGCalTypes::k_OffsetRotation)));
0865       }
0866     }
0867   }
0868   return mask;
0869 }
0870 
0871 int HGCalDDDConstants::maxCells(bool reco) const {
0872   int cells(0);
0873   for (unsigned int i = 0; i < layers(reco); ++i) {
0874     int lay = reco ? hgpar_->depth_[i] : hgpar_->layer_[i];
0875     if (cells < maxCells(lay, reco))
0876       cells = maxCells(lay, reco);
0877   }
0878   return cells;
0879 }
0880 
0881 int HGCalDDDConstants::maxCells(int lay, bool reco) const {
0882   const auto& index = getIndex(lay, reco);
0883   if (index.first < 0)
0884     return 0;
0885   if (waferHexagon6()) {
0886     unsigned int cells(0);
0887     for (unsigned int k = 0; k < hgpar_->waferTypeT_.size(); ++k) {
0888       if (waferInLayerTest(k, index.first, hgpar_->defineFull_)) {
0889         unsigned int cell = (hgpar_->waferTypeT_[k] - 1 == HGCSiliconDetId::HGCalFine) ? (hgpar_->cellFineX_.size())
0890                                                                                        : (hgpar_->cellCoarseX_.size());
0891         if (cell > cells)
0892           cells = cell;
0893       }
0894     }
0895     return static_cast<int>(cells);
0896   } else if (waferHexagon8()) {
0897     int cells(0);
0898     for (unsigned int k = 0; k < hgpar_->waferCopy_.size(); ++k) {
0899       if (waferInLayerTest(k, index.first, hgpar_->defineFull_)) {
0900         auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(
0901             lay, HGCalWaferIndex::waferU(hgpar_->waferCopy_[k]), HGCalWaferIndex::waferV(hgpar_->waferCopy_[k])));
0902         int type = ((itr == hgpar_->typesInLayers_.end()) ? HGCSiliconDetId::HGCalCoarseThick
0903                                                           : hgpar_->waferTypeL_[itr->second]);
0904         int N = (type == HGCSiliconDetId::HGCalFine) ? hgpar_->nCellsFine_ : hgpar_->nCellsCoarse_;
0905         cells = std::max(cells, 3 * N * N);
0906       }
0907     }
0908     return cells;
0909   } else if (tileTrapezoid()) {
0910     return hgpar_->scintCells(index.first + hgpar_->firstLayer_);
0911   } else {
0912     return 0;
0913   }
0914 }
0915 
0916 int HGCalDDDConstants::maxRows(int lay, bool reco) const {
0917   int kymax(0);
0918   const auto& index = getIndex(lay, reco);
0919   int i = index.first;
0920   if (i < 0)
0921     return kymax;
0922   if (waferHexagon6()) {
0923     for (unsigned int k = 0; k < hgpar_->waferCopy_.size(); ++k) {
0924       if (waferInLayerTest(k, i, hgpar_->defineFull_)) {
0925         int ky = ((hgpar_->waferCopy_[k]) / 100) % 100;
0926         if (ky > kymax)
0927           kymax = ky;
0928       }
0929     }
0930   } else if (waferHexagon8()) {
0931     kymax = 1 + 2 * hgpar_->waferUVMaxLayer_[index.first];
0932   }
0933   return kymax;
0934 }
0935 
0936 int HGCalDDDConstants::modifyUV(int uv, int type1, int type2) const {
0937   // Modify u/v for transition of type1 to type2
0938   return (((type1 == type2) || (type1 * type2 != 0)) ? uv : ((type1 == 0) ? (2 * uv + 1) / 3 : (3 * uv) / 2));
0939 }
0940 
0941 int HGCalDDDConstants::modules(int lay, bool reco) const {
0942   if (getIndex(lay, reco).first < 0)
0943     return 0;
0944   else
0945     return max_modules_layer_[static_cast<int>(reco)][lay];
0946 }
0947 
0948 int HGCalDDDConstants::modulesInit(int lay, bool reco) const {
0949   int nmod(0);
0950   const auto& index = getIndex(lay, reco);
0951   if (index.first < 0)
0952     return nmod;
0953   if (!tileTrapezoid()) {
0954     for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
0955       if (waferInLayerTest(k, index.first, hgpar_->defineFull_))
0956         ++nmod;
0957     }
0958   } else {
0959     nmod = 1 + hgpar_->lastModule_[index.first] - hgpar_->firstModule_[index.first];
0960   }
0961   return nmod;
0962 }
0963 
0964 double HGCalDDDConstants::mouseBite(bool reco) const {
0965   return (reco ? hgpar_->mouseBite_ : HGCalParameters::k_ScaleToDDD * hgpar_->mouseBite_);
0966 }
0967 
0968 int HGCalDDDConstants::numberCells(bool reco) const {
0969   int cells = (tileTrapezoid() && ((hgpar_->waferMaskMode_ == HGCalGeomParameters::scintillatorFile) ||
0970                                    (hgpar_->waferMaskMode_ == HGCalGeomParameters::scintillatorCassette)))
0971                   ? tileCount(0, -1)
0972                   : 0;
0973   if (cells == 0) {
0974     unsigned int nlayer = (reco) ? hgpar_->depth_.size() : hgpar_->layer_.size();
0975     for (unsigned k = 0; k < nlayer; ++k) {
0976       std::vector<int> ncells = numberCells(((reco) ? hgpar_->depth_[k] : hgpar_->layer_[k]), reco);
0977       cells = std::accumulate(ncells.begin(), ncells.end(), cells);
0978     }
0979   }
0980   return cells;
0981 }
0982 
0983 std::vector<int> HGCalDDDConstants::numberCells(int lay, bool reco) const {
0984   const auto& index = getIndex(lay, reco);
0985   int i = index.first;
0986   std::vector<int> ncell;
0987   if (i >= 0) {
0988     if (waferHexagon6()) {
0989       for (unsigned int k = 0; k < hgpar_->waferTypeT_.size(); ++k) {
0990         if (waferInLayerTest(k, i, hgpar_->defineFull_)) {
0991           unsigned int cell = (hgpar_->waferTypeT_[k] - 1 == HGCSiliconDetId::HGCalFine)
0992                                   ? (hgpar_->cellFineX_.size())
0993                                   : (hgpar_->cellCoarseX_.size());
0994           ncell.emplace_back(static_cast<int>(cell));
0995         }
0996       }
0997     } else if (tileTrapezoid()) {
0998       int nphi = hgpar_->scintCells(lay);
0999       for (int k = hgpar_->firstModule_[i]; k <= hgpar_->lastModule_[i]; ++k)
1000         ncell.emplace_back(nphi);
1001     } else {
1002       for (unsigned int k = 0; k < hgpar_->waferCopy_.size(); ++k) {
1003         if (waferInLayerTest(k, index.first, hgpar_->defineFull_)) {
1004           int cell = numberCellsHexagon(lay,
1005                                         HGCalWaferIndex::waferU(hgpar_->waferCopy_[k]),
1006                                         HGCalWaferIndex::waferV(hgpar_->waferCopy_[k]),
1007                                         true);
1008           ncell.emplace_back(cell);
1009         }
1010       }
1011     }
1012   }
1013   return ncell;
1014 }
1015 
1016 int HGCalDDDConstants::numberCellsHexagon(int wafer) const {
1017   if (wafer >= 0 && wafer < static_cast<int>(hgpar_->waferTypeT_.size())) {
1018     if (hgpar_->waferTypeT_[wafer] - 1 == HGCSiliconDetId::HGCalFine)
1019       return static_cast<int>(hgpar_->cellFineX_.size());
1020     else
1021       return static_cast<int>(hgpar_->cellCoarseX_.size());
1022   } else {
1023     return 0;
1024   }
1025 }
1026 
1027 int HGCalDDDConstants::numberCellsHexagon(int lay, int waferU, int waferV, bool flag) const {
1028   auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(lay, waferU, waferV));
1029   int type =
1030       ((itr == hgpar_->typesInLayers_.end()) ? HGCSiliconDetId::HGCalCoarseThick : hgpar_->waferTypeL_[itr->second]);
1031   int N = (type == 0) ? hgpar_->nCellsFine_ : hgpar_->nCellsCoarse_;
1032   if (flag)
1033     return (3 * N * N);
1034   else
1035     return N;
1036 }
1037 
1038 std::pair<double, double> HGCalDDDConstants::rangeR(double z, bool reco) const {
1039   double rmin(0), rmax(0), zz(0);
1040   if (hgpar_->detectorType_ > 0) {
1041     zz = (reco ? std::abs(z) : HGCalParameters::k_ScaleFromDDD * std::abs(z));
1042     if (hgpar_->detectorType_ <= 2) {
1043       rmin = HGCalGeomTools::radius(zz, hgpar_->zFrontMin_, hgpar_->rMinFront_, hgpar_->slopeMin_);
1044     } else {
1045       rmin = HGCalGeomTools::radius(
1046           zz, hgpar_->firstLayer_, hgpar_->firstMixedLayer_, hgpar_->zLayerHex_, hgpar_->radiusMixBoundary_);
1047     }
1048     if ((hgpar_->detectorType_ == 2) && (zz >= hgpar_->zLayerHex_[hgpar_->firstMixedLayer_ - 1])) {
1049       rmax = HGCalGeomTools::radius(
1050           zz, hgpar_->firstLayer_, hgpar_->firstMixedLayer_, hgpar_->zLayerHex_, hgpar_->radiusMixBoundary_);
1051     } else {
1052       rmax = HGCalGeomTools::radius(zz, hgpar_->zFrontTop_, hgpar_->rMaxFront_, hgpar_->slopeTop_);
1053     }
1054   }
1055   if (!reco) {
1056     rmin *= HGCalParameters::k_ScaleToDDD;
1057     rmax *= HGCalParameters::k_ScaleToDDD;
1058   }
1059 #ifdef EDM_ML_DEBUG
1060   edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:rangeR: " << z << ":" << zz << " R " << rmin << ":" << rmax;
1061 #endif
1062   return std::make_pair(rmin, rmax);
1063 }
1064 
1065 std::pair<double, double> HGCalDDDConstants::rangeRLayer(int lay, bool reco) const {
1066   double rmin(0), rmax(0);
1067   const auto& index = getIndex(lay, reco);
1068   if (index.first >= 0 && index.first < static_cast<int>(hgpar_->rMinLayHex_.size())) {
1069     rmin = hgpar_->rMinLayHex_[index.first];
1070     rmax = hgpar_->rMaxLayHex_[index.first];
1071   }
1072   if (!reco) {
1073     rmin *= HGCalParameters::k_ScaleToDDD;
1074     rmax *= HGCalParameters::k_ScaleToDDD;
1075   }
1076 #ifdef EDM_ML_DEBUG
1077   edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:rangeR: " << lay << ":" << index.first << " R " << rmin << ":"
1078                                 << rmax;
1079 #endif
1080   return std::make_pair(rmin, rmax);
1081 }
1082 
1083 std::pair<double, double> HGCalDDDConstants::rangeZ(bool reco) const {
1084   double zmin = (hgpar_->zLayerHex_[0] - hgpar_->waferThick_);
1085   double zmax = (hgpar_->zLayerHex_[hgpar_->zLayerHex_.size() - 1] + hgpar_->waferThick_);
1086 #ifdef EDM_ML_DEBUG
1087   edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:rangeZ: " << zmin << ":" << zmax << ":" << hgpar_->waferThick_;
1088 #endif
1089   if (!reco) {
1090     zmin *= HGCalParameters::k_ScaleToDDD;
1091     zmax *= HGCalParameters::k_ScaleToDDD;
1092   }
1093   return std::make_pair(zmin, zmax);
1094 }
1095 
1096 std::pair<int, int> HGCalDDDConstants::rowColumnWafer(int wafer) const {
1097   int row(0), col(0);
1098   if (wafer < static_cast<int>(hgpar_->waferCopy_.size())) {
1099     int copy = hgpar_->waferCopy_[wafer];
1100     col = HGCalTypes::getUnpackedU(copy);
1101     row = HGCalTypes::getUnpackedV(copy);
1102     ;
1103   }
1104   return std::make_pair(row, col);
1105 }
1106 
1107 std::pair<int, int> HGCalDDDConstants::simToReco(int cell, int lay, int mod, bool half) const {
1108   if (!waferHexagon6()) {
1109     return std::make_pair(cell, lay);
1110   } else {
1111     const auto& index = getIndex(lay, false);
1112     int i = index.first;
1113     if (i < 0) {
1114       edm::LogWarning("HGCalGeom") << "Wrong Layer # " << lay << " not in the list ***** ERROR *****";
1115       return std::make_pair(-1, -1);
1116     }
1117     if (mod >= static_cast<int>(hgpar_->waferTypeL_.size())) {
1118       edm::LogWarning("HGCalGeom") << "Invalid Wafer # " << mod << "should be < " << (hgpar_->waferTypeL_).size()
1119                                    << " ***** ERROR *****";
1120       return std::make_pair(-1, -1);
1121     }
1122     int depth(-1);
1123     int kx = cell;
1124     int type = hgpar_->waferTypeL_[mod];
1125     if (type == 1) {
1126       depth = hgpar_->layerGroup_[i];
1127     } else if (type == 2) {
1128       depth = hgpar_->layerGroupM_[i];
1129     } else {
1130       depth = hgpar_->layerGroupO_[i];
1131     }
1132     return std::make_pair(kx, depth);
1133   }
1134 }
1135 
1136 int HGCalDDDConstants::tileCount(int layer, int ring) const {
1137   int laymin(layer), laymax(layer), ringmin(ring), ringmax(ring), kount(0);
1138   if (layer == 0) {
1139     laymin = hgpar_->firstLayer_;
1140     laymax = lastLayer(true);
1141   }
1142 #ifdef EDM_ML_DEBUG
1143   edm::LogVerbatim("HGCalGeom") << "tileCount: layer " << layer << " ring " << ring << " layerMin/Max " << laymin << ":"
1144                                 << laymax;
1145 #endif
1146   for (int lay = laymin; lay <= laymax; ++lay) {
1147     if (ring < 0) {
1148       int ll = lay - hgpar_->firstLayer_;
1149       ringmin = hgpar_->tileRingRange_[ll].first;
1150       ringmax = hgpar_->tileRingRange_[ll].second;
1151     }
1152 #ifdef EDM_ML_DEBUG
1153     edm::LogVerbatim("HGCalGeom") << "tileCount: lay " << lay << ":" << (lay - hgpar_->firstLayer_) << " rings "
1154                                   << ringmin << ":" << ringmax;
1155 #endif
1156     for (int rin = ringmin; rin <= ringmax; ++rin) {
1157       int indx = HGCalTileIndex::tileIndex(lay, rin + 1, 0);
1158       auto itr = hgpar_->tileInfoMap_.find(indx);
1159 #ifdef EDM_ML_DEBUG
1160       edm::LogVerbatim("HGCalGeom") << "tileCount: rin  " << rin << " indx " << indx << " itr "
1161                                     << (itr != hgpar_->tileInfoMap_.end());
1162 #endif
1163       if (itr != hgpar_->tileInfoMap_.end()) {
1164         for (int k = 0; k < 4; ++k) {
1165           std::bitset<24> b(itr->second.hex[k]);
1166           kount += b.count();
1167         }
1168       }
1169 #ifdef EDM_ML_DEBUG
1170       edm::LogVerbatim("HGCalGeom") << "tileCount: lay|rin  " << lay << ":" << rin << " kount " << kount;
1171 #endif
1172     }
1173   }
1174   return (3 * kount);
1175 }
1176 
1177 bool HGCalDDDConstants::tileExist(int zside, int layer, int ring, int phi) const {
1178   int indx = HGCalTileIndex::tileIndex(layer, ring, 0);
1179   auto itr = hgpar_->tileInfoMap_.find(indx);
1180   bool ok = (itr == hgpar_->tileInfoMap_.end()) ? false : HGCalTileIndex::tileExist(itr->second.hex, zside, phi);
1181   return ok;
1182 }
1183 
1184 HGCalParameters::tileInfo HGCalDDDConstants::tileInfo(int zside, int layer, int ring) const {
1185   int indx = HGCalTileIndex::tileIndex(layer, ring, 0);
1186   auto itr = hgpar_->tileInfoMap_.find(indx);
1187   return ((itr == hgpar_->tileInfoMap_.end()) ? HGCalParameters::tileInfo() : itr->second);
1188 }
1189 
1190 std::pair<int, int> HGCalDDDConstants::tileRings(int layer) const {
1191   if ((mode_ == HGCalGeometryMode::TrapezoidFile) || (mode_ == HGCalGeometryMode::TrapezoidModule) ||
1192       (mode_ == HGCalGeometryMode::TrapezoidCassette)) {
1193     int ll = layer - hgpar_->firstLayer_;
1194     if (ll >= 0 && ll < static_cast<int>(hgpar_->tileRingRange_.size()))
1195       return hgpar_->tileRingRange_[ll];
1196   }
1197   return std::make_pair(0, 0);
1198 }
1199 
1200 std::pair<int, int> HGCalDDDConstants::tileType(int layer, int ring, int phi) const {
1201   int indx = HGCalTileIndex::tileIndex(layer, ring, phi);
1202   int type(-1), sipm(-1);
1203   auto itr = hgpar_->tileInfoMap_.find(indx);
1204   if (itr != hgpar_->tileInfoMap_.end()) {
1205     type = 1 + (itr->second).type;
1206     sipm = ((itr->second).sipm == HGCalTypes::SiPMLarge) ? 0 : 1;
1207   }
1208   return std::make_pair(type, sipm);
1209 }
1210 
1211 int HGCalDDDConstants::waferFromCopy(int copy) const {
1212   const int ncopies = hgpar_->waferCopy_.size();
1213   int wafer(ncopies);
1214   bool result(false);
1215   for (int k = 0; k < ncopies; ++k) {
1216     if (copy == hgpar_->waferCopy_[k]) {
1217       wafer = k;
1218       result = true;
1219       break;
1220     }
1221   }
1222   if (!result) {
1223     wafer = -1;
1224 #ifdef EDM_ML_DEBUG
1225     edm::LogVerbatim("HGCalGeom") << "Cannot find " << copy << " in a list of " << ncopies << " members";
1226     for (int k = 0; k < ncopies; ++k)
1227       edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << hgpar_->waferCopy_[k];
1228 #endif
1229   }
1230 #ifdef EDM_ML_DEBUG
1231   edm::LogVerbatim("HGCalGeom") << "WaferFromCopy " << copy << ":" << wafer << ":" << result;
1232 #endif
1233   return wafer;
1234 }
1235 
1236 void HGCalDDDConstants::waferFromPosition(const double x, const double y, int& wafer, int& icell, int& celltyp) const {
1237   // Input x, y in Geant4 unit and transformed to CMSSW standard
1238   double xx = HGCalParameters::k_ScaleFromDDD * x;
1239   double yy = HGCalParameters::k_ScaleFromDDD * y;
1240   int size_ = static_cast<int>(hgpar_->waferCopy_.size());
1241   wafer = size_;
1242   for (int k = 0; k < size_; ++k) {
1243     double dx = std::abs(xx - hgpar_->waferPosX_[k]);
1244     double dy = std::abs(yy - hgpar_->waferPosY_[k]);
1245     if (dx <= rmax_ && dy <= hexside_) {
1246       if ((dy <= 0.5 * hexside_) || (dx * tan30deg_ <= (hexside_ - dy))) {
1247         wafer = k;
1248         celltyp = hgpar_->waferTypeT_[k];
1249         xx -= hgpar_->waferPosX_[k];
1250         yy -= hgpar_->waferPosY_[k];
1251         break;
1252       }
1253     }
1254   }
1255   if (wafer < size_) {
1256     if (celltyp - 1 == HGCSiliconDetId::HGCalFine)
1257       icell = cellHex(
1258           xx, yy, 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[0], hgpar_->cellFineX_, hgpar_->cellFineY_);
1259     else
1260       icell = cellHex(xx,
1261                       yy,
1262                       0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[1],
1263                       hgpar_->cellCoarseX_,
1264                       hgpar_->cellCoarseY_);
1265   } else {
1266     wafer = -1;
1267 #ifdef EDM_ML_DEBUG
1268     edm::LogWarning("HGCalGeom") << "Cannot get wafer type corresponding to " << x << ":" << y << "    " << xx << ":"
1269                                  << yy;
1270 #endif
1271   }
1272 #ifdef EDM_ML_DEBUG
1273   edm::LogVerbatim("HGCalGeom") << "Position " << x << ":" << y << " Wafer " << wafer << ":" << size_ << " XX " << xx
1274                                 << ":" << yy << " Cell " << icell << " Type " << celltyp;
1275 #endif
1276 }
1277 
1278 void HGCalDDDConstants::waferFromPosition(const double x,
1279                                           const double y,
1280                                           const int layer,
1281                                           int& waferU,
1282                                           int& waferV,
1283                                           int& cellU,
1284                                           int& cellV,
1285                                           int& celltype,
1286                                           double& wt,
1287                                           bool extend,
1288                                           bool debug) const {
1289   // Expect x, y as in SIM step
1290   waferU = waferV = 1 + hgpar_->waferUVMax_;
1291   cellU = cellV = celltype = 0;
1292   if ((hgpar_->xLayerHex_.empty()) || (hgpar_->yLayerHex_.empty()))
1293     return;
1294   int ll = layer - hgpar_->firstLayer_;
1295   int layertype = layerType(layer);
1296   bool rotx = ((!hgpar_->layerType_.empty()) && (layertype == HGCalTypes::WaferCenterR));
1297   double xx(0), yy(0);
1298   if (rotx) {
1299     std::pair<double, double> xy =
1300         getXY(layer, HGCalParameters::k_ScaleFromDDD * x, HGCalParameters::k_ScaleFromDDD * y, true);
1301     xx = xy.first - hgpar_->xLayerHex_[ll];
1302     yy = xy.second - hgpar_->yLayerHex_[ll];
1303   } else {
1304     xx = HGCalParameters::k_ScaleFromDDD * x - hgpar_->xLayerHex_[ll];
1305     yy = HGCalParameters::k_ScaleFromDDD * y - hgpar_->yLayerHex_[ll];
1306   }
1307 #ifdef EDM_ML_DEBUG
1308   if (debug)
1309     edm::LogVerbatim("HGCalGeom") << "waferFromPosition:: Layer " << layer << ":" << ll << " Rot " << rotx << " X " << x
1310                                   << ":" << xx << " Y " << y << ":" << yy;
1311 #endif
1312   double rmax = extend ? rmaxT_ : rmax_;
1313   double hexside = extend ? hexsideT_ : hexside_;
1314   for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
1315     double dx0(0), dy0(0);
1316     waferU = HGCalWaferIndex::waferU(hgpar_->waferCopy_[k]);
1317     waferV = HGCalWaferIndex::waferV(hgpar_->waferCopy_[k]);
1318     if (mode_ == HGCalGeometryMode::Hexagon8Cassette) {
1319       int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1320       auto ktr = hgpar_->waferInfoMap_.find(indx);
1321       if (ktr != hgpar_->waferInfoMap_.end()) {
1322         auto cshift = hgcassette_.getShift(layer, 1, (ktr->second).cassette);
1323 #ifdef EDM_ML_DEBUG
1324         if (debug)
1325           edm::LogVerbatim("HGCalGeom") << "Cassette " << (ktr->second).cassette << " Shift " << cshift.first << ":"
1326                                         << cshift.second;
1327 #endif
1328         dx0 = cshift.first;
1329         dy0 = cshift.second;
1330       }
1331     }
1332     double dx = std::abs(xx - dx0 - hgpar_->waferPosX_[k]);
1333     double dy = std::abs(yy - dy0 - hgpar_->waferPosY_[k]);
1334     if (dx <= rmax && dy <= hexside) {
1335       if ((dy <= 0.5 * hexside) || (dx * tan30deg_ <= (hexside - dy))) {
1336         if (waferHexagon8File()) {
1337           int index = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1338           celltype = HGCalWaferType::getType(index, hgpar_->waferInfoMap_);
1339         } else {
1340           auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1341           celltype = ((itr == hgpar_->typesInLayers_.end()) ? HGCSiliconDetId::HGCalCoarseThick
1342                                                             : hgpar_->waferTypeL_[itr->second]);
1343         }
1344 #ifdef EDM_ML_DEBUG
1345         if (debug)
1346           edm::LogVerbatim("HGCalGeom") << "WaferFromPosition:: Input " << layer << ":" << ll << ":"
1347                                         << hgpar_->firstLayer_ << ":" << rotx << ":" << x << ":" << y << ":"
1348                                         << hgpar_->xLayerHex_[ll] << ":" << hgpar_->yLayerHex_[ll] << ":" << xx << ":"
1349                                         << yy << " compared with " << hgpar_->waferPosX_[k] << ":"
1350                                         << hgpar_->waferPosY_[k] << " difference " << dx << ":" << dy << ":"
1351                                         << dx * tan30deg_ << ":" << (hexside_ - dy) << " comparator " << rmax_ << ":"
1352                                         << rmaxT_ << ":" << hexside_ << ":" << hexsideT_ << " wafer " << waferU << ":"
1353                                         << waferV << ":" << celltype;
1354 #endif
1355         xx -= (dx0 + hgpar_->waferPosX_[k]);
1356         yy -= (dy0 + hgpar_->waferPosY_[k]);
1357         break;
1358       }
1359     }
1360   }
1361   if ((std::abs(waferU) <= hgpar_->waferUVMax_) && (celltype >= 0)) {
1362     int place(HGCalCell::cellPlacementOld), part(HGCalTypes::WaferFull);
1363     if (mode_ == HGCalGeometryMode::Hexagon8Cassette) {
1364       int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1365       auto ktr = hgpar_->waferInfoMap_.find(indx);
1366       if (ktr != hgpar_->waferInfoMap_.end()) {
1367         place = HGCalCell::cellPlacementIndex(1, HGCalTypes::layerFrontBack(layertype), (ktr->second).orient);
1368         part = (ktr->second).part;
1369 #ifdef EDM_ML_DEBUG
1370         if (debug)
1371           edm::LogVerbatim("HGCalGeom") << "waferFromPosition: frontback " << layertype << ":"
1372                                         << HGCalTypes::layerFrontBack(layertype) << " Orient " << (ktr->second).orient
1373                                         << " place " << place << " part " << part;
1374 #endif
1375       }
1376     }
1377     cellHex(xx, yy, celltype, place, part, cellU, cellV, extend, debug);
1378     wt = (((celltype < 2) && (hgpar_->useSimWt_ > 0)) ? (hgpar_->cellThickness_[celltype] / hgpar_->waferThick_) : 1.0);
1379   } else {
1380     cellU = cellV = 2 * hgpar_->nCellsFine_;
1381     wt = 1.0;
1382     celltype = -1;
1383   }
1384   if ((celltype < 0) && debug) {
1385     double x1(xx);
1386     double y1(yy);
1387     edm::LogVerbatim("HGCalGeom") << "waferfFromPosition: Bad type for X " << x << ":" << x1 << ":" << xx << " Y " << y
1388                                   << ":" << y1 << ":" << yy << " Wafer " << waferU << ":" << waferV << " Cell " << cellU
1389                                   << ":" << cellV;
1390     for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
1391       double dx = std::abs(x1 - hgpar_->waferPosX_[k]);
1392       double dy = std::abs(y1 - hgpar_->waferPosY_[k]);
1393       edm::LogVerbatim("HGCalGeom") << "Wafer [" << k << "] Position (" << hgpar_->waferPosX_[k] << ", "
1394                                     << hgpar_->waferPosY_[k] << ") difference " << dx << ":" << dy << ":"
1395                                     << dx * tan30deg_ << ":" << hexside - dy << " Paramerers " << rmax << ":"
1396                                     << hexside;
1397     }
1398   }
1399   edm::LogVerbatim("HGCalGeomX") << "Input x:y:layer " << x << ":" << y << ":" << layer << " Wafer " << waferU << ":"
1400                                  << waferV << " Cell " << cellU << ":" << cellV << ":" << celltype << " wt " << wt;
1401 }
1402 
1403 bool HGCalDDDConstants::waferInLayer(int wafer, int lay, bool reco) const {
1404   const auto& indx = getIndex(lay, reco);
1405   if (indx.first < 0)
1406     return false;
1407   return waferInLayerTest(wafer, indx.first, hgpar_->defineFull_);
1408 }
1409 
1410 bool HGCalDDDConstants::waferFullInLayer(int wafer, int lay, bool reco) const {
1411   const auto& indx = getIndex(lay, reco);
1412   if (indx.first < 0)
1413     return false;
1414   return waferInLayerTest(wafer, indx.first, false);
1415 }
1416 
1417 HGCalParameters::waferInfo HGCalDDDConstants::waferInfo(int lay, int waferU, int waferV) const {
1418   int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
1419   auto itr = hgpar_->waferInfoMap_.find(indx);
1420   return ((itr == hgpar_->waferInfoMap_.end()) ? HGCalParameters::waferInfo() : itr->second);
1421 }
1422 
1423 std::pair<double, double> HGCalDDDConstants::waferParameters(bool reco) const {
1424   if (reco)
1425     return std::make_pair(rmax_, hexside_);
1426   else
1427     return std::make_pair(HGCalParameters::k_ScaleToDDD * rmax_, HGCalParameters::k_ScaleToDDD * hexside_);
1428 }
1429 
1430 std::pair<double, double> HGCalDDDConstants::waferPosition(int wafer, bool reco) const {
1431   double xx(0), yy(0);
1432   if (wafer >= 0 && wafer < static_cast<int>(hgpar_->waferPosX_.size())) {
1433     xx = hgpar_->waferPosX_[wafer];
1434     yy = hgpar_->waferPosY_[wafer];
1435   }
1436   if (!reco) {
1437     xx *= HGCalParameters::k_ScaleToDDD;
1438     yy *= HGCalParameters::k_ScaleToDDD;
1439   }
1440   return std::make_pair(xx, yy);
1441 }
1442 
1443 std::pair<double, double> HGCalDDDConstants::waferPosition(
1444     int lay, int waferU, int waferV, bool reco, bool debug) const {
1445   int ll = lay - hgpar_->firstLayer_;
1446   bool rotx = ((!hgpar_->layerType_.empty()) && (hgpar_->layerType_[ll] == HGCalTypes::WaferCenterR));
1447 #ifdef EDM_ML_DEBUG
1448   if (debug)
1449     edm::LogVerbatim("HGCalGeom") << "Layer " << lay << ":" << ll << " Rotation " << rotx << " U:V " << waferU << ":"
1450                                   << waferV;
1451 #endif
1452   auto xy = waferPositionNoRot(lay, waferU, waferV, reco, debug);
1453   std::pair<double, double> xy0 = (rotx) ? getXY(lay, xy.first, xy.second, false) : xy;
1454 #ifdef EDM_ML_DEBUG
1455   if (debug)
1456     edm::LogVerbatim("HGCalGeom") << "Without and with rotation " << xy.first << ":" << xy.second << ":" << xy0.first
1457                                   << ":" << xy0.second;
1458 #endif
1459   return xy0;
1460 }
1461 
1462 int HGCalDDDConstants::waferFileIndex(unsigned int kk) const {
1463   if (kk < hgpar_->waferInfoMap_.size()) {
1464     auto itr = hgpar_->waferInfoMap_.begin();
1465     std::advance(itr, kk);
1466     return itr->first;
1467   } else
1468     return 0;
1469 }
1470 
1471 std::tuple<int, int, int> HGCalDDDConstants::waferFileInfo(unsigned int kk) const {
1472   if (kk < hgpar_->waferInfoMap_.size()) {
1473     auto itr = hgpar_->waferInfoMap_.begin();
1474     std::advance(itr, kk);
1475     return std::make_tuple(itr->second.type, itr->second.part, itr->second.orient);
1476   } else
1477     return std::make_tuple(0, 0, 0);
1478 }
1479 
1480 std::tuple<int, int, int> HGCalDDDConstants::waferFileInfoFromIndex(int kk) const {
1481   auto itr = hgpar_->waferInfoMap_.find(kk);
1482   if (itr != hgpar_->waferInfoMap_.end()) {
1483     return std::make_tuple(itr->second.type, itr->second.part, itr->second.orient);
1484   } else
1485     return std::make_tuple(0, 0, 0);
1486 }
1487 
1488 GlobalPoint HGCalDDDConstants::waferLocal2Global(
1489     HepGeom::Point3D<float>& loc, const DetId& id, bool useWafer, bool reco, bool debug) const {
1490   HGCSiliconDetId detid(id);
1491   double x(0), y(0);
1492   if (useWafer) {
1493     auto xyw = waferPositionNoRot(detid.layer(), detid.waferU(), detid.waferV(), reco, debug);
1494     x = xyw.first;
1495     y = xyw.second;
1496   }
1497   auto xy = getXY(detid.layer(), (x + loc.x()), (y + loc.y()), false);
1498   double zz = (detid.zside() < 0) ? -(loc.z() + waferZ(detid.layer(), reco)) : (loc.z() + waferZ(detid.layer(), reco));
1499   double xx = (detid.zside() < 0) ? -xy.first : xy.first;
1500   return GlobalPoint(xx, xy.second, zz);
1501 }
1502 
1503 int HGCalDDDConstants::wafers() const {
1504   int wafer(0);
1505   if (!tileTrapezoid()) {
1506     for (unsigned int i = 0; i < layers(true); ++i) {
1507       int lay = hgpar_->depth_[i];
1508       wafer += modules(lay, true);
1509     }
1510   } else {
1511     wafer = static_cast<int>(hgpar_->moduleLayR_.size());
1512   }
1513   return wafer;
1514 }
1515 
1516 int HGCalDDDConstants::wafers(int layer, int type) const {
1517   int wafer(0);
1518   if (!tileTrapezoid()) {
1519     auto itr = waferLayer_.find(layer);
1520     if (itr != waferLayer_.end()) {
1521       unsigned ity = (type > 0 && type <= 2) ? type : 0;
1522       wafer = (itr->second)[ity];
1523     }
1524   } else {
1525     const auto& index = getIndex(layer, true);
1526     wafer = 1 + hgpar_->lastModule_[index.first] - hgpar_->firstModule_[index.first];
1527   }
1528   return wafer;
1529 }
1530 
1531 int HGCalDDDConstants::waferType(DetId const& id, bool fromFile) const {
1532   int type(1);
1533   if (waferHexagon8()) {
1534     if (fromFile && (waferFileSize() > 0)) {
1535       int layer(0), waferU(0), waferV(0);
1536       if (id.det() != DetId::Forward) {
1537         HGCSiliconDetId hid(id);
1538         layer = hid.layer();
1539         waferU = hid.waferU();
1540         waferV = hid.waferV();
1541       } else {
1542         HFNoseDetId hid(id);
1543         layer = hid.layer();
1544         waferU = hid.waferU();
1545         waferV = hid.waferV();
1546       }
1547       auto itr = hgpar_->waferInfoMap_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1548       if (itr != hgpar_->waferInfoMap_.end())
1549         type = (itr->second).type;
1550     } else {
1551       type = ((id.det() != DetId::Forward) ? HGCSiliconDetId(id).type() : HFNoseDetId(id).type());
1552     }
1553   } else if (waferHexagon6()) {
1554     type = waferTypeL(HGCalDetId(id).wafer()) - 1;
1555   }
1556   return type;
1557 }
1558 
1559 int HGCalDDDConstants::waferType(int layer, int waferU, int waferV, bool fromFile) const {
1560   int type(HGCSiliconDetId::HGCalCoarseThick);
1561   if (waferHexagon8()) {
1562     if (fromFile && (waferFileSize() > 0)) {
1563       auto itr = hgpar_->waferInfoMap_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1564       if (itr != hgpar_->waferInfoMap_.end())
1565         type = (itr->second).type;
1566     } else {
1567       auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1568       if (itr != hgpar_->typesInLayers_.end())
1569         type = hgpar_->waferTypeL_[itr->second];
1570     }
1571   } else if (waferHexagon6()) {
1572     if ((waferU >= 0) && (waferU < static_cast<int>(hgpar_->waferTypeL_.size())))
1573       type = (hgpar_->waferTypeL_[waferU] - 1);
1574   }
1575   return type;
1576 }
1577 
1578 std::tuple<int, int, int> HGCalDDDConstants::waferType(HGCSiliconDetId const& id, bool fromFile) const {
1579   const auto& index = HGCalWaferIndex::waferIndex(id.layer(), id.waferU(), id.waferV());
1580   int type(-1), part(-1), orient(-1);
1581   if (fromFile && (waferFileSize() > 0)) {
1582     auto itr = hgpar_->waferInfoMap_.find(index);
1583     if (itr != hgpar_->waferInfoMap_.end()) {
1584       type = (itr->second).type;
1585       part = (itr->second).part;
1586       orient = (itr->second).orient;
1587     }
1588   } else {
1589     auto ktr = hgpar_->typesInLayers_.find(index);
1590     if (ktr != hgpar_->typesInLayers_.end())
1591       type = hgpar_->waferTypeL_[ktr->second];
1592     auto itr = hgpar_->waferTypes_.find(index);
1593     if (itr != hgpar_->waferTypes_.end()) {
1594       if ((itr->second).second < HGCalTypes::k_OffsetRotation) {
1595         orient = (itr->second).second;
1596         if ((itr->second).first == HGCalGeomTools::k_allCorners) {
1597           part = HGCalTypes::WaferFull;
1598         } else if ((itr->second).first == HGCalGeomTools::k_fiveCorners) {
1599           part = HGCalTypes::WaferFive;
1600         } else if ((itr->second).first == HGCalGeomTools::k_fourCorners) {
1601           part = HGCalTypes::WaferHalf;
1602         } else if ((itr->second).first == HGCalGeomTools::k_threeCorners) {
1603           part = HGCalTypes::WaferThree;
1604         }
1605       } else {
1606         part = (itr->second).first;
1607         orient = ((itr->second).second - HGCalTypes::k_OffsetRotation);
1608       }
1609     } else {
1610       part = HGCalTypes::WaferFull;
1611       orient = 0;
1612     }
1613   }
1614   return std::make_tuple(type, part, orient);
1615 }
1616 
1617 std::pair<int, int> HGCalDDDConstants::waferTypeRotation(
1618     int layer, int waferU, int waferV, bool fromFile, bool debug) const {
1619   int type(HGCalTypes::WaferOut), rotn(0);
1620   int wl = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1621   bool withinList(true);
1622   if (fromFile && (waferFileSize() > 0)) {
1623     auto itr = hgpar_->waferInfoMap_.find(wl);
1624     withinList = (itr != hgpar_->waferInfoMap_.end());
1625     if (withinList) {
1626       type = (itr->second).part;
1627       rotn = (itr->second).orient;
1628     }
1629   } else {
1630     auto itr = hgpar_->waferTypes_.find(wl);
1631     if (waferHexagon8()) {
1632       withinList = (itr != hgpar_->waferTypes_.end());
1633       if (withinList) {
1634         if ((itr->second).second < HGCalTypes::k_OffsetRotation) {
1635           rotn = (itr->second).second;
1636           if ((itr->second).first == HGCalGeomTools::k_allCorners) {
1637             type = HGCalTypes::WaferFull;
1638           } else if ((itr->second).first == HGCalGeomTools::k_fiveCorners) {
1639             type = HGCalTypes::WaferFive;
1640           } else if ((itr->second).first == HGCalGeomTools::k_fourCorners) {
1641             type = HGCalTypes::WaferHalf;
1642           } else if ((itr->second).first == HGCalGeomTools::k_threeCorners) {
1643             type = HGCalTypes::WaferThree;
1644           }
1645         } else {
1646           type = (itr->second).first;
1647           rotn = ((itr->second).second - HGCalTypes::k_OffsetRotation);
1648         }
1649       } else {
1650         type = HGCalTypes::WaferFull;
1651         rotn = HGCalTypes::WaferCorner0;
1652       }
1653     }
1654   }
1655 #ifdef EDM_ML_DEBUG
1656   if (debug)
1657     edm::LogVerbatim("HGCalGeom") << "waferTypeRotation: Layer " << layer << " Wafer " << waferU << ":" << waferV
1658                                   << " Index " << std::hex << wl << std::dec << ":" << withinList << " Type " << type
1659                                   << " Rotation " << rotn;
1660 #endif
1661   return std::make_pair(type, rotn);
1662 }
1663 
1664 bool HGCalDDDConstants::waferVirtual(int layer, int waferU, int waferV) const {
1665   bool type(false);
1666   if (waferHexagon8()) {
1667     int wl = HGCalWaferIndex::waferIndex(layer, waferU, waferV, false);
1668     type = (hgpar_->waferTypes_.find(wl) != hgpar_->waferTypes_.end());
1669   } else if (waferHexagon6()) {
1670     int wl = HGCalWaferIndex::waferIndex(layer, waferU, 0, true);
1671     type = (hgpar_->waferTypes_.find(wl) != hgpar_->waferTypes_.end());
1672   }
1673   return type;
1674 }
1675 
1676 double HGCalDDDConstants::waferZ(int lay, bool reco) const {
1677   const auto& index = getIndex(lay, reco);
1678   if (index.first < 0)
1679     return 0;
1680   else
1681     return (reco ? hgpar_->zLayerHex_[index.first] : HGCalParameters::k_ScaleToDDD * hgpar_->zLayerHex_[index.first]);
1682 }
1683 
1684 int HGCalDDDConstants::cellHex(
1685     double xx, double yy, const double& cellR, const std::vector<double>& posX, const std::vector<double>& posY) const {
1686   int num(0);
1687   const double tol(0.00001);
1688   double cellY = 2.0 * cellR * tan30deg_;
1689   for (unsigned int k = 0; k < posX.size(); ++k) {
1690     double dx = std::abs(xx - posX[k]);
1691     double dy = std::abs(yy - posY[k]);
1692     if (dx <= (cellR + tol) && dy <= (cellY + tol)) {
1693       double xmax = (dy <= 0.5 * cellY) ? cellR : (cellR - (dy - 0.5 * cellY) / tan30deg_);
1694       if (dx <= (xmax + tol)) {
1695         num = k;
1696         break;
1697       }
1698     }
1699   }
1700   return num;
1701 }
1702 
1703 void HGCalDDDConstants::cellHex(
1704     double xloc, double yloc, int cellType, int place, int part, int& cellU, int& cellV, bool extend, bool debug) const {
1705   if (mode_ == HGCalGeometryMode::Hexagon8Cassette) {
1706     auto uv = (part == HGCalTypes::WaferFull)
1707                   ? hgcellUV_->cellUVFromXY3(xloc, yloc, place, cellType, true, debug)
1708                   : hgcellUV_->cellUVFromXY1(xloc, yloc, place, cellType, part, true, debug);
1709     cellU = uv.first;
1710     cellV = uv.second;
1711   } else if (waferHexagon8File()) {
1712     auto uv = hgcellUV_->cellUVFromXY3(xloc, yloc, place, cellType, extend, debug);
1713     cellU = uv.first;
1714     cellV = uv.second;
1715   } else {
1716     int ncell = (cellType == 0) ? hgpar_->nCellsFine_ : hgpar_->nCellsCoarse_;
1717     double delY = 2 * rmax_ / (3 * ncell);
1718     double delX = 0.5 * delY * sqrt3_;
1719     double delYT = (extend) ? (2 * rmaxT_ / (3 * ncell)) : delY;
1720     double delXT = 0.5 * delYT * sqrt3_;
1721     double v0 = ((xloc / delY - 1.0) / 1.5);
1722     int cv0 = (v0 > 0) ? (ncell + static_cast<int>(v0 + 0.5)) : (ncell - static_cast<int>(-v0 + 0.5));
1723     double u0 = (0.5 * yloc / delX + 0.5 * cv0);
1724     int cu0 = (u0 > 0) ? (ncell / 2 + static_cast<int>(u0 + 0.5)) : (ncell / 2 - static_cast<int>(-u0 + 0.5));
1725     cu0 = std::max(0, std::min(cu0, 2 * ncell - 1));
1726     cv0 = std::max(0, std::min(cv0, 2 * ncell - 1));
1727     if (cv0 - cu0 >= ncell)
1728       cv0 = cu0 + ncell - 1;
1729     if (debug)
1730       edm::LogVerbatim("HGCalGeom") << "cellHex: input " << xloc << ":" << yloc << ":" << cellType << " parameter "
1731                                     << delX << ":" << delY << " u0 " << u0 << ":" << cu0 << " v0 " << v0 << ":" << cv0;
1732     bool found(false);
1733     static constexpr int shift[3] = {0, 1, -1};
1734     for (int i1 = 0; i1 < 3; ++i1) {
1735       cellU = cu0 + shift[i1];
1736       for (int i2 = 0; i2 < 3; ++i2) {
1737         cellV = cv0 + shift[i2];
1738         if (((cellV - cellU) < ncell) && ((cellU - cellV) <= ncell) && (cellU >= 0) && (cellV >= 0) &&
1739             (cellU < 2 * ncell) && (cellV < 2 * ncell)) {
1740           double xc = (1.5 * (cellV - ncell) + 1.0) * delY;
1741           double yc = (2 * cellU - cellV - ncell) * delX;
1742           if ((std::abs(yloc - yc) <= delXT) && (std::abs(xloc - xc) <= delYT) &&
1743               ((std::abs(xloc - xc) <= 0.5 * delYT) ||
1744                (std::abs(yloc - yc) <= sqrt3_ * (delYT - std::abs(xloc - xc))))) {
1745             if (debug)
1746               edm::LogVerbatim("HGCalGeom")
1747                   << "cellHex: local " << xc << ":" << yc << " difference " << std::abs(xloc - xc) << ":"
1748                   << std::abs(yloc - yc) << ":" << sqrt3_ * (delY - std::abs(yloc - yc)) << " comparator " << delX
1749                   << ":" << delY << " (u,v) = (" << cellU << "," << cellV << ")";
1750             found = true;
1751             break;
1752           }
1753         }
1754       }
1755       if (found)
1756         break;
1757     }
1758     if (!found) {
1759       cellU = cu0;
1760       cellV = cv0;
1761     }
1762   }
1763 }
1764 
1765 std::pair<int, float> HGCalDDDConstants::getIndex(int lay, bool reco) const {
1766   int indx = layerIndex(lay, reco);
1767   if (indx < 0)
1768     return std::make_pair(-1, 0);
1769   float cell(0);
1770   if (waferHexagon6()) {
1771     cell = (reco ? hgpar_->moduleCellR_[0] : hgpar_->moduleCellS_[0]);
1772   } else {
1773     if (waferHexagon8()) {
1774       cell = (reco ? hgpar_->moduleCellR_[0] : hgpar_->moduleCellS_[0]);
1775     } else {
1776       cell = hgpar_->scintCellSize(lay);
1777     }
1778   }
1779   return std::make_pair(indx, cell);
1780 }
1781 
1782 int HGCalDDDConstants::layerFromIndex(int index, bool reco) const {
1783   int ll(-1);
1784   if (waferHexagon6() && reco) {
1785     ll = static_cast<int>(std::find(hgpar_->depthLayerF_.begin(), hgpar_->depthLayerF_.end(), index) -
1786                           hgpar_->depthLayerF_.begin());
1787     if (ll == static_cast<int>(hgpar_->depthLayerF_.size()))
1788       ll = -1;
1789   } else {
1790     ll = static_cast<int>(std::find(hgpar_->layerIndex_.begin(), hgpar_->layerIndex_.end(), index) -
1791                           hgpar_->layerIndex_.begin());
1792     if (ll == static_cast<int>(hgpar_->layerIndex_.size()))
1793       ll = -1;
1794   }
1795 #ifdef EDM_ML_DEBUG
1796   edm::LogVerbatim("HGCalGeom") << "LayerFromIndex for " << index << ":" << reco << ":" << waferHexagon6() << " is"
1797                                 << ll << ":" << (ll + hgpar_->firstLayer_);
1798 #endif
1799   return ((ll < 0) ? ll : (ll + hgpar_->firstLayer_));
1800 }
1801 
1802 bool HGCalDDDConstants::isValidCell(int lay, int wafer, int cell) const {
1803   // Calculate the position of the cell
1804   // Works for options HGCalHexagon/HGCalHexagonFull
1805   double x = hgpar_->waferPosX_[wafer];
1806   double y = hgpar_->waferPosY_[wafer];
1807   if (hgpar_->waferTypeT_[wafer] - 1 == HGCSiliconDetId::HGCalFine) {
1808     x += hgpar_->cellFineX_[cell];
1809     y += hgpar_->cellFineY_[cell];
1810   } else {
1811     x += hgpar_->cellCoarseX_[cell];
1812     y += hgpar_->cellCoarseY_[cell];
1813   }
1814   double rr = sqrt(x * x + y * y);
1815   bool result = ((rr >= hgpar_->rMinLayHex_[lay - 1]) && (rr <= hgpar_->rMaxLayHex_[lay - 1]) &&
1816                  (wafer < static_cast<int>(hgpar_->waferPosX_.size())));
1817 #ifdef EDM_ML_DEBUG
1818   if (!result)
1819     edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << wafer << ":" << cell << " Position " << x << ":" << y
1820                                   << ":" << rr << " Compare Limits " << hgpar_->rMinLayHex_[lay - 1] << ":"
1821                                   << hgpar_->rMaxLayHex_[lay - 1] << " Flag " << result;
1822 #endif
1823   return result;
1824 }
1825 
1826 bool HGCalDDDConstants::isValidCell8(int lay, int waferU, int waferV, int cellU, int cellV, int type) const {
1827   float x(0), y(0);
1828   int kndx = cellV * 100 + cellU;
1829   if (type == 0) {
1830     auto ktr = hgpar_->cellFineIndex_.find(kndx);
1831     if (ktr != hgpar_->cellFineIndex_.end()) {
1832       x = hgpar_->cellFineX_[ktr->second];
1833       y = hgpar_->cellFineY_[ktr->second];
1834     }
1835 #ifdef EDM_ML_DEBUG
1836     edm::LogVerbatim("HGCalGeom") << "Fine " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
1837                                   << (ktr != hgpar_->cellFineIndex_.end());
1838 #endif
1839   } else {
1840     auto ktr = hgpar_->cellCoarseIndex_.find(kndx);
1841     if (ktr != hgpar_->cellCoarseIndex_.end()) {
1842       x = hgpar_->cellCoarseX_[ktr->second];
1843       y = hgpar_->cellCoarseY_[ktr->second];
1844     }
1845 #ifdef EDM_ML_DEBUG
1846     edm::LogVerbatim("HGCalGeom") << "Coarse " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
1847                                   << (ktr != hgpar_->cellCoarseIndex_.end());
1848 #endif
1849   }
1850   const auto& xy = waferPositionNoRot(lay, waferU, waferV, true, false);
1851   x += xy.first;
1852   y += xy.second;
1853 #ifdef EDM_ML_DEBUG
1854   edm::LogVerbatim("HGCalGeom") << "With wafer (" << waferU << "," << waferV << ") " << x << ":" << y;
1855 #endif
1856   double rr = sqrt(x * x + y * y);
1857   int ll = lay - hgpar_->firstLayer_;
1858   bool result = ((rr >= hgpar_->rMinLayHex_[ll]) && (rr <= hgpar_->rMaxLayHex_[ll]));
1859 #ifdef EDM_ML_DEBUG
1860   edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << ll << ":" << waferU << ":" << waferV << ":" << cellU << ":"
1861                                 << cellV << " Position " << x << ":" << y << ":" << rr << " Compare Limits "
1862                                 << hgpar_->rMinLayHex_[ll] << ":" << hgpar_->rMaxLayHex_[ll] << " Flag " << result;
1863 #endif
1864   if (result && waferHexagon8File()) {
1865     int N = (type == 0) ? hgpar_->nCellsFine_ : hgpar_->nCellsCoarse_;
1866     auto partn = waferTypeRotation(lay, waferU, waferV, false, false);
1867     result = HGCalWaferMask::goodCell(cellU, cellV, N, partn.first, partn.second);
1868 #ifdef EDM_ML_DEBUG
1869     edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << waferU << ":" << waferV << ":" << cellU << ":" << cellV
1870                                   << " N " << N << " part " << partn.first << ":" << partn.second << " Result "
1871                                   << result;
1872 #endif
1873   }
1874   return result;
1875 }
1876 
1877 int32_t HGCalDDDConstants::waferIndex(int wafer, int index) const {
1878   int layer = layerFromIndex(index, true);
1879   int waferU = HGCalWaferIndex::waferU(hgpar_->waferCopy_[wafer]);
1880   int waferV = HGCalWaferIndex::waferV(hgpar_->waferCopy_[wafer]);
1881   int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1882 #ifdef EDM_ML_DEBUG
1883   edm::LogVerbatim("HGCalGeom") << "WaferIndex for " << wafer << ":" << index << " (" << layer << ":" << waferU << ":"
1884                                 << waferV << ") " << indx;
1885 #endif
1886   return indx;
1887 }
1888 
1889 bool HGCalDDDConstants::waferInLayerTest(int wafer, int lay, bool full) const {
1890   bool in = (waferHexagon6()) ? true : false;
1891   if (!in) {
1892     double xpos = hgpar_->waferPosX_[wafer] + hgpar_->xLayerHex_[lay];
1893     double ypos = hgpar_->waferPosY_[wafer] + hgpar_->yLayerHex_[lay];
1894     std::pair<int, int> corner = HGCalGeomTools::waferCorner(
1895         xpos, ypos, rmax_, hexside_, hgpar_->rMinLayHex_[lay], hgpar_->rMaxLayHex_[lay], in);
1896     in = (full ? (corner.first > 0) : (corner.first == static_cast<int>(HGCalParameters::k_CornerSize)));
1897     if (in && fullAndPart_) {
1898       int indx = waferIndex(wafer, lay);
1899       in = (hgpar_->waferInfoMap_.find(indx) != hgpar_->waferInfoMap_.end());
1900 #ifdef EDM_ML_DEBUG
1901       if (!in)
1902         edm::LogVerbatim("HGCalGeom") << "WaferInLayerTest: Layer " << lay << " wafer " << wafer << " index " << indx
1903                                       << "( " << HGCalWaferIndex::waferLayer(indx) << ", "
1904                                       << HGCalWaferIndex::waferU(indx) << ", " << HGCalWaferIndex::waferV(indx)
1905                                       << ") in " << in;
1906 #endif
1907     }
1908 #ifdef EDM_ML_DEBUG
1909     edm::LogVerbatim("HGCalGeom") << "WaferInLayerTest: Layer " << lay << " wafer " << wafer << " R-limits "
1910                                   << hgpar_->rMinLayHex_[lay] << ":" << hgpar_->rMaxLayHex_[lay] << " Corners "
1911                                   << corner.first << ":" << corner.second << " In " << in;
1912 #endif
1913   }
1914   return in;
1915 }
1916 
1917 std::pair<double, double> HGCalDDDConstants::waferPositionNoRot(
1918     int lay, int waferU, int waferV, bool reco, bool debug) const {
1919   int ll = lay - hgpar_->firstLayer_;
1920   double x = hgpar_->xLayerHex_[ll];
1921   double y = hgpar_->yLayerHex_[ll];
1922 #ifdef EDM_ML_DEBUG
1923   if (debug)
1924     edm::LogVerbatim("HGCalGeom") << "Layer " << lay << ":" << ll << " Shift " << hgpar_->xLayerHex_[ll] << ":"
1925                                   << hgpar_->yLayerHex_[ll] << " U:V " << waferU << ":" << waferV;
1926 #endif
1927   if (!reco) {
1928     x *= HGCalParameters::k_ScaleToDDD;
1929     y *= HGCalParameters::k_ScaleToDDD;
1930   }
1931   const auto& xy = waferPosition(waferU, waferV, reco);
1932   x += xy.first;
1933   y += xy.second;
1934 #ifdef EDM_ML_DEBUG
1935   if (debug)
1936     edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << ":" << xy.first << ":" << xy.second;
1937 #endif
1938   return std::make_pair(x, y);
1939 }
1940 
1941 std::pair<double, double> HGCalDDDConstants::waferPosition(int waferU, int waferV, bool reco) const {
1942   double xx(0), yy(0);
1943   int indx = HGCalWaferIndex::waferIndex(0, waferU, waferV);
1944   auto itr = hgpar_->wafersInLayers_.find(indx);
1945   if (itr != hgpar_->wafersInLayers_.end()) {
1946     xx = hgpar_->waferPosX_[itr->second];
1947     yy = hgpar_->waferPosY_[itr->second];
1948   }
1949   if (!reco) {
1950     xx *= HGCalParameters::k_ScaleToDDD;
1951     yy *= HGCalParameters::k_ScaleToDDD;
1952   }
1953   return std::make_pair(xx, yy);
1954 }
1955 
1956 #include "FWCore/Utilities/interface/typelookup.h"
1957 
1958 TYPELOOKUP_DATA_REG(HGCalDDDConstants);