Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-01-25 02:56:24

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