File indexing completed on 2025-04-26 02:01:03
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
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") << "HGCalDDDConstants::Mode " << mode_ << " FullAndPart " << fullAndPart_
0029 << " waferHex6:waverHex8 " << waferHexagon6() << ":" << waferHexagon8() << " cassettte "
0030 << cassetteMode();
0031 #endif
0032 if (waferHexagon6() || waferHexagon8()) {
0033 rmax_ = (HGCalParameters::k_ScaleFromDDD * (hgpar_->waferR_) * std::cos(30._deg));
0034 rmaxT_ = rmax_ + 0.5 * hgpar_->sensorSeparation_;
0035 hexside_ = 2.0 * rmax_ * tan30deg_;
0036 hexsideT_ = 2.0 * rmaxT_ * tan30deg_;
0037 hgcell_ = std::make_unique<HGCalCell>(2.0 * rmaxT_, hgpar_->nCellsFine_, hgpar_->nCellsCoarse_);
0038 hgcellUV_ = std::make_unique<HGCalCellUV>(
0039 2.0 * rmax_, hgpar_->sensorSeparation_, hgpar_->nCellsFine_, hgpar_->nCellsCoarse_);
0040 cellOffset_ = std::make_unique<HGCalCellOffset>(hgpar_->waferSize_,
0041 hgpar_->nCellsFine_,
0042 hgpar_->nCellsCoarse_,
0043 hgpar_->guardRingOffset_,
0044 hgpar_->mouseBite_,
0045 hgpar_->sensorSizeOffset_);
0046 #ifdef EDM_ML_DEBUG
0047 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::rmax_ " << rmax_ << ":" << rmaxT_ << ":" << hexside_ << ":"
0048 << hexsideT_ << " CellSize "
0049 << 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[0] << ":"
0050 << 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[1];
0051 #endif
0052 } else {
0053 hgcell_.reset();
0054 hgcellUV_.reset();
0055 cellOffset_.reset();
0056 }
0057 if (cassetteMode()) {
0058 if (mode_ == HGCalGeometryMode::TrapezoidFineCell)
0059 hgcassette_.setParameter(hgpar_->cassettes_, hgpar_->cassetteShiftTile_, true);
0060 else
0061 hgcassette_.setParameter(hgpar_->cassettes_, hgpar_->cassetteShift_, true);
0062 #ifdef EDM_ML_DEBUG
0063 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::Setup HGCalCassette for " << hgpar_->cassettes_
0064 << " cassettes";
0065 #endif
0066 }
0067
0068 modHalf_ = 0;
0069 maxWafersPerLayer_ = 0;
0070 for (int simreco = 0; simreco < 2; ++simreco) {
0071 tot_layers_[simreco] = layersInit((bool)simreco);
0072 max_modules_layer_[simreco].resize(tot_layers_[simreco] + 1);
0073 for (unsigned int layer = 1; layer <= tot_layers_[simreco]; ++layer) {
0074 max_modules_layer_[simreco][layer] = modulesInit(layer, (bool)simreco);
0075 if (simreco == 1) {
0076 modHalf_ += max_modules_layer_[simreco][layer];
0077 maxWafersPerLayer_ = std::max(maxWafersPerLayer_, max_modules_layer_[simreco][layer]);
0078 #ifdef EDM_ML_DEBUG
0079 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::Layer " << layer << " with "
0080 << max_modules_layer_[simreco][layer] << ":" << modHalf_ << " modules in RECO";
0081 } else {
0082 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::Layer " << layer << " with "
0083 << max_modules_layer_[simreco][layer] << " modules in SIM";
0084 #endif
0085 }
0086 }
0087 #ifdef EDM_ML_DEBUG
0088 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::SimReco " << simreco << " with " << tot_layers_[simreco]
0089 << " Layers";
0090 #endif
0091 }
0092 tot_wafers_ = wafers();
0093
0094 #ifdef EDM_ML_DEBUG
0095 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants initialized for " << name << " with " << layers(false) << ":"
0096 << layers(true) << " layers, " << wafers() << ":" << 2 * modHalf_
0097 << " wafers with maximum " << maxWafersPerLayer_ << " per layer and "
0098 << "maximum of " << maxCells(false) << ":" << maxCells(true) << " cells";
0099 #endif
0100 if (waferHexagon6() || waferHexagon8()) {
0101 int wminT(9999999), wmaxT(-9999999), kount1(0), kount2(0);
0102 for (unsigned int i = 0; i < getTrFormN(); ++i) {
0103 int lay0 = getTrForm(i).lay;
0104 int wmin(9999999), wmax(-9999999), kount(0);
0105 for (int wafer = 0; wafer < sectors(); ++wafer) {
0106 bool waferIn = waferInLayer(wafer, lay0, true);
0107 if (waferHexagon8()) {
0108 int kndx = HGCalWaferIndex::waferIndex(lay0,
0109 HGCalWaferIndex::waferU(hgpar_->waferCopy_[wafer]),
0110 HGCalWaferIndex::waferV(hgpar_->waferCopy_[wafer]));
0111 waferIn_[kndx] = waferIn;
0112 }
0113 if (waferIn) {
0114 int waferU = ((waferHexagon6()) ? wafer : HGCalWaferIndex::waferU(hgpar_->waferCopy_[wafer]));
0115 if (waferU < wmin)
0116 wmin = waferU;
0117 if (waferU > wmax)
0118 wmax = waferU;
0119 ++kount;
0120 }
0121 }
0122 if (wminT > wmin)
0123 wminT = wmin;
0124 if (wmaxT < wmax)
0125 wmaxT = wmax;
0126 if (kount1 < kount)
0127 kount1 = kount;
0128 kount2 += kount;
0129 #ifdef EDM_ML_DEBUG
0130 int lay1 = getIndex(lay0, true).first;
0131 edm::LogVerbatim("HGCalGeom") << "Index " << i << " Layer " << lay0 << ":" << lay1 << " Wafer " << wmin << ":"
0132 << wmax << ":" << kount;
0133 #endif
0134 HGCWaferParam a1{{wmin, wmax, kount}};
0135 waferLayer_[lay0] = a1;
0136 }
0137 waferMax_ = std::array<int, 4>{{wminT, wmaxT, kount1, kount2}};
0138 #ifdef EDM_ML_DEBUG
0139 edm::LogVerbatim("HGCalGeom") << "Overall wafer statistics: " << wminT << ":" << wmaxT << ":" << kount1 << ":"
0140 << kount2;
0141 #endif
0142 }
0143 }
0144
0145 std::pair<int, int> HGCalDDDConstants::assignCell(float x, float y, int lay, int subSec, bool reco) const {
0146 const auto& index = getIndex(lay, reco);
0147 if (index.first < 0)
0148 return std::make_pair(-1, -1);
0149 if (waferHexagon6()) {
0150 float xx = (reco) ? x : HGCalParameters::k_ScaleFromDDD * x;
0151 float yy = (reco) ? y : HGCalParameters::k_ScaleFromDDD * y;
0152
0153
0154 int wafer = cellHex(xx, yy, rmax_, hgpar_->waferPosX_, hgpar_->waferPosY_);
0155 if (wafer < 0 || wafer >= static_cast<int>(hgpar_->waferTypeT_.size())) {
0156 edm::LogWarning("HGCalGeom") << "Wafer no. out of bound for " << wafer << ":" << (hgpar_->waferTypeT_).size()
0157 << ":" << (hgpar_->waferPosX_).size() << ":" << (hgpar_->waferPosY_).size()
0158 << " ***** ERROR *****";
0159 return std::make_pair(-1, -1);
0160 } else {
0161
0162 xx -= hgpar_->waferPosX_[wafer];
0163 yy -= hgpar_->waferPosY_[wafer];
0164 if (hgpar_->waferTypeT_[wafer] == 1)
0165 return std::make_pair(wafer,
0166 cellHex(xx,
0167 yy,
0168 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[0],
0169 hgpar_->cellFineX_,
0170 hgpar_->cellFineY_));
0171 else
0172 return std::make_pair(wafer,
0173 cellHex(xx,
0174 yy,
0175 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[1],
0176 hgpar_->cellCoarseX_,
0177 hgpar_->cellCoarseY_));
0178 }
0179 } else {
0180 return std::make_pair(-1, -1);
0181 }
0182 }
0183
0184 std::array<int, 5> HGCalDDDConstants::assignCellHex(
0185 float x, float y, int zside, int lay, bool reco, bool extend, bool debug) const {
0186 int waferU(0), waferV(0), waferType(-1), cellU(0), cellV(0);
0187 if (waferHexagon8()) {
0188 double xx = (reco) ? HGCalParameters::k_ScaleToDDD * x : x;
0189 double yy = (reco) ? HGCalParameters::k_ScaleToDDD * y : y;
0190 double wt(1.0);
0191 #ifdef EDM_ML_DEBUG
0192 edm::LogVerbatim("HGCalGeom") << "assignCellHex x " << x << ":" << xx << " y " << y << ":" << yy << " Lay " << lay;
0193 #endif
0194 waferFromPosition(xx, yy, zside, lay, waferU, waferV, cellU, cellV, waferType, wt, extend, debug);
0195 }
0196 return std::array<int, 5>{{waferU, waferV, waferType, cellU, cellV}};
0197 }
0198
0199 std::array<int, 3> HGCalDDDConstants::assignCellTrap(float x, float y, float z, int layer, bool reco) const {
0200 int irad(-1), iphi(-1), type(-1);
0201 const auto& indx = getIndex(layer, reco);
0202 #ifdef EDM_ML_DEBUG
0203 edm::LogVerbatim("HGCalGeom") << "assignCellTrap: Layer " << layer << ":" << reco << " indx " << indx.first << ":"
0204 << indx.second;
0205 #endif
0206 if (indx.first < 0)
0207 return std::array<int, 3>{{irad, iphi, type}};
0208 int zside = (z > 0) ? 1 : -1;
0209 double xx = (reco) ? (zside * x) : (zside * HGCalParameters::k_ScaleFromDDD * x);
0210 double yy = (reco) ? y : HGCalParameters::k_ScaleFromDDD * y;
0211 int ll = layer - hgpar_->firstLayer_;
0212 xx -= hgpar_->xLayerHex_[ll];
0213 yy -= hgpar_->yLayerHex_[ll];
0214 double phi = (((yy == 0.0) && (xx == 0.0)) ? 0. : std::atan2(yy, xx));
0215 if (phi < 0)
0216 phi += (2.0 * M_PI);
0217 if (indx.second != 0)
0218 iphi = (1 + static_cast<int>(phi / indx.second)) % hgpar_->scintCells(layer);
0219 if (iphi == 0)
0220 iphi = hgpar_->scintCells(layer);
0221 if (cassetteMode()) {
0222 int cassette = HGCalTileIndex::tileCassette(iphi, hgpar_->phiOffset_, hgpar_->nphiCassette_, hgpar_->cassettes_);
0223 auto cshift = hgcassette_.getShift(layer, -1, cassette);
0224 #ifdef EDM_ML_DEBUG
0225 std::ostringstream st1;
0226 st1 << "Cassette " << cassette << " Shift " << cshift.first << ":" << cshift.second << " Original " << xx << ":"
0227 << yy;
0228 #endif
0229 xx += (zside * cshift.first);
0230 yy -= cshift.second;
0231 #ifdef EDM_ML_DEBUG
0232 st1 << " Shifted " << xx << ":" << yy;
0233 edm::LogVerbatim("HGCalGeomT") << st1.str();
0234 #endif
0235 }
0236 type = hgpar_->scintType(layer);
0237 double r = std::sqrt(xx * xx + yy * yy);
0238 auto ir = std::lower_bound(hgpar_->radiusLayer_[type].begin(), hgpar_->radiusLayer_[type].end(), r);
0239 irad = static_cast<int>(ir - hgpar_->radiusLayer_[type].begin());
0240 irad = hgpar_->scintFine(indx.first)
0241 ? (std::clamp(irad, hgpar_->iradMinBHFine_[indx.first], hgpar_->iradMaxBHFine_[indx.first]))
0242 : (std::clamp(irad, hgpar_->iradMinBH_[indx.first], hgpar_->iradMaxBH_[indx.first]));
0243 #ifdef EDM_ML_DEBUG
0244 edm::LogVerbatim("HGCalGeomT") << "assignCellTrap Input " << x << ":" << y << ":" << z << ":" << layer << ":" << reco
0245 << " x|y|r " << xx << ":" << yy << ":" << r << " phi " << phi << ":"
0246 << convertRadToDeg(phi) << " o/p " << irad << ":" << iphi << ":" << type;
0247 #endif
0248 if (!tileExist(zside, layer, irad, iphi)) {
0249 if (tileRingEdge(r, layer, irad)) {
0250 if (std::abs(r - hgpar_->radiusLayer_[type][irad - 1]) < tol_) {
0251 --irad;
0252 if (hgpar_->scintFine(indx.first)) {
0253 if (irad <= hgpar_->iradMinBHFine_[indx.first])
0254 irad = hgpar_->iradMinBHFine_[indx.first];
0255 } else {
0256 if (irad <= hgpar_->iradMinBH_[indx.first])
0257 irad = hgpar_->iradMinBH_[indx.first];
0258 }
0259 } else {
0260 ++irad;
0261 if (hgpar_->scintFine(indx.first)) {
0262 if (irad > hgpar_->iradMaxBHFine_[indx.first])
0263 irad = hgpar_->iradMaxBHFine_[indx.first];
0264 } else {
0265 if (irad > hgpar_->iradMaxBH_[indx.first])
0266 irad = hgpar_->iradMaxBH_[indx.first];
0267 }
0268 }
0269 #ifdef EDM_ML_DEBUG
0270 std::ostringstream st1;
0271 st1 << "assignCellTrap: ring # in ring type " << hgpar_->scintFine(indx.first) << "modified to " << irad;
0272 if (irad > hgpar_->iradMaxBHFine_[indx.first])
0273 st1 << ":" << hgpar_->iradMinBHFine_[indx.first] << ":" << hgpar_->iradMaxBHFine_[indx.first];
0274 else
0275 st1 << ":" << hgpar_->iradMinBH_[indx.first] << ":" << hgpar_->iradMaxBH_[indx.first];
0276 edm::LogVerbatim("HGCalGeomT") << st1.str();
0277 #endif
0278 } else if (tilePhiEdge(phi, layer, iphi)) {
0279 if (std::abs(phi - hgpar_->scintCellSize(layer) * (iphi - 1)) < tol_) {
0280 --iphi;
0281 if (iphi <= 0)
0282 iphi = 1;
0283 } else {
0284 ++iphi;
0285 if (iphi > hgpar_->scintCells(layer))
0286 iphi = 1;
0287 }
0288 #ifdef EDM_ML_DEBUG
0289 edm::LogVerbatim("HGCalGeomT") << "assignCellTrap: iphi # modified to " << iphi << ":"
0290 << hgpar_->scintCells(layer);
0291 #endif
0292 }
0293 }
0294 return std::array<int, 3>{{irad, iphi, type}};
0295 }
0296
0297 bool HGCalDDDConstants::cassetteShiftSilicon(int zside, int layer, int waferU, int waferV) const {
0298 bool shift(false);
0299 if (cassetteMode()) {
0300 int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
0301 auto ktr = hgpar_->waferInfoMap_.find(indx);
0302 if (ktr != hgpar_->waferInfoMap_.end()) {
0303 auto cshift = hgcassette_.getShift(layer, zside, (ktr->second).cassette);
0304 if ((cshift.first != 0) || (cshift.second != 0))
0305 shift = true;
0306 }
0307 }
0308 return shift;
0309 }
0310
0311 bool HGCalDDDConstants::cassetteShiftScintillator(int zside, int layer, int iphi) const {
0312 bool shift(false);
0313 if (cassetteMode()) {
0314 auto cshift = hgcassette_.getShift(layer, zside, cassetteTile(iphi));
0315 if ((cshift.first != 0) || (cshift.second != 0))
0316 shift = true;
0317 }
0318 return shift;
0319 }
0320
0321 double HGCalDDDConstants::cellArea(const HGCSiliconDetId& id, bool reco) const {
0322 double area(0);
0323 int32_t indx = HGCalWaferIndex::waferIndex(id.layer(), id.waferU(), id.waferV());
0324 auto ktr = hgpar_->waferInfoMap_.find(indx);
0325 if (ktr != hgpar_->waferInfoMap_.end()) {
0326 if (ktr->second.part == HGCalTypes::WaferFull) {
0327 area = cellOffset_->cellAreaUV(id.cellU(), id.cellV(), placementIndex(id), ktr->second.type, reco);
0328 } else {
0329 area =
0330 cellOffset_->cellAreaUV(id.cellU(), id.cellV(), placementIndex(id), ktr->second.type, ktr->second.part, reco);
0331 }
0332 }
0333 #ifdef EDM_ML_DEBUG
0334 edm::LogVerbatim("HGCalGeom") << "CellArea: " << id << " Area " << area;
0335 #endif
0336 return area;
0337 }
0338
0339 std::pair<double, double> HGCalDDDConstants::cellEtaPhiTrap(int type, int irad) const {
0340 double dr(0), df(0);
0341 if (tileTrapezoid()) {
0342 double r = 0.5 * ((hgpar_->radiusLayer_[type][irad - 1] + hgpar_->radiusLayer_[type][irad]));
0343 dr = (hgpar_->radiusLayer_[type][irad] - hgpar_->radiusLayer_[type][irad - 1]);
0344 df = r * hgpar_->cellSize_[type];
0345 }
0346 return std::make_pair(dr, df);
0347 }
0348
0349 bool HGCalDDDConstants::cellInLayer(int waferU, int waferV, int cellU, int cellV, int lay, int zside, bool reco) const {
0350 const auto& indx = getIndex(lay, true);
0351 if (indx.first >= 0) {
0352 if (cassetteMode()) {
0353 int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0354 auto ktr = hgpar_->waferInfoMap_.find(indx);
0355 int part = (ktr != hgpar_->waferInfoMap_.end()) ? (ktr->second).part : HGCalTypes::WaferFull;
0356 return HGCalWaferMask::goodCell(cellU, cellV, part);
0357 } else if (mode_ == HGCalGeometryMode::Hexagon8Module) {
0358 int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0359 auto ktr = hgpar_->waferInfoMap_.find(indx);
0360 int thck(HGCalTypes::WaferHD120), part(HGCalTypes::WaferFull), rotn(HGCalTypes::WaferOrient0);
0361 if (ktr != hgpar_->waferInfoMap_.end()) {
0362 thck = (ktr->second).type;
0363 part = (ktr->second).part;
0364 rotn = (ktr->second).orient;
0365 }
0366 int ncell = ((thck == HGCalTypes::WaferHD120) || (thck == HGCalTypes::WaferHD200)) ? hgpar_->nCellsFine_
0367 : hgpar_->nCellsCoarse_;
0368 return HGCalWaferMask::goodCell(cellU, cellV, ncell, part, rotn);
0369 } else if (waferHexagon8() || waferHexagon6()) {
0370 const auto& xy =
0371 ((waferHexagon8()) ? locateCell(zside, lay, waferU, waferV, cellU, cellV, reco, true, false, false, false)
0372 : locateCell(cellU, lay, waferU, reco));
0373 double rpos = sqrt(xy.first * xy.first + xy.second * xy.second);
0374 return ((rpos >= hgpar_->rMinLayHex_[indx.first]) && (rpos <= hgpar_->rMaxLayHex_[indx.first]));
0375 } else {
0376 return true;
0377 }
0378 } else {
0379 return false;
0380 }
0381 }
0382
0383 double HGCalDDDConstants::cellThickness(int layer, int waferU, int waferV) const {
0384 double thick(-1);
0385 int type = waferType(layer, waferU, waferV, false);
0386 if (type >= 0) {
0387 if (waferHexagon8()) {
0388 thick = 10000.0 * hgpar_->cellThickness_[type];
0389 } else if (waferHexagon6()) {
0390 thick = 100.0 * (type + 1);
0391 }
0392 }
0393 return thick;
0394 }
0395
0396 double HGCalDDDConstants::cellSizeHex(int type) const {
0397 int indx = ((waferHexagon8()) ? ((type >= 1) ? 1 : 0) : ((type == 1) ? 1 : 0));
0398 double cell = (tileTrapezoid() ? 0.5 * hgpar_->cellSize_[indx]
0399 : 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[indx]);
0400 return cell;
0401 }
0402
0403 int32_t HGCalDDDConstants::cellType(int type, int cellU, int cellV, int iz, int fwdBack, int orient) const {
0404 int placement = (orient < 0) ? HGCalCell::cellPlacementOld : HGCalCell::cellPlacementIndex(iz, fwdBack, orient);
0405 int ncell = ((type == HGCSiliconDetId::HGCalHD120) || (type == HGCSiliconDetId::HGCalHD200)) ? hgpar_->nCellsFine_
0406 : hgpar_->nCellsCoarse_;
0407 auto cellType = HGCalCell::cellType(cellU, cellV, ncell, placement);
0408 return cellType.first;
0409 }
0410
0411 double HGCalDDDConstants::distFromEdgeHex(double x, double y, double z) const {
0412
0413
0414 if (z < 0)
0415 x = -x;
0416 double dist(0);
0417
0418 double xx = HGCalParameters::k_ScaleFromDDD * x;
0419 double yy = HGCalParameters::k_ScaleFromDDD * y;
0420 if (waferHexagon8()) {
0421 int ll = layerIndex(getLayer(z, false), false);
0422 xx -= hgpar_->xLayerHex_[ll];
0423 yy -= hgpar_->yLayerHex_[ll];
0424 }
0425 int sizew = static_cast<int>(hgpar_->waferPosX_.size());
0426 int wafer = sizew;
0427
0428 for (int k = 0; k < sizew; ++k) {
0429 double dx = std::abs(xx - hgpar_->waferPosX_[k]);
0430 double dy = std::abs(yy - hgpar_->waferPosY_[k]);
0431 if ((dx <= rmax_) && (dy <= hexside_) && ((dy <= 0.5 * hexside_) || (dx * tan30deg_ <= (hexside_ - dy)))) {
0432 wafer = k;
0433 xx -= hgpar_->waferPosX_[k];
0434 yy -= hgpar_->waferPosY_[k];
0435 break;
0436 }
0437 }
0438
0439 if (wafer < sizew) {
0440 if (std::abs(yy) < 0.5 * hexside_) {
0441 dist = rmax_ - std::abs(xx);
0442 } else {
0443 dist = 0.5 * ((rmax_ - std::abs(xx)) - sqrt3_ * (std::abs(yy) - 0.5 * hexside_));
0444 }
0445 } else {
0446 dist = 0;
0447 }
0448 dist *= HGCalParameters::k_ScaleToDDD;
0449 #ifdef EDM_ML_DEBUG
0450 edm::LogVerbatim("HGCalGeom") << "DistFromEdgeHex: Local " << xx << ":" << yy << " wafer " << wafer << " flag "
0451 << (wafer < sizew) << " Distance " << rmax_ << ":" << (rmax_ - std::abs(xx)) << ":"
0452 << (std::abs(yy) - 0.5 * hexside_) << ":" << 0.5 * hexside_ << ":" << dist;
0453 #endif
0454 return dist;
0455 }
0456
0457 double HGCalDDDConstants::distFromEdgeTrap(double x, double y, double z) const {
0458
0459
0460 int lay = getLayer(z, false);
0461 double xx = (z < 0) ? -x : x;
0462 int indx = layerIndex(lay, false);
0463 double r = HGCalParameters::k_ScaleFromDDD * std::sqrt(x * x + y * y);
0464 double phi = (r == 0. ? 0. : std::atan2(y, xx));
0465 if (phi < 0)
0466 phi += (2.0 * M_PI);
0467 int type = hgpar_->scintType(lay);
0468 double cell = hgpar_->scintCellSize(lay);
0469
0470
0471 auto ir = std::lower_bound(hgpar_->radiusLayer_[type].begin(), hgpar_->radiusLayer_[type].end(), r);
0472 int irad = static_cast<int>(ir - hgpar_->radiusLayer_[type].begin());
0473 if (hgpar_->scintFine(indx))
0474 irad = std::clamp(irad, hgpar_->iradMinBHFine_[indx], hgpar_->iradMaxBHFine_[indx]);
0475 else
0476 irad = std::clamp(irad, hgpar_->iradMinBH_[indx], hgpar_->iradMaxBH_[indx]);
0477 int iphi = 1 + static_cast<int>(phi / cell);
0478 double dphi = std::max(0.0, (0.5 * cell - std::abs(phi - (iphi - 0.5) * cell)));
0479 double dist = std::min((r - hgpar_->radiusLayer_[type][irad - 1]), (hgpar_->radiusLayer_[type][irad] - r));
0480 #ifdef EDM_ML_DEBUG
0481 edm::LogVerbatim("HGCalGeom") << "DistFromEdgeTrap: Global " << x << ":" << y << ":" << z << " Layer " << lay
0482 << " Index " << indx << ":" << type << " xx " << xx << " R " << r << ":" << irad << ":"
0483 << hgpar_->radiusLayer_[type][irad - 1] << ":" << hgpar_->radiusLayer_[type][irad]
0484 << " Phi " << phi << ":" << iphi << ":" << (iphi - 0.5) * cell << " cell " << cell
0485 << " Dphi " << dphi << " Dist " << dist << ":" << r * dphi;
0486 #endif
0487 return HGCalParameters::k_ScaleToDDD * std::min(r * dphi, dist);
0488 }
0489
0490 int HGCalDDDConstants::getLayer(double z, bool reco) const {
0491
0492 unsigned int k = 0;
0493 double zz = (reco ? std::abs(z) : HGCalParameters::k_ScaleFromDDD * std::abs(z));
0494 const auto& zLayerHex = hgpar_->zLayerHex_;
0495 auto itr = std::find_if(zLayerHex.begin() + 1, zLayerHex.end(), [&k, &zz, &zLayerHex](double zLayer) {
0496 ++k;
0497 return zz < 0.5 * (zLayerHex[k - 1] + zLayerHex[k]);
0498 });
0499 int lay = (itr == zLayerHex.end()) ? static_cast<int>(zLayerHex.size()) : k;
0500 if (waferHexagon6() && reco) {
0501 int indx = layerIndex(lay, false);
0502 if (indx >= 0)
0503 lay = hgpar_->layerGroupO_[indx];
0504 } else {
0505 lay += (hgpar_->firstLayer_ - 1);
0506 }
0507 return lay;
0508 }
0509
0510 HGCalParameters::hgtrap HGCalDDDConstants::getModule(unsigned int indx, bool hexType, bool reco) const {
0511 HGCalParameters::hgtrap mytr;
0512 if (hexType) {
0513 if (indx >= hgpar_->waferTypeL_.size())
0514 edm::LogWarning("HGCalGeom") << "Wafer no. out bound for index " << indx << ":" << (hgpar_->waferTypeL_).size()
0515 << ":" << (hgpar_->waferPosX_).size() << ":" << (hgpar_->waferPosY_).size()
0516 << " ***** ERROR *****";
0517 unsigned int type =
0518 ((indx < hgpar_->waferTypeL_.size()) ? hgpar_->waferTypeL_[indx] - 1 : HGCSiliconDetId::HGCalLD300);
0519 mytr = hgpar_->getModule(type, reco);
0520 } else {
0521 mytr = hgpar_->getModule(indx, reco);
0522 }
0523 return mytr;
0524 }
0525
0526 std::vector<HGCalParameters::hgtrap> HGCalDDDConstants::getModules() const {
0527 std::vector<HGCalParameters::hgtrap> mytrs;
0528 for (unsigned int k = 0; k < hgpar_->moduleLayR_.size(); ++k)
0529 mytrs.emplace_back(hgpar_->getModule(k, true));
0530 return mytrs;
0531 }
0532
0533 int HGCalDDDConstants::getPhiBins(int lay) const { return (tileTrapezoid() ? hgpar_->scintCells(lay) : 0); }
0534
0535 std::pair<double, double> HGCalDDDConstants::getRangeR(int lay, bool reco) const {
0536 int indx = layerIndex(lay, false);
0537 if ((indx >= 0) && (indx < static_cast<int>(hgpar_->rMinLayHex_.size())))
0538 return std::make_pair(hgpar_->rMinLayHex_[indx], hgpar_->rMaxLayHex_[indx]);
0539 else
0540 return std::make_pair(0, -1.);
0541 }
0542
0543 std::pair<int, int> HGCalDDDConstants::getREtaRange(int lay) const {
0544 int irmin(0), irmax(0);
0545 if (tileTrapezoid()) {
0546 int indx = layerIndex(lay, false);
0547 if (indx >= 0) {
0548 if (hgpar_->scintFine(indx)) {
0549 if (indx < static_cast<int>(hgpar_->iradMinBHFine_.size())) {
0550 irmin = hgpar_->iradMinBHFine_[indx];
0551 irmax = hgpar_->iradMaxBHFine_[indx];
0552 }
0553 } else {
0554 if (indx < static_cast<int>(hgpar_->iradMinBH_.size())) {
0555 irmin = hgpar_->iradMinBH_[indx];
0556 irmax = hgpar_->iradMaxBH_[indx];
0557 }
0558 }
0559 }
0560 }
0561 return std::make_pair(irmin, irmax);
0562 }
0563
0564 std::vector<HGCalParameters::hgtrform> HGCalDDDConstants::getTrForms() const {
0565 std::vector<HGCalParameters::hgtrform> mytrs;
0566 for (unsigned int k = 0; k < hgpar_->trformIndex_.size(); ++k)
0567 mytrs.emplace_back(hgpar_->getTrForm(k));
0568 return mytrs;
0569 }
0570
0571 int HGCalDDDConstants::getTypeTrap(int layer) const {
0572
0573 if (tileTrapezoid()) {
0574 return hgpar_->scintType(layer);
0575 } else {
0576 return -1;
0577 }
0578 }
0579
0580 int HGCalDDDConstants::getTypeHex(int layer, int waferU, int waferV) const {
0581
0582 if (waferHexagon8()) {
0583 auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
0584 return ((itr == hgpar_->typesInLayers_.end() ? 2 : hgpar_->waferTypeL_[itr->second]));
0585 } else {
0586 return -1;
0587 }
0588 }
0589
0590 std::pair<double, double> HGCalDDDConstants::getXY(int layer, double x, double y, bool forwd) const {
0591 int ll = layer - hgpar_->firstLayer_;
0592 double x0(x), y0(y);
0593 if ((!hgpar_->layerType_.empty()) && (ll < static_cast<int>(hgpar_->layerRotV_.size()))) {
0594 if (forwd) {
0595 x0 = x * hgpar_->layerRotV_[ll].first + y * hgpar_->layerRotV_[ll].second;
0596 y0 = y * hgpar_->layerRotV_[ll].first - x * hgpar_->layerRotV_[ll].second;
0597 } else {
0598 x0 = x * hgpar_->layerRotV_[ll].first - y * hgpar_->layerRotV_[ll].second;
0599 y0 = y * hgpar_->layerRotV_[ll].first + x * hgpar_->layerRotV_[ll].second;
0600 }
0601 }
0602 #ifdef EDM_ML_DEBUG
0603 double x1(x0), y1(y0);
0604 if (ll < static_cast<int>(hgpar_->layerRotV_.size())) {
0605 if (forwd) {
0606 x1 = x0 * hgpar_->layerRotV_[ll].first - y0 * hgpar_->layerRotV_[ll].second;
0607 y1 = y0 * hgpar_->layerRotV_[ll].first + x0 * hgpar_->layerRotV_[ll].second;
0608 } else {
0609 x1 = x0 * hgpar_->layerRotV_[ll].first + y0 * hgpar_->layerRotV_[ll].second;
0610 y1 = y0 * hgpar_->layerRotV_[ll].first - x0 * hgpar_->layerRotV_[ll].second;
0611 }
0612 }
0613 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants: Layer " << layer << ":" << ll << " mode " << forwd << " x " << x
0614 << ":" << x0 << ":" << x1 << " y " << y << ":" << y0 << ":" << y1;
0615 #endif
0616 return std::make_pair(x0, y0);
0617 }
0618
0619 double HGCalDDDConstants::guardRingOffset(bool reco) const {
0620 return (reco ? hgpar_->guardRingOffset_ : HGCalParameters::k_ScaleToDDD * hgpar_->guardRingOffset_);
0621 }
0622
0623 bool HGCalDDDConstants::isHalfCell(int waferType, int cell) const {
0624 if (waferType < 1 || cell < 0)
0625 return false;
0626 return ((waferType == HGCalTypes::WaferLD200) || (waferType == HGCalTypes::WaferLD300))
0627 ? hgpar_->cellCoarseHalf_[cell]
0628 : hgpar_->cellFineHalf_[cell];
0629 }
0630
0631 bool HGCalDDDConstants::isValidHex(int lay, int mod, int cell, bool reco) const {
0632
0633 bool result(false), resultMod(false);
0634 int cellmax(0);
0635 if (waferHexagon6()) {
0636 int32_t copyNumber = hgpar_->waferCopy_[mod];
0637 result = ((lay > 0 && lay <= static_cast<int>(layers(reco))));
0638 if (result) {
0639 const int32_t lay_idx = reco ? (lay - 1) * 3 + 1 : lay;
0640 const auto& the_modules = hgpar_->copiesInLayers_[lay_idx];
0641 auto moditr = the_modules.find(copyNumber);
0642 result = resultMod = (moditr != the_modules.end());
0643 #ifdef EDM_ML_DEBUG
0644 if (!result)
0645 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants: Layer " << lay << ":" << lay_idx << " Copy " << copyNumber
0646 << ":" << mod << " Flag " << result;
0647 #endif
0648 if (result) {
0649 if (moditr->second >= 0) {
0650 if (mod >= static_cast<int>(hgpar_->waferTypeT_.size()))
0651 edm::LogWarning("HGCalGeom") << "Module no. out of bound for " << mod << " to be compared with "
0652 << (hgpar_->waferTypeT_).size() << " ***** ERROR *****";
0653 cellmax = (((hgpar_->waferTypeT_[mod] - 1) == HGCSiliconDetId::HGCalHD120) ||
0654 ((hgpar_->waferTypeT_[mod] - 1) == HGCSiliconDetId::HGCalHD200))
0655 ? static_cast<int>(hgpar_->cellFineX_.size())
0656 : static_cast<int>(hgpar_->cellCoarseX_.size());
0657 result = (cell >= 0 && cell <= cellmax);
0658 } else {
0659 result = isValidCell(lay_idx, mod, cell);
0660 }
0661 }
0662 }
0663 }
0664
0665 #ifdef EDM_ML_DEBUG
0666 if (!result)
0667 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants: Layer " << lay << ":"
0668 << (lay > 0 && (lay <= static_cast<int>(layers(reco)))) << " Module " << mod << ":"
0669 << resultMod << " Cell " << cell << ":" << cellmax << ":"
0670 << (cell >= 0 && cell <= cellmax) << ":" << maxCells(reco);
0671 #endif
0672 return result;
0673 }
0674
0675 bool HGCalDDDConstants::isValidHex8(int layer, int modU, int modV, bool fullAndPart) const {
0676
0677 int indx = HGCalWaferIndex::waferIndex(layer, modU, modV);
0678 auto itr = hgpar_->typesInLayers_.find(indx);
0679 #ifdef EDM_ML_DEBUG
0680 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:WaferType " << layer << ":" << modU << ":" << modV
0681 << ":" << indx << " Test " << (itr != hgpar_->typesInLayers_.end());
0682 #endif
0683 if (itr == hgpar_->typesInLayers_.end()) {
0684 #ifdef EDM_ML_DEBUG
0685 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: Cannot find " << layer << ":" << modU << ":" << modV
0686 << " in wadferIndex";
0687 #endif
0688 return false;
0689 }
0690
0691 if (fullAndPart_) {
0692 auto ktr = hgpar_->waferInfoMap_.find(indx);
0693 #ifdef EDM_ML_DEBUG
0694 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:WaferInfoMap " << layer << ":" << modU << ":"
0695 << modV << ":" << indx << " Test " << (ktr != hgpar_->waferInfoMap_.end());
0696 #endif
0697 if (ktr == hgpar_->waferInfoMap_.end()) {
0698 #ifdef EDM_ML_DEBUG
0699 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: Cannot find " << layer << ":" << modU << ":" << modV
0700 << " in wadferInfoMap";
0701 #endif
0702 return false;
0703 }
0704 } else {
0705 auto jtr = waferIn_.find(indx);
0706 #ifdef EDM_ML_DEBUG
0707 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:WaferIn " << jtr->first << ":" << jtr->second;
0708 #endif
0709 if (!(jtr->second)) {
0710 #ifdef EDM_ML_DEBUG
0711 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: Cannot find " << layer << ":" << modU << ":" << modV
0712 << " in wadferIn";
0713 #endif
0714 return false;
0715 }
0716 }
0717
0718 if (fullAndPart || fullAndPart_) {
0719 auto ktr = hgpar_->waferTypes_.find(indx);
0720 if (ktr != hgpar_->waferTypes_.end()) {
0721 if (hgpar_->waferMaskMode_ > 0) {
0722 if (ktr->second.first == HGCalTypes::WaferOut) {
0723 #ifdef EDM_ML_DEBUG
0724 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: Cannot find " << layer << ":" << modU << ":" << modV
0725 << " due to WaferOut";
0726 #endif
0727 return false;
0728 }
0729 } else {
0730 if (ktr->second.first < HGCalTypes::WaferCornerMin) {
0731 #ifdef EDM_ML_DEBUG
0732 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: Cannot find " << layer << ":" << modU << ":" << modV
0733 << " due to WaferCornerMin";
0734 #endif
0735 return false;
0736 }
0737 }
0738 }
0739 }
0740 return true;
0741 }
0742
0743 bool HGCalDDDConstants::isValidHex8(int layer, int modU, int modV, int cellU, int cellV, bool fullAndPart) const {
0744
0745 #ifdef EDM_ML_DEBUG
0746 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: layer|wafer " << layer << ":" << modU << ":" << modV << ":"
0747 << fullAndPart << " Valid " << isValidHex8(layer, modU, modV, fullAndPart);
0748 #endif
0749 if (!isValidHex8(layer, modU, modV, fullAndPart))
0750 return false;
0751 int indx = HGCalWaferIndex::waferIndex(layer, modU, modV);
0752 auto itr = hgpar_->typesInLayers_.find(indx);
0753 int type = hgpar_->waferTypeL_[itr->second];
0754 int N = (((hgpar_->waferTypeL_[itr->second] == HGCalTypes::WaferHD120) ||
0755 (hgpar_->waferTypeL_[itr->second] == HGCalTypes::WaferHD200))
0756 ? hgpar_->nCellsFine_
0757 : hgpar_->nCellsCoarse_);
0758 #ifdef EDM_ML_DEBUG
0759 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants::isValidHex8:Cell " << cellU << ":" << cellV << ":" << N
0760 << " Tests " << (cellU >= 0) << ":" << (cellU < 2 * N) << ":" << (cellV >= 0) << ":"
0761 << (cellV < 2 * N) << ":" << ((cellV - cellU) < N) << ":" << ((cellU - cellV) <= N);
0762 #endif
0763 if ((cellU < 0) || (cellU >= 2 * N) || (cellV < 0) || (cellV >= 2 * N)) {
0764 #ifdef EDM_ML_DEBUG
0765 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: Cannot statisfy Cell 1 condition " << cellU << ":" << cellV
0766 << ":" << N;
0767 #endif
0768 return false;
0769 }
0770 if (((cellV - cellU) >= N) || ((cellU - cellV) > N)) {
0771 #ifdef EDM_ML_DEBUG
0772 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:: Cannot statisfy Cell 2 condition " << cellU << ":" << cellV
0773 << ":" << N;
0774 #endif
0775 return false;
0776 }
0777 return isValidCell8(layer, modU, modV, cellU, cellV, type);
0778 }
0779
0780 bool HGCalDDDConstants::isValidTrap(int zside, int layer, int irad, int iphi) const {
0781
0782 const auto& indx = getIndex(layer, true);
0783 if (indx.first < 0)
0784 return false;
0785 bool ok = ((hgpar_->scintValidRing(indx.first, irad)) && (iphi > 0) && (iphi <= hgpar_->scintCells(layer)));
0786 bool valid = ((ok && trapezoidFile()) ? tileExist(zside, layer, irad, iphi) : ok);
0787 #ifdef EDM_ML_DEBUG
0788 bool tileEx = trapezoidFile() ? tileExist(zside, layer, irad, iphi) : true;
0789 if (!valid)
0790 edm::LogWarning("HGCalGeomT") << "HGCalDDDConstants::isValidityTrap: Input " << zside << ":" << layer << ":" << irad
0791 << ":" << iphi << " Range on Ring " << hgpar_->scintValidRing(indx.first, irad)
0792 << " Range on phi 0:" << hgpar_->scintCells(layer) << " tileExist " << tileEx
0793 << " Valid " << ok << ":" << tileExist(zside, layer, irad, iphi) << ":" << valid;
0794 else
0795 edm::LogVerbatim("HGCalGeomT") << "HGCalDDDConstants::isValidityTrap: Input " << zside << ":" << layer << ":"
0796 << irad << ":" << iphi << " Range on Ring "
0797 << hgpar_->scintValidRing(indx.first, irad)
0798 << " Range on phi 0:" << hgpar_->scintCells(layer) << " tileExist " << tileEx
0799 << " Valid " << ok << ":" << tileExist(zside, layer, irad, iphi) << ":" << valid;
0800 #endif
0801 return valid;
0802 }
0803
0804 int HGCalDDDConstants::lastLayer(bool reco) const {
0805 return (hgpar_->firstLayer_ + tot_layers_[static_cast<int>(reco)] - 1);
0806 }
0807
0808 unsigned int HGCalDDDConstants::layers(bool reco) const { return tot_layers_[static_cast<int>(reco)]; }
0809
0810 int HGCalDDDConstants::layerIndex(int lay, bool reco) const {
0811 int ll = lay - hgpar_->firstLayer_;
0812 if (ll < 0 || ll >= static_cast<int>(hgpar_->layerIndex_.size()))
0813 return -1;
0814 if (waferHexagon6()) {
0815 if (reco && ll >= static_cast<int>(hgpar_->depthIndex_.size()))
0816 return -1;
0817 return (reco ? hgpar_->depthLayerF_[ll] : hgpar_->layerIndex_[ll]);
0818 } else {
0819 return (hgpar_->layerIndex_[ll]);
0820 }
0821 }
0822
0823 unsigned int HGCalDDDConstants::layersInit(bool reco) const {
0824 return (reco ? hgpar_->depthIndex_.size() : hgpar_->layerIndex_.size());
0825 }
0826
0827 std::pair<float, float> HGCalDDDConstants::localToGlobal8(
0828 int zside, int lay, int waferU, int waferV, double localX, double localY, bool reco, bool debug) const {
0829 double x(localX), y(localY);
0830 bool rotx =
0831 ((!hgpar_->layerType_.empty()) && (hgpar_->layerType_[lay - hgpar_->firstLayer_] == HGCalTypes::WaferCenterR));
0832 if (debug)
0833 edm::LogVerbatim("HGCalGeom") << "LocalToGlobal8 " << lay << ":" << (lay - hgpar_->firstLayer_) << ":" << rotx
0834 << " Local (" << x << ":" << y << ") Reco " << reco;
0835 if (!reco) {
0836 x *= HGCalParameters::k_ScaleToDDD;
0837 y *= HGCalParameters::k_ScaleToDDD;
0838 }
0839 const auto& xy = waferPositionNoRot(lay, waferU, waferV, reco, debug);
0840 x += xy.first;
0841 y += xy.second;
0842 int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0843 auto ktr = hgpar_->waferInfoMap_.find(indx);
0844 if (cassetteMode() && (ktr != hgpar_->waferInfoMap_.end())) {
0845 auto cshift = hgcassette_.getShift(lay, -1, (ktr->second).cassette);
0846 std::ostringstream st1;
0847 if (debug)
0848 st1 << "Cassette " << (ktr->second).cassette << " Shift " << cshift.first << ":" << cshift.second << " Original "
0849 << x << ":" << y;
0850 if (!reco) {
0851 x -= ((HGCalParameters::k_ScaleToDDD)*zside * cshift.first);
0852 y += ((HGCalParameters::k_ScaleToDDD)*cshift.second);
0853 } else {
0854 x -= (zside * cshift.first);
0855 y += cshift.second;
0856 }
0857 if (debug) {
0858 st1 << " Final " << x << ":" << y;
0859 edm::LogVerbatim("HGCalGeom") << st1.str();
0860 }
0861 }
0862 if (debug)
0863 edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by adding " << xy.first << ":" << xy.second;
0864 return (rotx ? getXY(lay, x, y, false) : std::make_pair(x, y));
0865 }
0866
0867 std::pair<float, float> HGCalDDDConstants::locateCell(int cell, int lay, int type, bool reco) const {
0868
0869 float x(999999.), y(999999.);
0870 const auto& index = getIndex(lay, reco);
0871 int i = index.first;
0872 if (i < 0)
0873 return std::make_pair(x, y);
0874 if (waferHexagon6()) {
0875 x = hgpar_->waferPosX_[type];
0876 y = hgpar_->waferPosY_[type];
0877 #ifdef EDM_ML_DEBUG
0878 float x0(x), y0(y);
0879 #endif
0880 if ((hgpar_->waferTypeT_[type] - 1 == HGCSiliconDetId::HGCalHD120) ||
0881 (hgpar_->waferTypeT_[type] - 1 == HGCSiliconDetId::HGCalHD200)) {
0882 x += hgpar_->cellFineX_[cell];
0883 y += hgpar_->cellFineY_[cell];
0884 } else {
0885 x += hgpar_->cellCoarseX_[cell];
0886 y += hgpar_->cellCoarseY_[cell];
0887 }
0888 #ifdef EDM_ML_DEBUG
0889 edm::LogVerbatim("HGCalGeom") << "LocateCell (Wafer) " << x0 << ":" << y0 << " Final " << x << ":" << y;
0890 #endif
0891 if (!reco) {
0892 x *= HGCalParameters::k_ScaleToDDD;
0893 y *= HGCalParameters::k_ScaleToDDD;
0894 }
0895 }
0896 return std::make_pair(x, y);
0897 }
0898
0899 std::pair<float, float> HGCalDDDConstants::locateCell(int zside,
0900 int lay,
0901 int waferU,
0902 int waferV,
0903 int cellU,
0904 int cellV,
0905 bool reco,
0906 bool all,
0907 bool norot,
0908 bool cog,
0909 bool debug) const {
0910 double x(0), y(0);
0911 int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
0912 auto itr = hgpar_->typesInLayers_.find(indx);
0913 int type = ((itr == hgpar_->typesInLayers_.end()) ? 2 : hgpar_->waferTypeL_[itr->second]);
0914 int layertype = layerType(lay);
0915 bool rotx = (norot) ? false : (layertype == HGCalTypes::WaferCenterR);
0916 if (debug) {
0917 edm::LogVerbatim("HGCalGeom") << "LocateCell " << lay << ":" << (lay - hgpar_->firstLayer_) << ":" << layertype
0918 << ":" << rotx << ":" << waferU << ":" << waferV << ":" << indx << ":"
0919 << (itr == hgpar_->typesInLayers_.end()) << ":" << type << " Flags " << reco << ":"
0920 << all;
0921 }
0922 auto ktr = hgpar_->waferInfoMap_.end();
0923 int place(HGCalCell::cellPlacementOld);
0924 if (waferHexagon8File()) {
0925 if (cassetteMode()) {
0926 ktr = hgpar_->waferInfoMap_.find(indx);
0927 if (ktr != hgpar_->waferInfoMap_.end())
0928 place = HGCalCell::cellPlacementIndex(1, HGCalTypes::layerFrontBack(layertype), (ktr->second).orient);
0929 }
0930 int part = partialWaferType(lay, waferU, waferV);
0931 auto xy = (waferHexagon8Fine() || cog) ? cellOffset_->cellOffsetUV2XY1(cellU, cellV, place, type, part)
0932 : hgcell_->cellUV2XY2(cellU, cellV, place, type);
0933 x = xy.first;
0934 y = xy.second;
0935 if (debug)
0936 edm::LogVerbatim("HGCalGeom") << "Type " << type << " Place " << place << " Cell " << cellU << ":" << cellV
0937 << " Position " << x << ":" << y;
0938 } else {
0939 int kndx = cellV * 100 + cellU;
0940 if ((type == HGCSiliconDetId::HGCalHD120) || (type == HGCSiliconDetId::HGCalHD200)) {
0941 auto jtr = hgpar_->cellFineIndex_.find(kndx);
0942 if (jtr != hgpar_->cellFineIndex_.end()) {
0943 x = hgpar_->cellFineX_[jtr->second];
0944 y = hgpar_->cellFineY_[jtr->second];
0945 }
0946 if (debug)
0947 edm::LogVerbatim("HGCalGeom") << "Fine " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
0948 << (jtr != hgpar_->cellFineIndex_.end());
0949 } else {
0950 auto jtr = hgpar_->cellCoarseIndex_.find(kndx);
0951 if (jtr != hgpar_->cellCoarseIndex_.end()) {
0952 x = hgpar_->cellCoarseX_[jtr->second];
0953 y = hgpar_->cellCoarseY_[jtr->second];
0954 }
0955 if (debug)
0956 edm::LogVerbatim("HGCalGeom") << "Coarse " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y
0957 << ":" << (jtr != hgpar_->cellCoarseIndex_.end());
0958 }
0959 }
0960 if (!reco) {
0961 x *= HGCalParameters::k_ScaleToDDD;
0962 y *= HGCalParameters::k_ScaleToDDD;
0963 }
0964 if (all) {
0965 const auto& xy = waferPositionNoRot(lay, waferU, waferV, reco, debug);
0966 x += xy.first;
0967 y += xy.second;
0968 if (cassetteMode() && (ktr != hgpar_->waferInfoMap_.end())) {
0969 auto cshift = hgcassette_.getShift(lay, -1, (ktr->second).cassette);
0970 std::ostringstream st1;
0971 if (debug)
0972 st1 << "Cassette " << (ktr->second).cassette << " Shift " << cshift.first << ":" << cshift.second
0973 << " Original " << x << ":" << y;
0974 if (!reco) {
0975 x -= ((HGCalParameters::k_ScaleToDDD)*cshift.first);
0976 y += ((HGCalParameters::k_ScaleToDDD)*cshift.second);
0977 } else {
0978 x -= cshift.first;
0979 y += cshift.second;
0980 }
0981 if (debug) {
0982 st1 << " Final " << x << ":" << y;
0983 edm::LogVerbatim("HGCalGeom") << st1.str();
0984 }
0985 }
0986 if (debug)
0987 edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by adding " << xy.first << ":" << xy.second;
0988 }
0989 return (rotx ? getXY(lay, x, y, false) : std::make_pair(x, y));
0990 }
0991
0992 std::pair<float, float> HGCalDDDConstants::locateCell(const HGCSiliconDetId& id, bool cog, bool debug) const {
0993 return locateCell(
0994 id.zside(), id.layer(), id.waferU(), id.waferV(), id.cellU(), id.cellV(), true, true, false, cog, debug);
0995 }
0996
0997 std::pair<float, float> HGCalDDDConstants::locateCell(const HGCScintillatorDetId& id, bool debug) const {
0998 return locateCellTrap(id.zside(), id.layer(), id.iradius(), id.iphi(), true, debug);
0999 }
1000
1001 std::pair<float, float> HGCalDDDConstants::locateCellHex(int cell, int wafer, bool reco) const {
1002 float x(0), y(0);
1003 if ((hgpar_->waferTypeT_[wafer] - 1 == HGCSiliconDetId::HGCalHD120) ||
1004 (hgpar_->waferTypeT_[wafer] - 1 == HGCSiliconDetId::HGCalHD200)) {
1005 x = hgpar_->cellFineX_[cell];
1006 y = hgpar_->cellFineY_[cell];
1007 } else {
1008 x = hgpar_->cellCoarseX_[cell];
1009 y = hgpar_->cellCoarseY_[cell];
1010 }
1011 if (!reco) {
1012 x *= HGCalParameters::k_ScaleToDDD;
1013 y *= HGCalParameters::k_ScaleToDDD;
1014 }
1015 return std::make_pair(x, y);
1016 }
1017
1018 std::pair<float, float> HGCalDDDConstants::locateCellTrap(
1019 int zside, int lay, int irad, int iphi, bool reco, bool debug) const {
1020 float x(0), y(0);
1021 const auto& indx = getIndex(lay, reco);
1022 #ifdef EDM_ML_DEBUG
1023 edm::LogVerbatim("HGCalGeom") << "locateCellTrap:: Input " << lay << ":" << irad << ":" << iphi << ":" << zside << ":"
1024 << reco << ":" << indx.first;
1025 #endif
1026 debug = true;
1027 if (indx.first >= 0) {
1028 int ir = std::abs(irad);
1029 int type = hgpar_->scintType(lay);
1030 double phi = (iphi - 0.5) * indx.second;
1031 double z = hgpar_->zLayerHex_[indx.first];
1032 double r = 0.5 * (hgpar_->radiusLayer_[type][ir - 1] + hgpar_->radiusLayer_[type][ir]);
1033 std::pair<double, double> range = rangeR(z, true);
1034 if (debug) {
1035 std::ostringstream st1;
1036 st1 << "locateCellTrap:: Input " << lay << ":" << irad << ":" << iphi << ":" << reco << " indx " << indx.first
1037 << " IR " << ir << ":";
1038 if (hgpar_->scintFine(indx.first))
1039 st1 << hgpar_->iradMinBHFine_[indx.first] << ":" << hgpar_->iradMaxBHFine_[indx.first];
1040 else
1041 st1 << hgpar_->iradMinBH_[indx.first] << ":" << hgpar_->iradMaxBH_[indx.first];
1042 edm::LogVerbatim("HGCalGeom") << st1.str() << " Type " << type << " Z " << indx.first << ":" << z << " phi "
1043 << phi << ":" << convertRadToDeg(phi) << " R " << r << ":" << range.first << ":"
1044 << range.second << " file " << (!trapezoidFile()) << " CassetteMode "
1045 << cassetteMode();
1046 }
1047 if (!trapezoidFile())
1048 r = std::max(range.first, std::min(r, range.second));
1049 x = r * std::cos(phi);
1050 y = r * std::sin(phi);
1051 int ll = lay - hgpar_->firstLayer_;
1052 x += hgpar_->xLayerHex_[ll];
1053 y += hgpar_->yLayerHex_[ll];
1054 if (irad < 0)
1055 x = -x;
1056 if (cassetteMode()) {
1057 int cassette = HGCalTileIndex::tileCassette(iphi, hgpar_->phiOffset_, hgpar_->nphiCassette_, hgpar_->cassettes_);
1058 auto cshift = hgcassette_.getShift(lay, -1, cassette);
1059 std::ostringstream st1;
1060 if (debug)
1061 st1 << "Cassette " << cassette << " Shift " << cshift.first << ":" << cshift.second << " Original " << x << ":"
1062 << y;
1063 x -= cshift.first;
1064 y += cshift.second;
1065 if (debug) {
1066 st1 << " Final " << x << ":" << y;
1067 edm::LogVerbatim("HGCalGeom") << st1.str();
1068 }
1069 }
1070 }
1071 if (!reco) {
1072 x *= HGCalParameters::k_ScaleToDDD;
1073 y *= HGCalParameters::k_ScaleToDDD;
1074 }
1075 return std::make_pair(x, y);
1076 }
1077
1078 bool HGCalDDDConstants::maskCell(const DetId& detId, int corners) const {
1079 bool mask(false);
1080 if (corners > 2 && corners <= static_cast<int>(HGCalParameters::k_CornerSize)) {
1081 if (waferHexagon8()) {
1082 int N(0), layer(0), waferU(0), waferV(0), u(0), v(0);
1083 if (detId.det() == DetId::Forward) {
1084 HFNoseDetId id(detId);
1085 N = getUVMax(id.type());
1086 layer = id.layer();
1087 waferU = id.waferU();
1088 waferV = id.waferV();
1089 u = id.cellU();
1090 v = id.cellV();
1091 } else {
1092 HGCSiliconDetId id(detId);
1093 N = getUVMax(id.type());
1094 layer = id.layer();
1095 waferU = id.waferU();
1096 waferV = id.waferV();
1097 u = id.cellU();
1098 v = id.cellV();
1099 }
1100 int wl = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1101 auto itr = hgpar_->waferTypes_.find(wl);
1102 auto ktr = hgpar_->waferInfoMap_.find(wl);
1103 #ifdef EDM_ML_DEBUG
1104 edm::LogVerbatim("HGCalGeom") << "MaskCell: Layer " << layer << " Wafer " << waferU << ":" << waferV << " Index "
1105 << wl << ":" << (itr != hgpar_->waferTypes_.end()) << ":"
1106 << (ktr != hgpar_->waferInfoMap_.end());
1107 #endif
1108 if (cassetteMode()) {
1109 int part = (ktr != hgpar_->waferInfoMap_.end()) ? (ktr->second).part : HGCalTypes::WaferFull;
1110 mask = !(HGCalWaferMask::goodCell(u, v, part));
1111 } else if (itr != hgpar_->waferTypes_.end()) {
1112 if ((itr->second).second <= HGCalTypes::k_OffsetRotation)
1113 mask = HGCalWaferMask::maskCell(u, v, N, (itr->second).first, (itr->second).second, corners);
1114 else
1115 mask = !(HGCalWaferMask::goodCell(
1116 u, v, N, (itr->second).first, ((itr->second).second - HGCalTypes::k_OffsetRotation)));
1117 }
1118 }
1119 }
1120 return mask;
1121 }
1122
1123 int HGCalDDDConstants::maxCells(bool reco) const {
1124 int cells(0);
1125 for (unsigned int i = 0; i < layers(reco); ++i) {
1126 int lay = reco ? hgpar_->depth_[i] : hgpar_->layer_[i];
1127 if (cells < maxCells(lay, reco))
1128 cells = maxCells(lay, reco);
1129 }
1130 return cells;
1131 }
1132
1133 int HGCalDDDConstants::maxCells(int lay, bool reco) const {
1134 const auto& index = getIndex(lay, reco);
1135 if (index.first < 0)
1136 return 0;
1137 if (waferHexagon6()) {
1138 unsigned int cells(0);
1139 for (unsigned int k = 0; k < hgpar_->waferTypeT_.size(); ++k) {
1140 if (waferInLayerTest(k, index.first, hgpar_->defineFull_)) {
1141 unsigned int cell = (((hgpar_->waferTypeT_[k] - 1) == HGCSiliconDetId::HGCalHD120) ||
1142 ((hgpar_->waferTypeT_[k] - 1) == HGCSiliconDetId::HGCalHD200))
1143 ? (hgpar_->cellFineX_.size())
1144 : (hgpar_->cellCoarseX_.size());
1145 if (cell > cells)
1146 cells = cell;
1147 }
1148 }
1149 return static_cast<int>(cells);
1150 } else if (waferHexagon8()) {
1151 int cells(0);
1152 for (unsigned int k = 0; k < hgpar_->waferCopy_.size(); ++k) {
1153 if (waferInLayerTest(k, index.first, hgpar_->defineFull_)) {
1154 auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(
1155 lay, HGCalWaferIndex::waferU(hgpar_->waferCopy_[k]), HGCalWaferIndex::waferV(hgpar_->waferCopy_[k])));
1156 int type =
1157 ((itr == hgpar_->typesInLayers_.end()) ? HGCSiliconDetId::HGCalLD300 : hgpar_->waferTypeL_[itr->second]);
1158 int N = ((type == HGCSiliconDetId::HGCalHD120) || (type == HGCSiliconDetId::HGCalHD200))
1159 ? hgpar_->nCellsFine_
1160 : hgpar_->nCellsCoarse_;
1161 cells = std::max(cells, 3 * N * N);
1162 }
1163 }
1164 return cells;
1165 } else if (tileTrapezoid()) {
1166 return hgpar_->scintCells(index.first + hgpar_->firstLayer_);
1167 } else {
1168 return 0;
1169 }
1170 }
1171
1172 int HGCalDDDConstants::maxRows(int lay, bool reco) const {
1173 int kymax(0);
1174 const auto& index = getIndex(lay, reco);
1175 int i = index.first;
1176 if (i < 0)
1177 return kymax;
1178 if (waferHexagon6()) {
1179 for (unsigned int k = 0; k < hgpar_->waferCopy_.size(); ++k) {
1180 if (waferInLayerTest(k, i, hgpar_->defineFull_)) {
1181 int ky = ((hgpar_->waferCopy_[k]) / 100) % 100;
1182 if (ky > kymax)
1183 kymax = ky;
1184 }
1185 }
1186 } else if (waferHexagon8()) {
1187 kymax = 1 + 2 * hgpar_->waferUVMaxLayer_[index.first];
1188 }
1189 return kymax;
1190 }
1191
1192 int HGCalDDDConstants::modifyUV(int uv, int type1, int type2) const {
1193
1194 int uvx(uv);
1195 if (type1 != type2) {
1196 if ((type1 == HGCSiliconDetId::HGCalHD120) || (type1 == HGCSiliconDetId::HGCalHD200)) {
1197 if ((type2 == HGCSiliconDetId::HGCalLD200) || (type2 == HGCSiliconDetId::HGCalLD300))
1198 uvx = (2 * uv + 1) / 3;
1199 } else {
1200 if ((type2 == HGCSiliconDetId::HGCalHD120) || (type2 == HGCSiliconDetId::HGCalHD200))
1201 uvx = (3 * uv) / 2;
1202 }
1203 }
1204 return uvx;
1205 }
1206
1207 int HGCalDDDConstants::modules(int lay, bool reco) const {
1208 if (getIndex(lay, reco).first < 0)
1209 return 0;
1210 else
1211 return max_modules_layer_[static_cast<int>(reco)][lay];
1212 }
1213
1214 int HGCalDDDConstants::modulesInit(int lay, bool reco) const {
1215 int nmod(0);
1216 const auto& index = getIndex(lay, reco);
1217 if (index.first < 0)
1218 return nmod;
1219 if (!tileTrapezoid()) {
1220 for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
1221 if (waferInLayerTest(k, index.first, hgpar_->defineFull_))
1222 ++nmod;
1223 }
1224 } else {
1225 nmod = 1 + hgpar_->lastModule_[index.first] - hgpar_->firstModule_[index.first];
1226 }
1227 return nmod;
1228 }
1229
1230 double HGCalDDDConstants::mouseBite(bool reco) const {
1231 return (reco ? hgpar_->mouseBite_ : HGCalParameters::k_ScaleToDDD * hgpar_->mouseBite_);
1232 }
1233
1234 int HGCalDDDConstants::numberCells(bool reco) const {
1235 int cells = (tileTrapezoid() && ((hgpar_->waferMaskMode_ == HGCalGeomParameters::scintillatorFile) ||
1236 (hgpar_->waferMaskMode_ == HGCalGeomParameters::scintillatorCassette)))
1237 ? tileCount(0, -1)
1238 : 0;
1239 if (cells == 0) {
1240 unsigned int nlayer = (reco) ? hgpar_->depth_.size() : hgpar_->layer_.size();
1241 for (unsigned k = 0; k < nlayer; ++k) {
1242 std::vector<int> ncells = numberCells(((reco) ? hgpar_->depth_[k] : hgpar_->layer_[k]), reco);
1243 cells = std::accumulate(ncells.begin(), ncells.end(), cells);
1244 }
1245 }
1246 return cells;
1247 }
1248
1249 std::vector<int> HGCalDDDConstants::numberCells(int lay, bool reco) const {
1250 const auto& index = getIndex(lay, reco);
1251 int i = index.first;
1252 std::vector<int> ncell;
1253 if (i >= 0) {
1254 if (waferHexagon6()) {
1255 for (unsigned int k = 0; k < hgpar_->waferTypeT_.size(); ++k) {
1256 if (waferInLayerTest(k, i, hgpar_->defineFull_)) {
1257 unsigned int cell = (((hgpar_->waferTypeT_[k] - 1) == HGCSiliconDetId::HGCalHD120) ||
1258 ((hgpar_->waferTypeT_[k] - 1) == HGCSiliconDetId::HGCalHD200))
1259 ? (hgpar_->cellFineX_.size())
1260 : (hgpar_->cellCoarseX_.size());
1261 ncell.emplace_back(static_cast<int>(cell));
1262 }
1263 }
1264 } else if (tileTrapezoid()) {
1265 int nphi = hgpar_->scintCells(lay);
1266 for (int k = hgpar_->firstModule_[i]; k <= hgpar_->lastModule_[i]; ++k)
1267 ncell.emplace_back(nphi);
1268 } else {
1269 for (unsigned int k = 0; k < hgpar_->waferCopy_.size(); ++k) {
1270 if (waferInLayerTest(k, index.first, hgpar_->defineFull_)) {
1271 int cell = numberCellsHexagon(lay,
1272 HGCalWaferIndex::waferU(hgpar_->waferCopy_[k]),
1273 HGCalWaferIndex::waferV(hgpar_->waferCopy_[k]),
1274 true);
1275 ncell.emplace_back(cell);
1276 }
1277 }
1278 }
1279 }
1280 return ncell;
1281 }
1282
1283 int HGCalDDDConstants::numberCellsHexagon(int wafer) const {
1284 if (wafer >= 0 && wafer < static_cast<int>(hgpar_->waferTypeT_.size())) {
1285 if (((hgpar_->waferTypeT_[wafer] - 1) == HGCSiliconDetId::HGCalHD120) ||
1286 ((hgpar_->waferTypeT_[wafer] - 1) == HGCSiliconDetId::HGCalHD200))
1287 return static_cast<int>(hgpar_->cellFineX_.size());
1288 else
1289 return static_cast<int>(hgpar_->cellCoarseX_.size());
1290 } else {
1291 return 0;
1292 }
1293 }
1294
1295 int HGCalDDDConstants::numberCellsHexagon(int lay, int waferU, int waferV, bool flag) const {
1296 auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(lay, waferU, waferV));
1297 int type = ((itr == hgpar_->typesInLayers_.end()) ? HGCSiliconDetId::HGCalLD300 : hgpar_->waferTypeL_[itr->second]);
1298 int N = ((type == HGCSiliconDetId::HGCalHD120) || (type == HGCSiliconDetId::HGCalHD200)) ? hgpar_->nCellsFine_
1299 : hgpar_->nCellsCoarse_;
1300 if (flag)
1301 return (3 * N * N);
1302 else
1303 return N;
1304 }
1305
1306 int32_t HGCalDDDConstants::placementIndex(const HGCSiliconDetId& id) const {
1307 int32_t place(0);
1308 int32_t layer = id.layer();
1309 int32_t layertype = layerType(layer);
1310 int32_t indx = HGCalWaferIndex::waferIndex(layer, id.waferU(), id.waferV());
1311 auto ktr = hgpar_->waferInfoMap_.find(indx);
1312 if (ktr != hgpar_->waferInfoMap_.end()) {
1313 place = HGCalCell::cellPlacementIndex(id.zside(), layertype, (ktr->second).orient);
1314 }
1315 #ifdef EDM_ML_DEBUG
1316 edm::LogVerbatim("HGCalGeom") << "ID: " << id << " Layer " << layer << ":" << layertype << " Index " << indx << ":"
1317 << (ktr != hgpar_->waferInfoMap_.end()) << " Place " << place;
1318 #endif
1319 return place;
1320 }
1321
1322 std::pair<double, double> HGCalDDDConstants::rangeR(double z, bool reco) const {
1323 double rmin(0), rmax(0), zz(0);
1324 if (hgpar_->detectorType_ > 0) {
1325 zz = (reco ? std::abs(z) : HGCalParameters::k_ScaleFromDDD * std::abs(z));
1326 if (hgpar_->detectorType_ <= 2) {
1327 rmin = HGCalGeomTools::radius(zz, hgpar_->zFrontMin_, hgpar_->rMinFront_, hgpar_->slopeMin_);
1328 } else {
1329 rmin = HGCalGeomTools::radius(
1330 zz, hgpar_->firstLayer_, hgpar_->firstMixedLayer_, hgpar_->zLayerHex_, hgpar_->radiusMixBoundary_);
1331 }
1332 if ((hgpar_->detectorType_ == 2) && (zz >= hgpar_->zLayerHex_[hgpar_->firstMixedLayer_ - 1])) {
1333 rmax = HGCalGeomTools::radius(
1334 zz, hgpar_->firstLayer_, hgpar_->firstMixedLayer_, hgpar_->zLayerHex_, hgpar_->radiusMixBoundary_);
1335 } else {
1336 rmax = HGCalGeomTools::radius(zz, hgpar_->zFrontTop_, hgpar_->rMaxFront_, hgpar_->slopeTop_);
1337 }
1338 }
1339 if (!reco) {
1340 rmin *= HGCalParameters::k_ScaleToDDD;
1341 rmax *= HGCalParameters::k_ScaleToDDD;
1342 }
1343 #ifdef EDM_ML_DEBUG
1344 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:rangeR: " << z << ":" << zz << " R " << rmin << ":" << rmax;
1345 #endif
1346 return std::make_pair(rmin, rmax);
1347 }
1348
1349 std::pair<double, double> HGCalDDDConstants::rangeRLayer(int lay, bool reco) const {
1350 double rmin(0), rmax(0);
1351 const auto& index = getIndex(lay, reco);
1352 if (index.first >= 0 && index.first < static_cast<int>(hgpar_->rMinLayHex_.size())) {
1353 rmin = hgpar_->rMinLayHex_[index.first];
1354 rmax = hgpar_->rMaxLayHex_[index.first];
1355 }
1356 if (!reco) {
1357 rmin *= HGCalParameters::k_ScaleToDDD;
1358 rmax *= HGCalParameters::k_ScaleToDDD;
1359 }
1360 #ifdef EDM_ML_DEBUG
1361 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:rangeR: " << lay << ":" << index.first << " R " << rmin << ":"
1362 << rmax;
1363 #endif
1364 return std::make_pair(rmin, rmax);
1365 }
1366
1367 std::pair<double, double> HGCalDDDConstants::rangeZ(bool reco) const {
1368 double zmin = (hgpar_->zLayerHex_[0] - hgpar_->waferThick_);
1369 double zmax = (hgpar_->zLayerHex_[hgpar_->zLayerHex_.size() - 1] + hgpar_->waferThick_);
1370 #ifdef EDM_ML_DEBUG
1371 edm::LogVerbatim("HGCalGeom") << "HGCalDDDConstants:rangeZ: " << zmin << ":" << zmax << ":" << hgpar_->waferThick_;
1372 #endif
1373 if (!reco) {
1374 zmin *= HGCalParameters::k_ScaleToDDD;
1375 zmax *= HGCalParameters::k_ScaleToDDD;
1376 }
1377 return std::make_pair(zmin, zmax);
1378 }
1379
1380 std::pair<int, int> HGCalDDDConstants::rowColumnWafer(int wafer) const {
1381 int row(0), col(0);
1382 if (wafer < static_cast<int>(hgpar_->waferCopy_.size())) {
1383 int copy = hgpar_->waferCopy_[wafer];
1384 col = HGCalTypes::getUnpackedU(copy);
1385 row = HGCalTypes::getUnpackedV(copy);
1386 ;
1387 }
1388 return std::make_pair(row, col);
1389 }
1390
1391 double HGCalDDDConstants::sensorSizeOffset(bool reco) const {
1392 return (reco ? hgpar_->sensorSizeOffset_ : HGCalParameters::k_ScaleToDDD * hgpar_->sensorSizeOffset_);
1393 }
1394
1395 std::pair<int, int> HGCalDDDConstants::simToReco(int cell, int lay, int mod, bool half) const {
1396 if (!waferHexagon6()) {
1397 return std::make_pair(cell, lay);
1398 } else {
1399 const auto& index = getIndex(lay, false);
1400 int i = index.first;
1401 if (i < 0) {
1402 edm::LogWarning("HGCalGeom") << "Wrong Layer # " << lay << " not in the list ***** ERROR *****";
1403 return std::make_pair(-1, -1);
1404 }
1405 if (mod >= static_cast<int>(hgpar_->waferTypeL_.size())) {
1406 edm::LogWarning("HGCalGeom") << "Invalid Wafer # " << mod << "should be < " << (hgpar_->waferTypeL_).size()
1407 << " ***** ERROR *****";
1408 return std::make_pair(-1, -1);
1409 }
1410 int depth(-1);
1411 int kx = cell;
1412 int type = hgpar_->waferTypeL_[mod];
1413 if (type == 1) {
1414 depth = hgpar_->layerGroup_[i];
1415 } else if (type == 2) {
1416 depth = hgpar_->layerGroupM_[i];
1417 } else {
1418 depth = hgpar_->layerGroupO_[i];
1419 }
1420 return std::make_pair(kx, depth);
1421 }
1422 }
1423
1424 int HGCalDDDConstants::tileCount(int layer, int ring) const {
1425 int laymin(layer), laymax(layer), ringmin(ring), ringmax(ring), kount(0);
1426 if (layer == 0) {
1427 laymin = hgpar_->firstLayer_;
1428 laymax = lastLayer(true);
1429 }
1430 #ifdef EDM_ML_DEBUG
1431 edm::LogVerbatim("HGCalGeom") << "tileCount: layer " << layer << " ring " << ring << " layerMin/Max " << laymin << ":"
1432 << laymax;
1433 #endif
1434 for (int lay = laymin; lay <= laymax; ++lay) {
1435 if (ring < 0) {
1436 int ll = lay - hgpar_->firstLayer_;
1437 ringmin = hgpar_->tileRingRange_[ll].first;
1438 ringmax = hgpar_->tileRingRange_[ll].second;
1439 }
1440 #ifdef EDM_ML_DEBUG
1441 edm::LogVerbatim("HGCalGeom") << "tileCount: lay " << lay << ":" << (lay - hgpar_->firstLayer_) << " rings "
1442 << ringmin << ":" << ringmax;
1443 #endif
1444 for (int rin = ringmin; rin <= ringmax; ++rin) {
1445 int indx = HGCalTileIndex::tileIndex(lay, rin + 1, 0);
1446 auto itr = hgpar_->tileInfoMap_.find(indx);
1447 #ifdef EDM_ML_DEBUG
1448 edm::LogVerbatim("HGCalGeom") << "tileCount: rin " << rin << " indx " << indx << " itr "
1449 << (itr != hgpar_->tileInfoMap_.end());
1450 #endif
1451 if (itr != hgpar_->tileInfoMap_.end()) {
1452 for (int k = 0; k < 4; ++k) {
1453 std::bitset<24> b(itr->second.hex[k]);
1454 kount += b.count();
1455 }
1456 }
1457 #ifdef EDM_ML_DEBUG
1458 edm::LogVerbatim("HGCalGeom") << "tileCount: lay|rin " << lay << ":" << rin << " kount " << kount;
1459 #endif
1460 }
1461 }
1462 return (3 * kount);
1463 }
1464
1465 bool HGCalDDDConstants::tileExist(int zside, int layer, int ring, int phi) const {
1466 const auto& index = getIndex(layer, true);
1467 bool fine = hgpar_->scintFine(index.first);
1468 #ifdef EDM_ML_DEBUG
1469 edm::LogVerbatim("HGCalGeomT") << "TileExist: Layer : " << layer << ":" << index.first << " Fine " << fine;
1470 #endif
1471 bool ok = ((ring >= hgpar_->iradMinBH_[index.first]) && (ring <= hgpar_->iradMaxBH_[index.first]));
1472 if (ok) {
1473 if (fine) {
1474 int indx = HGCalTileIndex::tileIndex(layer, ring, 1);
1475 auto itr = hgpar_->tileInfoMap_.find(indx);
1476 ok = (itr == hgpar_->tileInfoMap_.end()) ? false : HGCalTileIndex::tileFineExist(itr->second.hex, zside, phi);
1477 #ifdef EDM_ML_DEBUG
1478 if (!ok)
1479 edm::LogWarning("HGCalGeomT") << "TileExist:input " << zside << ":" << layer << ":" << ring << ":" << phi
1480 << " Index flag " << indx << ":" << (itr != hgpar_->tileInfoMap_.end()) << " ok "
1481 << ok;
1482 if (HGCalTileIndex::tileFineExist(itr->second.hex, zside, phi) !=
1483 HGCalTileIndex::tileExist(itr->second.hex, zside, phi))
1484 edm::LogVerbatim("HGCalGeom") << "Zside:Layer:Ring:Phi " << zside << ":" << layer << ":" << ring << ":" << phi
1485 << " hex " << std::hex << itr->second.hex[0] << ":" << itr->second.hex[1] << ":"
1486 << itr->second.hex[2] << ":" << itr->second.hex[3] << ":" << itr->second.hex[4]
1487 << ":" << itr->second.hex[5] << std::dec << " OK " << ok << ":"
1488 << HGCalTileIndex::tileExist(itr->second.hex, zside, phi) << " CHECK";
1489 #endif
1490 return ok;
1491 } else {
1492 int indx = HGCalTileIndex::tileIndex(layer, ring, 0);
1493 auto itr = hgpar_->tileInfoMap_.find(indx);
1494 ok = (itr == hgpar_->tileInfoMap_.end()) ? false : HGCalTileIndex::tileExist(itr->second.hex, zside, phi);
1495 #ifdef EDM_ML_DEBUG
1496 if (!ok)
1497 edm::LogWarning("HGCalGeomT") << "TileExist:input " << zside << ":" << layer << ":" << ring << ":" << phi
1498 << " Index " << index.first << ":" << (itr != hgpar_->tileInfoMap_.end())
1499 << " ok " << ok;
1500 if (HGCalTileIndex::tileFineExist(itr->second.hex, zside, phi) !=
1501 HGCalTileIndex::tileExist(itr->second.hex, zside, phi))
1502 edm::LogVerbatim("HGCalGeom") << "Zside:Layer:Ring:Phi " << zside << ":" << layer << ":" << ring << ":" << phi
1503 << " hex " << std::hex << itr->second.hex[0] << ":" << itr->second.hex[1] << ":"
1504 << itr->second.hex[2] << ":" << itr->second.hex[3] << ":" << itr->second.hex[4]
1505 << ":" << itr->second.hex[5] << std::dec << " OK " << ok << ":"
1506 << HGCalTileIndex::tileFineExist(itr->second.hex, zside, phi) << " CHECK";
1507 #endif
1508 return ok;
1509 }
1510 #ifdef EDM_ML_DEBUG
1511 } else {
1512 edm::LogWarning("HGCalGeomT") << "TileExist:input " << zside << ":" << layer << ":" << ring << ":" << phi
1513 << " Index " << index.first << " Ring limits " << hgpar_->iradMinBH_[index.first]
1514 << ":" << hgpar_->iradMaxBH_[index.first] << " ok " << ok;
1515 #endif
1516 }
1517 return ok;
1518 }
1519
1520 HGCalParameters::tileInfo HGCalDDDConstants::tileInfo(int zside, int layer, int ring) const {
1521 const auto& index = getIndex(layer, true);
1522 bool fine = hgpar_->scintFine(index.first);
1523 HGCalParameters::tileInfo info;
1524 int indx = (fine) ? HGCalTileIndex::tileIndex(layer, ring, 1) : HGCalTileIndex::tileIndex(layer, ring, 0);
1525 ;
1526 auto itr = hgpar_->tileInfoMap_.find(indx);
1527 if (itr != hgpar_->tileInfoMap_.end())
1528 info = itr->second;
1529 #ifdef EDM_ML_DEBUG
1530 edm::LogVerbatim("HGCalGeomT") << "TileInfo:input " << zside << ":" << layer << ":" << ring << " Fine " << fine
1531 << " Index flag " << indx << ":" << (itr != hgpar_->tileInfoMap_.end());
1532 #endif
1533 return info;
1534 }
1535
1536 bool HGCalDDDConstants::tilePhiEdge(double phi, int layer, int iphi) const {
1537 double dif1 = std::abs(phi - hgpar_->scintCellSize(layer) * (iphi - 1));
1538 double dif2 = std::abs(phi - hgpar_->scintCellSize(layer) * iphi);
1539 #ifdef EDM_ML_DEBUG
1540 edm::LogVerbatim("HGCalGeomT") << "HGCalDDDConstants::tilePhiEdge:: input: " << phi << ":" << layer << ":" << iphi
1541 << " Differences " << dif1 << ":" << dif2;
1542 #endif
1543 return ((dif1 < tol_) || (dif2 < tol_));
1544 }
1545
1546 bool HGCalDDDConstants::tileRingEdge(double r, int layer, int ring) const {
1547 int type = hgpar_->scintType(layer);
1548 double dif1 = std::abs(r - hgpar_->radiusLayer_[type][ring - 1]);
1549 double dif2 = std::abs(r - hgpar_->radiusLayer_[type][ring]);
1550 #ifdef EDM_ML_DEBUG
1551 edm::LogVerbatim("HGCalGeomT") << "HGCalDDDConstants::tileRingEdge:: input: " << r << ":" << layer << ":" << ring
1552 << " Differences " << dif1 << ":" << dif2;
1553 #endif
1554 return ((dif1 < tol_) || (dif2 < tol_));
1555 }
1556 std::pair<int, int> HGCalDDDConstants::tileRings(int layer) const {
1557 if (trapezoidFile()) {
1558 int ll = layer - hgpar_->firstLayer_;
1559 if (ll >= 0 && ll < static_cast<int>(hgpar_->tileRingRange_.size()))
1560 return hgpar_->tileRingRange_[ll];
1561 }
1562 return std::make_pair(0, 0);
1563 }
1564
1565 std::pair<int, int> HGCalDDDConstants::tileType(int layer, int ring, int phi) const {
1566 int indx = HGCalTileIndex::tileIndex(layer, ring, phi);
1567 int type(-1), sipm(-1);
1568 auto itr = hgpar_->tileInfoMap_.find(indx);
1569 if (itr != hgpar_->tileInfoMap_.end()) {
1570 type = 1 + (itr->second).type;
1571 sipm = ((itr->second).sipm == HGCalTypes::SiPMLarge) ? 0 : 1;
1572 }
1573 #ifdef EDM_ML_DEBUG
1574 edm::LogVerbatim("HGCalGeom") << "tileType::Input layet:ring:phi " << layer << ":" << ring << ":" << phi
1575 << " Output Type:SiPM " << type << ":" << sipm;
1576 #endif
1577 return std::make_pair(type, sipm);
1578 }
1579
1580 int HGCalDDDConstants::waferFromCopy(int copy) const {
1581 const int ncopies = hgpar_->waferCopy_.size();
1582 int wafer(ncopies);
1583 bool result(false);
1584 for (int k = 0; k < ncopies; ++k) {
1585 if (copy == hgpar_->waferCopy_[k]) {
1586 wafer = k;
1587 result = true;
1588 break;
1589 }
1590 }
1591 if (!result) {
1592 wafer = -1;
1593 #ifdef EDM_ML_DEBUG
1594 edm::LogVerbatim("HGCalGeom") << "Cannot find " << copy << " in a list of " << ncopies << " members";
1595 for (int k = 0; k < ncopies; ++k)
1596 edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << hgpar_->waferCopy_[k];
1597 #endif
1598 }
1599 #ifdef EDM_ML_DEBUG
1600 edm::LogVerbatim("HGCalGeom") << "WaferFromCopy " << copy << ":" << wafer << ":" << result;
1601 #endif
1602 return wafer;
1603 }
1604
1605 void HGCalDDDConstants::waferFromPosition(const double x, const double y, int& wafer, int& icell, int& celltyp) const {
1606
1607 double xx = HGCalParameters::k_ScaleFromDDD * x;
1608 double yy = HGCalParameters::k_ScaleFromDDD * y;
1609 int size_ = static_cast<int>(hgpar_->waferCopy_.size());
1610 wafer = size_;
1611 for (int k = 0; k < size_; ++k) {
1612 double dx = std::abs(xx - hgpar_->waferPosX_[k]);
1613 double dy = std::abs(yy - hgpar_->waferPosY_[k]);
1614 if (dx <= rmax_ && dy <= hexside_) {
1615 if ((dy <= 0.5 * hexside_) || (dx * tan30deg_ <= (hexside_ - dy))) {
1616 wafer = k;
1617 celltyp = hgpar_->waferTypeT_[k];
1618 xx -= hgpar_->waferPosX_[k];
1619 yy -= hgpar_->waferPosY_[k];
1620 break;
1621 }
1622 }
1623 }
1624 if (wafer < size_) {
1625 if ((celltyp - 1 == HGCSiliconDetId::HGCalHD120) || (celltyp - 1 == HGCSiliconDetId::HGCalHD200))
1626 icell = cellHex(
1627 xx, yy, 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[0], hgpar_->cellFineX_, hgpar_->cellFineY_);
1628 else
1629 icell = cellHex(xx,
1630 yy,
1631 0.5 * HGCalParameters::k_ScaleFromDDD * hgpar_->cellSize_[1],
1632 hgpar_->cellCoarseX_,
1633 hgpar_->cellCoarseY_);
1634 } else {
1635 wafer = -1;
1636 #ifdef EDM_ML_DEBUG
1637 edm::LogVerbatim("HGCalGeom") << "Cannot get wafer type corresponding to " << x << ":" << y << " " << xx << ":"
1638 << yy;
1639 #endif
1640 }
1641 #ifdef EDM_ML_DEBUG
1642 edm::LogVerbatim("HGCalGeom") << "Position " << x << ":" << y << " Wafer " << wafer << ":" << size_ << " XX " << xx
1643 << ":" << yy << " Cell " << icell << " Type " << celltyp;
1644 #endif
1645 }
1646
1647 void HGCalDDDConstants::waferFromPosition(const double x,
1648 const double y,
1649 const int zside,
1650 const int layer,
1651 int& waferU,
1652 int& waferV,
1653 int& cellU,
1654 int& cellV,
1655 int& celltype,
1656 double& wt,
1657 bool extend,
1658 bool debug) const {
1659
1660 bool waferin = ((waferU == 0) && (waferV == 0));
1661 if (waferin)
1662 waferU = waferV = 1 + hgpar_->waferUVMax_;
1663 cellU = cellV = celltype = 0;
1664 if ((hgpar_->xLayerHex_.empty()) || (hgpar_->yLayerHex_.empty()))
1665 return;
1666 int ll = layer - hgpar_->firstLayer_;
1667 int layertype = layerType(layer);
1668 bool rotx = ((!hgpar_->layerType_.empty()) && (layertype == HGCalTypes::WaferCenterR));
1669 double xx(0), yy(0);
1670 if (rotx) {
1671 std::pair<double, double> xy =
1672 getXY(layer, HGCalParameters::k_ScaleFromDDD * x, HGCalParameters::k_ScaleFromDDD * y, true);
1673 xx = xy.first - hgpar_->xLayerHex_[ll];
1674 yy = xy.second - hgpar_->yLayerHex_[ll];
1675 } else {
1676 xx = HGCalParameters::k_ScaleFromDDD * x - hgpar_->xLayerHex_[ll];
1677 yy = HGCalParameters::k_ScaleFromDDD * y - hgpar_->yLayerHex_[ll];
1678 }
1679 if (debug)
1680 edm::LogVerbatim("HGCalGeom") << "waferFromPosition:: Layer " << layer << ":" << ll << " Rot " << rotx << " X " << x
1681 << ":" << xx << " Y " << y << ":" << yy << " side " << zside << " extend " << extend
1682 << " initial wafer index " << waferU << ":" << waferV;
1683 ;
1684 double rmax = extend ? rmaxT_ : rmax_;
1685 double hexside = extend ? hexsideT_ : hexside_;
1686 if (waferin) {
1687 double tolmin(100.0);
1688 for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
1689 double dx0(0), dy0(0);
1690 waferU = HGCalWaferIndex::waferU(hgpar_->waferCopy_[k]);
1691 waferV = HGCalWaferIndex::waferV(hgpar_->waferCopy_[k]);
1692 if (cassetteMode()) {
1693 int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1694 auto ktr = hgpar_->waferInfoMap_.find(indx);
1695 if (ktr != hgpar_->waferInfoMap_.end()) {
1696 auto cshift = hgcassette_.getShift(layer, -1, (ktr->second).cassette);
1697 if (debug)
1698 edm::LogVerbatim("HGCalGeom")
1699 << "Cassette " << (ktr->second).cassette << " Shift " << cshift.first << ":" << cshift.second;
1700 dx0 = -cshift.first;
1701 dy0 = cshift.second;
1702 }
1703 }
1704 double dx = std::abs(xx - dx0 - hgpar_->waferPosX_[k]);
1705 double dy = std::abs(yy - dy0 - hgpar_->waferPosY_[k]);
1706 constexpr double tolc = 0.01;
1707 if (debug) {
1708 edm::LogVerbatim("HGCalGeom") << "Wafer " << waferU << ":" << waferV << " position " << xx << ":" << yy
1709 << " Distance " << dx << ":" << dy << " diff0 " << (dx - rmax) << ":"
1710 << (dy - hexside) << " diff1 " << (dy - 0.5 * hexside) << ":"
1711 << (dx * tan30deg_ - (hexside - dy));
1712 if ((dx - rmax) <= tolc && (dy - hexside) <= tolc) {
1713 tolmin = std::min(tolmin, (dy - 0.5 * hexside));
1714 tolmin = std::min(tolmin, (dx * tan30deg_ - (hexside - dy)));
1715 }
1716 }
1717 if ((dx - rmax) <= tolc && (dy - hexside) <= tolc) {
1718 if (((dy - 0.5 * hexside) <= tolc) || ((dx * tan30deg_ - (hexside - dy)) <= tolc)) {
1719 if (waferHexagon8File()) {
1720 int index = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1721 celltype = HGCalWaferType::getType(index, hgpar_->waferInfoMap_);
1722 if (debug)
1723 edm::LogVerbatim("HGCalGeom")
1724 << "Position (" << x << ", " << y << ") Wafer type:partial:orient:cassette " << celltype << ":"
1725 << HGCalWaferType::getPartial(index, hgpar_->waferInfoMap_) << ":"
1726 << HGCalWaferType::getOrient(index, hgpar_->waferInfoMap_) << ":"
1727 << HGCalWaferType::getCassette(index, hgpar_->waferInfoMap_);
1728 } else {
1729 auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1730 celltype = ((itr == hgpar_->typesInLayers_.end()) ? HGCSiliconDetId::HGCalLD300
1731 : hgpar_->waferTypeL_[itr->second]);
1732 }
1733 if (debug)
1734 edm::LogVerbatim("HGCalGeom")
1735 << "WaferFromPosition:: Input " << layer << ":" << ll << ":" << hgpar_->firstLayer_ << ":" << rotx
1736 << ":" << x << ":" << y << ":" << hgpar_->xLayerHex_[ll] << ":" << hgpar_->yLayerHex_[ll] << ":" << xx
1737 << ":" << yy << " compared with " << hgpar_->waferPosX_[k] << ":" << hgpar_->waferPosY_[k]
1738 << " difference " << dx << ":" << dy << ":" << dx * tan30deg_ << ":" << (hexside_ - dy)
1739 << " comparator " << rmax_ << ":" << rmaxT_ << ":" << hexside_ << ":" << hexsideT_ << " wafer "
1740 << waferU << ":" << waferV << ":" << celltype;
1741 xx -= (dx0 + hgpar_->waferPosX_[k]);
1742 yy -= (dy0 + hgpar_->waferPosY_[k]);
1743 break;
1744 }
1745 }
1746 }
1747 if (debug)
1748 edm::LogVerbatim("HGCalGeom") << "Tolmin " << tolmin;
1749 } else {
1750 for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
1751 double dx0(0), dy0(0);
1752 if ((waferU == HGCalWaferIndex::waferU(hgpar_->waferCopy_[k])) &&
1753 (waferV == HGCalWaferIndex::waferV(hgpar_->waferCopy_[k]))) {
1754 if (cassetteMode()) {
1755 int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1756 auto ktr = hgpar_->waferInfoMap_.find(indx);
1757 if (ktr != hgpar_->waferInfoMap_.end()) {
1758 auto cshift = hgcassette_.getShift(layer, -1, (ktr->second).cassette);
1759 if (debug)
1760 edm::LogVerbatim("HGCalGeom")
1761 << "Cassette " << (ktr->second).cassette << " Shift " << cshift.first << ":" << cshift.second;
1762 dx0 = -cshift.first;
1763 dy0 = cshift.second;
1764 }
1765 }
1766 if (waferHexagon8File()) {
1767 int index = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1768 celltype = HGCalWaferType::getType(index, hgpar_->waferInfoMap_);
1769 if (debug)
1770 edm::LogVerbatim("HGCalGeom") << "Position (" << x << ", " << y << ") Wafer type:partial:orient:cassette "
1771 << celltype << ":" << HGCalWaferType::getPartial(index, hgpar_->waferInfoMap_)
1772 << ":" << HGCalWaferType::getOrient(index, hgpar_->waferInfoMap_) << ":"
1773 << HGCalWaferType::getCassette(index, hgpar_->waferInfoMap_);
1774 } else {
1775 auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1776 celltype =
1777 ((itr == hgpar_->typesInLayers_.end()) ? HGCSiliconDetId::HGCalLD300 : hgpar_->waferTypeL_[itr->second]);
1778 }
1779 xx -= (dx0 + hgpar_->waferPosX_[k]);
1780 yy -= (dy0 + hgpar_->waferPosY_[k]);
1781 break;
1782 }
1783 }
1784 }
1785 if ((std::abs(waferU) <= hgpar_->waferUVMax_) && (celltype >= 0)) {
1786 int place(HGCalCell::cellPlacementOld), part(HGCalTypes::WaferFull);
1787 if (cassetteMode()) {
1788 int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
1789 auto ktr = hgpar_->waferInfoMap_.find(indx);
1790 if (ktr != hgpar_->waferInfoMap_.end()) {
1791 place = HGCalCell::cellPlacementIndex(1, HGCalTypes::layerFrontBack(layertype), (ktr->second).orient);
1792 part = (ktr->second).part;
1793 if (debug)
1794 edm::LogVerbatim("HGCalGeom") << "waferFromPosition: frontback " << layertype << ":"
1795 << HGCalTypes::layerFrontBack(layertype) << " Orient " << (ktr->second).orient
1796 << " place " << place << " part " << part;
1797 }
1798 }
1799 cellHex(xx, yy, celltype, place, part, cellU, cellV, extend, debug);
1800 wt = ((((celltype == HGCSiliconDetId::HGCalHD120) || (celltype == HGCSiliconDetId::HGCalHD200)) &&
1801 (hgpar_->useSimWt_ > 0))
1802 ? (hgpar_->cellThickness_[celltype] / hgpar_->waferThick_)
1803 : 1.0);
1804 } else {
1805 cellU = cellV = 2 * hgpar_->nCellsFine_;
1806 wt = 1.0;
1807 celltype = -1;
1808 }
1809 if ((celltype < 0) && debug) {
1810 double x1(xx);
1811 double y1(yy);
1812 edm::LogVerbatim("HGCalGeom") << "waferfFromPosition: Bad type for X " << x << ":" << x1 << ":" << xx << " Y " << y
1813 << ":" << y1 << ":" << yy << " Wafer " << waferU << ":" << waferV << " Cell " << cellU
1814 << ":" << cellV;
1815 for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
1816 double dx = std::abs(x1 - hgpar_->waferPosX_[k]);
1817 double dy = std::abs(y1 - hgpar_->waferPosY_[k]);
1818 edm::LogVerbatim("HGCalGeom") << "Wafer [" << k << "] Position (" << hgpar_->waferPosX_[k] << ", "
1819 << hgpar_->waferPosY_[k] << ") difference " << dx << ":" << dy << ":"
1820 << dx * tan30deg_ << ":" << hexside - dy << " Paramerers " << rmax << ":"
1821 << hexside;
1822 }
1823 }
1824 edm::LogVerbatim("HGCalGeomX") << "Input x:y:layer " << x << ":" << y << ":" << layer << " Wafer " << waferU << ":"
1825 << waferV << " Cell " << cellU << ":" << cellV << ":" << celltype << " wt " << wt;
1826 }
1827
1828 bool HGCalDDDConstants::waferInLayer(int wafer, int lay, bool reco) const {
1829 const auto& indx = getIndex(lay, reco);
1830 if (indx.first < 0)
1831 return false;
1832 return waferInLayerTest(wafer, indx.first, hgpar_->defineFull_);
1833 }
1834
1835 bool HGCalDDDConstants::waferFullInLayer(int wafer, int lay, bool reco) const {
1836 const auto& indx = getIndex(lay, reco);
1837 if (indx.first < 0)
1838 return false;
1839 return waferInLayerTest(wafer, indx.first, false);
1840 }
1841
1842 HGCalParameters::waferInfo HGCalDDDConstants::waferInfo(int lay, int waferU, int waferV) const {
1843 int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
1844 auto itr = hgpar_->waferInfoMap_.find(indx);
1845 return ((itr == hgpar_->waferInfoMap_.end()) ? HGCalParameters::waferInfo() : itr->second);
1846 }
1847
1848 std::pair<double, double> HGCalDDDConstants::waferParameters(bool reco) const {
1849 if (reco)
1850 return std::make_pair(rmax_, hexside_);
1851 else
1852 return std::make_pair(HGCalParameters::k_ScaleToDDD * rmax_, HGCalParameters::k_ScaleToDDD * hexside_);
1853 }
1854
1855 std::pair<double, double> HGCalDDDConstants::waferPosition(int wafer, bool reco) const {
1856 double xx(0), yy(0);
1857 if (wafer >= 0 && wafer < static_cast<int>(hgpar_->waferPosX_.size())) {
1858 xx = hgpar_->waferPosX_[wafer];
1859 yy = hgpar_->waferPosY_[wafer];
1860 }
1861 if (!reco) {
1862 xx *= HGCalParameters::k_ScaleToDDD;
1863 yy *= HGCalParameters::k_ScaleToDDD;
1864 }
1865 return std::make_pair(xx, yy);
1866 }
1867
1868 std::pair<double, double> HGCalDDDConstants::waferPosition(
1869 int lay, int waferU, int waferV, bool reco, bool debug) const {
1870 int ll = lay - hgpar_->firstLayer_;
1871 bool rotx = ((!hgpar_->layerType_.empty()) && (hgpar_->layerType_[ll] == HGCalTypes::WaferCenterR));
1872 #ifdef EDM_ML_DEBUG
1873 if (debug)
1874 edm::LogVerbatim("HGCalGeom") << "Layer " << lay << ":" << ll << " Rotation " << rotx << " U:V " << waferU << ":"
1875 << waferV;
1876 #endif
1877 auto xy = waferPositionNoRot(lay, waferU, waferV, reco, debug);
1878 std::pair<double, double> xy0 = (rotx) ? getXY(lay, xy.first, xy.second, false) : xy;
1879 #ifdef EDM_ML_DEBUG
1880 if (debug)
1881 edm::LogVerbatim("HGCalGeom") << "Without and with rotation " << xy.first << ":" << xy.second << ":" << xy0.first
1882 << ":" << xy0.second;
1883 #endif
1884 return xy0;
1885 }
1886
1887 int HGCalDDDConstants::waferFileIndex(unsigned int kk) const {
1888 if (kk < hgpar_->waferInfoMap_.size()) {
1889 auto itr = hgpar_->waferInfoMap_.begin();
1890 std::advance(itr, kk);
1891 return itr->first;
1892 } else
1893 return 0;
1894 }
1895
1896 std::tuple<int, int, int, int> HGCalDDDConstants::waferFileInfo(unsigned int kk) const {
1897 if (kk < hgpar_->waferInfoMap_.size()) {
1898 auto itr = hgpar_->waferInfoMap_.begin();
1899 std::advance(itr, kk);
1900 return std::make_tuple(itr->second.type, itr->second.part, itr->second.orient, itr->second.cassette);
1901 } else
1902 return std::make_tuple(0, 0, 0, 0);
1903 }
1904
1905 std::tuple<int, int, int, int> HGCalDDDConstants::waferFileInfoFromIndex(int kk) const {
1906 auto itr = hgpar_->waferInfoMap_.find(kk);
1907 if (itr != hgpar_->waferInfoMap_.end()) {
1908 return std::make_tuple(itr->second.type, itr->second.part, itr->second.orient, itr->second.cassette);
1909 } else
1910 return std::make_tuple(0, 0, 0, 0);
1911 }
1912
1913 GlobalPoint HGCalDDDConstants::waferLocal2Global(
1914 HepGeom::Point3D<float>& loc, const DetId& id, bool useWafer, bool reco, bool debug) const {
1915 HGCSiliconDetId detid(id);
1916 double x(0), y(0);
1917 if (useWafer) {
1918 auto xyw = waferPositionNoRot(detid.layer(), detid.waferU(), detid.waferV(), reco, debug);
1919 x = xyw.first;
1920 y = xyw.second;
1921 }
1922 auto xy = getXY(detid.layer(), (x + loc.x()), (y + loc.y()), false);
1923 double zz = (detid.zside() < 0) ? -(loc.z() + waferZ(detid.layer(), reco)) : (loc.z() + waferZ(detid.layer(), reco));
1924 double xx = (detid.zside() < 0) ? -xy.first : xy.first;
1925 return GlobalPoint(xx, xy.second, zz);
1926 }
1927
1928 int HGCalDDDConstants::wafers() const {
1929 int wafer(0);
1930 if (!tileTrapezoid()) {
1931 for (unsigned int i = 0; i < layers(true); ++i) {
1932 int lay = hgpar_->depth_[i];
1933 wafer += modules(lay, true);
1934 }
1935 } else {
1936 wafer = static_cast<int>(hgpar_->moduleLayR_.size());
1937 }
1938 return wafer;
1939 }
1940
1941 int HGCalDDDConstants::wafers(int layer, int type) const {
1942 int wafer(0);
1943 if (!tileTrapezoid()) {
1944 auto itr = waferLayer_.find(layer);
1945 if (itr != waferLayer_.end()) {
1946 unsigned ity = (type > 0 && type <= 2) ? type : 0;
1947 wafer = (itr->second)[ity];
1948 }
1949 } else {
1950 const auto& index = getIndex(layer, true);
1951 wafer = 1 + hgpar_->lastModule_[index.first] - hgpar_->firstModule_[index.first];
1952 }
1953 return wafer;
1954 }
1955
1956 int HGCalDDDConstants::waferType(DetId const& id, bool fromFile) const {
1957 int type(1);
1958 if (waferHexagon8()) {
1959 if (fromFile && (waferFileSize() > 0)) {
1960 int layer(0), waferU(0), waferV(0);
1961 if (id.det() != DetId::Forward) {
1962 HGCSiliconDetId hid(id);
1963 layer = hid.layer();
1964 waferU = hid.waferU();
1965 waferV = hid.waferV();
1966 } else {
1967 HFNoseDetId hid(id);
1968 layer = hid.layer();
1969 waferU = hid.waferU();
1970 waferV = hid.waferV();
1971 }
1972 auto itr = hgpar_->waferInfoMap_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1973 if (itr != hgpar_->waferInfoMap_.end())
1974 type = (itr->second).type;
1975 } else {
1976 type = ((id.det() != DetId::Forward) ? HGCSiliconDetId(id).type() : HFNoseDetId(id).type());
1977 }
1978 } else if (waferHexagon6()) {
1979 type = waferTypeL(HGCalDetId(id).wafer()) - 1;
1980 }
1981 return type;
1982 }
1983
1984 int HGCalDDDConstants::waferType(int layer, int waferU, int waferV, bool fromFile) const {
1985 int type(HGCSiliconDetId::HGCalLD300);
1986 if (waferHexagon8()) {
1987 if (fromFile && (waferFileSize() > 0)) {
1988 auto itr = hgpar_->waferInfoMap_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1989 if (itr != hgpar_->waferInfoMap_.end())
1990 type = (itr->second).type;
1991 } else {
1992 auto itr = hgpar_->typesInLayers_.find(HGCalWaferIndex::waferIndex(layer, waferU, waferV));
1993 if (itr != hgpar_->typesInLayers_.end())
1994 type = hgpar_->waferTypeL_[itr->second];
1995 }
1996 } else if (waferHexagon6()) {
1997 if ((waferU >= 0) && (waferU < static_cast<int>(hgpar_->waferTypeL_.size())))
1998 type = (hgpar_->waferTypeL_[waferU] - 1);
1999 }
2000 return type;
2001 }
2002
2003 std::tuple<int, int, int> HGCalDDDConstants::waferType(HGCSiliconDetId const& id, bool fromFile) const {
2004 const auto& index = HGCalWaferIndex::waferIndex(id.layer(), id.waferU(), id.waferV());
2005 int type(-1), part(-1), orient(-1);
2006 if (fromFile && (waferFileSize() > 0)) {
2007 auto itr = hgpar_->waferInfoMap_.find(index);
2008 if (itr != hgpar_->waferInfoMap_.end()) {
2009 type = (itr->second).type;
2010 part = (itr->second).part;
2011 orient = (itr->second).orient;
2012 }
2013 } else {
2014 auto ktr = hgpar_->typesInLayers_.find(index);
2015 if (ktr != hgpar_->typesInLayers_.end())
2016 type = hgpar_->waferTypeL_[ktr->second];
2017 auto itr = hgpar_->waferTypes_.find(index);
2018 if (itr != hgpar_->waferTypes_.end()) {
2019 if ((itr->second).second < HGCalTypes::k_OffsetRotation) {
2020 orient = (itr->second).second;
2021 if ((itr->second).first == HGCalGeomTools::k_allCorners) {
2022 part = HGCalTypes::WaferFull;
2023 } else if ((itr->second).first == HGCalGeomTools::k_fiveCorners) {
2024 part = HGCalTypes::WaferFive;
2025 } else if ((itr->second).first == HGCalGeomTools::k_fourCorners) {
2026 part = HGCalTypes::WaferHalf;
2027 } else if ((itr->second).first == HGCalGeomTools::k_threeCorners) {
2028 part = HGCalTypes::WaferThree;
2029 }
2030 } else {
2031 part = (itr->second).first;
2032 orient = ((itr->second).second - HGCalTypes::k_OffsetRotation);
2033 }
2034 } else {
2035 part = HGCalTypes::WaferFull;
2036 orient = 0;
2037 }
2038 }
2039 return std::make_tuple(type, part, orient);
2040 }
2041
2042 std::pair<int, int> HGCalDDDConstants::waferTypeRotation(
2043 int layer, int waferU, int waferV, bool fromFile, bool debug) const {
2044 int type(HGCalTypes::WaferOut), rotn(0);
2045 int wl = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
2046 bool withinList(true);
2047 if (fromFile && (waferFileSize() > 0)) {
2048 auto itr = hgpar_->waferInfoMap_.find(wl);
2049 withinList = (itr != hgpar_->waferInfoMap_.end());
2050 if (withinList) {
2051 type = (itr->second).part;
2052 rotn = (itr->second).orient;
2053 }
2054 } else {
2055 auto itr = hgpar_->waferTypes_.find(wl);
2056 if (waferHexagon8()) {
2057 withinList = (itr != hgpar_->waferTypes_.end());
2058 if (withinList) {
2059 if ((itr->second).second < HGCalTypes::k_OffsetRotation) {
2060 rotn = (itr->second).second;
2061 if ((itr->second).first == HGCalGeomTools::k_allCorners) {
2062 type = HGCalTypes::WaferFull;
2063 } else if ((itr->second).first == HGCalGeomTools::k_fiveCorners) {
2064 type = HGCalTypes::WaferFive;
2065 } else if ((itr->second).first == HGCalGeomTools::k_fourCorners) {
2066 type = HGCalTypes::WaferHalf;
2067 } else if ((itr->second).first == HGCalGeomTools::k_threeCorners) {
2068 type = HGCalTypes::WaferThree;
2069 }
2070 } else {
2071 type = (itr->second).first;
2072 rotn = ((itr->second).second - HGCalTypes::k_OffsetRotation);
2073 }
2074 } else {
2075 type = HGCalTypes::WaferFull;
2076 rotn = HGCalTypes::WaferCorner0;
2077 }
2078 }
2079 }
2080 #ifdef EDM_ML_DEBUG
2081 if (debug)
2082 edm::LogVerbatim("HGCalGeom") << "waferTypeRotation: Layer " << layer << " Wafer " << waferU << ":" << waferV
2083 << " Index " << std::hex << wl << std::dec << ":" << withinList << " Type " << type
2084 << " Rotation " << rotn;
2085 #endif
2086 return std::make_pair(type, rotn);
2087 }
2088
2089 bool HGCalDDDConstants::waferVirtual(int layer, int waferU, int waferV) const {
2090 bool type(false);
2091 if (waferHexagon8()) {
2092 int wl = HGCalWaferIndex::waferIndex(layer, waferU, waferV, false);
2093 type = (hgpar_->waferTypes_.find(wl) != hgpar_->waferTypes_.end());
2094 } else if (waferHexagon6()) {
2095 int wl = HGCalWaferIndex::waferIndex(layer, waferU, 0, true);
2096 type = (hgpar_->waferTypes_.find(wl) != hgpar_->waferTypes_.end());
2097 }
2098 return type;
2099 }
2100
2101 double HGCalDDDConstants::waferZ(int lay, bool reco) const {
2102 const auto& index = getIndex(lay, reco);
2103 if (index.first < 0)
2104 return 0;
2105 else
2106 return (reco ? hgpar_->zLayerHex_[index.first] : HGCalParameters::k_ScaleToDDD * hgpar_->zLayerHex_[index.first]);
2107 }
2108
2109 int HGCalDDDConstants::cellHex(
2110 double xx, double yy, const double& cellR, const std::vector<double>& posX, const std::vector<double>& posY) const {
2111 int num(0);
2112 const double tol(0.00001);
2113 double cellY = 2.0 * cellR * tan30deg_;
2114 for (unsigned int k = 0; k < posX.size(); ++k) {
2115 double dx = std::abs(xx - posX[k]);
2116 double dy = std::abs(yy - posY[k]);
2117 if (dx <= (cellR + tol) && dy <= (cellY + tol)) {
2118 double xmax = (dy <= 0.5 * cellY) ? cellR : (cellR - (dy - 0.5 * cellY) / tan30deg_);
2119 if (dx <= (xmax + tol)) {
2120 num = k;
2121 break;
2122 }
2123 }
2124 }
2125 return num;
2126 }
2127
2128 void HGCalDDDConstants::cellHex(
2129 double xloc, double yloc, int cellType, int place, int part, int& cellU, int& cellV, bool extend, bool debug) const {
2130 if (cassetteMode()) {
2131 auto uv = (part == HGCalTypes::WaferFull)
2132 ? hgcellUV_->cellUVFromXY3(xloc, yloc, place, cellType, true, debug)
2133 : hgcellUV_->cellUVFromXY1(xloc, yloc, place, cellType, part, true, debug);
2134 cellU = uv.first;
2135 cellV = uv.second;
2136 } else if (waferHexagon8File()) {
2137 auto uv = hgcellUV_->cellUVFromXY3(xloc, yloc, place, cellType, extend, debug);
2138 cellU = uv.first;
2139 cellV = uv.second;
2140 } else {
2141 int ncell = ((cellType == HGCSiliconDetId::HGCalHD120) || (cellType == HGCSiliconDetId::HGCalHD200))
2142 ? hgpar_->nCellsFine_
2143 : hgpar_->nCellsCoarse_;
2144 double delY = 2 * rmax_ / (3 * ncell);
2145 double delX = 0.5 * delY * sqrt3_;
2146 double delYT = (extend) ? (2 * rmaxT_ / (3 * ncell)) : delY;
2147 double delXT = 0.5 * delYT * sqrt3_;
2148 double v0 = ((xloc / delY - 1.0) / 1.5);
2149 int cv0 = (v0 > 0) ? (ncell + static_cast<int>(v0 + 0.5)) : (ncell - static_cast<int>(-v0 + 0.5));
2150 double u0 = (0.5 * yloc / delX + 0.5 * cv0);
2151 int cu0 = (u0 > 0) ? (ncell / 2 + static_cast<int>(u0 + 0.5)) : (ncell / 2 - static_cast<int>(-u0 + 0.5));
2152 cu0 = std::max(0, std::min(cu0, 2 * ncell - 1));
2153 cv0 = std::max(0, std::min(cv0, 2 * ncell - 1));
2154 if (cv0 - cu0 >= ncell)
2155 cv0 = cu0 + ncell - 1;
2156 if (debug)
2157 edm::LogVerbatim("HGCalGeom") << "cellHex: input " << xloc << ":" << yloc << ":" << cellType << " parameter "
2158 << delX << ":" << delY << " u0 " << u0 << ":" << cu0 << " v0 " << v0 << ":" << cv0;
2159 bool found(false);
2160 static constexpr int shift[3] = {0, 1, -1};
2161 for (int i1 = 0; i1 < 3; ++i1) {
2162 cellU = cu0 + shift[i1];
2163 for (int i2 = 0; i2 < 3; ++i2) {
2164 cellV = cv0 + shift[i2];
2165 if (((cellV - cellU) < ncell) && ((cellU - cellV) <= ncell) && (cellU >= 0) && (cellV >= 0) &&
2166 (cellU < 2 * ncell) && (cellV < 2 * ncell)) {
2167 double xc = (1.5 * (cellV - ncell) + 1.0) * delY;
2168 double yc = (2 * cellU - cellV - ncell) * delX;
2169 if ((std::abs(yloc - yc) <= delXT) && (std::abs(xloc - xc) <= delYT) &&
2170 ((std::abs(xloc - xc) <= 0.5 * delYT) ||
2171 (std::abs(yloc - yc) <= sqrt3_ * (delYT - std::abs(xloc - xc))))) {
2172 if (debug)
2173 edm::LogVerbatim("HGCalGeom")
2174 << "cellHex: local " << xc << ":" << yc << " difference " << std::abs(xloc - xc) << ":"
2175 << std::abs(yloc - yc) << ":" << sqrt3_ * (delY - std::abs(yloc - yc)) << " comparator " << delX
2176 << ":" << delY << " (u,v) = (" << cellU << "," << cellV << ")";
2177 found = true;
2178 break;
2179 }
2180 }
2181 }
2182 if (found)
2183 break;
2184 }
2185 if (!found) {
2186 cellU = cu0;
2187 cellV = cv0;
2188 }
2189 }
2190 }
2191
2192 std::pair<int, float> HGCalDDDConstants::getIndex(int lay, bool reco) const {
2193 int indx = layerIndex(lay, reco);
2194 if (indx < 0)
2195 return std::make_pair(-1, 0);
2196 float cell(0);
2197 if (waferHexagon6()) {
2198 cell = (reco ? hgpar_->moduleCellR_[0] : hgpar_->moduleCellS_[0]);
2199 } else {
2200 if (waferHexagon8()) {
2201 cell = (reco ? hgpar_->moduleCellR_[0] : hgpar_->moduleCellS_[0]);
2202 } else {
2203 cell = hgpar_->scintCellSize(lay);
2204 }
2205 }
2206 return std::make_pair(indx, cell);
2207 }
2208
2209 int HGCalDDDConstants::layerFromIndex(int index, bool reco) const {
2210 int ll(-1);
2211 if (waferHexagon6() && reco) {
2212 ll = static_cast<int>(std::find(hgpar_->depthLayerF_.begin(), hgpar_->depthLayerF_.end(), index) -
2213 hgpar_->depthLayerF_.begin());
2214 if (ll == static_cast<int>(hgpar_->depthLayerF_.size()))
2215 ll = -1;
2216 } else {
2217 ll = static_cast<int>(std::find(hgpar_->layerIndex_.begin(), hgpar_->layerIndex_.end(), index) -
2218 hgpar_->layerIndex_.begin());
2219 if (ll == static_cast<int>(hgpar_->layerIndex_.size()))
2220 ll = -1;
2221 }
2222 #ifdef EDM_ML_DEBUG
2223 edm::LogVerbatim("HGCalGeom") << "LayerFromIndex for " << index << ":" << reco << ":" << waferHexagon6() << " is"
2224 << ll << ":" << (ll + hgpar_->firstLayer_);
2225 #endif
2226 return ((ll < 0) ? ll : (ll + hgpar_->firstLayer_));
2227 }
2228
2229 bool HGCalDDDConstants::isValidCell(int lay, int wafer, int cell) const {
2230
2231
2232 double x = hgpar_->waferPosX_[wafer];
2233 double y = hgpar_->waferPosY_[wafer];
2234 if (((hgpar_->waferTypeT_[wafer] - 1) == HGCSiliconDetId::HGCalHD120) ||
2235 ((hgpar_->waferTypeT_[wafer] - 1) == HGCSiliconDetId::HGCalHD200)) {
2236 x += hgpar_->cellFineX_[cell];
2237 y += hgpar_->cellFineY_[cell];
2238 } else {
2239 x += hgpar_->cellCoarseX_[cell];
2240 y += hgpar_->cellCoarseY_[cell];
2241 }
2242 double rr = sqrt(x * x + y * y);
2243 bool result = ((rr >= hgpar_->rMinLayHex_[lay - 1]) && (rr <= hgpar_->rMaxLayHex_[lay - 1]) &&
2244 (wafer < static_cast<int>(hgpar_->waferPosX_.size())));
2245 #ifdef EDM_ML_DEBUG
2246 if (!result)
2247 edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << wafer << ":" << cell << " Position " << x << ":" << y
2248 << ":" << rr << " Compare Limits " << hgpar_->rMinLayHex_[lay - 1] << ":"
2249 << hgpar_->rMaxLayHex_[lay - 1] << " Flag " << result;
2250 #endif
2251 return result;
2252 }
2253
2254 bool HGCalDDDConstants::isValidCell8(int lay, int waferU, int waferV, int cellU, int cellV, int type) const {
2255 bool result(false);
2256 auto partn = waferTypeRotation(lay, waferU, waferV, false, false);
2257 #ifdef EDM_ML_DEBUG
2258 edm::LogVerbatim("HGCalGeom") << "waferHexagon8 " << waferHexagon8File() << ":" << mode_ << ":" << cassetteMode()
2259 << " part " << partn.first << ":" << partn.second;
2260 #endif
2261 if (cassetteMode()) {
2262 result = HGCalWaferMask::goodCell(cellU, cellV, partn.first);
2263 #ifdef EDM_ML_DEBUG
2264 edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << waferU << ":" << waferV << ":" << cellU << ":" << cellV
2265 << " Result " << result << " from goodCell";
2266 #endif
2267 } else {
2268 float x(0), y(0);
2269 int kndx = cellV * 100 + cellU;
2270 if ((type == HGCSiliconDetId::HGCalHD120) || (type == HGCSiliconDetId::HGCalHD200)) {
2271 auto ktr = hgpar_->cellFineIndex_.find(kndx);
2272 if (ktr != hgpar_->cellFineIndex_.end()) {
2273 x = hgpar_->cellFineX_[ktr->second];
2274 y = hgpar_->cellFineY_[ktr->second];
2275 }
2276 #ifdef EDM_ML_DEBUG
2277 edm::LogVerbatim("HGCalGeom") << "Fine " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
2278 << (ktr != hgpar_->cellFineIndex_.end());
2279 #endif
2280 } else {
2281 auto ktr = hgpar_->cellCoarseIndex_.find(kndx);
2282 if (ktr != hgpar_->cellCoarseIndex_.end()) {
2283 x = hgpar_->cellCoarseX_[ktr->second];
2284 y = hgpar_->cellCoarseY_[ktr->second];
2285 }
2286 #ifdef EDM_ML_DEBUG
2287 edm::LogVerbatim("HGCalGeom") << "Coarse " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
2288 << (ktr != hgpar_->cellCoarseIndex_.end());
2289 #endif
2290 }
2291 const auto& xy = waferPositionNoRot(lay, waferU, waferV, true, false);
2292 x += xy.first;
2293 y += xy.second;
2294 #ifdef EDM_ML_DEBUG
2295 edm::LogVerbatim("HGCalGeom") << "With wafer (" << waferU << "," << waferV << ") " << x << ":" << y;
2296 #endif
2297 double rr = sqrt(x * x + y * y);
2298 int ll = lay - hgpar_->firstLayer_;
2299 double tol = waferHexagon8File() ? 0.5 : 0.0;
2300 result = (((rr + tol) >= hgpar_->rMinLayHex_[ll]) && (rr <= hgpar_->rMaxLayHex_[ll]));
2301 #ifdef EDM_ML_DEBUG
2302 edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << ll << ":" << waferU << ":" << waferV << ":" << cellU
2303 << ":" << cellV << " Position " << x << ":" << y << ":" << rr << " Compare Limits "
2304 << hgpar_->rMinLayHex_[ll] << ":" << hgpar_->rMaxLayHex_[ll] << " Flag " << result
2305 << " from Radius Limits";
2306 #endif
2307 if (result && waferHexagon8File()) {
2308 int N = ((type == HGCSiliconDetId::HGCalHD120) || (type == HGCSiliconDetId::HGCalHD200)) ? hgpar_->nCellsFine_
2309 : hgpar_->nCellsCoarse_;
2310 result = HGCalWaferMask::goodCell(cellU, cellV, N, partn.first, partn.second);
2311 #ifdef EDM_ML_DEBUG
2312 edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << waferU << ":" << waferV << ":" << cellU << ":" << cellV
2313 << " N " << N << " part " << partn.first << ":" << partn.second << " Result "
2314 << result << " from goodCell";
2315 #endif
2316 }
2317 }
2318 return result;
2319 }
2320
2321 int32_t HGCalDDDConstants::waferIndex(int wafer, int index) const {
2322 int layer = layerFromIndex(index, true);
2323 int waferU = HGCalWaferIndex::waferU(hgpar_->waferCopy_[wafer]);
2324 int waferV = HGCalWaferIndex::waferV(hgpar_->waferCopy_[wafer]);
2325 int indx = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
2326 #ifdef EDM_ML_DEBUG
2327 edm::LogVerbatim("HGCalGeom") << "WaferIndex for " << wafer << ":" << index << " (" << layer << ":" << waferU << ":"
2328 << waferV << ") " << indx;
2329 #endif
2330 return indx;
2331 }
2332
2333 bool HGCalDDDConstants::waferInLayerTest(int wafer, int lay, bool full) const {
2334 bool in = (waferHexagon6()) ? true : false;
2335 if (!in) {
2336 double xpos = hgpar_->waferPosX_[wafer] + hgpar_->xLayerHex_[lay];
2337 double ypos = hgpar_->waferPosY_[wafer] + hgpar_->yLayerHex_[lay];
2338 std::pair<int, int> corner = HGCalGeomTools::waferCorner(
2339 xpos, ypos, rmax_, hexside_, hgpar_->rMinLayHex_[lay], hgpar_->rMaxLayHex_[lay], in);
2340 in = (full ? (corner.first > 0) : (corner.first == static_cast<int>(HGCalParameters::k_CornerSize)));
2341 if (in && fullAndPart_) {
2342 int indx = waferIndex(wafer, lay);
2343 in = (hgpar_->waferInfoMap_.find(indx) != hgpar_->waferInfoMap_.end());
2344 #ifdef EDM_ML_DEBUG
2345 if (!in)
2346 edm::LogVerbatim("HGCalGeom") << "WaferInLayerTest: Layer " << lay << " wafer " << wafer << " index " << indx
2347 << "( " << HGCalWaferIndex::waferLayer(indx) << ", "
2348 << HGCalWaferIndex::waferU(indx) << ", " << HGCalWaferIndex::waferV(indx)
2349 << ") in " << in;
2350 #endif
2351 }
2352 #ifdef EDM_ML_DEBUG
2353 edm::LogVerbatim("HGCalGeom") << "WaferInLayerTest: Layer " << lay << " wafer " << wafer << " R-limits "
2354 << hgpar_->rMinLayHex_[lay] << ":" << hgpar_->rMaxLayHex_[lay] << " Corners "
2355 << corner.first << ":" << corner.second << " In " << in;
2356 #endif
2357 }
2358 return in;
2359 }
2360
2361 std::pair<double, double> HGCalDDDConstants::waferPositionNoRot(
2362 int lay, int waferU, int waferV, bool reco, bool debug) const {
2363 int ll = lay - hgpar_->firstLayer_;
2364 double x = hgpar_->xLayerHex_[ll];
2365 double y = hgpar_->yLayerHex_[ll];
2366 #ifdef EDM_ML_DEBUG
2367 if (debug)
2368 edm::LogVerbatim("HGCalGeom") << "Layer " << lay << ":" << ll << " Shift " << hgpar_->xLayerHex_[ll] << ":"
2369 << hgpar_->yLayerHex_[ll] << " U:V " << waferU << ":" << waferV;
2370 #endif
2371 if (!reco) {
2372 x *= HGCalParameters::k_ScaleToDDD;
2373 y *= HGCalParameters::k_ScaleToDDD;
2374 }
2375 const auto& xy = waferPosition(waferU, waferV, reco);
2376 x += xy.first;
2377 y += xy.second;
2378 #ifdef EDM_ML_DEBUG
2379 if (debug)
2380 edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << ":" << xy.first << ":" << xy.second;
2381 #endif
2382 return std::make_pair(x, y);
2383 }
2384
2385 std::pair<double, double> HGCalDDDConstants::waferPosition(int waferU, int waferV, bool reco) const {
2386 double xx(0), yy(0);
2387 int indx = HGCalWaferIndex::waferIndex(0, waferU, waferV);
2388 auto itr = hgpar_->wafersInLayers_.find(indx);
2389 if (itr != hgpar_->wafersInLayers_.end()) {
2390 xx = hgpar_->waferPosX_[itr->second];
2391 yy = hgpar_->waferPosY_[itr->second];
2392 }
2393 if (!reco) {
2394 xx *= HGCalParameters::k_ScaleToDDD;
2395 yy *= HGCalParameters::k_ScaleToDDD;
2396 }
2397 return std::make_pair(xx, yy);
2398 }
2399
2400 #include "FWCore/Utilities/interface/typelookup.h"
2401
2402 TYPELOOKUP_DATA_REG(HGCalDDDConstants);
