TimingCalibration/plugins/PPSDiamondSampicTimingCalibrationPCLHarvester.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    CalibPPS/TimingCalibration/PPSDiamondSampicTimingCalibrationPCLHarvester
0004 // Class:      PPSDiamondSampicTimingCalibrationPCLHarvester
0005 //
0006 /**\class PPSDiamondSampicTimingCalibrationPCLHarvester PPSDiamondSampicTimingCalibrationPCLHarvester.cc CalibPPS/TimingCalibration/PPSDiamondSampicTimingCalibrationPCLHarvester/plugins/PPSDiamondSampicTimingCalibrationPCLHarvester.cc
0007
0008  Description: Harvester of the DiamondSampic calibration which produces sqlite file with a new channel alignment
0009
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013 //
0014 // Original Author:  Christopher Misan
0015 //         Created:  Mon, 26 Jul 2021 16:36:13 GMT
0016 //
0017 //
0018 #include "TAxis.h"
0019 #include "TH1.h"
0020 #include "TArrayD.h"
0021 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0022
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "FWCore/Framework/interface/ESHandle.h"
0026 #include "FWCore/Framework/interface/EventSetup.h"
0027 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0028 #include "FWCore/ServiceRegistry/interface/Service.h"
0029
0030 #include <boost/property_tree/json_parser.hpp>
0031 #include <boost/property_tree/ptree.hpp>
0032
0033 #include "Geometry/VeryForwardGeometryBuilder/interface/CTPPSGeometry.h"
0034 #include "Geometry/Records/interface/VeryForwardRealGeometryRecord.h"
0035
0036 #include "CalibPPS/TimingCalibration/interface/TimingCalibrationStruct.h"
0037 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
0038
0039 #include "DataFormats/CTPPSDetId/interface/CTPPSDiamondDetId.h"
0040 #include "CondFormats/PPSObjects/interface/PPSTimingCalibration.h"
0041 #include "CondFormats/DataRecord/interface/PPSTimingCalibrationRcd.h"
0042
0043 namespace pt = boost::property_tree;
0044
0045 class PPSDiamondSampicTimingCalibrationPCLHarvester : public DQMEDHarvester {
0046 public:
0047   PPSDiamondSampicTimingCalibrationPCLHarvester(const edm::ParameterSet&);
0048   void beginRun(const edm::Run&, const edm::EventSetup&) override;
0049   static void fillDescriptions(edm::ConfigurationDescriptions&);
0050
0051 private:
0052   void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override;
0053   void calibJson(DQMStore::IGetter& iGetter);
0054   void calibDb(DQMStore::IGetter& iGetter);
0055   bool getDbSampicChannel(
0056       DQMStore::IGetter& iGetter, int& db, int& sampic, int& channel, std::string ch_name, CTPPSDiamondDetId detid);
0057   edm::ESGetToken<CTPPSGeometry, VeryForwardRealGeometryRecord> geomEsToken_;
0058   edm::ESGetToken<PPSTimingCalibration, PPSTimingCalibrationRcd> timingCalibrationToken_;
0059   edm::ESHandle<PPSTimingCalibration> hTimingCalib_;
0060   std::vector<CTPPSDiamondDetId> detids_;
0061   const std::string dqmDir_;
0062   const unsigned int min_entries_;
0063   const std::string jsonCalibFile_, jsonOutputPath_;
0064 };
0065
0066 //------------------------------------------------------------------------------
0067
0068 PPSDiamondSampicTimingCalibrationPCLHarvester::PPSDiamondSampicTimingCalibrationPCLHarvester(
0069     const edm::ParameterSet& iConfig)
0070     : geomEsToken_(esConsumes<edm::Transition::BeginRun>()),
0071       dqmDir_(iConfig.getParameter<std::string>("dqmDir")),
0072       min_entries_(iConfig.getParameter<unsigned int>("minEntries")),
0073       jsonCalibFile_(iConfig.getParameter<std::string>("jsonCalibFile")),
0074       jsonOutputPath_(iConfig.getParameter<std::string>("jsonOutputPath")) {
0075   if (jsonCalibFile_.empty())
0076     timingCalibrationToken_ = esConsumes<edm::Transition::BeginRun>(
0077         edm::ESInputTag(iConfig.getParameter<std::string>("timingCalibrationTag")));
0078 }
0079
0080 //------------------------------------------------------------------------------
0081
0082 void PPSDiamondSampicTimingCalibrationPCLHarvester::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0083   if (jsonCalibFile_.empty())
0084     hTimingCalib_ = iSetup.getHandle(timingCalibrationToken_);
0085   const auto& geom = iSetup.getData(geomEsToken_);
0086   for (auto it = geom.beginSensor(); it != geom.endSensor(); ++it) {
0087     if (!CTPPSDiamondDetId::check(it->first))
0088       continue;
0089     const CTPPSDiamondDetId detid(it->first);
0090     detids_.emplace_back(detid);
0091   }
0092 }
0093
0094 //------------------------------------------------------------------------------
0095
0096 bool PPSDiamondSampicTimingCalibrationPCLHarvester::getDbSampicChannel(
0097     DQMStore::IGetter& iGetter, int& db, int& sampic, int& channel, std::string path, CTPPSDiamondDetId detid) {
0098   auto histDb = iGetter.get(path + "db");
0099   auto histSampic = iGetter.get(path + "sampic");
0100   auto histChannel = iGetter.get(path + "channel");
0101
0102   if (histDb == nullptr) {
0103     edm::LogWarning("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0104         << "Failed to retrieve db for detid: " << detid;
0105     return false;
0106   }
0107
0108   if (histSampic == nullptr) {
0109     edm::LogWarning("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0110         << "Failed to retrieve sampic for detid: " << detid;
0111     return false;
0112   }
0113
0114   if (histChannel == nullptr) {
0115     edm::LogWarning("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0116         << "Failed to retrieve channel hwId for detid: " << detid;
0117     return false;
0118   }
0119
0120   db = histDb->getIntValue();
0121   sampic = histSampic->getIntValue();
0122   channel = histChannel->getIntValue();
0123
0124   return true;
0125 }
0126 //------------------------------------------------------------------------------
0127
0128 void PPSDiamondSampicTimingCalibrationPCLHarvester::calibJson(DQMStore::IGetter& iGetter) {
0129   std::unordered_map<uint32_t, dqm::reco::MonitorElement*> timeHisto;
0130   std::string ch_name, path;
0131
0132   pt::ptree node;
0133   pt::read_json(jsonCalibFile_, node);
0134   const std::string formula = node.get<std::string>("formula");
0135
0136   for (const auto& detid : detids_) {
0137     detid.channelName(path, CTPPSDiamondDetId::nPath);
0138     detid.channelName(ch_name);
0139     path = dqmDir_ + "/" + path + "/" + ch_name;
0140     const auto chid = detid.rawId();
0141
0142     int db, sampic, channel;
0143     if (!getDbSampicChannel(iGetter, db, sampic, channel, path, detid))
0144       continue;
0145
0146     int ct = 0;
0147     for (pt::ptree::value_type& par : node.get_child("parameters." + std::to_string(db))) {
0148       double new_time_offset;
0149       if (ct == 16 * (1 - sampic) + channel) {  //flip the calibration - sampic 1 is first in json
0150         double old_time_offset = par.second.get<double>("time_offset");
0151
0152         timeHisto[chid] = iGetter.get(path);
0153
0154         if (timeHisto[chid] == nullptr) {
0155           edm::LogWarning("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0156               << "Failed to retrieve time monitor for detid" << detid;
0157           par.second.put<double>("time_offset", old_time_offset);
0158           continue;
0159         }
0160
0161         if (min_entries_ > 0 && timeHisto[chid]->getEntries() < min_entries_) {
0162           edm::LogWarning("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0163               << "Not enough entries for channel (" << detid << "): " << timeHisto[chid]->getEntries() << " < "
0164               << min_entries_ << ". Skipping calibration.";
0165           par.second.put<double>("time_offset", old_time_offset);
0166           continue;
0167         }
0168         new_time_offset = old_time_offset - timeHisto[chid]->getMean();
0169         //scale x axis of the plots by calculated offset
0170         timeHisto[chid]->getTH1F()->GetXaxis()->Set(
0171             timeHisto[chid]->getTH1F()->GetXaxis()->GetNbins(),
0172             timeHisto[chid]->getTH1F()->GetXaxis()->GetXmin() + new_time_offset,   // new Xmin
0173             timeHisto[chid]->getTH1F()->GetXaxis()->GetXmax() + new_time_offset);  // new Xmax
0174         timeHisto[chid]->getTH1F()->ResetStats();
0175
0176         par.second.put<double>("time_offset", new_time_offset);
0177         break;
0178       }
0179       ct++;
0180     }
0181   }
0182   pt::write_json(jsonOutputPath_, node);
0183 }
0184
0185 //------------------------------------------------------------------------------
0186
0187 void PPSDiamondSampicTimingCalibrationPCLHarvester::calibDb(DQMStore::IGetter& iGetter) {
0188   PPSTimingCalibration calib = *hTimingCalib_;
0189
0190   // book the parameters containers
0191   PPSTimingCalibration::ParametersMap params;
0192   PPSTimingCalibration::TimingMap time_info;
0193
0194   std::unordered_map<uint32_t, dqm::reco::MonitorElement*> timeHisto;
0195   std::string rp_name, plane_name, ch_name, path;
0196   const std::string& formula = calib.formula();
0197
0198   for (const auto& detid : detids_) {
0199     detid.channelName(path, CTPPSDiamondDetId::nPath);
0200     detid.channelName(ch_name);
0201     path = dqmDir_ + "/" + path + "/" + ch_name;
0202     const auto chid = detid.rawId();
0203
0204     int db, sampic, channel;
0205     if (!getDbSampicChannel(iGetter, db, sampic, channel, path, detid))
0206       continue;
0207
0208     PPSTimingCalibration::Key key;
0209     key.db = db;
0210     key.sampic = sampic;
0211     key.channel = channel;
0212
0213     double timeOffset = calib.timeOffset(db, sampic, channel);
0214     double timePrecision = calib.timePrecision(db, sampic, channel);
0215     if (timeOffset == 0 && timePrecision == 0) {
0216       edm::LogWarning("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0217           << "No calibration found for db: " << db << " sampic: " << sampic << " channel: " << channel;
0218       continue;
0219     }
0220
0221     int cell_ct = 0;
0222     while (!calib.parameters(db, sampic, channel, cell_ct).empty()) {
0223       auto parameters = calib.parameters(db, sampic, channel, cell_ct);
0224       key.cell = cell_ct;
0225       params[key] = parameters;
0226       cell_ct++;
0227     }
0228
0229     key.cell = -1;
0230
0231     time_info[key] = {timeOffset, timePrecision};
0232     timeHisto[chid] = iGetter.get(path);
0233     if (timeHisto[chid] == nullptr) {
0234       edm::LogWarning("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0235           << "Failed to retrieve time monitor for detid: " << detid;
0236       continue;
0237     }
0238
0239     if (min_entries_ > 0 && timeHisto[chid]->getEntries() < min_entries_) {
0240       edm::LogInfo("PPSDiamondSampicTimingCalibrationPCLHarvester:dqmEndJob")
0241           << "Not enough entries (" << detid << "): " << timeHisto[chid]->getEntries() << " < " << min_entries_
0242           << ". Skipping calibration.";
0243       continue;
0244     }
0245
0246     double new_time_offset = timeOffset - timeHisto[chid]->getMean();
0247     //scale x axis of the plots by calculated offset
0248     timeHisto[chid]->getTH1F()->GetXaxis()->Set(
0249         timeHisto[chid]->getTH1F()->GetXaxis()->GetNbins(),
0250         timeHisto[chid]->getTH1F()->GetXaxis()->GetXmin() + new_time_offset,   // new Xmin
0251         timeHisto[chid]->getTH1F()->GetXaxis()->GetXmax() + new_time_offset);  // new Xmax
0252     timeHisto[chid]->getTH1F()->ResetStats();
0253
0254     time_info[key] = {new_time_offset, timePrecision};
0255   }
0256
0257   auto calibPPS = PPSTimingCalibration(formula, params, time_info);
0258   // write the object
0259   edm::Service<cond::service::PoolDBOutputService> poolDbService;
0260   poolDbService->writeOneIOV(calibPPS, poolDbService->currentTime(), "PPSTimingCalibrationRcd_SAMPIC");
0261 }
0262
0263 //------------------------------------------------------------------------------
0264
0265 void PPSDiamondSampicTimingCalibrationPCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) {
0266   iBooker.cd();
0267   iBooker.setCurrentFolder(dqmDir_);
0268   if (jsonCalibFile_.empty())
0269     calibDb(iGetter);
0270   else
0271     calibJson(iGetter);
0272 }
0273
0274 //------------------------------------------------------------------------------
0275
0276 void PPSDiamondSampicTimingCalibrationPCLHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0277   edm::ParameterSetDescription desc;
0278   desc.add<std::string>("timingCalibrationTag", "GlobalTag:PPSDiamondSampicCalibration")
0279       ->setComment("input tag for timing calibration retrieval");
0280   desc.add<std::string>("dqmDir", "AlCaReco/PPSDiamondSampicTimingCalibrationPCL")
0281       ->setComment("input path for the various DQM plots");
0282   desc.add<unsigned int>("minEntries", 1)->setComment("minimal number of hits to extract calibration");
0283   desc.add<std::string>("jsonCalibFile", "")
0284       ->setComment(
0285           "input path for json file containing calibration, if none, calibration will be obtained from db instead");
0286   desc.add<std::string>("jsonOutputPath", "offset_cal.json")->setComment("output path for the new json calibration");
0287   descriptions.add("PPSDiamondSampicTimingCalibrationPCLHarvester", desc);
0288 }
0289
0290 DEFINE_FWK_MODULE(PPSDiamondSampicTimingCalibrationPCLHarvester);