Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-09-07 04:35:45

0001 #include <memory>
0002 #include "FWCore/Framework/interface/Frameworkfwd.h"
0003 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "CondFormats/JetMETObjects/interface/QGLikelihoodObject.h"
0010 #include "JetMETCorrections/Objects/interface/JetCorrectionsRecord.h"
0011 #include "CondFormats/DataRecord/interface/QGLikelihoodRcd.h"
0012 
0013 class QGLikelihoodDBReader : public edm::one::EDAnalyzer<> {
0014 public:
0015   explicit QGLikelihoodDBReader(const edm::ParameterSet&);
0016   ~QGLikelihoodDBReader() override {}
0017 
0018 private:
0019   void beginJob() override {}
0020   void analyze(const edm::Event&, const edm::EventSetup&) override;
0021   void endJob() override {}
0022 
0023   std::string mPayloadName;
0024   edm::ESGetToken<QGLikelihoodObject, QGLikelihoodRcd> mPayloadToken;
0025   bool mCreateTextFile, mPrintScreen;
0026 };
0027 
0028 QGLikelihoodDBReader::QGLikelihoodDBReader(const edm::ParameterSet& iConfig) {
0029   mPayloadName = iConfig.getUntrackedParameter<std::string>("payloadName");
0030   mPayloadToken = esConsumes(edm::ESInputTag("", mPayloadName));
0031   mPrintScreen = iConfig.getUntrackedParameter<bool>("printScreen");
0032   mCreateTextFile = iConfig.getUntrackedParameter<bool>("createTextFile");
0033 }
0034 
0035 void QGLikelihoodDBReader::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0036   edm::LogInfo("UserOutput") << "Inspecting QGLikelihood payload with label:" << mPayloadName << std::endl;
0037   edm::ESHandle<QGLikelihoodObject> QGLParamsColl = iSetup.getHandle(mPayloadToken);
0038 
0039   edm::LogInfo("UserOutput") << "Ranges in which the QGTagger could be applied:"
0040                              << "  pt: " << QGLParamsColl->qgValidRange.PtMin << " --> "
0041                              << QGLParamsColl->qgValidRange.PtMax << ", eta: " << QGLParamsColl->qgValidRange.EtaMin
0042                              << " --> " << QGLParamsColl->qgValidRange.EtaMax
0043                              << ", rho: " << QGLParamsColl->qgValidRange.RhoMin << " --> "
0044                              << QGLParamsColl->qgValidRange.RhoMax << std::endl;
0045 
0046   std::vector<QGLikelihoodObject::Entry> const& data = QGLParamsColl->data;
0047   edm::LogInfo("UserOutput") << "There are " << data.size()
0048                              << " entries (categories with associated PDF):" << std::endl;
0049   for (auto idata = data.begin(); idata != data.end(); ++idata) {
0050     int varIndex = idata->category.VarIndex;
0051     int qgBin = idata->category.QGIndex;
0052     double etaMin = idata->category.EtaMin;
0053     double etaMax = idata->category.EtaMax;
0054     double rhoMin = idata->category.RhoMin;
0055     double rhoMax = idata->category.RhoMax;
0056     double ptMin = idata->category.PtMin;
0057     double ptMax = idata->category.PtMax;
0058 
0059     char buff[1000];
0060     sprintf(buff,
0061             "var=%1d, qg=%1d, ptMin=%8.2f, ptMax=%8.2f, etaMin=%3.1f, etaMax=%3.1f, rhoMin=%6.2f, rhoMax=%6.2f",
0062             varIndex,
0063             qgBin,
0064             ptMin,
0065             ptMax,
0066             etaMin,
0067             etaMax,
0068             rhoMin,
0069             rhoMax);
0070     edm::LogVerbatim("UserOutput") << buff << std::endl;
0071   }
0072 }
0073 
0074 DEFINE_FWK_MODULE(QGLikelihoodDBReader);