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