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