File indexing completed on 2024-04-06 12:31:03
0001 #include <string>
0002 #include <vector>
0003 #include <iostream>
0004 #include <algorithm>
0005
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/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/Utilities/interface/InputTag.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013
0014 #include "DetectorDescription/DDCMS/interface/DDFilteredView.h"
0015 #include "DetectorDescription/DDCMS/interface/DDCompactView.h"
0016 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0017
0018 class DD4hep_ListIds : public edm::one::EDAnalyzer<> {
0019 public:
0020 DD4hep_ListIds(const edm::ParameterSet &);
0021 ~DD4hep_ListIds() override;
0022
0023 private:
0024 void analyze(const edm::Event &, const edm::EventSetup &) override;
0025 void beginJob() override {}
0026 void endJob() override;
0027
0028
0029
0030
0031
0032 bool printMaterial_;
0033 std::vector<std::string> materials_;
0034 edm::ESGetToken<cms::DDCompactView, IdealGeometryRecord> cpvToken_;
0035 };
0036
0037 DD4hep_ListIds::DD4hep_ListIds(const edm::ParameterSet &pset)
0038 : printMaterial_(pset.getUntrackedParameter<bool>("printMaterial")),
0039 materials_(pset.getUntrackedParameter<std::vector<std::string> >("materials")),
0040 cpvToken_(esConsumes(edm::ESInputTag("", "CMS"))) {}
0041
0042 DD4hep_ListIds::~DD4hep_ListIds() {}
0043
0044 void DD4hep_ListIds::analyze(const edm::Event &evt, const edm::EventSetup &setup) {
0045 auto cpv = setup.getTransientHandle(cpvToken_);
0046
0047 std::string attribute = "TkDDDStructure";
0048 cms::DDFilter filter(attribute, "");
0049 cms::DDFilteredView fv(*cpv, filter);
0050 fv.firstChild();
0051 std::set<std::string_view> tkdss;
0052
0053 for (const auto &t : fv.specpars()) {
0054 tkdss.insert(t.second->strValue(attribute));
0055 }
0056
0057 for (const auto &i : tkdss) {
0058 edm::LogVerbatim("ListIds") << "\nFiltering " << i;
0059 cms::DDFilter filter1(attribute, {i.data(), i.size()});
0060 cms::DDFilteredView fv1(*cpv, filter1);
0061 fv1.firstChild();
0062
0063 std::vector<const cms::Node *> nodes;
0064 while (fv1.firstChild()) {
0065 nodes = fv1.geoHistory();
0066 if (std::find(materials_.begin(),
0067 materials_.end(),
0068 nodes[nodes.size() - 1]->GetVolume()->GetMaterial()->GetName()) != materials_.end()) {
0069 for (const auto &n : nodes) {
0070 edm::LogVerbatim("ListIds") << n->GetVolume()->GetName() << "[" << n->GetNumber() << "]/";
0071 }
0072 edm::LogVerbatim("ListIds") << "Material:|" << nodes[nodes.size() - 1]->GetVolume()->GetMaterial()->GetName()
0073 << "|\n";
0074 }
0075 }
0076 }
0077 }
0078
0079 void DD4hep_ListIds::endJob() {}
0080
0081
0082
0083 #include "FWCore/PluginManager/interface/ModuleDef.h"
0084 #include "FWCore/Framework/interface/MakerMacros.h"
0085 DEFINE_FWK_MODULE(DD4hep_ListIds);