File indexing completed on 2024-11-22 02:45:25
0001 #include "SimG4CMS/Calo/interface/HGCGuardRing.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003 #include "Geometry/HGCalCommonData/interface/HGCalWaferMask.h"
0004 #include "Geometry/HGCalCommonData/interface/HGCalWaferType.h"
0005 #include <iostream>
0006
0007
0008
0009 HGCGuardRing::HGCGuardRing(const HGCalDDDConstants& hgc)
0010 : hgcons_(hgc),
0011 modeUV_(hgcons_.geomMode()),
0012 v17OrLess_(hgcons_.v17OrLess()),
0013 waferSize_(hgcons_.waferSize(false)),
0014 sensorSizeOffset_(hgcons_.sensorSizeOffset(false)),
0015 guardRingOffset_(hgcons_.guardRingOffset(false)) {
0016 offset_ = sensorSizeOffset_ + guardRingOffset_;
0017 xmax_ = 0.5 * waferSize_ - offset_;
0018 ymax_ = xmax_ / sqrt3_;
0019 #ifdef EDM_ML_DEBUG
0020 edm::LogVerbatim("HGCSim") << "Creating HGCGuardRing with wafer size " << waferSize_ << ", Offsets "
0021 << sensorSizeOffset_ << ":" << guardRingOffset_ << ":" << offset_ << ", and mode "
0022 << modeUV_ << " xmax|ymax " << xmax_ << ":" << ymax_;
0023 #endif
0024 }
0025
0026 bool HGCGuardRing::exclude(G4ThreeVector& point, int zside, int frontBack, int layer, int waferU, int waferV) {
0027 bool check(false);
0028 if (hgcons_.waferHexagon8Module()) {
0029 int index = HGCalWaferIndex::waferIndex(layer, waferU, waferV);
0030 int partial = HGCalWaferType::getPartial(index, hgcons_.getParameter()->waferInfoMap_);
0031 #ifdef EDM_ML_DEBUG
0032 edm::LogVerbatim("HGCSim") << "HGCGuardRing:: Layer " << layer << " wafer " << waferU << ":" << waferV << " index "
0033 << index << " partial " << partial;
0034 #endif
0035 if (partial == HGCalTypes::WaferFull) {
0036 double dx = std::abs(point.x());
0037 double dy = std::abs(point.y());
0038 if (dx > xmax_) {
0039 check = true;
0040 } else if (dy > (2 * ymax_)) {
0041 check = true;
0042 } else {
0043 check = (dx > (sqrt3_ * (2 * ymax_ - dy)));
0044 }
0045 #ifdef EDM_ML_DEBUG
0046 edm::LogVerbatim("HGCSim") << "HGCGuardRing:: Point " << point << " zside " << zside << " layer " << layer
0047 << " wafer " << waferU << ":" << waferV << " partial type " << partial << ":"
0048 << HGCalTypes::WaferFull << " x " << dx << ":" << xmax_ << " y " << dy << ":" << ymax_
0049 << " check " << check;
0050 #endif
0051 } else if (partial > 0) {
0052 int orient = HGCalWaferType::getOrient(index, hgcons_.getParameter()->waferInfoMap_);
0053 #ifdef EDM_ML_DEBUG
0054 edm::LogVerbatim("HGCSim") << "HGCGuardRing:: Orient " << orient << " Mode " << modeUV_;
0055 #endif
0056 if (modeUV_ == HGCalGeometryMode::Hexagon8Module) {
0057 std::vector<std::pair<double, double> > wxy =
0058 HGCalWaferMask::waferXY(partial, orient, zside, waferSize_, offset_, 0.0, 0.0, v17OrLess_);
0059 check = !(insidePolygon(point.x(), point.y(), wxy));
0060 #ifdef EDM_ML_DEBUG
0061 std::ostringstream st1;
0062 st1 << "HGCGuardRing:: Point " << point << " Partial/orient/zside/size/offset " << partial << ":" << orient
0063 << ":" << zside << ":" << waferSize_ << offset_ << " with " << wxy.size() << " points:";
0064 for (unsigned int k = 0; k < wxy.size(); ++k)
0065 st1 << " (" << wxy[k].first << ", " << wxy[k].second << ")";
0066 edm::LogVerbatim("HGCSim") << st1.str();
0067 #endif
0068 } else {
0069 int placement = HGCalCell::cellPlacementIndex(zside, frontBack, orient);
0070 std::vector<std::pair<double, double> > wxy =
0071 HGCalWaferMask::waferXY(partial, placement, waferSize_, offset_, 0.0, 0.0, v17OrLess_);
0072 check = !(insidePolygon(point.x(), point.y(), wxy));
0073 #ifdef EDM_ML_DEBUG
0074 std::ostringstream st1;
0075 st1 << "HGCGuardRing:: Point " << point << " Partial/frontback/orient/zside/placeemnt/size/offset " << partial
0076 << ":" << frontBack << ":" << orient << ":" << zside << ":" << placement << ":" << waferSize_ << offset_
0077 << " with " << wxy.size() << " points:";
0078 for (unsigned int k = 0; k < wxy.size(); ++k)
0079 st1 << " (" << wxy[k].first << ", " << wxy[k].second << ")";
0080 edm::LogVerbatim("HGCSim") << st1.str();
0081 #endif
0082 }
0083 } else {
0084 check = true;
0085 }
0086 }
0087 return check;
0088 }
0089
0090 bool HGCGuardRing::insidePolygon(double x, double y, const std::vector<std::pair<double, double> >& xyv) {
0091 int counter(0);
0092 double x1(xyv[0].first), y1(xyv[0].second);
0093 for (unsigned i1 = 1; i1 <= xyv.size(); i1++) {
0094 unsigned i2 = (i1 % xyv.size());
0095 double x2(xyv[i2].first), y2(xyv[i2].second);
0096 if (y > std::min(y1, y2)) {
0097 if (y <= std::max(y1, y2)) {
0098 if (x <= std::max(x1, x2)) {
0099 if (y1 != y2) {
0100 double xinter = (y - y1) * (x2 - x1) / (y2 - y1) + x1;
0101 if ((x1 == x2) || (x <= xinter))
0102 ++counter;
0103 }
0104 }
0105 }
0106 }
0107 x1 = x2;
0108 y1 = y2;
0109 }
0110
0111 if (counter % 2 == 0)
0112 return false;
0113 else
0114 return true;
0115 }