File indexing completed on 2023-03-17 11:23:04
0001 #include "PixelInactiveAreaFinder.h"
0002
0003 #include "FWCore/Utilities/interface/VecArray.h"
0004 #include "FWCore/Utilities/interface/transform.h"
0005 #include "CondFormats/SiPixelObjects/interface/SiPixelQuality.h"
0006 #include "DataFormats/Math/interface/deltaPhi.h"
0007 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0008 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0009
0010 #include "TrackingTools/TransientTrackingRecHit/interface/SeedingLayerSetsLooper.h"
0011
0012 #include <fstream>
0013 #include <queue>
0014 #include <algorithm>
0015
0016 std::ostream& operator<<(std::ostream& os, SeedingLayerSetsBuilder::SeedingLayerId layer) {
0017 if (std::get<0>(layer) == GeomDetEnumerators::PixelBarrel) {
0018 os << "BPix";
0019 } else {
0020 os << "FPix";
0021 }
0022 os << std::get<2>(layer);
0023 if (std::get<1>(layer) == TrackerDetSide::PosEndcap) {
0024 os << "_pos";
0025 } else if (std::get<1>(layer) == TrackerDetSide::NegEndcap) {
0026 os << "_neg";
0027 }
0028 return os;
0029 }
0030
0031 namespace {
0032 using LayerPair = std::pair<SeedingLayerSetsBuilder::SeedingLayerId, SeedingLayerSetsBuilder::SeedingLayerId>;
0033 using ActiveLayerSetToInactiveSetsMap = std::map<LayerPair, edm::VecArray<LayerPair, 5>>;
0034 using Stream = std::stringstream;
0035 using Span_t = std::pair<float, float>;
0036
0037 ActiveLayerSetToInactiveSetsMap createActiveToInactiveMap() {
0038 ActiveLayerSetToInactiveSetsMap map;
0039
0040 auto bpix = [](int layer) {
0041 return SeedingLayerSetsBuilder::SeedingLayerId(GeomDetEnumerators::PixelBarrel, TrackerDetSide::Barrel, layer);
0042 };
0043 auto fpix_pos = [](int disk) {
0044 return SeedingLayerSetsBuilder::SeedingLayerId(GeomDetEnumerators::PixelEndcap, TrackerDetSide::PosEndcap, disk);
0045 };
0046 auto fpix_neg = [](int disk) {
0047 return SeedingLayerSetsBuilder::SeedingLayerId(GeomDetEnumerators::PixelEndcap, TrackerDetSide::NegEndcap, disk);
0048 };
0049
0050 auto add_permutations = [&](std::array<SeedingLayerSetsBuilder::SeedingLayerId, 4> quads) {
0051 do {
0052
0053
0054 if (quads[0] > quads[1] || quads[2] > quads[3])
0055 continue;
0056
0057 map[std::make_pair(quads[0], quads[1])].emplace_back(quads[2], quads[3]);
0058 } while (std::next_permutation(quads.begin(), quads.end()));
0059 };
0060
0061
0062 add_permutations({{bpix(1), bpix(2), bpix(3), bpix(4)}});
0063
0064
0065 add_permutations({{bpix(1), bpix(2), bpix(3), fpix_pos(1)}});
0066 add_permutations({{bpix(1), bpix(2), bpix(3), fpix_neg(1)}});
0067
0068
0069 add_permutations({{bpix(1), bpix(2), fpix_pos(1), fpix_pos(2)}});
0070 add_permutations({{bpix(1), bpix(2), fpix_neg(1), fpix_neg(2)}});
0071
0072
0073 add_permutations({{bpix(1), fpix_pos(1), fpix_pos(2), fpix_pos(3)}});
0074 add_permutations({{bpix(1), fpix_neg(1), fpix_neg(2), fpix_neg(3)}});
0075
0076 #ifdef EDM_ML_DEBUG
0077 LogDebug("PixelInactiveAreaFinder") << "Active to inactive mapping";
0078 for (const auto& elem : map) {
0079 std::stringstream ss;
0080 for (const auto& layerPair : elem.second) {
0081 ss << layerPair.first << "+" << layerPair.second << ",";
0082 }
0083 LogTrace("PixelInactiveAreaFinder") << " " << elem.first.first << "+" << elem.first.second << " => " << ss.str();
0084 }
0085 #endif
0086
0087 return map;
0088 }
0089
0090 void detGroupSpanInfo(const PixelInactiveAreaFinder::DetGroupSpan& cspan, Stream& ss) {
0091 using std::fixed;
0092 using std::left;
0093 using std::noshowpos;
0094 using std::right;
0095 using std::setfill;
0096 using std::setprecision;
0097 using std::setw;
0098 using std::showpos;
0099 std::string deli = "; ";
0100 ss << "subdetid:[" << cspan.subdetId << "]" << deli;
0101 if (cspan.subdetId == PixelSubdetector::PixelBarrel) {
0102 ss << "layer:[" << cspan.layer << "]" << deli;
0103 } else {
0104 ss << "disk:[" << cspan.disk << "]" << deli;
0105 }
0106 ss
0107
0108 << setprecision(16) << showpos << "phi:<" << right << setw(12) << cspan.phiSpan.first << "," << left << setw(12)
0109 << cspan.phiSpan.second << ">" << deli << "z:<" << right << setw(7) << cspan.zSpan.first << "," << left
0110 << setw(7) << cspan.zSpan.second << ">" << deli << noshowpos << "r:<" << right << setw(10) << cspan.rSpan.first
0111 << "," << left << setw(10) << cspan.rSpan.second << ">" << deli;
0112 }
0113 void printOverlapSpans(const PixelInactiveAreaFinder::InactiveAreas& areasLayers) {
0114 auto spansLayerSets = areasLayers.spansAndLayerSets(GlobalPoint(0, 0, 0), std::numeric_limits<float>::infinity());
0115
0116 Stream ss;
0117 for (auto const& spansLayers : spansLayerSets) {
0118 ss << "Overlapping detGroups:\n";
0119 for (auto const& cspan : spansLayers.first) {
0120 detGroupSpanInfo(cspan, ss);
0121 ss << std::endl;
0122 }
0123 }
0124 edm::LogPrint("PixelInactiveAreaFinder") << ss.str();
0125 }
0126
0127
0128 bool phiRangesOverlap(const Span_t& phiSpanA, const Span_t& phiSpanB) {
0129 float x1, x2, y1, y2;
0130 std::tie(x1, x2) = phiSpanA;
0131 std::tie(y1, y2) = phiSpanB;
0132
0133 if (x1 <= x2 && y1 <= y2) {
0134 return x1 <= y2 && y1 <= x2;
0135 } else if ((x1 > x2 && y1 <= y2) || (y1 > y2 && x1 <= x2)) {
0136 return y1 <= x2 || x1 <= y2;
0137 } else if (x1 > x2 && y1 > y2) {
0138 return true;
0139 } else {
0140 return false;
0141 }
0142 }
0143
0144
0145 bool phiMoreClockwise(float phiA, float phiB) {
0146
0147 return reco::deltaPhi(phiA, phiB) <= 0.f;
0148 }
0149 bool phiMoreCounterclockwise(float phiA, float phiB) {
0150
0151 return reco::deltaPhi(phiA, phiB) >= 0.f;
0152 }
0153
0154
0155 float zAxisIntersection(const float zrPointA[2], const float zrPointB[2]) {
0156 return (zrPointB[0] - zrPointA[0]) / (zrPointB[1] - zrPointA[1]) * (-zrPointA[1]) + zrPointA[0];
0157 }
0158 bool getZAxisOverlapRangeBarrel(const PixelInactiveAreaFinder::DetGroupSpan& cspanA,
0159 const PixelInactiveAreaFinder::DetGroupSpan& cspanB,
0160 std::pair<float, float>& range) {
0161 PixelInactiveAreaFinder::DetGroupSpan cspanUpper;
0162 PixelInactiveAreaFinder::DetGroupSpan cspanLower;
0163 if (cspanA.rSpan.second < cspanB.rSpan.first) {
0164 cspanLower = cspanA;
0165 cspanUpper = cspanB;
0166 } else if (cspanA.rSpan.first > cspanB.rSpan.second) {
0167 cspanUpper = cspanA;
0168 cspanLower = cspanB;
0169 } else {
0170 return false;
0171 }
0172 float lower = 0;
0173 float upper = 0;
0174 if (cspanUpper.zSpan.second < cspanLower.zSpan.first) {
0175
0176 const float pointUpperDetGroupL[2] = {cspanUpper.zSpan.second, cspanUpper.rSpan.second};
0177 const float pointLowerDetGroupL[2] = {cspanLower.zSpan.first, cspanLower.rSpan.first};
0178 lower = zAxisIntersection(pointUpperDetGroupL, pointLowerDetGroupL);
0179
0180 const float pointUpperDetGroupU[2] = {cspanUpper.zSpan.first, cspanUpper.rSpan.first};
0181 const float pointLowerDetGroupU[2] = {cspanLower.zSpan.second, cspanLower.rSpan.second};
0182 upper = zAxisIntersection(pointUpperDetGroupU, pointLowerDetGroupU);
0183 } else if (cspanUpper.zSpan.first <= cspanLower.zSpan.second && cspanLower.zSpan.first <= cspanUpper.zSpan.second) {
0184
0185 const float pointUpperDetGroupL[2] = {cspanUpper.zSpan.second, cspanUpper.rSpan.first};
0186 const float pointLowerDetGroupL[2] = {cspanLower.zSpan.first, cspanLower.rSpan.second};
0187 lower = zAxisIntersection(pointUpperDetGroupL, pointLowerDetGroupL);
0188
0189 const float pointUpperDetGroupU[2] = {cspanUpper.zSpan.first, cspanUpper.rSpan.first};
0190 const float pointLowerDetGroupU[2] = {cspanLower.zSpan.second, cspanLower.rSpan.second};
0191 upper = zAxisIntersection(pointUpperDetGroupU, pointLowerDetGroupU);
0192 } else if (cspanUpper.zSpan.first > cspanLower.zSpan.second) {
0193
0194 const float pointUpperDetGroupL[2] = {cspanUpper.zSpan.second, cspanUpper.rSpan.first};
0195 const float pointLowerDetGroupL[2] = {cspanLower.zSpan.first, cspanLower.rSpan.second};
0196 lower = zAxisIntersection(pointUpperDetGroupL, pointLowerDetGroupL);
0197
0198 const float pointUpperDetGroupU[2] = {cspanUpper.zSpan.first, cspanUpper.rSpan.second};
0199 const float pointLowerDetGroupU[2] = {cspanLower.zSpan.second, cspanLower.rSpan.first};
0200 upper = zAxisIntersection(pointUpperDetGroupU, pointLowerDetGroupU);
0201 } else {
0202
0203 return false;
0204 }
0205 range = std::pair<float, float>(lower, upper);
0206 return true;
0207 }
0208 bool getZAxisOverlapRangeEndcap(const PixelInactiveAreaFinder::DetGroupSpan& cspanA,
0209 const PixelInactiveAreaFinder::DetGroupSpan& cspanB,
0210 std::pair<float, float>& range) {
0211
0212 PixelInactiveAreaFinder::DetGroupSpan cspanNearer;
0213 PixelInactiveAreaFinder::DetGroupSpan cspanFurther;
0214 float lower = 0;
0215 float upper = 0;
0216 if (cspanA.zSpan.first < 0 && cspanB.zSpan.first < 0) {
0217 if (cspanA.zSpan.second < cspanB.zSpan.first) {
0218 cspanFurther = cspanA;
0219 cspanNearer = cspanB;
0220 } else if (cspanB.zSpan.second < cspanA.zSpan.first) {
0221 cspanFurther = cspanB;
0222 cspanNearer = cspanA;
0223 } else {
0224 #ifdef EDM_ML_DEBUG
0225 LogTrace("PixelInactiveAreaFinder") << "No overlap, same disk propably. Spans:";
0226 Stream ss;
0227 detGroupSpanInfo(cspanA, ss);
0228 ss << std::endl;
0229 detGroupSpanInfo(cspanB, ss);
0230 ss << std::endl;
0231 LogTrace("PixelInactiveAreaFinder") << ss.str();
0232 ss.str(std::string());
0233 LogTrace("PixelInactiveAreaFinder") << "**";
0234 #endif
0235 return false;
0236 }
0237 if (cspanFurther.rSpan.second > cspanNearer.rSpan.first) {
0238 const float pointA[2] = {cspanFurther.zSpan.second, cspanFurther.rSpan.second};
0239 const float pointB[2] = {cspanNearer.zSpan.first, cspanNearer.rSpan.first};
0240 lower = zAxisIntersection(pointA, pointB);
0241 if (cspanFurther.rSpan.first > cspanNearer.rSpan.second) {
0242 const float pointC[2] = {cspanFurther.zSpan.first, cspanFurther.rSpan.first};
0243 const float pointD[2] = {cspanNearer.zSpan.second, cspanFurther.rSpan.second};
0244 upper = zAxisIntersection(pointC, pointD);
0245 } else {
0246 upper = std::numeric_limits<float>::infinity();
0247 }
0248 } else {
0249 #ifdef EDM_ML_DEBUG
0250 LogTrace("PixelInactiveAreaFinder") << "No overlap, further detGroup is lower. Spans:";
0251 Stream ss;
0252 detGroupSpanInfo(cspanA, ss);
0253 ss << std::endl;
0254 detGroupSpanInfo(cspanB, ss);
0255 ss << std::endl;
0256 LogTrace("PixelInactiveAreaFinder") << ss.str();
0257 ss.str(std::string());
0258 LogTrace("PixelInactiveAreaFinder") << "**";
0259 #endif
0260 return false;
0261 }
0262 } else if (cspanA.zSpan.first > 0 && cspanB.zSpan.first > 0) {
0263 if (cspanA.zSpan.first > cspanB.zSpan.second) {
0264 cspanFurther = cspanA;
0265 cspanNearer = cspanB;
0266 } else if (cspanB.zSpan.first > cspanA.zSpan.second) {
0267 cspanFurther = cspanB;
0268 cspanNearer = cspanA;
0269 } else {
0270 #ifdef EDM_ML_DEBUG
0271 LogTrace("PixelInactiveAreaFinder") << "No overlap, same disk propably. Spans:";
0272 Stream ss;
0273 detGroupSpanInfo(cspanA, ss);
0274 ss << std::endl;
0275 detGroupSpanInfo(cspanB, ss);
0276 ss << std::endl;
0277 LogTrace("PixelInactiveAreaFinder") << ss.str();
0278 ss.str(std::string());
0279 LogTrace("PixelInactiveAreaFinder") << "**";
0280 #endif
0281 return false;
0282 }
0283 if (cspanFurther.rSpan.second > cspanNearer.rSpan.first) {
0284 const float pointA[2] = {cspanFurther.zSpan.first, cspanFurther.rSpan.second};
0285 const float pointB[2] = {cspanNearer.zSpan.second, cspanNearer.rSpan.first};
0286 upper = zAxisIntersection(pointA, pointB);
0287 if (cspanFurther.rSpan.first > cspanNearer.rSpan.second) {
0288 const float pointC[2] = {cspanFurther.zSpan.second, cspanFurther.rSpan.first};
0289 const float pointD[2] = {cspanNearer.zSpan.first, cspanFurther.rSpan.second};
0290 lower = zAxisIntersection(pointC, pointD);
0291 } else {
0292 lower = -std::numeric_limits<float>::infinity();
0293 }
0294 } else {
0295 #ifdef EDM_ML_DEBUG
0296 LogTrace("PixelInactiveAreaFinder") << "No overlap, further detGroup lower. Spans:";
0297 Stream ss;
0298 detGroupSpanInfo(cspanA, ss);
0299 ss << std::endl;
0300 detGroupSpanInfo(cspanB, ss);
0301 ss << std::endl;
0302 LogTrace("PixelInactiveAreaFinder") << ss.str();
0303 ss.str(std::string());
0304 LogTrace("PixelInactiveAreaFinder") << "**";
0305 #endif
0306 return false;
0307 }
0308 } else {
0309 #ifdef EDM_ML_DEBUG
0310 LogTrace("PixelInactiveAreaFinder") << "No overlap, different sides of z axis. Spans:";
0311 Stream ss;
0312 detGroupSpanInfo(cspanA, ss);
0313 ss << std::endl;
0314 detGroupSpanInfo(cspanB, ss);
0315 ss << std::endl;
0316 LogTrace("PixelInactiveAreaFinder") << ss.str();
0317 ss.str(std::string());
0318 LogTrace("PixelInactiveAreaFinder") << "**";
0319 #endif
0320 return false;
0321 }
0322 range = std::pair<float, float>(lower, upper);
0323 return true;
0324 }
0325
0326 bool getZAxisOverlapRangeBarrelEndcap(const PixelInactiveAreaFinder::DetGroupSpan& cspanBar,
0327 const PixelInactiveAreaFinder::DetGroupSpan& cspanEnd,
0328 std::pair<float, float>& range) {
0329 float lower = 0;
0330 float upper = 0;
0331 if (cspanEnd.rSpan.second > cspanBar.rSpan.first) {
0332 if (cspanEnd.zSpan.second < cspanBar.zSpan.first) {
0333
0334 const float pointA[2] = {cspanEnd.zSpan.second, cspanEnd.rSpan.second};
0335 const float pointB[2] = {cspanBar.zSpan.first, cspanBar.rSpan.first};
0336 lower = zAxisIntersection(pointA, pointB);
0337 if (cspanEnd.rSpan.first > cspanBar.rSpan.second) {
0338
0339 const float pointC[2] = {cspanEnd.zSpan.first, cspanEnd.rSpan.first};
0340 const float pointD[2] = {cspanBar.zSpan.second, cspanBar.rSpan.second};
0341 upper = zAxisIntersection(pointC, pointD);
0342 } else {
0343 upper = std::numeric_limits<float>::infinity();
0344 }
0345 } else if (cspanEnd.zSpan.first > cspanBar.zSpan.second) {
0346
0347 const float pointA[2] = {cspanEnd.zSpan.first, cspanEnd.rSpan.second};
0348 const float pointB[2] = {cspanBar.zSpan.second, cspanBar.rSpan.first};
0349 upper = zAxisIntersection(pointA, pointB);
0350 if (cspanEnd.rSpan.first > cspanBar.rSpan.second) {
0351 const float pointC[2] = {cspanEnd.zSpan.second, cspanEnd.rSpan.first};
0352 const float pointD[2] = {cspanBar.zSpan.first, cspanBar.rSpan.second};
0353 lower = zAxisIntersection(pointC, pointD);
0354 } else {
0355 lower = -std::numeric_limits<float>::infinity();
0356 }
0357 } else {
0358 return false;
0359 }
0360 } else {
0361 return false;
0362 }
0363 range = std::pair<float, float>(lower, upper);
0364 return true;
0365 }
0366 }
0367
0368 std::vector<
0369 std::pair<edm::VecArray<PixelInactiveAreaFinder::Area, 2>, std::vector<PixelInactiveAreaFinder::LayerSetIndex>>>
0370 PixelInactiveAreaFinder::InactiveAreas::areasAndLayerSets(const GlobalPoint& point, float zwidth) const {
0371 auto spansLayerSets = spansAndLayerSets(point, zwidth);
0372
0373
0374 std::vector<std::pair<VecArray2<Area>, std::vector<LayerSetIndex>>> ret;
0375 for (auto& item : spansLayerSets) {
0376 auto& innerSpan = item.first[0];
0377 auto& outerSpan = item.first[1];
0378 VecArray2<Area> areas;
0379 areas.emplace_back(innerSpan.rSpan.first,
0380 innerSpan.rSpan.second,
0381 innerSpan.phiSpan.first,
0382 innerSpan.phiSpan.second,
0383 innerSpan.zSpan.first,
0384 innerSpan.zSpan.second);
0385 areas.emplace_back(outerSpan.rSpan.first,
0386 outerSpan.rSpan.second,
0387 outerSpan.phiSpan.first,
0388 outerSpan.phiSpan.second,
0389 outerSpan.zSpan.first,
0390 outerSpan.zSpan.second);
0391 ret.emplace_back(areas, std::move(item.second));
0392 }
0393
0394 return ret;
0395 }
0396
0397 std::vector<std::pair<edm::VecArray<PixelInactiveAreaFinder::DetGroupSpan, 2>,
0398 std::vector<PixelInactiveAreaFinder::LayerSetIndex>>>
0399 PixelInactiveAreaFinder::InactiveAreas::spansAndLayerSets(const GlobalPoint& point, float zwidth) const {
0400
0401 const float zmin = point.z() - zwidth;
0402 const float zmax = point.z() + zwidth;
0403
0404 std::vector<std::pair<VecArray2<DetGroupSpan>, std::vector<LayerSetIndex>>> ret;
0405
0406 LogDebug("PixelInactiveAreaFinder") << "Origin at " << point.x() << "," << point.y() << "," << point.z()
0407 << " z half width " << zwidth;
0408
0409 for (LayerSetIndex i = 0, end = inactiveLayerPairIndices_->size(); i < end; ++i) {
0410 const auto& layerIdxPair = (*inactiveLayerPairIndices_)[i];
0411 const auto& innerSpans = inactiveSpans_[layerIdxPair.first];
0412 const auto& outerSpans = inactiveSpans_[layerIdxPair.second];
0413
0414 for (const auto& innerSpan : innerSpans) {
0415 for (const auto& outerSpan : outerSpans) {
0416 if (phiRangesOverlap(innerSpan.phiSpan, outerSpan.phiSpan)) {
0417 std::pair<float, float> range(0, 0);
0418
0419 bool zOverlap = false;
0420 const auto innerDet = std::get<0>((*inactiveLayers_)[layerIdxPair.first]);
0421 const auto outerDet = std::get<0>((*inactiveLayers_)[layerIdxPair.second]);
0422 if (innerDet == GeomDetEnumerators::PixelBarrel) {
0423 if (outerDet == GeomDetEnumerators::PixelBarrel)
0424 zOverlap = getZAxisOverlapRangeBarrel(innerSpan, outerSpan, range);
0425 else
0426 zOverlap = getZAxisOverlapRangeBarrelEndcap(innerSpan, outerSpan, range);
0427 } else {
0428 if (outerDet == GeomDetEnumerators::PixelEndcap)
0429 zOverlap = getZAxisOverlapRangeEndcap(innerSpan, outerSpan, range);
0430 else
0431 throw cms::Exception("LogicError") << "Forward->barrel transition is not supported";
0432 }
0433
0434 if (zOverlap && zmin <= range.second && range.first <= zmax) {
0435 #ifdef EDM_ML_DEBUG
0436 Stream ss;
0437 for (auto ind : (*layerSetIndexInactiveToActive_)[i]) {
0438 ss << ind << ",";
0439 }
0440 ss << "\n ";
0441 detGroupSpanInfo(innerSpan, ss);
0442 ss << "\n ";
0443 detGroupSpanInfo(outerSpan, ss);
0444 LogTrace("PixelInactiveAreaFinder") << " adding areas for active layer sets " << ss.str();
0445 #endif
0446 VecArray2<DetGroupSpan> vec;
0447 vec.emplace_back(innerSpan);
0448 vec.emplace_back(outerSpan);
0449 ret.emplace_back(std::move(vec), (*layerSetIndexInactiveToActive_)[i]);
0450 }
0451 }
0452 }
0453 }
0454 }
0455
0456 return ret;
0457 }
0458
0459 PixelInactiveAreaFinder::PixelInactiveAreaFinder(
0460 const edm::ParameterSet& iConfig,
0461 const std::vector<SeedingLayerSetsBuilder::SeedingLayerId>& seedingLayers,
0462 const SeedingLayerSetsLooper& seedingLayerSetsLooper,
0463 edm::ConsumesCollector&& iC)
0464 : debug_(iConfig.getUntrackedParameter<bool>("debug")),
0465 createPlottingFiles_(iConfig.getUntrackedParameter<bool>("createPlottingFiles")),
0466 ignoreSingleFPixPanelModules_(iConfig.getParameter<bool>("ignoreSingleFPixPanelModules")),
0467 inactivePixelDetectorTokens_(
0468 edm::vector_transform(iConfig.getParameter<std::vector<edm::InputTag>>("inactivePixelDetectorLabels"),
0469 [&](const auto& tag) { return iC.consumes<DetIdCollection>(tag); })),
0470 badPixelFEDChannelsTokens_(
0471 edm::vector_transform(iConfig.getParameter<std::vector<edm::InputTag>>("badPixelFEDChannelCollectionLabels"),
0472 [&](const auto& tag) { return iC.consumes<PixelFEDChannelCollection>(tag); })),
0473 trackerGeometryToken_(iC.esConsumes()),
0474 trackerTopologyToken_(iC.esConsumes()),
0475 pixelQualityToken_(iC.esConsumes()) {
0476 #ifdef EDM_ML_DEBUG
0477 for (const auto& layer : seedingLayers) {
0478 LogTrace("PixelInactiveAreaFinder") << "Input layer subdet " << std::get<0>(layer) << " side "
0479 << static_cast<unsigned int>(std::get<1>(layer)) << " layer "
0480 << std::get<2>(layer);
0481 }
0482 #endif
0483
0484 auto findOrAdd = [&](SeedingLayerId layer) -> unsigned short {
0485 auto found = std::find(inactiveLayers_.cbegin(), inactiveLayers_.cend(), layer);
0486 if (found == inactiveLayers_.cend()) {
0487 auto ret = inactiveLayers_.size();
0488 inactiveLayers_.push_back(layer);
0489 return ret;
0490 }
0491 return std::distance(inactiveLayers_.cbegin(), found);
0492 };
0493
0494
0495 const auto activeToInactiveMap = createActiveToInactiveMap();
0496
0497
0498
0499 LayerSetIndex i = 0;
0500 for (const auto& layerSet : seedingLayerSetsLooper.makeRange(seedingLayers)) {
0501 assert(layerSet.size() == 2);
0502 auto found = activeToInactiveMap.find(std::make_pair(layerSet[0], layerSet[1]));
0503 if (found == activeToInactiveMap.end()) {
0504 throw cms::Exception("Configuration")
0505 << "Encountered layer pair " << layerSet[0] << "+" << layerSet[1]
0506 << " not found from the internal 'active layer pairs' to 'inactive layer pairs' mapping; either fix the "
0507 "input or the mapping (in PixelInactiveAreaFinder.cc)";
0508 }
0509
0510 LogTrace("PixelInactiveAreaFinder") << "Input layer set " << layerSet[0] << "+" << layerSet[1];
0511 for (const auto& inactiveLayerSet : found->second) {
0512 auto innerInd = findOrAdd(inactiveLayerSet.first);
0513 auto outerInd = findOrAdd(inactiveLayerSet.second);
0514
0515 auto found = std::find(
0516 inactiveLayerSetIndices_.cbegin(), inactiveLayerSetIndices_.cend(), std::make_pair(innerInd, outerInd));
0517 if (found == inactiveLayerSetIndices_.end()) {
0518 inactiveLayerSetIndices_.emplace_back(innerInd, outerInd);
0519 layerSetIndexInactiveToActive_.push_back(std::vector<LayerSetIndex>{i});
0520 } else {
0521 layerSetIndexInactiveToActive_.at(std::distance(inactiveLayerSetIndices_.cbegin(), found))
0522 .push_back(i);
0523 }
0524
0525 LogTrace("PixelInactiveAreaFinder")
0526 << " inactive layer set " << inactiveLayerSet.first << "+" << inactiveLayerSet.second;
0527 }
0528
0529 ++i;
0530 }
0531
0532 #ifdef EDM_ML_DEBUG
0533 LogDebug("PixelInactiveAreaFinder") << "All inactive layer sets";
0534 for (const auto& idxPair : inactiveLayerSetIndices_) {
0535 LogTrace("PixelInactiveAreaFinder") << " " << inactiveLayers_[idxPair.first] << "+"
0536 << inactiveLayers_[idxPair.second];
0537 }
0538 #endif
0539 }
0540
0541 void PixelInactiveAreaFinder::fillDescriptions(edm::ParameterSetDescription& desc) {
0542 desc.add<std::vector<edm::InputTag>>("inactivePixelDetectorLabels",
0543 std::vector<edm::InputTag>{{edm::InputTag("siPixelDigis")}})
0544 ->setComment("One or more DetIdCollections of modules to mask on the fly for a given event");
0545 desc.add<std::vector<edm::InputTag>>("badPixelFEDChannelCollectionLabels",
0546 std::vector<edm::InputTag>{{edm::InputTag("siPixelDigis")}})
0547 ->setComment("One or more PixelFEDChannelCollections of modules+ROCs to mask on the fly for a given event");
0548 desc.add<bool>("ignoreSingleFPixPanelModules", false);
0549
0550 desc.addUntracked<bool>("debug", false);
0551 desc.addUntracked<bool>("createPlottingFiles", false);
0552 }
0553
0554 PixelInactiveAreaFinder::InactiveAreas PixelInactiveAreaFinder::inactiveAreas(const edm::Event& iEvent,
0555 const edm::EventSetup& iSetup) {
0556
0557 trackerGeometry_ = &iSetup.getData(trackerGeometryToken_);
0558 trackerTopology_ = &iSetup.getData(trackerTopologyToken_);
0559
0560
0561 updatePixelDets(iSetup);
0562 getBadPixelDets(iEvent, iSetup);
0563
0564
0565 if (createPlottingFiles_) {
0566 createPlottingFiles();
0567 }
0568
0569
0570
0571 DetGroupSpanContainerPair cspans = detGroupSpans();
0572
0573
0574
0575 std::vector<DetGroupSpanContainer> spans(inactiveLayers_.size());
0576
0577 auto doWork = [&](const DetGroupSpanContainer& container) {
0578 for (const auto& span : container) {
0579 const auto subdet = span.subdetId == PixelSubdetector::PixelBarrel ? GeomDetEnumerators::PixelBarrel
0580 : GeomDetEnumerators::PixelEndcap;
0581 const auto side = (subdet == GeomDetEnumerators::PixelBarrel
0582 ? TrackerDetSide::Barrel
0583 : (span.zSpan.first < 0 ? TrackerDetSide::NegEndcap : TrackerDetSide::PosEndcap));
0584 const auto layer = subdet == GeomDetEnumerators::PixelBarrel ? span.layer : span.disk;
0585 auto found = std::find(inactiveLayers_.begin(), inactiveLayers_.end(), SeedingLayerId(subdet, side, layer));
0586 if (found != inactiveLayers_.end()) {
0587 spans[std::distance(inactiveLayers_.begin(), found)].push_back(span);
0588 }
0589 }
0590 };
0591 doWork(cspans.first);
0592 doWork(cspans.second);
0593
0594 auto ret =
0595 InactiveAreas(&inactiveLayers_, std::move(spans), &inactiveLayerSetIndices_, &layerSetIndexInactiveToActive_);
0596
0597 if (debug_) {
0598 printOverlapSpans(ret);
0599 }
0600
0601 return ret;
0602 }
0603
0604
0605 void PixelInactiveAreaFinder::updatePixelDets(const edm::EventSetup& iSetup) {
0606 if (!geometryWatcher_.check(iSetup))
0607 return;
0608
0609
0610 {
0611
0612 if (trackerGeometry_->numberOfLayers(PixelSubdetector::PixelBarrel) != 4) {
0613 throw cms::Exception("NotImplemented")
0614 << "This module supports only a detector with 4 pixel barrel layers, the current geometry has "
0615 << trackerGeometry_->numberOfLayers(PixelSubdetector::PixelBarrel);
0616 }
0617 if (trackerGeometry_->numberOfLayers(PixelSubdetector::PixelEndcap) != 3) {
0618 throw cms::Exception("NotImplemented")
0619 << "This module supports only a detector with 3 pixel forward disks, the current geometry has "
0620 << trackerGeometry_->numberOfLayers(PixelSubdetector::PixelEndcap);
0621 }
0622
0623 std::array<std::array<unsigned short, 100>, 4> counts = {};
0624 for (const auto& det : trackerGeometry_->detsPXB()) {
0625 const auto layer = trackerTopology_->layer(det->geographicalId());
0626 const auto ladder = trackerTopology_->pxbLadder(det->geographicalId());
0627 if (ladder >= 100) {
0628 throw cms::Exception("LogicError")
0629 << "Got a ladder with number " << ladder
0630 << " while the expected maximum was 100; either something is wrong or the maximum has to be increased.";
0631 }
0632 counts[layer - 1][ladder - 1] += 1;
0633 }
0634
0635
0636
0637 nModulesPerLadder = counts[0][0];
0638 if (nModulesPerLadder == 0) {
0639 throw cms::Exception("LogicError") << "Ladder 1 of layer 1 has 0 modules, something fishy is going on.";
0640 }
0641
0642 LogDebug("PixelInactiveAreaFinder") << "Number of modules per ladder " << nModulesPerLadder
0643 << "; below are number of ladders per layer";
0644
0645
0646 for (unsigned layer = 0; layer < 4; ++layer) {
0647 nBPixLadders[layer] =
0648 std::count_if(counts[layer].begin(), counts[layer].end(), [](unsigned short val) { return val > 0; });
0649 LogTrace("PixelInactiveAreaFinder")
0650 << "BPix layer " << (layer + 1) << " has " << nBPixLadders[layer] << " ladders";
0651
0652 auto fail = std::count_if(counts[layer].begin(), counts[layer].end(), [&](unsigned short val) {
0653 return val != nModulesPerLadder && val > 0;
0654 });
0655 if (fail != 0) {
0656 throw cms::Exception("LogicError")
0657 << "Layer " << (layer + 1) << " had " << fail
0658 << " ladders whose number of modules/ladder differed from the ladder 1 of layer 1 (" << nModulesPerLadder
0659 << "). Something fishy is going on.";
0660 }
0661 }
0662 }
0663
0664
0665 if (!createPlottingFiles_)
0666 return;
0667
0668 pixelDetsBarrel_.clear();
0669 pixelDetsEndcap_.clear();
0670
0671 for (auto const& geomDetPtr : trackerGeometry_->detsPXB()) {
0672 if (geomDetPtr->geographicalId().subdetId() == PixelSubdetector::PixelBarrel) {
0673 pixelDetsBarrel_.push_back(geomDetPtr->geographicalId().rawId());
0674 }
0675 }
0676 for (auto const& geomDetPtr : trackerGeometry_->detsPXF()) {
0677 if (geomDetPtr->geographicalId().subdetId() == PixelSubdetector::PixelEndcap) {
0678 pixelDetsEndcap_.push_back(geomDetPtr->geographicalId().rawId());
0679 }
0680 }
0681 std::sort(pixelDetsBarrel_.begin(), pixelDetsBarrel_.end());
0682 std::sort(pixelDetsEndcap_.begin(), pixelDetsEndcap_.end());
0683 }
0684 void PixelInactiveAreaFinder::getBadPixelDets(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0685 auto addDetId = [&](const auto id) {
0686 const auto detid = DetId(id);
0687 const auto subdet = detid.subdetId();
0688 if (subdet == PixelSubdetector::PixelBarrel) {
0689 badPixelDetsBarrel_.push_back(detid.rawId());
0690 } else if (subdet == PixelSubdetector::PixelEndcap) {
0691 badPixelDetsEndcap_.push_back(detid.rawId());
0692 }
0693 };
0694
0695
0696 auto const& pixelQuality = iSetup.getData(pixelQualityToken_);
0697
0698 for (auto const& disabledModule : pixelQuality.getBadComponentList()) {
0699 addDetId(disabledModule.DetID);
0700 }
0701
0702
0703 for (const auto& token : inactivePixelDetectorTokens_) {
0704 edm::Handle<DetIdCollection> detIds;
0705 iEvent.getByToken(token, detIds);
0706 for (const auto& id : *detIds) {
0707 addDetId(id);
0708 }
0709 }
0710
0711
0712
0713 for (const auto& token : badPixelFEDChannelsTokens_) {
0714 edm::Handle<PixelFEDChannelCollection> pixelFEDChannelCollectionHandle;
0715 iEvent.getByToken(token, pixelFEDChannelCollectionHandle);
0716 for (const auto& disabledChannels : *pixelFEDChannelCollectionHandle) {
0717 addDetId(disabledChannels.detId());
0718 }
0719 }
0720
0721
0722 std::sort(badPixelDetsBarrel_.begin(), badPixelDetsBarrel_.end());
0723 std::sort(badPixelDetsEndcap_.begin(), badPixelDetsEndcap_.end());
0724 badPixelDetsBarrel_.erase(std::unique(badPixelDetsBarrel_.begin(), badPixelDetsBarrel_.end()),
0725 badPixelDetsBarrel_.end());
0726 badPixelDetsEndcap_.erase(std::unique(badPixelDetsEndcap_.begin(), badPixelDetsEndcap_.end()),
0727 badPixelDetsEndcap_.end());
0728 }
0729
0730 void PixelInactiveAreaFinder::detInfo(const det_t& det, Stream& ss) {
0731 using std::fixed;
0732 using std::left;
0733 using std::noshowpos;
0734 using std::right;
0735 using std::setfill;
0736 using std::setprecision;
0737 using std::setw;
0738 using std::showpos;
0739 using std::tie;
0740 std::string deli = "; ";
0741 ss << "id:[" << det << "]" << deli;
0742 ss << "subdetid:[" << DetId(det).subdetId() << "]" << deli;
0743 if (DetId(det).subdetId() == PixelSubdetector::PixelBarrel) {
0744 unsigned int layer = trackerTopology_->pxbLayer(DetId(det));
0745 unsigned int ladder = trackerTopology_->pxbLadder(DetId(det));
0746 unsigned int module = trackerTopology_->pxbModule(DetId(det));
0747 ss << "layer:[" << layer << "]" << deli << "ladder:[" << right << setw(2) << ladder << "]" << deli << "module:["
0748 << module << "]" << deli;
0749 } else if (DetId(det).subdetId() == PixelSubdetector::PixelEndcap) {
0750 unsigned int disk = trackerTopology_->pxfDisk(DetId(det));
0751 unsigned int blade = trackerTopology_->pxfBlade(DetId(det));
0752 unsigned int panel = trackerTopology_->pxfPanel(DetId(det));
0753 ss << left << setw(6) << "disk:"
0754 << "[" << right << disk << "]" << deli << left << setw(7) << "blade:"
0755 << "[" << setw(2) << right << blade << "]" << deli << left << setw(7) << "panel:"
0756 << "[" << right << panel << "]" << deli;
0757 }
0758 float phiA, phiB, zA, zB, rA, rB;
0759 auto detSurface = trackerGeometry_->idToDet(DetId(det))->surface();
0760 tie(phiA, phiB) = detSurface.phiSpan();
0761 tie(zA, zB) = detSurface.zSpan();
0762 tie(rA, rB) = detSurface.rSpan();
0763 ss << setprecision(16) << fixed << showpos << setfill(' ') << "phi:[" << right << setw(12) << phiA << "," << left
0764 << setw(12) << phiB << "]" << deli << "z:[" << right << setw(7) << zA << "," << left << setw(7) << zB << "]"
0765 << deli << noshowpos << "r:[" << right << setw(10) << rA << "," << left << setw(10) << rB << "]" << deli;
0766 }
0767 void PixelInactiveAreaFinder::printPixelDets() {
0768 edm::LogPrint("PixelInactiveAreaFinder") << "Barrel detectors:";
0769 Stream ss;
0770 for (auto const& det : pixelDetsBarrel_) {
0771 detInfo(det, ss);
0772 edm::LogPrint("PixelInactiveAreaFinder") << ss.str();
0773 ss.str(std::string());
0774 }
0775 edm::LogPrint("PixelInactiveAreaFinder") << "Endcap detectors;";
0776 for (auto const& det : pixelDetsEndcap_) {
0777 detInfo(det, ss);
0778 edm::LogPrint("PixelInactiveAreaFinder") << ss.str();
0779 ss.str(std::string());
0780 }
0781 }
0782 void PixelInactiveAreaFinder::printBadPixelDets() {
0783 edm::LogPrint("PixelInactiveAreaFinder") << "Bad barrel detectors:";
0784 Stream ss;
0785 for (auto const& det : badPixelDetsBarrel_) {
0786 detInfo(det, ss);
0787 edm::LogPrint("PixelInactiveAreaFinder") << ss.str();
0788 ss.str(std::string());
0789 }
0790 edm::LogPrint("PixelInactiveAreaFinder") << "Endcap detectors;";
0791 for (auto const& det : badPixelDetsEndcap_) {
0792 detInfo(det, ss);
0793 edm::LogPrint("PixelInactiveAreaFinder") << ss.str();
0794 ss.str(std::string());
0795 }
0796 }
0797 void PixelInactiveAreaFinder::printBadDetGroups() {
0798 DetGroupContainer badDetGroupsBar = badDetGroupsBarrel();
0799 DetGroupContainer badDetGroupsEnd = badDetGroupsEndcap();
0800 Stream ss;
0801 for (auto const& detGroup : badDetGroupsBar) {
0802 ss << std::setfill(' ') << std::left << std::setw(16) << "DetGroup:";
0803 DetGroupSpan cspan;
0804 getPhiSpanBarrel(detGroup, cspan);
0805 getZSpan(detGroup, cspan);
0806 getRSpan(detGroup, cspan);
0807 detGroupSpanInfo(cspan, ss);
0808 ss << std::endl;
0809 for (auto const& det : detGroup) {
0810 detInfo(det, ss);
0811 ss << std::endl;
0812 }
0813 ss << std::endl;
0814 }
0815 for (auto const& detGroup : badDetGroupsEnd) {
0816 ss << std::setfill(' ') << std::left << std::setw(16) << "DetGroup:";
0817 DetGroupSpan cspan;
0818 getPhiSpanEndcap(detGroup, cspan);
0819 getZSpan(detGroup, cspan);
0820 getRSpan(detGroup, cspan);
0821 detGroupSpanInfo(cspan, ss);
0822 ss << std::endl;
0823 for (auto const& det : detGroup) {
0824 detInfo(det, ss);
0825 ss << std::endl;
0826 }
0827 ss << std::endl;
0828 }
0829 edm::LogPrint("PixelInactiveAreaFinder") << ss.str();
0830 }
0831 void PixelInactiveAreaFinder::printBadDetGroupSpans() {
0832 DetGroupSpanContainerPair cspans = detGroupSpans();
0833 Stream ss;
0834 for (auto const& cspan : cspans.first) {
0835 detGroupSpanInfo(cspan, ss);
0836 ss << std::endl;
0837 }
0838 for (auto const& cspan : cspans.second) {
0839 detGroupSpanInfo(cspan, ss);
0840 ss << std::endl;
0841 }
0842 edm::LogPrint("PixelInactiveAreaFinder") << ss.str();
0843 }
0844 void PixelInactiveAreaFinder::createPlottingFiles() {
0845
0846 Stream ss;
0847 std::ofstream fsDet("DETECTORS.txt");
0848 for (auto const& det : pixelDetsBarrel_) {
0849 detInfo(det, ss);
0850 ss << std::endl;
0851 }
0852 edm::LogPrint("PixelInactiveAreaFinder") << "Endcap detectors;";
0853 for (auto const& det : pixelDetsEndcap_) {
0854 detInfo(det, ss);
0855 ss << std::endl;
0856 }
0857 fsDet << ss.rdbuf();
0858 ss.str(std::string());
0859
0860 std::ofstream fsBadDet("BADDETECTORS.txt");
0861 for (auto const& det : badPixelDetsBarrel_) {
0862 detInfo(det, ss);
0863 ss << std::endl;
0864 }
0865 for (auto const& det : badPixelDetsEndcap_) {
0866 detInfo(det, ss);
0867 ss << std::endl;
0868 }
0869 fsBadDet << ss.rdbuf();
0870 ss.str(std::string());
0871
0872 std::ofstream fsSpans("DETGROUPSPANS.txt");
0873 DetGroupSpanContainerPair cspans = detGroupSpans();
0874 for (auto const& cspan : cspans.first) {
0875 detGroupSpanInfo(cspan, ss);
0876 ss << std::endl;
0877 }
0878 for (auto const& cspan : cspans.second) {
0879 detGroupSpanInfo(cspan, ss);
0880 ss << std::endl;
0881 }
0882 fsSpans << ss.rdbuf();
0883 ss.str(std::string());
0884 }
0885
0886 bool PixelInactiveAreaFinder::detWorks(det_t det) {
0887 return std::find(badPixelDetsBarrel_.begin(), badPixelDetsBarrel_.end(), det) == badPixelDetsBarrel_.end() &&
0888 std::find(badPixelDetsEndcap_.begin(), badPixelDetsEndcap_.end(), det) == badPixelDetsEndcap_.end();
0889 }
0890 PixelInactiveAreaFinder::DetGroup PixelInactiveAreaFinder::badAdjecentDetsBarrel(const det_t& det) {
0891 using std::remove_if;
0892
0893 DetGroup adj;
0894 auto const tTopo = trackerTopology_;
0895 auto const& detId = DetId(det);
0896 unsigned int layer = tTopo->pxbLayer(detId);
0897 unsigned int ladder = tTopo->pxbLadder(detId);
0898 unsigned int module = tTopo->pxbModule(detId);
0899 unsigned int nLads = nBPixLadders[layer];
0900
0901 adj.push_back(tTopo->pxbDetId(layer, ((ladder - 1) + 1) % nLads + 1, module)());
0902 adj.push_back(tTopo->pxbDetId(layer, ((ladder - 1) - 1 + nLads) % nLads + 1, module)());
0903
0904 if (module == 1) {
0905 adj.push_back(tTopo->pxbDetId(layer, ladder, module + 1)());
0906 } else if (module == nModulesPerLadder) {
0907 adj.push_back(tTopo->pxbDetId(layer, ladder, module - 1)());
0908 } else {
0909 adj.push_back(tTopo->pxbDetId(layer, ladder, module + 1)());
0910 adj.push_back(tTopo->pxbDetId(layer, ladder, module - 1)());
0911 }
0912
0913 adj.erase(remove_if(adj.begin(), adj.end(), [&](auto c) { return this->detWorks(c); }), adj.end());
0914 return adj;
0915 }
0916 PixelInactiveAreaFinder::DetGroup PixelInactiveAreaFinder::badAdjecentDetsEndcap(const det_t& det) {
0917
0918 using std::ignore;
0919 using std::tie;
0920 DetGroup adj;
0921 Span_t phiSpan, phiSpanComp;
0922 float z, zComp;
0923 unsigned int disk, diskComp;
0924 auto const& detSurf = trackerGeometry_->idToDet(DetId(det))->surface();
0925 phiSpan = detSurf.phiSpan();
0926 tie(z, ignore) = detSurf.zSpan();
0927 disk = trackerTopology_->pxfDisk(DetId(det));
0928
0929 for (auto const& detComp : badPixelDetsEndcap_) {
0930 auto const& detIdComp = DetId(detComp);
0931 auto const& detSurfComp = trackerGeometry_->idToDet(detIdComp)->surface();
0932 diskComp = trackerTopology_->pxfDisk(detIdComp);
0933 phiSpanComp = detSurfComp.phiSpan();
0934 tie(zComp, ignore) = detSurfComp.zSpan();
0935 if (det != detComp && disk == diskComp && z * zComp > 0 && phiRangesOverlap(phiSpan, phiSpanComp)) {
0936 adj.push_back(detComp);
0937 }
0938 }
0939 return adj;
0940 }
0941 PixelInactiveAreaFinder::DetGroup PixelInactiveAreaFinder::reachableDetGroup(const det_t& initDet,
0942 DetectorSet& foundDets) {
0943 DetGroup reachableDetGroup;
0944 std::queue<det_t> workQueue;
0945 det_t workDet;
0946 DetGroup badAdjDets;
0947 foundDets.insert(initDet);
0948 workQueue.push(initDet);
0949 reachableDetGroup.push_back(initDet);
0950 while (!workQueue.empty()) {
0951 workDet = workQueue.front();
0952 workQueue.pop();
0953 if (DetId(workDet).subdetId() == PixelSubdetector::PixelBarrel) {
0954 badAdjDets = this->badAdjecentDetsBarrel(workDet);
0955 } else if (DetId(workDet).subdetId() == PixelSubdetector::PixelEndcap) {
0956 badAdjDets = this->badAdjecentDetsEndcap(workDet);
0957 } else {
0958 badAdjDets = {};
0959 }
0960 for (auto const& badDet : badAdjDets) {
0961 if (foundDets.find(badDet) == foundDets.end()) {
0962 reachableDetGroup.push_back(badDet);
0963 foundDets.insert(badDet);
0964 workQueue.push(badDet);
0965 }
0966 }
0967 }
0968 return reachableDetGroup;
0969 }
0970 PixelInactiveAreaFinder::DetGroupContainer PixelInactiveAreaFinder::badDetGroupsBarrel() {
0971 DetGroupContainer detGroups;
0972 DetectorSet foundDets;
0973 for (auto const& badDet : badPixelDetsBarrel_) {
0974 if (foundDets.find(badDet) == foundDets.end()) {
0975 detGroups.push_back(this->reachableDetGroup(badDet, foundDets));
0976 }
0977 }
0978 return detGroups;
0979 }
0980 PixelInactiveAreaFinder::DetGroupContainer PixelInactiveAreaFinder::badDetGroupsEndcap() {
0981 DetGroupContainer detGroups;
0982 DetectorSet foundDets;
0983 for (auto const& badDet : badPixelDetsEndcap_) {
0984 if (foundDets.find(badDet) == foundDets.end()) {
0985 auto adjacentDets = this->reachableDetGroup(badDet, foundDets);
0986 if (ignoreSingleFPixPanelModules_ && adjacentDets.size() == 1) {
0987
0988
0989
0990 continue;
0991 }
0992 detGroups.push_back(adjacentDets);
0993 }
0994 }
0995 return detGroups;
0996 }
0997
0998 void PixelInactiveAreaFinder::getPhiSpanBarrel(const DetGroup& detGroup, DetGroupSpan& cspan) {
0999
1000 if (detGroup.empty()) {
1001 cspan = DetGroupSpan();
1002 return;
1003 } else {
1004 cspan.layer = trackerTopology_->pxbLayer(DetId(detGroup[0]));
1005 cspan.disk = 0;
1006 }
1007 using uint = unsigned int;
1008 using LadderSet = std::set<uint>;
1009 using LadVec = std::vector<uint>;
1010 LadderSet lads;
1011 for (auto const& det : detGroup) {
1012 lads.insert(trackerTopology_->pxbLadder(DetId(det)));
1013 }
1014 LadVec ladv(lads.begin(), lads.end());
1015 uint nLadders = nBPixLadders[cspan.layer];
1016
1017 uint i = 0;
1018 uint currentLadder = ladv[0];
1019 uint previousLadder = ladv[(ladv.size() + i - 1) % ladv.size()];
1020
1021 while ((nLadders + currentLadder - 1) % nLadders == previousLadder) {
1022 ++i;
1023 currentLadder = ladv[i % ladv.size()];
1024 previousLadder = ladv[(ladv.size() + i - 1) % ladv.size()];
1025 if (i == ladv.size()) {
1026 cspan.phiSpan.first = std::numeric_limits<float>::epsilon();
1027 cspan.phiSpan.second = -std::numeric_limits<float>::epsilon();
1028 return;
1029 }
1030 }
1031 uint startLadder = currentLadder;
1032 uint endLadder = previousLadder;
1033 auto detStart = trackerTopology_->pxbDetId(cspan.layer, startLadder, 1);
1034 auto detEnd = trackerTopology_->pxbDetId(cspan.layer, endLadder, 1);
1035 cspan.phiSpan.first = trackerGeometry_->idToDet(detStart)->surface().phiSpan().first;
1036 cspan.phiSpan.second = trackerGeometry_->idToDet(detEnd)->surface().phiSpan().second;
1037 }
1038 void PixelInactiveAreaFinder::getPhiSpanEndcap(const DetGroup& detGroup, DetGroupSpan& cspan) {
1039
1040
1041
1042
1043
1044
1045 Stream ss;
1046 bool found = false;
1047 auto const tGeom = trackerGeometry_;
1048 DetGroup::const_iterator startDetIter = detGroup.begin();
1049 Span_t phiSpan, phiSpanComp;
1050 unsigned int counter = 0;
1051 while (!found) {
1052 phiSpan = tGeom->idToDet(DetId(*startDetIter))->surface().phiSpan();
1053 for (DetGroup::const_iterator compDetIter = detGroup.begin(); compDetIter != detGroup.end(); ++compDetIter) {
1054 phiSpanComp = tGeom->idToDet(DetId(*compDetIter))->surface().phiSpan();
1055 if (phiRangesOverlap(phiSpan, phiSpanComp) && phiMoreClockwise(phiSpanComp.first, phiSpan.first) &&
1056 startDetIter != compDetIter) {
1057 ++counter;
1058 if (counter > detGroup.size()) {
1059 cspan.phiSpan.first = std::numeric_limits<float>::epsilon();
1060 cspan.phiSpan.second = -std::numeric_limits<float>::epsilon();
1061 return;
1062 }
1063 startDetIter = compDetIter;
1064 break;
1065 } else if (compDetIter == detGroup.end() - 1) {
1066 found = true;
1067 }
1068 }
1069 }
1070 cspan.phiSpan.first = phiSpan.first;
1071
1072 found = false;
1073 DetGroup::const_iterator endDetIter = detGroup.begin();
1074 counter = 0;
1075 while (!found) {
1076 phiSpan = tGeom->idToDet(DetId(*endDetIter))->surface().phiSpan();
1077 for (DetGroup::const_iterator compDetIter = detGroup.begin(); compDetIter != detGroup.end(); ++compDetIter) {
1078 phiSpanComp = tGeom->idToDet(DetId(*compDetIter))->surface().phiSpan();
1079 if (phiRangesOverlap(phiSpan, phiSpanComp) && phiMoreCounterclockwise(phiSpanComp.second, phiSpan.second) &&
1080 endDetIter != compDetIter) {
1081 ++counter;
1082 if (counter > detGroup.size()) {
1083 cspan.phiSpan.first = std::numeric_limits<float>::epsilon();
1084 cspan.phiSpan.second = -std::numeric_limits<float>::epsilon();
1085 return;
1086 }
1087 endDetIter = compDetIter;
1088 break;
1089 } else if (compDetIter == detGroup.end() - 1) {
1090 found = true;
1091 }
1092 }
1093 }
1094 cspan.phiSpan.second = phiSpan.second;
1095 }
1096 void PixelInactiveAreaFinder::getZSpan(const DetGroup& detGroup, DetGroupSpan& cspan) {
1097 auto cmpFun = [this](det_t detA, det_t detB) {
1098 return trackerGeometry_->idToDet(DetId(detA))->surface().zSpan().first <
1099 trackerGeometry_->idToDet(DetId(detB))->surface().zSpan().first;
1100 };
1101
1102 auto minmaxIters = std::minmax_element(detGroup.begin(), detGroup.end(), cmpFun);
1103 cspan.zSpan.first = trackerGeometry_->idToDet(DetId(*(minmaxIters.first)))->surface().zSpan().first;
1104 cspan.zSpan.second = trackerGeometry_->idToDet(DetId(*(minmaxIters.second)))->surface().zSpan().second;
1105 }
1106 void PixelInactiveAreaFinder::getRSpan(const DetGroup& detGroup, DetGroupSpan& cspan) {
1107 auto cmpFun = [this](det_t detA, det_t detB) {
1108 return trackerGeometry_->idToDet(DetId(detA))->surface().rSpan().first <
1109 trackerGeometry_->idToDet(DetId(detB))->surface().rSpan().first;
1110 };
1111
1112 auto minmaxIters = std::minmax_element(detGroup.begin(), detGroup.end(), cmpFun);
1113 cspan.rSpan.first = trackerGeometry_->idToDet(DetId(*(minmaxIters.first)))->surface().rSpan().first;
1114 cspan.rSpan.second = trackerGeometry_->idToDet(DetId(*(minmaxIters.second)))->surface().rSpan().second;
1115 }
1116 void PixelInactiveAreaFinder::getSpan(const DetGroup& detGroup, DetGroupSpan& cspan) {
1117 auto firstDetIt = detGroup.begin();
1118 if (firstDetIt != detGroup.end()) {
1119 cspan.subdetId = DetId(*firstDetIt).subdetId();
1120 if (cspan.subdetId == 1) {
1121 cspan.layer = trackerTopology_->pxbLayer(DetId(*firstDetIt));
1122 cspan.disk = 0;
1123 getPhiSpanBarrel(detGroup, cspan);
1124 } else if (cspan.subdetId == 2) {
1125 cspan.disk = trackerTopology_->pxfDisk(DetId(*firstDetIt));
1126 cspan.layer = 0;
1127 getPhiSpanEndcap(detGroup, cspan);
1128 }
1129 getZSpan(detGroup, cspan);
1130 getRSpan(detGroup, cspan);
1131 }
1132 }
1133 PixelInactiveAreaFinder::DetGroupSpanContainerPair PixelInactiveAreaFinder::detGroupSpans() {
1134 DetGroupSpanContainer cspansBarrel;
1135 DetGroupSpanContainer cspansEndcap;
1136 DetGroupContainer badDetGroupsBar = badDetGroupsBarrel();
1137 DetGroupContainer badDetGroupsEnd = badDetGroupsEndcap();
1138 for (auto const& detGroup : badDetGroupsBar) {
1139 DetGroupSpan cspan;
1140 getSpan(detGroup, cspan);
1141 cspansBarrel.push_back(cspan);
1142 }
1143 for (auto const& detGroup : badDetGroupsEnd) {
1144 DetGroupSpan cspan;
1145 getSpan(detGroup, cspan);
1146 cspansEndcap.push_back(cspan);
1147 }
1148 return DetGroupSpanContainerPair(cspansBarrel, cspansEndcap);
1149 }