File indexing completed on 2022-09-27 05:25:56
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0012
0013 #include "FWCore/Framework/interface/Frameworkfwd.h"
0014 #include "FWCore/Framework/interface/MakerMacros.h"
0015 #include "FWCore/Framework/interface/ESHandle.h"
0016 #include "FWCore/Framework/interface/EventSetup.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/ServiceRegistry/interface/Service.h"
0019
0020 #include "Geometry/VeryForwardGeometryBuilder/interface/CTPPSGeometry.h"
0021 #include "Geometry/Records/interface/VeryForwardRealGeometryRecord.h"
0022
0023 #include "CalibPPS/TimingCalibration/interface/TimingCalibrationStruct.h"
0024 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
0025
0026 #include "DataFormats/CTPPSDetId/interface/CTPPSDiamondDetId.h"
0027 #include "CondFormats/PPSObjects/interface/PPSTimingCalibration.h"
0028 #include "TFitResult.h"
0029
0030
0031 class PPSTimingCalibrationPCLHarvester : public DQMEDHarvester {
0032 public:
0033 PPSTimingCalibrationPCLHarvester(const edm::ParameterSet&);
0034 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0035
0036 static void fillDescriptions(edm::ConfigurationDescriptions&);
0037
0038 private:
0039 void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override;
0040 edm::ESGetToken<CTPPSGeometry, VeryForwardRealGeometryRecord> geomEsToken_;
0041 std::vector<CTPPSDiamondDetId> detids_;
0042 const std::string dqmDir_;
0043 const std::string formula_;
0044 const unsigned int min_entries_;
0045 static constexpr double resolution_ = 0.1;
0046 static constexpr double offset_ = 0.;
0047 TF1 interp_;
0048 };
0049
0050
0051
0052 PPSTimingCalibrationPCLHarvester::PPSTimingCalibrationPCLHarvester(const edm::ParameterSet& iConfig)
0053 : geomEsToken_(esConsumes<edm::Transition::BeginRun>()),
0054 dqmDir_(iConfig.getParameter<std::string>("dqmDir")),
0055 formula_(iConfig.getParameter<std::string>("formula")),
0056 min_entries_(iConfig.getParameter<unsigned int>("minEntries")),
0057 interp_("interp", formula_.c_str(), 10.5, 25.) {
0058
0059 edm::Service<cond::service::PoolDBOutputService> poolDbService;
0060 if (!poolDbService.isAvailable())
0061 throw cms::Exception("PPSTimingCalibrationPCLHarvester") << "PoolDBService required";
0062
0063
0064 interp_.SetParLimits(1, 9., 15.);
0065 interp_.SetParLimits(2, 0.2, 2.5);
0066 }
0067
0068
0069
0070 void PPSTimingCalibrationPCLHarvester::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0071 const auto& geom = iSetup.getData(geomEsToken_);
0072 for (auto it = geom.beginSensor(); it != geom.endSensor(); ++it) {
0073 if (!CTPPSDiamondDetId::check(it->first))
0074 continue;
0075 const CTPPSDiamondDetId detid(it->first);
0076 if (detid.station() == 1)
0077 detids_.emplace_back(detid);
0078 }
0079 }
0080
0081
0082
0083 void PPSTimingCalibrationPCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) {
0084
0085 PPSTimingCalibration::ParametersMap calib_params;
0086 PPSTimingCalibration::TimingMap calib_time;
0087
0088 iGetter.cd();
0089 iGetter.setCurrentFolder(dqmDir_);
0090
0091
0092 TimingCalibrationHistograms hists;
0093 std::string ch_name;
0094 for (const auto& detid : detids_) {
0095 detid.channelName(ch_name);
0096 const auto chid = detid.rawId();
0097 const PPSTimingCalibration::Key key{
0098 (int)detid.arm(), (int)detid.station(), (int)detid.plane(), (int)detid.channel()};
0099
0100 calib_params[key] = {0, 0, 0, 0};
0101 calib_time[key] = std::make_pair(offset_, resolution_);
0102
0103 hists.leadingTime[chid] = iGetter.get(dqmDir_ + "/t_" + ch_name);
0104 if (hists.leadingTime[chid] == nullptr) {
0105 edm::LogInfo("PPSTimingCalibrationPCLHarvester:dqmEndJob")
0106 << "Failed to retrieve leading time monitor for channel (" << detid << ").";
0107 continue;
0108 }
0109 hists.toT[chid] = iGetter.get(dqmDir_ + "/tot_" + ch_name);
0110 if (hists.toT[chid] == nullptr) {
0111 edm::LogInfo("PPSTimingCalibrationPCLHarvester:dqmEndJob")
0112 << "Failed to retrieve time over threshold monitor for channel (" << detid << ").";
0113 continue;
0114 }
0115 hists.leadingTimeVsToT[chid] = iGetter.get(dqmDir_ + "/tvstot_" + ch_name);
0116 if (hists.leadingTimeVsToT[chid] == nullptr) {
0117 edm::LogInfo("PPSTimingCalibrationPCLHarvester:dqmEndJob")
0118 << "Failed to retrieve leading time vs. time over threshold monitor for channel (" << detid << ").";
0119 continue;
0120 }
0121 if (min_entries_ > 0 && hists.leadingTimeVsToT[chid]->getEntries() < min_entries_) {
0122 edm::LogWarning("PPSTimingCalibrationPCLHarvester:dqmEndJob")
0123 << "Not enough entries for channel (" << detid << "): " << hists.leadingTimeVsToT[chid]->getEntries() << " < "
0124 << min_entries_ << ". Skipping calibration.";
0125 continue;
0126 }
0127 const double upper_tot_range = hists.toT[chid]->getMean() + 2.5;
0128 {
0129
0130 std::string ch_name;
0131 detid.channelName(ch_name);
0132 auto profile = iBooker.bookProfile(ch_name + "_prof_x", ch_name + "_prof_x", 240, 0., 60., 450, -20., 25.);
0133
0134 std::unique_ptr<TProfile> prof(hists.leadingTimeVsToT[chid]->getTH2F()->ProfileX("_prof_x", 1, -1));
0135 *(profile->getTProfile()) = *((TProfile*)prof->Clone());
0136 profile->getTProfile()->SetTitle(ch_name.c_str());
0137 profile->getTProfile()->SetName(ch_name.c_str());
0138
0139 interp_.SetParameters(hists.leadingTime[chid]->getRMS(),
0140 hists.toT[chid]->getMean(),
0141 0.8,
0142 hists.leadingTime[chid]->getMean() - hists.leadingTime[chid]->getRMS());
0143 const auto& res = profile->getTProfile()->Fit(&interp_, "B+", "", 10.4, upper_tot_range);
0144 if (!(bool)res) {
0145 calib_params[key] = {
0146 interp_.GetParameter(0), interp_.GetParameter(1), interp_.GetParameter(2), interp_.GetParameter(3)};
0147 calib_time[key] =
0148 std::make_pair(offset_, resolution_);
0149
0150
0151 } else
0152 edm::LogWarning("PPSTimingCalibrationPCLHarvester:dqmEndJob")
0153 << "Fit did not converge for channel (" << detid << ").";
0154 }
0155 }
0156
0157
0158 PPSTimingCalibration calib(formula_, calib_params, calib_time);
0159
0160
0161 edm::Service<cond::service::PoolDBOutputService> poolDbService;
0162 poolDbService->writeOneIOV(calib, poolDbService->currentTime(), "PPSTimingCalibrationRcd_HPTDC");
0163 }
0164
0165
0166
0167 void PPSTimingCalibrationPCLHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0168 edm::ParameterSetDescription desc;
0169 desc.add<std::string>("dqmDir", "AlCaReco/PPSTimingCalibrationPCL")
0170 ->setComment("input path for the various DQM plots");
0171 desc.add<std::string>("formula", "[0]/(exp((x-[1])/[2])+1)+[3]")
0172 ->setComment("interpolation formula for the time walk component");
0173 desc.add<unsigned int>("minEntries", 100)->setComment("minimal number of hits to extract calibration");
0174 descriptions.addWithDefaultLabel(desc);
0175 }
0176
0177 DEFINE_FWK_MODULE(PPSTimingCalibrationPCLHarvester);