File indexing completed on 2024-04-06 12:14:52
0001 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0002 #include "FWCore/Framework/interface/EventSetup.h"
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "Geometry/HcalCommonData/interface/HcalDDDRecConstants.h"
0006 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0007 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0008 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0009 #include "Geometry/HcalTowerAlgo/interface/HcalGeometry.h"
0010 #include "Geometry/HcalTowerAlgo/interface/HcalFlexiHardcodeGeometryLoader.h"
0011 #include "Geometry/HcalTowerAlgo/interface/HcalHardcodeGeometryLoader.h"
0012 #include "Geometry/HcalTowerAlgo/interface/HcalTrigTowerGeometry.h"
0013 #include "DataFormats/HcalDetId/interface/HcalTrigTowerDetId.h"
0014 #include <iostream>
0015 #include <string>
0016
0017 class HcalGeometryTester : public edm::one::EDAnalyzer<> {
0018 public:
0019 explicit HcalGeometryTester(const edm::ParameterSet&);
0020 ~HcalGeometryTester(void) override {}
0021
0022 void beginJob() override {}
0023 void analyze(edm::Event const&, edm::EventSetup const&) override;
0024 void endJob() override {}
0025
0026 private:
0027 void testValidDetIds(CaloSubdetectorGeometry* geom,
0028 const HcalTopology& topo,
0029 DetId::Detector det,
0030 int subdet,
0031 const std::string& label);
0032 void testClosestCells(CaloSubdetectorGeometry* geom, const HcalTopology& top);
0033 void testClosestCell(const HcalDetId& detId, CaloSubdetectorGeometry* geom);
0034 void testTriggerGeometry(const HcalTopology& topology);
0035 void testFlexiValidDetIds(CaloSubdetectorGeometry* geom,
0036 const HcalTopology& topology,
0037 DetId::Detector det,
0038 int subdet,
0039 const std::string& label,
0040 std::vector<int>& dins);
0041 void testFlexiGeomHF(CaloSubdetectorGeometry* geom);
0042
0043 bool useOld_;
0044 edm::ESGetToken<HcalDDDRecConstants, HcalRecNumberingRecord> tok_ddrec_;
0045 edm::ESGetToken<HcalTopology, HcalRecNumberingRecord> tok_htopo_;
0046 };
0047
0048 HcalGeometryTester::HcalGeometryTester(const edm::ParameterSet& iConfig) {
0049 useOld_ = iConfig.getParameter<bool>("UseOldLoader");
0050 tok_ddrec_ = esConsumes<HcalDDDRecConstants, HcalRecNumberingRecord>();
0051 tok_htopo_ = esConsumes<HcalTopology, HcalRecNumberingRecord>();
0052 }
0053
0054 void HcalGeometryTester::analyze(const edm::Event& , const edm::EventSetup& iSetup) {
0055 const HcalDDDRecConstants hcons = iSetup.getData(tok_ddrec_);
0056 const HcalTopology topology = iSetup.getData(tok_htopo_);
0057 CaloSubdetectorGeometry* geom(nullptr);
0058 if (useOld_) {
0059 HcalHardcodeGeometryLoader m_loader;
0060 geom = m_loader.load(topology);
0061 } else {
0062 HcalFlexiHardcodeGeometryLoader m_loader;
0063 geom = m_loader.load(topology, hcons);
0064 }
0065
0066 testValidDetIds(geom, topology, DetId::Hcal, HcalBarrel, " BARREL ");
0067 testValidDetIds(geom, topology, DetId::Hcal, HcalEndcap, " ENDCAP ");
0068 testValidDetIds(geom, topology, DetId::Hcal, HcalOuter, " OUTER ");
0069 testValidDetIds(geom, topology, DetId::Hcal, HcalForward, " FORWARD ");
0070
0071 testTriggerGeometry(topology);
0072
0073 testClosestCells(geom, topology);
0074
0075 edm::LogVerbatim("HCalGeom") << "HcalGeometryTester::Test SLHC Hcal Geometry";
0076 std::vector<int> dins;
0077
0078 testFlexiValidDetIds(geom, topology, DetId::Hcal, HcalBarrel, " BARREL ", dins);
0079 testFlexiValidDetIds(geom, topology, DetId::Hcal, HcalEndcap, " ENDCAP ", dins);
0080 testFlexiValidDetIds(geom, topology, DetId::Hcal, HcalOuter, " OUTER ", dins);
0081 testFlexiValidDetIds(geom, topology, DetId::Hcal, HcalForward, " FORWARD ", dins);
0082
0083 testFlexiGeomHF(geom);
0084 }
0085
0086 void HcalGeometryTester::testValidDetIds(CaloSubdetectorGeometry* caloGeom,
0087 const HcalTopology& topology,
0088 DetId::Detector det,
0089 int subdet,
0090 const std::string& label) {
0091 std::stringstream s;
0092 s << label << " : \n";
0093 const std::vector<DetId>& idshb = caloGeom->getValidDetIds(det, subdet);
0094
0095 int counter = 0;
0096 for (std::vector<DetId>::const_iterator i = idshb.begin(); i != idshb.end(); i++, ++counter) {
0097 HcalDetId hid = (*i);
0098 s << counter << ": din " << topology.detId2denseId(*i) << ":" << hid;
0099 auto cell = caloGeom->getGeometry(*i);
0100 s << *cell << std::endl;
0101 }
0102
0103 s << "=== Total " << counter << " cells in " << label << std::endl;
0104 edm::LogVerbatim("HCalGeom") << s.str();
0105 }
0106
0107 void HcalGeometryTester::testClosestCells(CaloSubdetectorGeometry* g, const HcalTopology& topology) {
0108
0109 HcalDetId barrelDet1(HcalBarrel, 1, 1, 1);
0110 HcalDetId barrelDet2(HcalBarrel, 16, 50, 1);
0111 HcalDetId barrelDet3(HcalBarrel, -6, 20, 3);
0112 HcalDetId barrelDet4(HcalBarrel, -12, 33, 4);
0113 HcalDetId endcapDet1(HcalEndcap, -17, 72, 1);
0114 HcalDetId endcapDet2(HcalEndcap, 29, 35, 1);
0115 HcalDetId endcapDet3(HcalEndcap, 20, 51, 3);
0116 HcalDetId endcapDet4(HcalEndcap, -24, 63, 5);
0117 HcalDetId forwardDet1(HcalForward, 30, 71, 1);
0118 HcalDetId forwardDet3(HcalForward, -40, 71, 1);
0119
0120 if (topology.valid(barrelDet1))
0121 testClosestCell(barrelDet1, g);
0122 if (topology.valid(barrelDet2))
0123 testClosestCell(barrelDet2, g);
0124 if (topology.valid(barrelDet3))
0125 testClosestCell(barrelDet3, g);
0126 if (topology.valid(barrelDet4))
0127 testClosestCell(barrelDet4, g);
0128 if (topology.valid(endcapDet1))
0129 testClosestCell(endcapDet1, g);
0130 if (topology.valid(endcapDet2))
0131 testClosestCell(endcapDet2, g);
0132 if (topology.valid(endcapDet3))
0133 testClosestCell(endcapDet3, g);
0134 if (topology.valid(endcapDet4))
0135 testClosestCell(endcapDet4, g);
0136 if (topology.valid(forwardDet1))
0137 testClosestCell(forwardDet1, g);
0138 if (topology.valid(forwardDet3))
0139 testClosestCell(forwardDet3, g);
0140
0141 const std::vector<DetId>& idsb = g->getValidDetIds(DetId::Hcal, HcalBarrel);
0142 for (auto id : idsb) {
0143 testClosestCell(HcalDetId(id), g);
0144 }
0145
0146 const std::vector<DetId>& idse = g->getValidDetIds(DetId::Hcal, HcalEndcap);
0147 for (auto id : idse) {
0148 testClosestCell(HcalDetId(id), g);
0149 }
0150 }
0151
0152 void HcalGeometryTester::testClosestCell(const HcalDetId& detId, CaloSubdetectorGeometry* geom) {
0153 auto cell = geom->getGeometry(detId);
0154 HcalDetId closest = geom->getClosestCell(cell->getPosition());
0155 edm::LogVerbatim("HCalGeom") << "i/p " << detId << " position " << cell->getPosition() << " closest " << closest;
0156
0157 if (closest != detId) {
0158 edm::LogVerbatim("HCalGeom") << "HcalGeometryTester::Mismatch. Original HCAL cell is " << detId
0159 << " while nearest is " << closest << " ***ERROR***";
0160 }
0161 }
0162
0163 void HcalGeometryTester::testTriggerGeometry(const HcalTopology& topology) {
0164 HcalTrigTowerGeometry trigTowers(&topology);
0165 edm::LogVerbatim("HCalGeom") << "HCAL trigger tower eta bounds:";
0166 for (int ieta = 1; ieta <= 32; ++ieta) {
0167 double eta1, eta2;
0168 trigTowers.towerEtaBounds(ieta, 0, eta1, eta2);
0169 edm::LogVerbatim("HCalGeom") << "[" << ieta << "] " << eta1 << " " << eta2;
0170 }
0171
0172
0173 HcalDetId barrelDet(HcalBarrel, 1, 1, 1);
0174 HcalDetId endcapDet(HcalEndcap, 29, 1, 1);
0175 HcalDetId forwardDet1(HcalForward, 29, 71, 1);
0176 HcalDetId forwardDet2(HcalForward, 29, 71, 2);
0177 HcalDetId forwardDet3(HcalForward, 40, 71, 1);
0178
0179 using TowerDets = std::vector<HcalTrigTowerDetId>;
0180 if (topology.valid(barrelDet)) {
0181 TowerDets barrelTowers = trigTowers.towerIds(barrelDet);
0182 edm::LogVerbatim("HCalGeom") << "Trigger Tower Size: Barrel " << barrelTowers.size();
0183 assert(barrelTowers.size() == 1);
0184 edm::LogVerbatim("HCalGeom") << "Tower[0] " << barrelTowers[0];
0185 }
0186 if (topology.valid(endcapDet)) {
0187 TowerDets endcapTowers = trigTowers.towerIds(endcapDet);
0188 edm::LogVerbatim("HCalGeom") << "Trigger Tower Size: Endcap " << endcapTowers.size();
0189 assert(!endcapTowers.empty());
0190 for (unsigned int k = 0; k < endcapTowers.size(); ++k)
0191 edm::LogVerbatim("HCalGeom") << "Tower[" << k << "] " << endcapTowers[k];
0192 }
0193 if (topology.valid(forwardDet1)) {
0194 TowerDets forwardTowers1 = trigTowers.towerIds(forwardDet1);
0195 edm::LogVerbatim("HCalGeom") << "Trigger Tower Size: Forward1 " << forwardTowers1.size();
0196 assert(!forwardTowers1.empty());
0197 for (unsigned int k = 0; k < forwardTowers1.size(); ++k)
0198 edm::LogVerbatim("HCalGeom") << "Tower[" << k << "] " << forwardTowers1[k];
0199 }
0200 if (topology.valid(forwardDet2)) {
0201 TowerDets forwardTowers2 = trigTowers.towerIds(forwardDet2);
0202 edm::LogVerbatim("HCalGeom") << "Trigger Tower Size: Forward2 " << forwardTowers2.size();
0203 assert(!forwardTowers2.empty());
0204 for (unsigned int k = 0; k < forwardTowers2.size(); ++k)
0205 edm::LogVerbatim("HCalGeom") << "Tower[" << k << "] " << forwardTowers2[k];
0206 }
0207 if (topology.valid(forwardDet3)) {
0208 TowerDets forwardTowers3 = trigTowers.towerIds(forwardDet3);
0209 edm::LogVerbatim("HCalGeom") << "Trigger Tower Size: Forward3 " << forwardTowers3.size();
0210 assert(!forwardTowers3.empty());
0211 for (unsigned int k = 0; k < forwardTowers3.size(); ++k)
0212 edm::LogVerbatim("HCalGeom") << "Tower[" << k << "] " << forwardTowers3[k];
0213 }
0214 }
0215
0216 void HcalGeometryTester::testFlexiValidDetIds(CaloSubdetectorGeometry* caloGeom,
0217 const HcalTopology& topology,
0218 DetId::Detector det,
0219 int subdet,
0220 const std::string& label,
0221 std::vector<int>& dins) {
0222 std::stringstream s;
0223 s << label << " : " << std::endl;
0224 const std::vector<DetId>& idshb = caloGeom->getValidDetIds(det, subdet);
0225
0226 int counter = 0;
0227 for (std::vector<DetId>::const_iterator i = idshb.begin(); i != idshb.end(); i++, ++counter) {
0228 HcalDetId hid = (*i);
0229 s << counter << ": din " << topology.detId2denseId(*i) << ":" << hid;
0230 dins.emplace_back(topology.detId2denseId(*i));
0231
0232 auto cell = caloGeom->getGeometry(*i);
0233 s << *cell << std::endl;
0234 }
0235
0236 std::sort(dins.begin(), dins.end());
0237 s << "=== Total " << counter << " cells in " << label << std::endl;
0238
0239 counter = 0;
0240 for (std::vector<int>::const_iterator i = dins.begin(); i != dins.end(); ++i, ++counter) {
0241 HcalDetId hid = (topology.denseId2detId(*i));
0242 HcalDetId ihid = (topology.denseId2detId(dins[counter]));
0243 s << counter << ": din " << (*i) << " :" << hid << " == " << ihid << std::endl;
0244 }
0245 edm::LogVerbatim("HCalGeom") << s.str();
0246 }
0247
0248 void HcalGeometryTester::testFlexiGeomHF(CaloSubdetectorGeometry* caloGeom) {
0249 std::stringstream s;
0250 s << "Test HF Geometry : " << std::endl;
0251 for (int ieta = 29; ieta <= 41; ++ieta) {
0252 HcalDetId cell3(HcalForward, ieta, 3, 1);
0253 auto cellGeometry3 = caloGeom->getGeometry(cell3);
0254 if (cellGeometry3) {
0255 s << "cell geometry iphi=3 -> ieta=" << ieta << " eta " << cellGeometry3->getPosition().eta() << "+-"
0256 << std::abs(cellGeometry3->getCorners()[0].eta() - cellGeometry3->getCorners()[2].eta()) / 2 << " phi "
0257 << cellGeometry3->getPosition().phi() / 3.1415 * 180 << "+-"
0258 << std::abs(cellGeometry3->getCorners()[0].phi() - cellGeometry3->getCorners()[2].phi()) / 3.1415 * 180 / 2.
0259 << std::endl;
0260 }
0261 }
0262 edm::LogVerbatim("HCalGeom") << s.str();
0263 }
0264
0265 DEFINE_FWK_MODULE(HcalGeometryTester);
