UETableProducer.cc

Go to the documentation of this file.

// system include files
#include <iostream>
#include <fstream>
#include <vector>
#include <string>
#include <algorithm>
 
// user include files
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/one/EDAnalyzer.h"
 
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/Framework/interface/ESHandle.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
 
#include "FWCore/ServiceRegistry/interface/Service.h"
 
#include "CondFormats/DataRecord/interface/HeavyIonUERcd.h"
#include "CondFormats/HIObjects/interface/UETable.h"
#include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
 
// For DB entry using JetCorrector to store the vector of float
#include "CondFormats/JetMETObjects/interface/JetCorrectorParameters.h"
 
#include "FWCore/ParameterSet/interface/FileInPath.h"
 
using namespace std;
//
// class decleration
//
 
namespace {
  class UETableProducer : public edm::one::EDAnalyzer<> {
  public:
    explicit UETableProducer(const edm::ParameterSet&);
    ~UETableProducer() override;
 
    static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
 
  private:
    void beginJob() override {}
    void analyze(const edm::Event&, const edm::EventSetup&) override;
    void endJob() override;
 
    // ----------member data ---------------------------
 
    bool debug_;
    bool jetCorrectorFormat_;
 
    string calibrationFile_;
    unsigned int runnum_;
 
    unsigned int index = 0, np[5], ni0[2], ni1[2], ni2[2];
 
    //    ue_interpolation_pf0[15][344],
    //    ue_interpolation_pf1[15][344],
    //    ue_interpolation_pf2[15][82];
  };
}  // namespace
//
// constants, enums and typedefs
//
 
//
// static data member definitions
//
 
//
// constructors and destructor
//
UETableProducer::UETableProducer(const edm::ParameterSet& iConfig) : runnum_(0) {
  //now do what ever initialization is needed
  calibrationFile_ = iConfig.getParameter<std::string>("txtFile");
  //calibrationFile_ = "RecoHI/HiJetAlgos/data/ue_calibrations_pf_data.txt";
 
  debug_ = iConfig.getUntrackedParameter<bool>("debug", false);
  jetCorrectorFormat_ = iConfig.getUntrackedParameter<bool>("jetCorrectorFormat", false);
 
  np[0] = 3;   // Number of reduced PF ID (track, ECAL, HCAL)
  np[1] = 15;  // Number of pseudorapidity block
  np[2] = 5;   // Fourier series order
  np[3] = 2;   // Re or Im
  np[4] = 82;  // Number of feature parameter
 
  ni0[0] = np[1];
  ni0[1] = 344;  // Number of track binning (kept = ECAL)
 
  ni1[0] = np[1];
  ni1[1] = 344;  // Number of ECAL binning
 
  ni2[0] = np[1];
  ni2[1] = 82;  // Number of HCAL binning
}
 
UETableProducer::~UETableProducer() {
  // do anything here that needs to be done at desctruction time
  // (e.g. close files, deallocate resources etc.)
}
 
//
// member functions
//
 
// ------------ method called to for each event  ------------
void UETableProducer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
  // nothing
}
 
// ------------ method called once each job just after ending the event loop  ------------
void UETableProducer::endJob() {
  std::string qpDataName = calibrationFile_;
  std::ifstream textTable_(qpDataName.c_str());
 
  std::vector<float> ue_vec;
  std::unique_ptr<UETable> ue_predictor_pf;
 
  if (!jetCorrectorFormat_) {
    ue_predictor_pf = std::make_unique<UETable>();
  }
  //  unsigned int Nnp_full = np[0] * np[1] * np[2] * np[3] * np[4];
  unsigned int Nnp = np[0] * np[1] * (1 + (np[2] - 1) * np[3]) * np[4];
  unsigned int Nni0 = ni0[0] * ni0[1];
  unsigned int Nni1 = ni1[0] * ni1[1];
  unsigned int Nni2 = ni2[0] * ni2[1];
 
  if (!jetCorrectorFormat_) {
    ue_predictor_pf->np.resize(5);
    ue_predictor_pf->ni0.resize(2);
    ue_predictor_pf->ni1.resize(2);
    ue_predictor_pf->ni2.resize(2);
 
    std::copy(np, np + 5, ue_predictor_pf->np.begin());
    std::copy(ni0, ni0 + 2, ue_predictor_pf->ni0.begin());
    std::copy(ni1, ni1 + 2, ue_predictor_pf->ni1.begin());
    std::copy(ni2, ni2 + 2, ue_predictor_pf->ni2.begin());
 
    static const float edge_pseudorapidity[16] = {-5.191,
                                                  -2.650,
                                                  -2.043,
                                                  -1.740,
                                                  -1.479,
                                                  -1.131,
                                                  -0.783,
                                                  -0.522,
                                                  0.522,
                                                  0.783,
                                                  1.131,
                                                  1.479,
                                                  1.740,
                                                  2.043,
                                                  2.650,
                                                  5.191};
 
    ue_predictor_pf->edgeEta.resize(16);
 
    std::copy(edge_pseudorapidity, edge_pseudorapidity + 16, ue_predictor_pf->edgeEta.begin());
  }
 
  std::string line;
 
  while (std::getline(textTable_, line)) {
    if (line.empty() || line[0] == '#') {
      std::cout << " continue " << std::endl;
      continue;
    }
    std::istringstream linestream(line);
    float val;
    int bin0, bin1, bin2, bin3, bin4;
    if (index < Nnp) {
      linestream >> bin0 >> bin1 >> bin2 >> bin3 >> bin4 >> val;
      ue_vec.push_back(val);
    } else if (index < Nnp + Nni0) {
      linestream >> bin0 >> bin1 >> val;
      ue_vec.push_back(val);
    } else if (index < Nnp + Nni0 + Nni1) {
      linestream >> bin0 >> bin1 >> val;
      ue_vec.push_back(val);
    } else if (index < Nnp + Nni0 + Nni1 + Nni2) {
      linestream >> bin0 >> bin1 >> val;
      ue_vec.push_back(val);
    }
    ++index;
  }
 
  edm::Service<cond::service::PoolDBOutputService> pool;
 
  if (pool.isAvailable()) {
    if (jetCorrectorFormat_) {
      // A minimal dummy line that satisfies the JME # token >= 6 requirement, and has the correct key type
      JetCorrectorParameters::Definitions definition("1 0 0 0 Correction L1Offset");
      std::vector<JetCorrectorParameters::Record> record(
          1,
          JetCorrectorParameters::Record(
              ue_vec.size(), std::vector<float>(ue_vec.size(), 0), std::vector<float>(ue_vec.size(), 0), ue_vec));
      JetCorrectorParameters parameter(definition, record);
 
      std::unique_ptr<JetCorrectorParametersCollection> jme_payload =
          std::make_unique<JetCorrectorParametersCollection>();
 
      jme_payload->push_back(JetCorrectorParametersCollection::L1Offset, parameter);
 
      if (pool->isNewTagRequest("JetCorrectionsRecord")) {
        pool->createNewIOV<JetCorrectorParametersCollection>(
            jme_payload.get(), pool->beginOfTime(), pool->endOfTime(), "JetCorrectionsRecord");
      } else {
        pool->appendSinceTime<JetCorrectorParametersCollection>(
            jme_payload.get(), pool->currentTime(), "JetCorrectionsRecord");
      }
    } else {
      ue_predictor_pf->values = ue_vec;
 
      if (pool->isNewTagRequest("HeavyIonUERcd")) {
        pool->createNewIOV<UETable>(ue_predictor_pf.get(), pool->beginOfTime(), pool->endOfTime(), "HeavyIonUERcd");
      } else {
        pool->appendSinceTime<UETable>(ue_predictor_pf.get(), pool->currentTime(), "HeavyIonUERcd");
      }
    }
  }
}
 
void UETableProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
  edm::ParameterSetDescription desc;
  desc.add<std::string>("txtFile", "example");
  desc.addUntracked<bool>("debug", false);
  desc.addUntracked<bool>("jetCorrectorFormat", false);
  descriptions.add("produceUETable", desc);
}
 
//define this as a plug-in
DEFINE_FWK_MODULE(UETableProducer);