File indexing completed on 2023-03-17 11:25:47
0001 #include <string>
0002 #include <vector>
0003 #include <map>
0004 #include <iostream>
0005 #include <iomanip>
0006
0007 #include <fmt/printf.h>
0008
0009 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/Framework/interface/ESTransientHandle.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/Utilities/interface/InputTag.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016
0017 #include "DetectorDescription/DDCMS/interface/DDFilteredView.h"
0018 #include "DetectorDescription/DDCMS/interface/DDCompactView.h"
0019 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0020
0021 #include <TROOT.h>
0022 #include <TProfile2D.h>
0023 #include <TH2F.h>
0024 #include <TLine.h>
0025 #include <TText.h>
0026 #include <TColor.h>
0027 #include <TCanvas.h>
0028 #include <TStyle.h>
0029 #include <TFrame.h>
0030 #include <TLegend.h>
0031 #include <TLegendEntry.h>
0032
0033 #include "DD4hep_MaterialAccountingGroup.h"
0034 #include "DD4hep_TrackingMaterialPlotter.h"
0035
0036 class DD4hep_ListGroups : public edm::one::EDAnalyzer<> {
0037 public:
0038 DD4hep_ListGroups(const edm::ParameterSet &iConfig);
0039 ~DD4hep_ListGroups() override;
0040
0041 private:
0042 void analyze(const edm::Event &, const edm::EventSetup &) override;
0043 void beginJob() override {}
0044 void endJob() override;
0045
0046 const edm::ESInputTag m_tag;
0047 const bool m_saveSummaryPlot;
0048 const edm::ESGetToken<cms::DDCompactView, IdealGeometryRecord> ddToken_;
0049
0050 std::vector<TH2F *> m_plots;
0051 std::set<std::string_view> m_group_names;
0052 std::vector<unsigned int> m_color;
0053 std::vector<int> m_gradient;
0054 std::vector<DD4hep_MaterialAccountingGroup *> m_groups;
0055 void fillColor();
0056 void fillGradient();
0057 void fillMaterialDifferences();
0058 void produceAndSaveSummaryPlot(const cms::DDCompactView &cpv);
0059 std::vector<std::pair<std::shared_ptr<TLine>, std::shared_ptr<TText>>> overlayEtaReferences();
0060 std::map<std::string, std::pair<float, float>> m_diff;
0061 std::map<std::string, std::pair<float, float>> m_values;
0062 };
0063
0064 #include "DD4hep_ListGroupsMaterialDifference.h"
0065
0066 DD4hep_ListGroups::DD4hep_ListGroups(const edm::ParameterSet &iConfig)
0067 : m_tag(iConfig.getParameter<edm::ESInputTag>("DDDetector")),
0068 m_saveSummaryPlot(iConfig.getUntrackedParameter<bool>("SaveSummaryPlot")),
0069 ddToken_(esConsumes<cms::DDCompactView, IdealGeometryRecord>(m_tag)) {
0070 m_plots.clear();
0071 m_groups.clear();
0072 TColor::InitializeColors();
0073 fillColor();
0074 fillMaterialDifferences();
0075 fillGradient();
0076 }
0077
0078 DD4hep_ListGroups::~DD4hep_ListGroups() {}
0079
0080 void DD4hep_ListGroups::produceAndSaveSummaryPlot(const cms::DDCompactView &cpv) {
0081 const double scale = 10.;
0082
0083 static int markerStyles[10] = {kFullCircle,
0084 kFullSquare,
0085 kFullTriangleUp,
0086 kFullTriangleDown,
0087 kOpenCircle,
0088 kOpenSquare,
0089 kOpenTriangleUp,
0090 kOpenDiamond,
0091 kOpenCross,
0092 kFullStar};
0093
0094 for (auto n : m_group_names) {
0095 m_groups.emplace_back(new DD4hep_MaterialAccountingGroup(n.data(), cpv));
0096 }
0097
0098 auto canvas = std::make_unique<TCanvas>(
0099 "Grouping_rz", "Grouping - RZ view", (int)(600 * scale * 1.25), (int)(120 * scale * 1.50));
0100 canvas->GetFrame()->SetFillColor(kWhite);
0101 gStyle->SetOptStat(0);
0102
0103 unsigned int color_index = 1;
0104
0105 auto leg = std::make_unique<TLegend>(0.1, 0.1, 0.23, 0.34);
0106 leg->SetHeader("Tracker Material Grouping");
0107 leg->SetTextFont(42);
0108 leg->SetTextSize(0.008);
0109 leg->SetNColumns(3);
0110 auto radlen = std::make_unique<TProfile2D>("OverallRadLen", "OverallRadLen", 600., -300., 300, 120., 0., 120.);
0111 auto eneloss =
0112 std::make_unique<TProfile2D>("OverallEnergyLoss", "OverallEnergyLoss", 600., -300., 300, 120., 0., 120.);
0113 auto radlen_diff = std::make_unique<TProfile2D>(
0114 "OverallDifferencesRadLen", "OverallDifferencesRadLen", 600., -300., 300, 120., 0., 120.);
0115 auto eneloss_diff = std::make_unique<TProfile2D>(
0116 "OverallDifferencesEnergyLoss", "OverallDifferencesEnergyLoss", 600., -300., 300, 120., 0., 120.);
0117
0118 for (auto g : m_groups) {
0119 m_plots.push_back(
0120 new TH2F(g->name().c_str(), g->name().c_str(), 6000., -300., 300, 1200., 0., 120.));
0121 TH2F ¤t = *m_plots.back();
0122 current.SetMarkerColor(m_color[color_index]);
0123 current.SetMarkerStyle(markerStyles[color_index % 10]);
0124 current.SetMarkerSize(0.8);
0125 current.SetLineWidth(1);
0126 for (const auto &element : g->elements()) {
0127 current.Fill(element.z(), element.perp());
0128 radlen->Fill(element.z(), element.perp(), m_values[g->name()].first);
0129 eneloss->Fill(element.z(), element.perp(), m_values[g->name()].second);
0130 radlen_diff->Fill(element.z(), element.perp(), m_diff[g->name()].first);
0131 eneloss_diff->Fill(element.z(), element.perp(), m_diff[g->name()].second);
0132 }
0133
0134 if (color_index == 1)
0135 current.Draw();
0136 else
0137 current.Draw("SAME");
0138
0139 leg->AddEntry(¤t, g->name().c_str(), "lp")->SetTextColor(m_color[color_index]);
0140 color_index++;
0141
0142 color_index = color_index % m_color.size();
0143 }
0144 leg->Draw();
0145 canvas->SaveAs("Grouping.png");
0146
0147 std::vector<std::pair<std::shared_ptr<TLine>, std::shared_ptr<TText>>> lines = overlayEtaReferences();
0148
0149 canvas->Clear();
0150 radlen->SetMinimum(0);
0151 radlen->SetMaximum(0.25);
0152 radlen->Draw("COLZ");
0153 for (const auto &line : lines) {
0154 line.first->SetLineWidth(5);
0155 line.first->Draw();
0156 line.second->Draw();
0157 }
0158 canvas->SaveAs("RadLenValues.png");
0159
0160 canvas->Clear();
0161 eneloss->SetMinimum(0.00001);
0162 eneloss->SetMaximum(0.0005);
0163 eneloss->Draw("COLZ");
0164 for (const auto &line : lines) {
0165 line.first->SetLineWidth(5);
0166 line.first->Draw();
0167 line.second->Draw();
0168 }
0169 canvas->SaveAs("EnergyLossValues.png");
0170
0171 canvas->Clear();
0172 gStyle->SetPalette(m_gradient.size(), &m_gradient.front());
0173 gStyle->SetNumberContours(m_gradient.size());
0174 radlen_diff->SetMinimum(-100);
0175 radlen_diff->SetMaximum(100);
0176 radlen_diff->Draw("COLZ");
0177 for (const auto &line : lines) {
0178 line.first->SetLineWidth(5);
0179 line.first->Draw();
0180 line.second->Draw();
0181 }
0182 canvas->SaveAs("RadLenChanges.png");
0183
0184 canvas->Clear();
0185 eneloss_diff->SetMinimum(-100);
0186 eneloss_diff->SetMaximum(100);
0187 eneloss_diff->Draw("COLZ");
0188 for (const auto &line : lines) {
0189 line.first->SetLineWidth(5);
0190 line.first->Draw();
0191 line.second->Draw();
0192 }
0193 canvas->SaveAs("EnergyLossChanges.png");
0194 }
0195
0196 void DD4hep_ListGroups::fillColor(void) {
0197
0198
0199
0200
0201
0202
0203 m_color.emplace_back(kBlack);
0204
0205 m_color.emplace_back(kAzure);
0206 m_color.emplace_back(kAzure - 1);
0207 m_color.emplace_back(kAzure + 1);
0208 m_color.emplace_back(kAzure + 2);
0209
0210 m_color.emplace_back(kGreen);
0211 m_color.emplace_back(kGreen + 2);
0212 m_color.emplace_back(kGreen + 4);
0213 m_color.emplace_back(kSpring + 9);
0214 m_color.emplace_back(kSpring + 4);
0215 m_color.emplace_back(kSpring);
0216
0217 m_color.emplace_back(kRed);
0218 m_color.emplace_back(kRed + 2);
0219 m_color.emplace_back(kRed - 7);
0220 m_color.emplace_back(kRed - 5);
0221 m_color.emplace_back(kRed - 10);
0222 m_color.emplace_back(kRed - 1);
0223 m_color.emplace_back(kRed - 2);
0224 m_color.emplace_back(kRed - 3);
0225 m_color.emplace_back(kPink - 2);
0226 m_color.emplace_back(kPink - 3);
0227 m_color.emplace_back(kPink - 4);
0228 m_color.emplace_back(kPink + 9);
0229 m_color.emplace_back(kPink + 8);
0230 m_color.emplace_back(kPink + 7);
0231 m_color.emplace_back(kMagenta - 2);
0232 m_color.emplace_back(kMagenta - 3);
0233 m_color.emplace_back(kMagenta - 4);
0234 m_color.emplace_back(kMagenta - 5);
0235 m_color.emplace_back(kMagenta - 6);
0236 m_color.emplace_back(kMagenta - 7);
0237 m_color.emplace_back(kRed);
0238 m_color.emplace_back(kMagenta - 9);
0239 m_color.emplace_back(kViolet);
0240
0241 m_color.emplace_back(kOrange + 9);
0242 m_color.emplace_back(kOrange + 7);
0243 m_color.emplace_back(kOrange + 5);
0244 m_color.emplace_back(kOrange - 2);
0245 m_color.emplace_back(kOrange - 3);
0246 m_color.emplace_back(kOrange - 6);
0247 m_color.emplace_back(kOrange + 4);
0248 m_color.emplace_back(kOrange - 7);
0249 m_color.emplace_back(kOrange);
0250 m_color.emplace_back(kOrange + 10);
0251
0252 m_color.emplace_back(kViolet + 10);
0253 m_color.emplace_back(kViolet + 6);
0254 m_color.emplace_back(kViolet + 3);
0255 m_color.emplace_back(kViolet - 7);
0256 m_color.emplace_back(kViolet - 1);
0257 m_color.emplace_back(kViolet + 9);
0258 m_color.emplace_back(kViolet - 5);
0259 m_color.emplace_back(kViolet - 3);
0260 m_color.emplace_back(kViolet);
0261
0262 m_color.emplace_back(kAzure);
0263 m_color.emplace_back(kAzure + 8);
0264 m_color.emplace_back(kAzure + 2);
0265 m_color.emplace_back(kAzure + 4);
0266 m_color.emplace_back(kCyan + 1);
0267 m_color.emplace_back(kCyan - 9);
0268 m_color.emplace_back(kCyan + 3);
0269 m_color.emplace_back(kCyan + 4);
0270 m_color.emplace_back(kAzure);
0271 m_color.emplace_back(kAzure + 8);
0272 m_color.emplace_back(kAzure + 2);
0273 m_color.emplace_back(kAzure + 5);
0274 m_color.emplace_back(kCyan + 1);
0275 m_color.emplace_back(kCyan - 9);
0276 m_color.emplace_back(kCyan + 3);
0277 m_color.emplace_back(kCyan + 4);
0278 m_color.emplace_back(kAzure);
0279 m_color.emplace_back(kAzure + 8);
0280 m_color.emplace_back(kAzure + 2);
0281 m_color.emplace_back(kAzure + 5);
0282 m_color.emplace_back(kCyan + 1);
0283 m_color.emplace_back(kCyan - 9);
0284 m_color.emplace_back(kCyan + 3);
0285 m_color.emplace_back(kCyan + 4);
0286 }
0287
0288 void DD4hep_ListGroups::fillGradient() {
0289 m_gradient.reserve(200);
0290 unsigned int steps = 100;
0291
0292 unsigned int index = ((TObjArray *)gROOT->GetListOfColors())->GetLast() + 1;
0293
0294 float r1, g1, b1, r2, g2, b2;
0295 static_cast<TColor *>(gROOT->GetListOfColors()->At(kBlue + 1))->GetRGB(r1, g1, b1);
0296 static_cast<TColor *>(gROOT->GetListOfColors()->At(kAzure + 10))->GetRGB(r2, g2, b2);
0297 float delta_r = (r2 - r1) / (steps - 1);
0298 float delta_g = (g2 - g1) / (steps - 1);
0299 float delta_b = (b2 - b1) / (steps - 1);
0300
0301 m_gradient.emplace_back(kBlue + 4);
0302 unsigned int ii = 0;
0303 for (unsigned int i = 0; i < steps; ++i, ++ii) {
0304 new TColor(static_cast<Int_t>(index + ii), r1 + delta_r * i, g1 + delta_g * i, b1 + delta_b * i);
0305 m_gradient.emplace_back(index + ii);
0306 }
0307
0308 m_gradient.emplace_back(kWhite);
0309
0310 static_cast<TColor *>(gROOT->GetListOfColors()->At(kOrange))->GetRGB(r1, g1, b1);
0311 static_cast<TColor *>(gROOT->GetListOfColors()->At(kOrange + 7))->GetRGB(r2, g2, b2);
0312 delta_r = (r2 - r1) / (steps - 1);
0313 delta_g = (g2 - g1) / (steps - 1);
0314 delta_b = (b2 - b1) / (steps - 1);
0315 for (unsigned int i = 0; i < steps; ++i, ++ii) {
0316 new TColor(static_cast<Int_t>(index + ii), r1 + delta_r * i, g1 + delta_g * i, b1 + delta_b * i);
0317 m_gradient.emplace_back(index + ii);
0318 }
0319 m_gradient.emplace_back(kRed);
0320 }
0321
0322 std::vector<std::pair<std::shared_ptr<TLine>, std::shared_ptr<TText>>> DD4hep_ListGroups::overlayEtaReferences() {
0323 std::vector<std::pair<std::shared_ptr<TLine>, std::shared_ptr<TText>>> lines;
0324
0325 lines.reserve(40);
0326 std::pair<float, float> deltaZ(293, 298);
0327 std::pair<float, float> deltaR(115, 118);
0328 float text_size = 0.033;
0329
0330 for (float eta = 0.; eta <= 3.8; eta += 0.2) {
0331 float theta = 2. * atan(exp(-eta));
0332 if (eta >= 1.8) {
0333 lines.emplace_back(
0334 std::make_shared<TLine>(deltaZ.first, deltaZ.first * tan(theta), deltaZ.second, deltaZ.second * tan(theta)),
0335 std::make_shared<TText>(deltaZ.first, deltaZ.first * tan(theta), fmt::sprintf("%2.1f", eta).c_str()));
0336 lines.back().second->SetTextFont(42);
0337 lines.back().second->SetTextSize(text_size);
0338 lines.back().second->SetTextAlign(33);
0339 lines.emplace_back(
0340 std::make_shared<TLine>(-deltaZ.first, deltaZ.first * tan(theta), -deltaZ.second, deltaZ.second * tan(theta)),
0341 std::make_shared<TText>(-deltaZ.first, deltaZ.first * tan(theta), fmt::sprintf("-%2.1f", eta).c_str()));
0342 lines.back().second->SetTextFont(42);
0343 lines.back().second->SetTextSize(text_size);
0344 lines.back().second->SetTextAlign(13);
0345 } else {
0346 lines.emplace_back(
0347 std::make_shared<TLine>(deltaR.first / tan(theta), deltaR.first, deltaR.second / tan(theta), deltaR.second),
0348 std::make_shared<TText>(deltaR.first / tan(theta), deltaR.first, fmt::sprintf("%2.1f", eta).c_str()));
0349 lines.back().second->SetTextFont(42);
0350 lines.back().second->SetTextSize(text_size);
0351 lines.back().second->SetTextAlign(23);
0352 if (eta != 0) {
0353 lines.emplace_back(
0354 std::make_shared<TLine>(
0355 -deltaR.first / tan(theta), deltaR.first, -deltaR.second / tan(theta), deltaR.second),
0356 std::make_shared<TText>(-deltaR.first / tan(theta), deltaR.first, fmt::sprintf("-%2.1f", eta).c_str()));
0357 lines.back().second->SetTextFont(42);
0358 lines.back().second->SetTextSize(text_size);
0359 lines.back().second->SetTextAlign(23);
0360 }
0361 }
0362 }
0363 return lines;
0364 }
0365
0366 void DD4hep_ListGroups::analyze(const edm::Event &evt, const edm::EventSetup &setup) {
0367 edm::ESTransientHandle<cms::DDCompactView> cpv = setup.getTransientHandle(ddToken_);
0368 cms::DDFilter filter("TrackingMaterialGroup", "");
0369 cms::DDFilteredView fv(*cpv, filter);
0370
0371 for (const auto &t : fv.specpars()) {
0372 m_group_names.insert(t.second->strValue("TrackingMaterialGroup"));
0373 }
0374
0375 for (const auto &i : m_group_names) {
0376 cms::DDFilter filter1("TrackingMaterialGroup", {i.data(), i.size()});
0377 cms::DDFilteredView fv1(*cpv, filter1);
0378 bool firstChild = fv1.firstChild();
0379
0380 for (const auto &j : fv1.specpars()) {
0381 for (const auto &k : j.second->paths) {
0382 if (firstChild) {
0383 std::vector<std::vector<cms::Node *>> children = fv1.children(k);
0384 for (auto const &path : children) {
0385 for (auto const &node : path) {
0386 edm::LogVerbatim("TrackingMaterialGroup") << node->GetName() << "/";
0387 }
0388 cms::Translation trans = fv1.translation(path);
0389 edm::LogVerbatim("TrackingMaterialGroup")
0390 << "(" << trans.x() << ", " << trans.y() << ", " << trans.z() << ")\n";
0391 }
0392 }
0393 }
0394 }
0395 }
0396
0397 if (m_saveSummaryPlot)
0398 produceAndSaveSummaryPlot(*cpv);
0399 }
0400
0401 void DD4hep_ListGroups::endJob() {}
0402
0403
0404
0405 #include "FWCore/PluginManager/interface/ModuleDef.h"
0406 #include "FWCore/Framework/interface/MakerMacros.h"
0407 DEFINE_FWK_MODULE(DD4hep_ListGroups);