File indexing completed on 2024-06-04 04:35:02
0001 #include "Geometry/HGCalCommonData/interface/HGCalCellUV.h"
0002 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0003 #include "Geometry/HGCalCommonData/interface/HGCalWaferMask.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include <iostream>
0006 #include <array>
0007 #include <algorithm>
0008 #include <cassert>
0009
0010 HGCalCellUV::HGCalCellUV(double waferSize, double separation, int32_t nFine, int32_t nCoarse) : waferSize_(waferSize) {
0011 hgcalcell_ = std::make_unique<HGCalCell>(waferSize, nFine, nCoarse);
0012 assert(nFine > 0 && nCoarse > 0);
0013 ncell_[0] = nFine;
0014 ncell_[1] = nCoarse;
0015 for (int k = 0; k < 2; ++k) {
0016 cellX_[k] = waferSize / (3 * ncell_[k]);
0017 cellY_[k] = 0.5 * sqrt3_ * cellX_[k];
0018 cellXTotal_[k] = (waferSize + separation) / (3 * ncell_[k]);
0019 cellYTotal_[k] = 0.5 * sqrt3_ * cellXTotal_[k];
0020 }
0021
0022 for (int placement = 0; placement < HGCalCell::cellPlacementTotal; ++placement) {
0023
0024 for (int iu = 0; iu < 2 * ncell_[0]; ++iu) {
0025 for (int iv = 0; iv < 2 * ncell_[0]; ++iv) {
0026 int u = (placement < HGCalCell::cellPlacementExtra) ? iv : iu;
0027 int v = (placement < HGCalCell::cellPlacementExtra) ? iu : iv;
0028 if (((v - u) < ncell_[0]) && (u - v) <= ncell_[0]) {
0029 cellPosFine_[placement][std::pair<int, int>(u, v)] = hgcalcell_->cellUV2XY1(u, v, placement, 0);
0030 }
0031 }
0032 }
0033
0034 for (int iu = 0; iu < 2 * ncell_[1]; ++iu) {
0035 for (int iv = 0; iv < 2 * ncell_[1]; ++iv) {
0036 int u = (placement < HGCalCell::cellPlacementExtra) ? iv : iu;
0037 int v = (placement < HGCalCell::cellPlacementExtra) ? iu : iv;
0038 if (((v - u) < ncell_[1]) && (u - v) <= ncell_[1]) {
0039 cellPosCoarse_[placement][std::pair<int, int>(u, v)] = hgcalcell_->cellUV2XY1(u, v, placement, 1);
0040 }
0041 }
0042 }
0043 }
0044 }
0045
0046 std::pair<int32_t, int32_t> HGCalCellUV::cellUVFromXY1(
0047 double xloc, double yloc, int32_t placement, int32_t type, bool extend, bool debug) const {
0048 if (type != 0)
0049 type = 1;
0050
0051 double xloc1 = (placement >= HGCalCell::cellPlacementExtra) ? xloc : -xloc;
0052 double cellY = (extend) ? cellYTotal_[type] : cellY_[type];
0053 int rot = placement % HGCalCell::cellPlacementExtra;
0054 static constexpr std::array<double, 6> fcos = {{1.0, cos60_, -cos60_, -1.0, -cos60_, cos60_}};
0055 static constexpr std::array<double, 6> fsin = {{0.0, sin60_, sin60_, 0.0, -sin60_, -sin60_}};
0056 double x = xloc1 * fcos[rot] - yloc * fsin[rot];
0057 double y = xloc1 * fsin[rot] + yloc * fcos[rot];
0058
0059
0060 double u = x * sin60_ + y * cos60_;
0061 double v = -x * sin60_ + y * cos60_;
0062 double w = y;
0063
0064
0065 int iu = std::floor(u / cellY) + 3 * (ncell_[type]) - 1;
0066 int iv = std::floor(v / cellY) + 3 * (ncell_[type]);
0067 int iw = std::floor(w / cellY) + 1;
0068
0069 int isv = (iu + iw) / 3;
0070 int isu = (iv + iw) / 3;
0071
0072
0073 if ((iu + iw) < 0) {
0074 isu = (iv + iw + 1) / 3;
0075 isv = 0;
0076 } else if (isv - isu > ncell_[type] - 1) {
0077 isu = (iv + iw + 1) / 3;
0078 isv = (iu + iw - 1) / 3;
0079 } else if (isu > 2 * ncell_[type] - 1) {
0080 isu = 2 * ncell_[type] - 1;
0081 isv = (iu + iw - 1) / 3;
0082 }
0083 if (debug)
0084 edm::LogVerbatim("HGCalGeom") << "cellUVFromXY1: Input " << xloc << ":" << yloc << ":" << extend << " Output "
0085 << isu << ":" << isv;
0086 return std::make_pair(isu, isv);
0087 }
0088
0089 std::pair<int32_t, int32_t> HGCalCellUV::cellUVFromXY2(
0090 double xloc, double yloc, int32_t placement, int32_t type, bool extend, bool debug) const {
0091
0092
0093 if (type != 0)
0094 type = 1;
0095 double cellY = (extend) ? cellYTotal_[type] : cellY_[type];
0096 double xloc1 = (placement >= HGCalCell::cellPlacementExtra) ? xloc : -1 * xloc;
0097 int rot = placement % HGCalCell::cellPlacementExtra;
0098 static constexpr std::array<double, 6> fcos = {{cos60_, 1.0, cos60_, -cos60_, -1.0, -cos60_}};
0099 static constexpr std::array<double, 6> fsin = {{-sin60_, 0.0, sin60_, sin60_, 0.0, -sin60_}};
0100 double x = xloc1 * fcos[rot] - yloc * fsin[rot];
0101 double y = xloc1 * fsin[rot] + yloc * fcos[rot];
0102
0103 int32_t u(-100), v(-100);
0104 int ncell = ncell_[type];
0105 double r = cellY;
0106 double l1 = (y / r) + ncell - 1.0;
0107 int l2 = std::floor((0.5 * y + 0.5 * x / sqrt3_) / r + ncell - 4.0 / 3.0);
0108 int l3 = std::floor((x / sqrt3_) / r + ncell - 4.0 / 3.0);
0109 double l4 = (y + sqrt3_ * x) / (2 * r) + 2 * ncell - 2;
0110 double l5 = (y - sqrt3_ * x) / (2 * r) - ncell;
0111 double u1 = (y / r) + ncell + 1.0;
0112 int u2 = std::ceil((0.5 * y + 0.5 * x / sqrt3_) / r + ncell + 2.0 / 3.0);
0113 int u3 = std::ceil((x / sqrt3_) / r + ncell);
0114 double u4 = l4 + 2;
0115 double u5 = l5 + 2;
0116
0117 for (int ui = l2 + 1; ui < u2; ui++) {
0118 for (int vi = l3 + 1; vi < u3; vi++) {
0119 int c1 = 2 * ui - vi;
0120 int c4 = ui + vi;
0121 int c5 = ui - 2 * vi;
0122 if ((c1 < u1) && (c1 > l1) && (c4 < u4) && (c4 > l4) && (c5 < u5) && (c5 > l5)) {
0123 u = ui;
0124 v = vi;
0125 }
0126 }
0127 }
0128
0129
0130 if (v == -1) {
0131 if (y < (2 * u - v - ncell) * r) {
0132 ++v;
0133 } else {
0134 ++u;
0135 ++v;
0136 }
0137 }
0138 if (v - u == ncell) {
0139 if ((y + sqrt3_ * x) < ((u + v - 2 * ncell + 1) * 2 * r)) {
0140 --v;
0141 } else {
0142 ++u;
0143 }
0144 }
0145 if (u == 2 * ncell) {
0146 if ((y - sqrt3_ * x) < ((u - 2 * v + ncell - 1) * 2 * r)) {
0147 --u;
0148 } else {
0149 --u;
0150 --v;
0151 }
0152 }
0153 if (debug)
0154 edm::LogVerbatim("HGCalGeom") << "cellUVFromXY2: Input " << xloc << ":" << yloc << ":" << extend << " Output " << u
0155 << ":" << v;
0156 return std::make_pair(u, v);
0157 }
0158
0159 std::pair<int32_t, int32_t> HGCalCellUV::cellUVFromXY3(
0160 double xloc, double yloc, int32_t placement, int32_t type, bool extend, bool debug) const {
0161
0162
0163 if (type != 0)
0164 type = 1;
0165 double cellX = (extend) ? cellXTotal_[type] : cellX_[type];
0166 double cellY = (extend) ? cellYTotal_[type] : cellY_[type];
0167 double xloc1 = (placement >= HGCalCell::cellPlacementExtra) ? xloc : -xloc;
0168 int rot = placement % HGCalCell::cellPlacementExtra;
0169 static constexpr std::array<double, 6> fcos = {{1.0, cos60_, -cos60_, -1.0, -cos60_, cos60_}};
0170 static constexpr std::array<double, 6> fsin = {{0.0, sin60_, sin60_, 0.0, -sin60_, -sin60_}};
0171 double xprime = xloc1 * fcos[rot] - yloc * fsin[rot];
0172 double yprime = xloc1 * fsin[rot] + yloc * fcos[rot];
0173 double x = xprime + cellX;
0174 double y = yprime;
0175
0176 x = x / cellX;
0177 y = y / cellY;
0178
0179 double cu = 2 * x / 3;
0180 double cv = -x / 3 + y / 2;
0181 double cw = -x / 3 - y / 2;
0182
0183 int iu = std::round(cu);
0184 int iv = std::round(cv);
0185 int iw = std::round(cw);
0186
0187 if (iu + iv + iw != 0) {
0188 double arr[] = {std::abs(cu - iu), std::abs(cv - iv), std::abs(cw - iw)};
0189 int i = std::distance(arr, std::max_element(arr, arr + 3));
0190
0191 if (i == 1)
0192 iv = (std::round(cv) == std::floor(cv)) ? std::ceil(cv) : std::floor(cv);
0193 else if (i == 2)
0194 iw = (std::round(cw) == std::floor(cw)) ? std::ceil(cw) : std::floor(cw);
0195 }
0196
0197
0198 int u(ncell_[type] + iv), v(ncell_[type] - 1 - iw);
0199 double xcell = (1.5 * (v - u) + 0.5) * cellX_[type];
0200 double ycell = (v + u - 2 * ncell_[type] + 1) * cellY;
0201 if (v == -1) {
0202 if ((yprime - sqrt3_ * xprime) < (ycell - sqrt3_ * xcell)) {
0203 ++v;
0204 } else {
0205 ++u;
0206 ++v;
0207 }
0208 }
0209 if (v - u == ncell_[type]) {
0210 if (yprime < ycell) {
0211 --v;
0212 } else {
0213 ++u;
0214 }
0215 }
0216 if (u == 2 * ncell_[type]) {
0217 if ((yprime + sqrt3_ * xprime) > (ycell + sqrt3_ * xcell)) {
0218 --u;
0219 } else {
0220 --u;
0221 --v;
0222 }
0223 }
0224
0225 if (debug)
0226 edm::LogVerbatim("HGCalGeom") << "cellUVFromXY3: Input " << xloc << ":" << yloc << ":" << extend << " Output " << u
0227 << ":" << v;
0228 return std::make_pair(u, v);
0229 }
0230
0231 std::pair<int, int> HGCalCellUV::cellUVFromXY4(
0232 double xloc, double yloc, int32_t placement, int32_t type, bool extend, bool debug) {
0233 if (type != 0)
0234 return cellUVFromXY4(
0235 xloc, yloc, ncell_[1], cellX_[1], cellY_[1], cellXTotal_[1], cellY_[1], cellPosCoarse_[placement], extend, debug);
0236 else
0237 return cellUVFromXY4(
0238 xloc, yloc, ncell_[0], cellX_[0], cellY_[0], cellXTotal_[0], cellY_[0], cellPosFine_[placement], extend, debug);
0239 }
0240
0241 std::pair<int, int> HGCalCellUV::cellUVFromXY4(double xloc,
0242 double yloc,
0243 int n,
0244 double cellX,
0245 double cellY,
0246 double cellXTotal,
0247 double cellYTotal,
0248 std::map<std::pair<int, int>, std::pair<double, double> >& cellPos,
0249 bool extend,
0250 bool debug) {
0251 std::pair<int, int> uv = std::make_pair(-1, -1);
0252 std::map<std::pair<int, int>, std::pair<double, double> >::const_iterator itr;
0253 for (itr = cellPos.begin(); itr != cellPos.end(); ++itr) {
0254 double delX = std::abs(xloc - (itr->second).first);
0255 double delY = std::abs(yloc - (itr->second).second);
0256 if ((delX < cellX) && (delY < cellY)) {
0257 if ((delX < (0.5 * cellX)) || (delY < (2.0 * cellY - sqrt3_ * delX))) {
0258 uv = itr->first;
0259 break;
0260 }
0261 }
0262 }
0263 if ((uv.first < 0) && extend) {
0264 for (itr = cellPos.begin(); itr != cellPos.end(); ++itr) {
0265 double delX = std::abs(xloc - (itr->second).first);
0266 double delY = std::abs(yloc - (itr->second).second);
0267 if ((delX < cellXTotal) && (delY < cellYTotal)) {
0268 if ((delX < (0.5 * cellXTotal)) || (delY < (2.0 * cellYTotal - sqrt3_ * delX))) {
0269 uv = itr->first;
0270 break;
0271 }
0272 }
0273 }
0274 }
0275 if (debug)
0276 edm::LogVerbatim("HGCalGeom") << "cellUVFromXY4: Input " << xloc << ":" << yloc << ":" << extend << " Output "
0277 << uv.first << ":" << uv.second;
0278 return uv;
0279 }
0280
0281 std::pair<int32_t, int32_t> HGCalCellUV::cellUVFromXY1(
0282 double xloc,
0283 double yloc,
0284 int32_t placement,
0285 int32_t type,
0286 int32_t partial,
0287 bool extend,
0288 bool debug) const {
0289 if (type != 0)
0290 type = 1;
0291 double cellX = (extend) ? cellXTotal_[type] : cellX_[type];
0292 double cellY = (extend) ? cellYTotal_[type] : cellY_[type];
0293 std::pair<int, int> uv = HGCalCellUV::cellUVFromXY1(xloc, yloc, placement, type, extend, debug);
0294 int u = uv.first;
0295 int v = uv.second;
0296 if (partial == HGCalTypes::WaferLDTop) {
0297 if (u * HGCalTypes::edgeWaferLDTop[0] + v * HGCalTypes::edgeWaferLDTop[1] == HGCalTypes::edgeWaferLDTop[2] + 1) {
0298 double xloc1 = (placement >= HGCalCell::cellPlacementExtra) ? xloc : -xloc;
0299 int rot = placement % HGCalCell::cellPlacementExtra;
0300 static constexpr std::array<double, 6> fcos = {{1.0, cos60_, -cos60_, -1.0, -cos60_, cos60_}};
0301 static constexpr std::array<double, 6> fsin = {{0.0, sin60_, sin60_, 0.0, -sin60_, -sin60_}};
0302 double xprime = -1 * (xloc1 * fcos[rot] - yloc * fsin[rot]);
0303 double yprime = xloc1 * fsin[rot] + yloc * fcos[rot];
0304 double xcell = -1 * (1.5 * (v - u) + 0.5) * cellX;
0305 double ycell = (v + u - 2 * ncell_[type] + 1) * cellY;
0306 if ((yprime - sqrt3_ * xprime) > (ycell - sqrt3_ * xcell)) {
0307 --u;
0308 if ((v - u) >= ncell_[1])
0309 --v;
0310 } else {
0311 --u;
0312 --v;
0313 v = std::max(v, 0);
0314 }
0315 }
0316 } else if (partial == HGCalTypes::WaferHDBottom) {
0317 if (u * HGCalTypes::edgeWaferHDBottom[0] + v * HGCalTypes::edgeWaferHDBottom[1] ==
0318 HGCalTypes::edgeWaferHDBottom[2] + 1) {
0319 double xloc1 = (placement >= HGCalCell::cellPlacementExtra) ? xloc : -xloc;
0320 int rot = placement % HGCalCell::cellPlacementExtra;
0321 static constexpr std::array<double, 6> fcos = {{1.0, cos60_, -cos60_, -1.0, -cos60_, cos60_}};
0322 static constexpr std::array<double, 6> fsin = {{0.0, sin60_, sin60_, 0.0, -sin60_, -sin60_}};
0323 double xprime = -1 * (xloc1 * fcos[rot] - yloc * fsin[rot]);
0324 double yprime = xloc1 * fsin[rot] + yloc * fcos[rot];
0325 double xcell = -1 * (1.5 * (v - u) + 0.5) * cellX;
0326 double ycell = (v + u - 2 * ncell_[type] + 1) * cellY;
0327 if ((yprime - sqrt3_ * xprime) > (ycell - sqrt3_ * xcell)) {
0328 ++u;
0329 ++v;
0330 } else {
0331 ++u;
0332 }
0333 }
0334 }
0335 if (debug)
0336 edm::LogVerbatim("HGCalGeom") << "cellUVFromXY5: Input " << xloc << ":" << yloc << ":" << extend << " Output " << u
0337 << ":" << v;
0338 return std::make_pair(u, v);
0339 }
0340
0341 std::pair<int32_t, int32_t> HGCalCellUV::cellUVFromXY2(
0342 double xloc,
0343 double yloc,
0344 int32_t placement,
0345 int32_t type,
0346 int32_t partial,
0347 bool extend,
0348 bool debug) const {
0349 if (type != 0)
0350 type = 1;
0351 std::pair<int, int> uv = HGCalCellUV::cellUVFromXY1(xloc, yloc, placement, type, extend, debug);
0352 int u = uv.first;
0353 int v = uv.second;
0354 if ((partial != HGCalTypes::WaferFull) && (type == 1)) {
0355 if (u == 1 && v == 8) {
0356 std::array<double, 4> criterion =
0357 HGCalWaferMask::maskCut(HGCalTypes::WaferLDThree, placement, waferSize_, 0.0, false);
0358 if ((criterion[0] * yloc) + (criterion[1] * xloc) < -criterion[2]) {
0359 ++u;
0360 ++v;
0361 }
0362 }
0363 if (u == 15 && v == 15) {
0364 std::array<double, 4> criterion =
0365 HGCalWaferMask::maskCut(HGCalTypes::WaferLDThree, placement, waferSize_, 0.0, false);
0366 if ((criterion[0] * yloc) + (criterion[1] * xloc) < -criterion[2]) {
0367 --u;
0368 }
0369 }
0370 if (u * HGCalTypes::edgeWaferLDTop[0] + v * HGCalTypes::edgeWaferLDTop[1] == HGCalTypes::edgeWaferLDTop[2] + 1) {
0371 std::array<double, 4> criterion =
0372 HGCalWaferMask::maskCut(HGCalTypes::WaferLDTop, placement, waferSize_, 0.0, false);
0373 if ((criterion[0] * yloc) + (criterion[1] * xloc) < -criterion[2]) {
0374 std::pair<double, double> xy1 = hgcalcell_->cellUV2XY1(u, v, placement, 1);
0375 std::array<double, 4> criterion2 =
0376 HGCalWaferMask::maskCut(HGCalTypes::WaferLDThree, placement, waferSize_, 0.0, false);
0377 if (((criterion2[0] * yloc) + (criterion2[1] * xloc) - (criterion2[0] * xy1.second) -
0378 (criterion2[1] * xy1.first)) < 0.0) {
0379 --u;
0380 if ((v - u) >= ncell_[1])
0381 --v;
0382 } else {
0383 --u;
0384 --v;
0385 v = std::max(v, 0);
0386 }
0387 }
0388 }
0389 } else if ((partial != HGCalTypes::WaferFull) && (type == 0)) {
0390 if (u == 10 && v == 0) {
0391 std::array<double, 4> criterion =
0392 HGCalWaferMask::maskCut(HGCalTypes::WaferHDTop, placement, waferSize_, 0.0, false);
0393 if ((criterion[0] * yloc) + (criterion[1] * xloc) < -criterion[2]) {
0394 --u;
0395 }
0396 }
0397 if (u == 10 && v == 21) {
0398 std::array<double, 4> criterion =
0399 HGCalWaferMask::maskCut(HGCalTypes::WaferHDTop, placement, waferSize_, 0.0, false);
0400 if ((criterion[0] * yloc) + (criterion[1] * xloc) < -criterion[2]) {
0401 --u;
0402 --v;
0403 }
0404 }
0405 if (u * HGCalTypes::edgeWaferHDBottom[0] + v * HGCalTypes::edgeWaferHDBottom[1] ==
0406 HGCalTypes::edgeWaferHDBottom[2] + 1) {
0407 std::array<double, 4> criterion =
0408 HGCalWaferMask::maskCut(HGCalTypes::WaferHDBottom, placement, waferSize_, 0.0, false);
0409 if ((criterion[0] * yloc) + (criterion[1] * xloc) < -criterion[2]) {
0410 std::pair<double, double> xy1 = hgcalcell_->cellUV2XY1(u, v, placement, 0);
0411 std::array<double, 4> criterion2 =
0412 HGCalWaferMask::maskCut(HGCalTypes::WaferHDRight, placement, waferSize_, 0.0, false);
0413 if (((criterion2[0] * yloc) + (criterion2[1] * xloc) - (criterion2[0] * xy1.second) -
0414 (criterion2[1] * xy1.first)) < 0.0) {
0415 ++u;
0416 ++v;
0417 } else {
0418 ++u;
0419 }
0420 }
0421 }
0422 }
0423 if (debug)
0424 edm::LogVerbatim("HGCalGeom") << "cellUVFromXY5: Input " << xloc << ":" << yloc << ":" << extend << " Output " << u
0425 << ":" << v;
0426 return std::make_pair(u, v);
0427 }