File indexing completed on 2024-08-30 02:10:35
0001 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0002 #include "DQMServices/Core/interface/DQMStore.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/Frameworkfwd.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/ServiceRegistry/interface/Service.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010
0011
0012 #include <fmt/printf.h>
0013
0014 class SiPixelTrackComparisonHarvester : public DQMEDHarvester {
0015 public:
0016 explicit SiPixelTrackComparisonHarvester(const edm::ParameterSet&);
0017 ~SiPixelTrackComparisonHarvester() override = default;
0018 void dqmEndJob(DQMStore::IBooker& ibooker, DQMStore::IGetter& igetter) override;
0019 void project2DalongDiagonal(MonitorElement* input2D, DQMStore::IBooker& ibooker);
0020 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0021
0022 private:
0023
0024 const std::string topFolder_;
0025 };
0026
0027 SiPixelTrackComparisonHarvester::SiPixelTrackComparisonHarvester(const edm::ParameterSet& iConfig)
0028 : topFolder_(iConfig.getParameter<std::string>("topFolderName")) {}
0029
0030 void SiPixelTrackComparisonHarvester::dqmEndJob(DQMStore::IBooker& ibooker, DQMStore::IGetter& igetter) {
0031 MonitorElement* hpt_eta_tkAllReference = igetter.get(topFolder_ + "/ptetatrkAllReference");
0032 if (hpt_eta_tkAllReference == nullptr) {
0033 edm::LogError("SiPixelTrackComparisonHarvester")
0034 << "MonitorElement not found: " << topFolder_ << "/ptetatrkAllReference. Skipping.";
0035 return;
0036 }
0037
0038 MonitorElement* hpt_eta_tkAllReferencematched = igetter.get(topFolder_ + "/ptetatrkAllReferencematched");
0039 if (hpt_eta_tkAllReferencematched == nullptr) {
0040 edm::LogError("SiPixelTrackComparisonHarvester")
0041 << "MonitorElement not found: " << topFolder_ << "/ptetatrkAllReferencematched. Skipping.";
0042 return;
0043 }
0044
0045 MonitorElement* hphi_z_tkAllReference = igetter.get(topFolder_ + "/phiztrkAllReference");
0046 if (hphi_z_tkAllReference == nullptr) {
0047 edm::LogError("SiPixelTrackComparisonHarvester")
0048 << "MonitorElement not found: " << topFolder_ << "/phiztrkAllReference. Skipping.";
0049 return;
0050 }
0051
0052 MonitorElement* hphi_z_tkAllReferencematched = igetter.get(topFolder_ + "/phiztrkAllReferencematched");
0053 if (hphi_z_tkAllReferencematched == nullptr) {
0054 edm::LogError("SiPixelTrackComparisonHarvester")
0055 << "MonitorElement not found: " << topFolder_ << "/phiztrkAllReferencematched. Skipping.";
0056 return;
0057 }
0058
0059 ibooker.cd();
0060 ibooker.setCurrentFolder(topFolder_);
0061 MonitorElement* hpt_eta_matchRatio = ibooker.book2D(
0062 "matchingeff_pt_eta", "Efficiency of track matching; #eta; p_{T} [GeV];", 30, -M_PI, M_PI, 200, 0., 200.);
0063 MonitorElement* hphi_z_matchRatio = ibooker.book2D(
0064 "matchingeff_phi_z", "Efficiency of track matching; #phi; z [cm];", 30, -M_PI, M_PI, 30, -30., 30.);
0065
0066 hpt_eta_matchRatio->divide(hpt_eta_tkAllReferencematched, hpt_eta_tkAllReference, 1., 1., "B");
0067 hphi_z_matchRatio->divide(hphi_z_tkAllReferencematched, hphi_z_tkAllReference, 1., 1., "B");
0068
0069
0070 std::vector<std::string> listOfMEsToProject = {"nTracks",
0071 "nLooseAndAboveTracks",
0072 "nLooseAndAboveTracks_matched",
0073 "nRecHits",
0074 "nLayers",
0075 "nChi2ndof",
0076 "charge",
0077 "pt",
0078 "curvature",
0079 "eta",
0080 "phi",
0081 "z",
0082 "tip"};
0083 for (const auto& me : listOfMEsToProject) {
0084 MonitorElement* input2D = igetter.get(topFolder_ + "/" + me);
0085 edm::LogPrint("SiPixelTrackComparisonHarvester") << "processing " << topFolder_ + "/" + me;
0086 this->project2DalongDiagonal(input2D, ibooker);
0087 }
0088 }
0089
0090 void SiPixelTrackComparisonHarvester::project2DalongDiagonal(MonitorElement* input2D, DQMStore::IBooker& ibooker) {
0091 if (input2D == nullptr) {
0092 edm::LogError("SiPixelTrackComparisonHarvester")
0093 << "ME needed for diagonal projection is not found in the input file at" << topFolder_ << ". Skipping.";
0094 return;
0095 }
0096
0097 ibooker.cd();
0098 ibooker.setCurrentFolder(topFolder_ + "/projectedDifferences");
0099 const auto& h_name = fmt::sprintf("%s_proj", input2D->getName());
0100 const auto& h_title = fmt::sprintf(";%s CPU -GPU difference", input2D->getTitle());
0101 const auto& span = (input2D->getAxisMax() - input2D->getAxisMin());
0102 const auto& b_w = span / input2D->getNbinsX();
0103 const auto& nbins = ((input2D->getNbinsX() % 2) == 0) ? input2D->getNbinsX() + 1 : input2D->getNbinsX();
0104
0105 MonitorElement* diagonalized = ibooker.book1D(h_name, h_title, nbins, -span / 2., span / 2.);
0106
0107
0108 for (int i = 1; i <= input2D->getNbinsX(); i++) {
0109 for (int j = 1; j <= input2D->getNbinsY(); j++) {
0110 diagonalized->Fill((i - j) * b_w, input2D->getBinContent(i, j));
0111 }
0112 }
0113
0114
0115
0116 for (int bin = 1; bin <= diagonalized->getNbinsX(); bin++) {
0117 diagonalized->setBinError(bin, 0.f);
0118 }
0119 }
0120
0121 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0122 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0123 void SiPixelTrackComparisonHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0124 edm::ParameterSetDescription desc;
0125 desc.add<std::string>("topFolderName", "SiPixelHeterogeneous/PixelTrackCompareGPUvsCPU/");
0126 descriptions.addWithDefaultLabel(desc);
0127 }
0128
0129
0130 DEFINE_FWK_MODULE(SiPixelTrackComparisonHarvester);