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->createOneIOV<JetCorrectorParametersCollection>(*jme_payload, pool->beginOfTime(), "JetCorrectionsRecord");
      } else {
        pool->appendOneIOV<JetCorrectorParametersCollection>(*jme_payload, pool->currentTime(), "JetCorrectionsRecord");
      }
    } else {
      ue_predictor_pf->values = ue_vec;

      if (pool->isNewTagRequest("HeavyIonUERcd")) {
        pool->createOneIOV<UETable>(*ue_predictor_pf, pool->beginOfTime(), "HeavyIonUERcd");
      } else {
        pool->appendOneIOV<UETable>(*ue_predictor_pf, 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);