File indexing completed on 2024-04-06 12:14:20
0001 #include <iostream>
0002 #include <string>
0003 #include <vector>
0004
0005 #include "FWCore/Framework/interface/Frameworkfwd.h"
0006 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0007 #include "FWCore/Framework/interface/Event.h"
0008 #include "FWCore/Framework/interface/EventSetup.h"
0009 #include "FWCore/Framework/interface/ESTransientHandle.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0014
0015 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0016 #include "Geometry/HGCalTBCommonData/interface/HGCalTBDDDConstants.h"
0017 #include "Geometry/CaloTopology/interface/HGCalTBTopology.h"
0018 #include "DataFormats/ForwardDetId/interface/ForwardSubdetector.h"
0019 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0020
0021 class HGCalTBTopologyTester : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
0022 public:
0023 explicit HGCalTBTopologyTester(const edm::ParameterSet&);
0024
0025 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0026
0027 private:
0028 void analyze(edm::Event const&, edm::EventSetup const&) override;
0029 void beginJob() override {}
0030 void beginRun(edm::Run const&, edm::EventSetup const&) override {}
0031 void endRun(edm::Run const&, edm::EventSetup const&) override {}
0032 void doTest(const HGCalTBTopology& topology);
0033
0034
0035 const std::string detectorName_;
0036 const std::vector<int> type_, layer_, sector_, cells_;
0037 const edm::ESGetToken<HGCalTBTopology, IdealGeometryRecord> tokTopo_;
0038 };
0039
0040 HGCalTBTopologyTester::HGCalTBTopologyTester(const edm::ParameterSet& iC)
0041 : detectorName_(iC.getParameter<std::string>("detectorName")),
0042 type_(iC.getParameter<std::vector<int> >("types")),
0043 layer_(iC.getParameter<std::vector<int> >("layers")),
0044 sector_(iC.getParameter<std::vector<int> >("sector")),
0045 cells_(iC.getParameter<std::vector<int> >("cells")),
0046 tokTopo_{esConsumes<HGCalTBTopology, IdealGeometryRecord>(edm::ESInputTag{"", detectorName_})} {}
0047
0048 void HGCalTBTopologyTester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0049 edm::ParameterSetDescription desc;
0050 std::vector<int> types = {0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2};
0051 std::vector<int> layer = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18};
0052 std::vector<int> sect = {1, 1, 2, 2, 3, 3, 5, 5, 6, 6, 7, 7, 8, 8, 9, 9, 10, 10};
0053 std::vector<int> cells = {0, 4, 12, 14, 18, 23, 1, 4, 7, 10, 13, 16, 0, 3, 6, 9, 12, 15};
0054 desc.add<std::string>("detectorName", "HGCalEESensitive");
0055 desc.add<std::vector<int> >("types", types);
0056 desc.add<std::vector<int> >("layers", layer);
0057 desc.add<std::vector<int> >("sector", sect);
0058 desc.add<std::vector<int> >("cells", cells);
0059 descriptions.add("hgcalTBTopologyTesterEE", desc);
0060 }
0061
0062 void HGCalTBTopologyTester::analyze(edm::Event const&, edm::EventSetup const& iSetup) {
0063 doTest(iSetup.getData(tokTopo_));
0064 }
0065
0066 void HGCalTBTopologyTester::doTest(const HGCalTBTopology& topology) {
0067 for (unsigned int i = 0; i < type_.size(); ++i) {
0068 DetId id;
0069 if (detectorName_ == "HGCalEESensitive") {
0070 id = HGCalDetId(ForwardSubdetector::HGCEE, 1, layer_[i], type_[i], sector_[i], cells_[i]);
0071 } else if (detectorName_ == "HGCalHESiliconSensitive") {
0072 id = HGCalDetId(ForwardSubdetector::HGCHEF, 1, layer_[i], type_[i], sector_[i], cells_[i]);
0073 } else {
0074 break;
0075 }
0076 if (topology.valid(id)) {
0077 std::vector<DetId> ids = topology.neighbors(id);
0078 unsigned int k(0);
0079 edm::LogVerbatim("HGCalGeom") << static_cast<HGCalDetId>(id) << " has " << ids.size() << " neighbours:";
0080 for (const auto& idn : ids) {
0081 edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << (HGCSiliconDetId)(idn);
0082 ++k;
0083 }
0084 }
0085 }
0086 }
0087
0088
0089 DEFINE_FWK_MODULE(HGCalTBTopologyTester);