Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-02-12 09:07:47

0001 #include "FWCore/Framework/interface/Frameworkfwd.h"
0002 #include "FWCore/Framework/interface/MakerMacros.h"
0003 #include "FWCore/Utilities/interface/InputTag.h"
0004 #include "DQMServices/Core/interface/DQMStore.h"
0005 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "DataFormats/NanoAOD/interface/FlatTable.h"
0008 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0009 
0010 #include <memory>
0011 
0012 #include <numeric>
0013 #include <regex>
0014 #include <sstream>
0015 
0016 namespace {
0017   std::string replaceStringsToColumGets(const std::string &expr, const nanoaod::FlatTable &table) {
0018     std::regex token("\\w+");
0019     std::sregex_iterator tbegin(expr.begin(), expr.end(), token), tend;
0020     if (tbegin == tend)
0021       return expr;
0022     std::stringstream out;
0023     std::sregex_iterator last;
0024     for (std::sregex_iterator i = tbegin; i != tend; last = i, ++i) {
0025       std::smatch match = *i;
0026       out << match.prefix().str();
0027       if (table.columnIndex(match.str()) != -1) {
0028         out << "getAnyValue(\"" << match.str() << "\")";
0029       } else {
0030         out << match.str();
0031       }
0032     }
0033     out << last->suffix().str();
0034     return out.str();
0035   };
0036 }  // namespace
0037 
0038 class NanoAODDQM : public DQMEDAnalyzer {
0039 public:
0040   typedef nanoaod::FlatTable FlatTable;
0041 
0042   NanoAODDQM(const edm::ParameterSet &);
0043   void analyze(const edm::Event &, const edm::EventSetup &) override;
0044 
0045 protected:
0046   //Book histograms
0047   void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override;
0048 
0049 private:
0050   class Plot {
0051   public:
0052     Plot(MonitorElement *me) : plot_(me) {}
0053     virtual ~Plot() {}
0054     virtual void fill(const FlatTable &table, const std::vector<bool> &rowsel) = 0;
0055     const std::string &name() const { return plot_->getName(); }
0056 
0057   protected:
0058     MonitorElement *plot_;
0059   };
0060   class Count1D : public Plot {
0061   public:
0062     Count1D(DQMStore::IBooker &booker, const edm::ParameterSet &cfg)
0063         : Plot(booker.book1D(cfg.getParameter<std::string>("name"),
0064                              cfg.getParameter<std::string>("title"),
0065                              cfg.getParameter<uint32_t>("nbins"),
0066                              cfg.getParameter<double>("min"),
0067                              cfg.getParameter<double>("max"))) {}
0068     ~Count1D() override {}
0069     void fill(const FlatTable &table, const std::vector<bool> &rowsel) override {
0070       plot_->Fill(std::accumulate(rowsel.begin(), rowsel.end(), 0u));
0071     }
0072   };
0073 
0074   class Plot1D : public Plot {
0075   public:
0076     Plot1D(DQMStore::IBooker &booker, const edm::ParameterSet &cfg)
0077         : Plot(booker.book1D(cfg.getParameter<std::string>("name"),
0078                              cfg.getParameter<std::string>("title"),
0079                              cfg.getParameter<uint32_t>("nbins"),
0080                              cfg.getParameter<double>("min"),
0081                              cfg.getParameter<double>("max"))),
0082           col_(cfg.getParameter<std::string>("column")) {}
0083     ~Plot1D() override {}
0084     void fill(const FlatTable &table, const std::vector<bool> &rowsel) override {
0085       int icol = table.columnIndex(col_);
0086       if (icol == -1)
0087         return;  // columns may be missing (e.g. mc-only)
0088       switch (table.columnType(icol)) {
0089         case FlatTable::ColumnType::Float:
0090           vfill<float>(table, icol, rowsel);
0091           break;
0092         case FlatTable::ColumnType::Int:
0093           vfill<int>(table, icol, rowsel);
0094           break;
0095         case FlatTable::ColumnType::Int8:
0096           vfill<int8_t>(table, icol, rowsel);
0097           break;
0098         case FlatTable::ColumnType::UInt8:
0099           vfill<uint8_t>(table, icol, rowsel);
0100           break;
0101         case FlatTable::ColumnType::Bool:
0102           vfill<bool>(table, icol, rowsel);
0103           break;
0104         default:
0105           throw cms::Exception("LogicError", "Unsupported type");
0106       }
0107     }
0108 
0109   protected:
0110     std::string col_;
0111     template <typename T>
0112     void vfill(const FlatTable &table, int icol, const std::vector<bool> &rowsel) {
0113       const auto &data = table.columnData<T>(icol);
0114       for (unsigned int i = 0, n = data.size(); i < n; ++i) {
0115         if (rowsel[i])
0116           plot_->Fill(data[i]);
0117       }
0118     }
0119   };
0120   class Profile1D : public Plot {
0121   public:
0122     Profile1D(DQMStore::IBooker &booker, const edm::ParameterSet &cfg)
0123         : Plot(booker.bookProfile(cfg.getParameter<std::string>("name"),
0124                                   cfg.getParameter<std::string>("title"),
0125                                   cfg.getParameter<uint32_t>("nbins"),
0126                                   cfg.getParameter<double>("min"),
0127                                   cfg.getParameter<double>("max"),
0128                                   0.,
0129                                   0.,
0130                                   "")),
0131           ycol_(cfg.getParameter<std::string>("ycolumn")),
0132           xcol_(cfg.getParameter<std::string>("xcolumn")) {}
0133     ~Profile1D() override {}
0134     void fill(const FlatTable &table, const std::vector<bool> &rowsel) override {
0135       int icolx = table.columnIndex(xcol_);
0136       int icoly = table.columnIndex(ycol_);
0137       if (icolx == -1)
0138         throw cms::Exception("LogicError", "Missing " + xcol_);
0139       if (icoly == -1)
0140         throw cms::Exception("LogicError", "Missing " + ycol_);
0141       for (unsigned int irow = 0, n = table.size(); irow < n; ++irow) {
0142         if (rowsel[irow])
0143           plot_->Fill(table.getAnyValue(irow, icolx), table.getAnyValue(irow, icoly));
0144       }
0145     }
0146 
0147   protected:
0148     std::string ycol_, xcol_;
0149   };
0150 
0151   static std::unique_ptr<Plot> makePlot(DQMStore::IBooker &booker, const edm::ParameterSet &cfg) {
0152     const std::string &kind = cfg.getParameter<std::string>("kind");
0153     if (kind == "none")
0154       return nullptr;
0155     if (kind == "count1d")
0156       return std::make_unique<Count1D>(booker, cfg);
0157     if (kind == "hist1d")
0158       return std::make_unique<Plot1D>(booker, cfg);
0159     if (kind == "prof1d")
0160       return std::make_unique<Profile1D>(booker, cfg);
0161     throw cms::Exception("Configuration", "Unsupported plot kind '" + kind + "'");
0162   }
0163 
0164   struct SelGroupConfig {
0165     typedef StringCutObjectSelector<FlatTable::RowView> Selector;
0166     std::string name;
0167     std::string cutstr;
0168     std::unique_ptr<StringCutObjectSelector<FlatTable::RowView>> cutptr;
0169     std::vector<std::unique_ptr<Plot>> plots;
0170     SelGroupConfig() : name(), cutstr(), cutptr(), plots() {}
0171     SelGroupConfig(const std::string &nam, const std::string &cut) : name(nam), cutstr(cut), cutptr(), plots() {}
0172     bool nullCut() const { return cutstr.empty(); }
0173     void fillSel(const FlatTable &table, std::vector<bool> &out) {
0174       out.resize(table.size());
0175       if (nullCut()) {
0176         std::fill(out.begin(), out.end(), true);
0177       } else {
0178         if (!cutptr) {
0179           cutptr = std::make_unique<Selector>(replaceStringsToColumGets(cutstr, table));
0180         }
0181         for (unsigned int i = 0, n = table.size(); i < n; ++i) {
0182           out[i] = (*cutptr)(table.row(i));
0183         }
0184       }
0185     }
0186   };
0187   struct GroupConfig {
0188     std::vector<edm::ParameterSet> plotPSets;
0189     std::vector<SelGroupConfig> selGroups;
0190   };
0191   std::map<std::string, GroupConfig> groups_;
0192 };
0193 
0194 NanoAODDQM::NanoAODDQM(const edm::ParameterSet &iConfig) {
0195   const edm::ParameterSet &vplots = iConfig.getParameter<edm::ParameterSet>("vplots");
0196   for (const std::string &name : vplots.getParameterNamesForType<edm::ParameterSet>()) {
0197     auto &group = groups_[name];
0198     const auto &pset = vplots.getParameter<edm::ParameterSet>(name);
0199     group.plotPSets = pset.getParameter<std::vector<edm::ParameterSet>>("plots");
0200     group.selGroups.emplace_back();  // no selection (all entries)
0201     const auto &cuts = pset.getParameter<edm::ParameterSet>("sels");
0202     for (const std::string &cname : cuts.getParameterNamesForType<std::string>()) {
0203       group.selGroups.emplace_back(cname, cuts.getParameter<std::string>(cname));
0204     }
0205   }
0206   consumesMany<FlatTable>();
0207 }
0208 
0209 void NanoAODDQM::bookHistograms(DQMStore::IBooker &booker, edm::Run const &, edm::EventSetup const &) {
0210   booker.setCurrentFolder("Physics/NanoAODDQM");
0211 
0212   for (auto &pair : groups_) {
0213     booker.setCurrentFolder("Physics/NanoAODDQM/" + pair.first);
0214     for (auto &sels : pair.second.selGroups) {
0215       std::string dir("Physics/NanoAODDQM/" + pair.first);
0216       if (!sels.nullCut())
0217         dir += "/" + sels.name;
0218       booker.setCurrentFolder(dir);
0219       auto &plots = sels.plots;
0220       plots.clear();
0221       plots.reserve(pair.second.plotPSets.size());
0222       for (const auto &cfg : pair.second.plotPSets) {
0223         auto plot = makePlot(booker, cfg);
0224         if (plot)
0225           plots.push_back(std::move(plot));
0226       }
0227     }
0228   }
0229 }
0230 
0231 void NanoAODDQM::analyze(const edm::Event &iEvent, const edm::EventSetup &) {
0232   std::vector<edm::Handle<FlatTable>> alltables;
0233   iEvent.getManyByType(alltables);
0234   std::map<std::string, std::pair<const FlatTable *, std::vector<const FlatTable *>>> maintables;
0235 
0236   for (const auto &htab : alltables) {
0237     if (htab->extension())
0238       continue;
0239     maintables[htab->name()] = std::make_pair(htab.product(), std::vector<const FlatTable *>());
0240   }
0241   for (const auto &htab : alltables) {
0242     if (htab->extension()) {
0243       if (maintables.find(htab->name()) == maintables.end())
0244         throw cms::Exception("LogicError", "Missing main table for " + htab->name());
0245       maintables[htab->name()].second.push_back(htab.product());
0246     }
0247   }
0248 
0249   FlatTable merged;
0250   for (auto &pair : groups_) {
0251     const std::string &name = pair.first;
0252     if (maintables.find(name) == maintables.end())
0253       continue;  // may happen for missing collections
0254     auto &tables = maintables[name];
0255     const FlatTable *table = tables.first;
0256     if (!tables.second.empty()) {
0257       merged = *tables.first;
0258       for (auto *other : tables.second) {
0259         merged.addExtension(*other);
0260       }
0261       table = &merged;
0262     }
0263     std::vector<bool> selbits;
0264     for (auto &sel : pair.second.selGroups) {
0265       sel.fillSel(*table, selbits);
0266 
0267       for (auto &plot : sel.plots) {
0268         plot->fill(*table, selbits);
0269       }
0270     }
0271   }
0272 }
0273 
0274 #include "FWCore/Framework/interface/MakerMacros.h"
0275 DEFINE_FWK_MODULE(NanoAODDQM);