File indexing completed on 2024-04-06 12:07:33
0001 #include "DQM/HcalTasks/interface/HcalOnlineHarvesting.h"
0002
0003 using namespace hcaldqm;
0004 using namespace hcaldqm::constants;
0005
0006 HcalOnlineHarvesting::HcalOnlineHarvesting(edm::ParameterSet const& ps)
0007 : DQHarvester(ps), _nBad(0), _nTotal(0), _reportSummaryMap(nullptr) {
0008
0009
0010 _vsumgen.resize(nSummary);
0011 _vnames.resize(nSummary);
0012 _vmarks.resize(nSummary);
0013 for (uint32_t i = 0; i < _vmarks.size(); i++)
0014 _vmarks[i] = false;
0015 _vnames[fRaw] = "RawTask";
0016 _vnames[fDigi] = "DigiTask";
0017 _vnames[fReco] = "RecHitTask";
0018 _vnames[fTP] = "TPTask";
0019 _vnames[fPedestal] = "PedestalTask";
0020
0021 auto iC = consumesCollector();
0022 _vsumgen[fRaw] = new hcaldqm::RawRunSummary("RawRunHarvesting", _vnames[fRaw], ps, iC);
0023 _vsumgen[fDigi] = new hcaldqm::DigiRunSummary("DigiRunHarvesting", _vnames[fDigi], ps, iC);
0024 _vsumgen[fReco] = new hcaldqm::RecoRunSummary("RecoRunHarvesting", _vnames[fReco], ps, iC);
0025 _vsumgen[fTP] = new hcaldqm::TPRunSummary("TPRunHarvesting", _vnames[fTP], ps, iC);
0026 _vsumgen[fPedestal] = new hcaldqm::PedestalRunSummary("PedestalRunHarvesting", _vnames[fPedestal], ps, iC);
0027
0028 _thresh_bad_bad = ps.getUntrackedParameter("thresh_bad_bad", 0.05);
0029 }
0030
0031 void HcalOnlineHarvesting::beginRun(edm::Run const& r, edm::EventSetup const& es) {
0032 DQHarvester::beginRun(r, es);
0033 for (std::vector<DQClient*>::const_iterator it = _vsumgen.begin(); it != _vsumgen.end(); ++it)
0034 (*it)->beginRun(r, es);
0035 }
0036
0037 void HcalOnlineHarvesting::_dqmEndLuminosityBlock(DQMStore::IBooker& ib,
0038 DQMStore::IGetter& ig,
0039 edm::LuminosityBlock const&,
0040 edm::EventSetup const&) {
0041
0042 if (ig.get(_subsystem + "/" + _vnames[fRaw] + "/EventsTotal") != nullptr)
0043 _vmarks[fRaw] = true;
0044 if (ig.get(_subsystem + "/" + _vnames[fDigi] + "/EventsTotal") != nullptr)
0045 _vmarks[fDigi] = true;
0046 if (ig.get(_subsystem + "/" + _vnames[fTP] + "/EventsTotal") != nullptr)
0047 _vmarks[fTP] = true;
0048 if (ig.get(_subsystem + "/" + _vnames[fReco] + "/EventsTotal") != nullptr)
0049 _vmarks[fReco] = true;
0050 if (ig.get(_subsystem + "/" + _vnames[fPedestal] + "/EventsTotal") != nullptr)
0051 _vmarks[fPedestal] = true;
0052
0053
0054
0055 if (!_reportSummaryMap) {
0056 ig.setCurrentFolder(_subsystem + "/EventInfo");
0057 _reportSummaryMap =
0058 ib.book2D("reportSummaryMap", "reportSummaryMap", _maxLS, 1, _maxLS + 1, _vFEDs.size(), 0, _vFEDs.size());
0059 for (uint32_t i = 0; i < _vFEDs.size(); i++) {
0060 char name[5];
0061 sprintf(name, "%d", _vFEDs[i]);
0062 _reportSummaryMap->setBinLabel(i + 1, name, 2);
0063 }
0064
0065 _reportSummaryMap->getTH1()->SetBit(BIT(BIT_OFFSET + BIT_AXIS_LS));
0066
0067
0068 for (uint32_t i = 0; i < _vnames.size(); i++)
0069 _vcSummaryvsLS.push_back(ContainerSingle2D(_vnames[i],
0070 "SummaryvsLS",
0071 new hcaldqm::quantity::LumiSection(_maxLS),
0072 new hcaldqm::quantity::FEDQuantity(_vFEDs),
0073 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fState)));
0074
0075
0076 for (uint32_t i = 0; i < _vmarks.size(); i++) {
0077 if (_vmarks[i])
0078 _vcSummaryvsLS[i].load(ig, _subsystem);
0079 }
0080
0081
0082 _cKnownBadChannels_depth.initialize("RunInfo",
0083 "KnownBadChannels",
0084 hcaldqm::hashfunctions::fdepth,
0085 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fieta),
0086 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fiphi),
0087 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN),
0088 0);
0089 _cKnownBadChannels_depth.book(ib, _emap, _subsystem);
0090 for (uintCompactMap::const_iterator it = _xQuality.begin(); it != _xQuality.end(); ++it)
0091 _cKnownBadChannels_depth.fill(HcalDetId(it->first));
0092
0093 ig.setCurrentFolder(_subsystem + "/EventInfo");
0094 _runSummary = ib.book2D("runSummary", "runSummary", 1, 0, 1, 1, 0, 1);
0095 }
0096
0097 int ifed = 0;
0098 hcaldqm::flag::Flag fTotal("Status", hcaldqm::flag::fNCDAQ);
0099 if (_ptype != fOffline) {
0100 for (std::vector<uint32_t>::const_iterator it = _vhashFEDs.begin(); it != _vhashFEDs.end(); ++it) {
0101 HcalElectronicsId eid(*it);
0102 hcaldqm::flag::Flag fSum("Status", hcaldqm::flag::fNCDAQ);
0103 for (uint32_t im = 0; im < _vmarks.size(); im++)
0104 if (_vmarks[im]) {
0105 int x = _vcSummaryvsLS[im].getBinContent(eid, _currentLS);
0106 hcaldqm::flag::Flag flag("Status", (hcaldqm::flag::State)x);
0107 fSum += flag;
0108 }
0109 _reportSummaryMap->setBinContent(_currentLS, ifed + 1, int(fSum._state));
0110 ifed++;
0111 fTotal += fSum;
0112 }
0113 }
0114
0115
0116
0117 if (fTotal._state == hcaldqm::flag::fBAD)
0118 _nBad++;
0119 _nTotal++;
0120 if (double(_nBad) / double(_nTotal) >= _thresh_bad_bad)
0121 _runSummary->setBinContent(1, 1, int(hcaldqm::flag::fBAD));
0122 else if (fTotal._state == hcaldqm::flag::fNCDAQ)
0123 _runSummary->setBinContent(1, 1, int(hcaldqm::flag::fNCDAQ));
0124 else
0125 _runSummary->setBinContent(1, 1, int(hcaldqm::flag::fGOOD));
0126
0127
0128 if (_vmarks[fTP]) {
0129 MonitorElement* meOccupancy_HF_depth = ig.get("Hcal/TPTask/OccupancyDataHF_depth/OccupancyDataHF_depth");
0130 MonitorElement* meOccupancyNoTDC_HF_depth =
0131 ig.get("Hcal/TPTask/OccupancyEmulHFNoTDC_depth/OccupancyEmulHFNoTDC_depth");
0132 MonitorElement* meOccupancy_HF_ieta = ig.get("Hcal/TPTask/OccupancyDataHF_ieta/OccupancyDataHF_ieta");
0133 MonitorElement* meOccupancyNoTDC_HF_ieta =
0134 ig.get("Hcal/TPTask/OccupancyEmulHFNoTDC_ieta/OccupancyEmulHFNoTDC_ieta");
0135
0136 if (meOccupancy_HF_depth && meOccupancyNoTDC_HF_depth && meOccupancy_HF_ieta && meOccupancyNoTDC_HF_ieta) {
0137 TH2F* hOccupancy_HF_depth = meOccupancy_HF_depth->getTH2F();
0138 TH2F* hOccupancyNoTDC_HF_depth = meOccupancyNoTDC_HF_depth->getTH2F();
0139 TH1D* hOccupancy_HF_ieta = meOccupancy_HF_ieta->getTH1D();
0140 TH1D* hOccupancyNoTDC_HF_ieta = meOccupancyNoTDC_HF_ieta->getTH1D();
0141
0142 TH2F* hEfficiency_HF_depth = (TH2F*)hOccupancy_HF_depth->Clone();
0143 hEfficiency_HF_depth->Divide(hOccupancyNoTDC_HF_depth);
0144 TH1D* hEfficiency_HF_ieta = (TH1D*)hOccupancy_HF_ieta->Clone();
0145 hEfficiency_HF_ieta->Divide(hOccupancyNoTDC_HF_ieta);
0146
0147 ib.setCurrentFolder("Hcal/TPTask");
0148
0149 MonitorElement* meEfficiency_HF_depth = ib.book2D("TDCCutEfficiency_depth", hEfficiency_HF_depth);
0150 meEfficiency_HF_depth->setEfficiencyFlag();
0151 MonitorElement* meEfficiency_HF_ieta = ib.book1DD("TDCCutEfficiency_ieta", hEfficiency_HF_ieta);
0152 meEfficiency_HF_ieta->setEfficiencyFlag();
0153
0154 delete hEfficiency_HF_depth;
0155 delete hEfficiency_HF_ieta;
0156 }
0157 }
0158 }
0159
0160
0161
0162
0163 void HcalOnlineHarvesting::_dqmEndJob(DQMStore::IBooker& ib, DQMStore::IGetter& ig) {}
0164
0165 DEFINE_FWK_MODULE(HcalOnlineHarvesting);