File indexing completed on 2025-07-15 22:06:19
0001 #include "DQM/HcalTasks/interface/LaserTask.h"
0002
0003 using namespace hcaldqm;
0004 using namespace hcaldqm::constants;
0005 using namespace hcaldqm::filter;
0006 LaserTask::LaserTask(edm::ParameterSet const& ps)
0007 : DQTask(ps), hcalDbServiceToken_(esConsumes<HcalDbService, HcalDbRecord, edm::Transition::BeginRun>()) {
0008 _nevents = ps.getUntrackedParameter<int>("nevents", 2000);
0009
0010
0011 _tagQIE11 = ps.getUntrackedParameter<edm::InputTag>("tagHE", edm::InputTag("hcalDigis"));
0012 _tagHO = ps.getUntrackedParameter<edm::InputTag>("tagHO", edm::InputTag("hcalDigis"));
0013 _tagQIE10 = ps.getUntrackedParameter<edm::InputTag>("tagHF", edm::InputTag("hcalDigis"));
0014 _taguMN = ps.getUntrackedParameter<edm::InputTag>("taguMN", edm::InputTag("hcalDigis"));
0015 _tagLaserMon = ps.getUntrackedParameter<edm::InputTag>("tagLaserMon", edm::InputTag("LASERMON"));
0016 _tokQIE11 = consumes<QIE11DigiCollection>(_tagQIE11);
0017 _tokHO = consumes<HODigiCollection>(_tagHO);
0018 _tokQIE10 = consumes<QIE10DigiCollection>(_tagQIE10);
0019 _tokuMN = consumes<HcalUMNioDigi>(_taguMN);
0020 _tokLaserMon = consumes<QIE10DigiCollection>(_tagLaserMon);
0021 _tagFEDs = ps.getUntrackedParameter<edm::InputTag>("tagFEDs", edm::InputTag("hltHcalCalibrationRaw"));
0022 _tokFEDs = consumes<FEDRawDataCollection>(_tagFEDs);
0023
0024 _vflags.resize(nLaserFlag);
0025 _vflags[fBadTiming] = hcaldqm::flag::Flag("BadTiming");
0026 _vflags[fMissingLaserMon] = hcaldqm::flag::Flag("MissingLaserMon");
0027
0028
0029 _lowHBHE = ps.getUntrackedParameter<double>("lowHBHE", 50);
0030 _lowHO = ps.getUntrackedParameter<double>("lowHO", 100);
0031 _lowHF = ps.getUntrackedParameter<double>("lowHF", 50);
0032 _laserType = (uint32_t)ps.getUntrackedParameter<uint32_t>("laserType");
0033
0034
0035 _vLaserMonIPhi = ps.getUntrackedParameter<std::vector<int>>("vLaserMonIPhi");
0036 _laserMonIEta = ps.getUntrackedParameter<int>("laserMonIEta");
0037 _laserMonCBox = ps.getUntrackedParameter<int>("laserMonCBox");
0038 _laserMonDigiOverlap = ps.getUntrackedParameter<int>("laserMonDigiOverlap");
0039 _laserMonTS0 = ps.getUntrackedParameter<int>("laserMonTS0");
0040 _laserMonThreshold = ps.getUntrackedParameter<double>("laserMonThreshold", 1.e5);
0041 _thresh_frac_timingreflm = ps.getUntrackedParameter<double>("_thresh_frac_timingreflm", 0.01);
0042 _thresh_min_lmsumq = ps.getUntrackedParameter<double>("thresh_min_lmsumq", 50000.);
0043
0044 std::vector<double> vTimingRangeHB =
0045 ps.getUntrackedParameter<std::vector<double>>("thresh_timingreflm_HB", std::vector<double>({-70, -10.}));
0046 std::vector<double> vTimingRangeHE =
0047 ps.getUntrackedParameter<std::vector<double>>("thresh_timingreflm_HE", std::vector<double>({-60., 0.}));
0048 std::vector<double> vTimingRangeHO =
0049 ps.getUntrackedParameter<std::vector<double>>("thresh_timingreflm_HO", std::vector<double>({-50., 20.}));
0050 std::vector<double> vTimingRangeHF =
0051 ps.getUntrackedParameter<std::vector<double>>("thresh_timingreflm_HF", std::vector<double>({-50., 20.}));
0052 _thresh_timingreflm[HcalBarrel] = std::make_pair(vTimingRangeHB[0], vTimingRangeHB[1]);
0053 _thresh_timingreflm[HcalEndcap] = std::make_pair(vTimingRangeHE[0], vTimingRangeHE[1]);
0054 _thresh_timingreflm[HcalOuter] = std::make_pair(vTimingRangeHO[0], vTimingRangeHO[1]);
0055 _thresh_timingreflm[HcalForward] = std::make_pair(vTimingRangeHF[0], vTimingRangeHF[1]);
0056 }
0057
0058 void LaserTask::bookHistograms(DQMStore::IBooker& ib, edm::Run const& r, edm::EventSetup const& es) {
0059 if (_ptype == fLocal)
0060 if (r.runAuxiliary().run() == 1)
0061 return;
0062
0063 DQTask::bookHistograms(ib, r, es);
0064
0065 edm::ESHandle<HcalDbService> dbService = es.getHandle(hcalDbServiceToken_);
0066 _emap = dbService->getHcalMapping();
0067
0068 std::vector<uint32_t> vhashVME;
0069 std::vector<uint32_t> vhashuTCA;
0070 std::vector<uint32_t> vhashC36;
0071 vhashVME.push_back(
0072 HcalElectronicsId(constants::FIBERCH_MIN, constants::FIBER_VME_MIN, SPIGOT_MIN, CRATE_VME_MIN).rawId());
0073 vhashuTCA.push_back(HcalElectronicsId(CRATE_uTCA_MIN, SLOT_uTCA_MIN, FIBER_uTCA_MIN1, FIBERCH_MIN, false).rawId());
0074 _filter_VME.initialize(filter::fFilter, hcaldqm::hashfunctions::fElectronics, vhashVME);
0075 _filter_uTCA.initialize(filter::fFilter, hcaldqm::hashfunctions::fElectronics, vhashuTCA);
0076
0077
0078 _cSignalMean_Subdet.initialize(_name,
0079 "SignalMean",
0080 hcaldqm::hashfunctions::fSubdet,
0081 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10fC_100000Coarse),
0082 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0083 0);
0084 _cSignalRMS_Subdet.initialize(_name,
0085 "SignalRMS",
0086 hcaldqm::hashfunctions::fSubdet,
0087 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_1000),
0088 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0089 0);
0090 _cTimingMean_Subdet.initialize(_name,
0091 "TimingMean",
0092 hcaldqm::hashfunctions::fSubdet,
0093 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0094 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0095 0);
0096 _cTimingRMS_Subdet.initialize(_name,
0097 "TimingRMS",
0098 hcaldqm::hashfunctions::fSubdet,
0099 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0100 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0101 0);
0102
0103 _cADC_SubdetPM.initialize(_name,
0104 "ADC",
0105 hcaldqm::hashfunctions::fSubdetPM,
0106 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fADC_128),
0107 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0108 0);
0109
0110 if (_ptype == fLocal) {
0111 _cSignalMean_FEDuTCA.initialize(_name,
0112 "SignalMean",
0113 hcaldqm::hashfunctions::fFED,
0114 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fSlotuTCA),
0115 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fFiberuTCAFiberCh),
0116 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10fC_100000Coarse),
0117 0);
0118 _cSignalRMS_FEDuTCA.initialize(_name,
0119 "SignalRMS",
0120 hcaldqm::hashfunctions::fFED,
0121 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fSlotuTCA),
0122 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fFiberuTCAFiberCh),
0123 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_3000),
0124 0);
0125 _cTimingMean_FEDuTCA.initialize(_name,
0126 "TimingMean",
0127 hcaldqm::hashfunctions::fFED,
0128 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fSlotuTCA),
0129 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fFiberuTCAFiberCh),
0130 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0131 0);
0132 _cTimingRMS_FEDuTCA.initialize(_name,
0133 "TimingRMS",
0134 hcaldqm::hashfunctions::fFED,
0135 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fSlotuTCA),
0136 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fFiberuTCAFiberCh),
0137 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0138 0);
0139 _cMissing_FEDuTCA.initialize(_name,
0140 "Missing",
0141 hcaldqm::hashfunctions::fFED,
0142 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fSlotuTCA),
0143 new hcaldqm::quantity::ElectronicsQuantity(hcaldqm::quantity::fFiberuTCAFiberCh),
0144 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN),
0145 0);
0146
0147 _cShapeCut_FEDSlot.initialize(_name,
0148 "Shape",
0149 hcaldqm::hashfunctions::fFEDSlot,
0150 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0151 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_3000),
0152 0);
0153 }
0154 _cTimingvsEvent_SubdetPM.initialize(_name,
0155 "TimingvsEvent",
0156 hcaldqm::hashfunctions::fSubdetPM,
0157 new hcaldqm::quantity::EventNumber(_nevents),
0158 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0159 0);
0160 _cSignalvsEvent_SubdetPM.initialize(_name,
0161 "SignalvsEvent",
0162 hcaldqm::hashfunctions::fSubdetPM,
0163 new hcaldqm::quantity::EventNumber(_nevents),
0164 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_3000),
0165 0);
0166 _cTimingvsLS_SubdetPM.initialize(_name,
0167 "TimingvsLS",
0168 hcaldqm::hashfunctions::fSubdetPM,
0169 new hcaldqm::quantity::LumiSection(_maxLS),
0170 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0171 0);
0172 _cSignalvsLS_SubdetPM.initialize(_name,
0173 "SignalvsLS",
0174 hcaldqm::hashfunctions::fSubdetPM,
0175 new hcaldqm::quantity::LumiSection(_maxLS),
0176 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_3000),
0177 0);
0178 _cTimingvsBX_SubdetPM.initialize(_name,
0179 "TimingvsBX",
0180 hcaldqm::hashfunctions::fSubdetPM,
0181 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fBX),
0182 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0183 0);
0184 _cSignalvsBX_SubdetPM.initialize(_name,
0185 "SignalvsBX",
0186 hcaldqm::hashfunctions::fSubdetPM,
0187 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fBX),
0188 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_3000),
0189 0);
0190
0191 _cSignalMean_depth.initialize(_name,
0192 "SignalMean",
0193 hcaldqm::hashfunctions::fdepth,
0194 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fieta),
0195 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fiphi),
0196 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10fC_100000Coarse),
0197 0);
0198 _cSignalRMS_depth.initialize(_name,
0199 "SignalRMS",
0200 hcaldqm::hashfunctions::fdepth,
0201 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fieta),
0202 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fiphi),
0203 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_1000),
0204 0);
0205 _cTimingMean_depth.initialize(_name,
0206 "TimingMean",
0207 hcaldqm::hashfunctions::fdepth,
0208 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fieta),
0209 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fiphi),
0210 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0211 0);
0212 _cTimingRMS_depth.initialize(_name,
0213 "TimingRMS",
0214 hcaldqm::hashfunctions::fdepth,
0215 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fieta),
0216 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fiphi),
0217 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS200),
0218 0);
0219
0220 _cMissing_depth.initialize(_name,
0221 "Missing",
0222 hcaldqm::hashfunctions::fdepth,
0223 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fieta),
0224 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fiphi),
0225 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN),
0226 0);
0227
0228 if (_ptype == fOnline || _ptype == fLocal) {
0229 _cADCvsTS_SubdetPM.initialize(_name,
0230 "ADCvsTS",
0231 hcaldqm::hashfunctions::fSubdetPM,
0232 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0233 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0234 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0235 0);
0236 }
0237
0238
0239 _xSignalSum.initialize(hcaldqm::hashfunctions::fDChannel);
0240 _xSignalSum2.initialize(hcaldqm::hashfunctions::fDChannel);
0241 _xTimingSum.initialize(hcaldqm::hashfunctions::fDChannel);
0242 _xTimingSum2.initialize(hcaldqm::hashfunctions::fDChannel);
0243 _xEntries.initialize(hcaldqm::hashfunctions::fDChannel);
0244
0245
0246 if (_ptype == fOnline || _ptype == fLocal) {
0247 _cLaserMonSumQ.initialize(_name,
0248 "LaserMonSumQ",
0249 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_1000000),
0250 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN),
0251 0);
0252 _cLaserMonSumQ.showOverflowX(true);
0253
0254 _cLaserMonTiming.initialize(_name,
0255 "LaserMonTiming",
0256 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTime_ns_250),
0257 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN),
0258 0);
0259 _cLaserMonTiming.showOverflowX(true);
0260
0261 if (_ptype == fOnline) {
0262 _cLaserMonSumQ_LS.initialize(_name,
0263 "LaserMonSumQ_LS",
0264 new hcaldqm::quantity::LumiSectionCoarse(_maxLS, 10),
0265 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_1000000));
0266 _cLaserMonTiming_LS.initialize(_name,
0267 "LaserMonTiming_LS",
0268 new hcaldqm::quantity::LumiSectionCoarse(_maxLS, 10),
0269 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_100TS));
0270 _cTimingDiffLS_SubdetPM.initialize(_name,
0271 "TimingDiff_DigiMinusLaserMon",
0272 hcaldqm::hashfunctions::fSubdetPM,
0273 new hcaldqm::quantity::LumiSectionCoarse(_maxLS, 10),
0274 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fRBX),
0275 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTimingDiff_ns),
0276 0);
0277 _cTimingDiffLS_SubdetPM.showOverflowY(true);
0278 } else if (_ptype == fLocal) {
0279 _cLaserMonSumQ_Event.initialize(_name,
0280 "LaserMonSumQ_Event",
0281 new hcaldqm::quantity::EventNumber(_nevents),
0282 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::ffC_1000000));
0283 _cLaserMonTiming_Event.initialize(_name,
0284 "LaserMonTiming_Event",
0285 new hcaldqm::quantity::EventNumber(_nevents),
0286 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_100TS));
0287 _cTimingDiffEvent_SubdetPM.initialize(_name,
0288 "TimingDiff_DigiMinusLaserMon",
0289 hcaldqm::hashfunctions::fSubdetPM,
0290 new hcaldqm::quantity::EventNumber(_nevents),
0291 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fRBX),
0292 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTimingDiff_ns),
0293 0);
0294 _cTimingDiffEvent_SubdetPM.showOverflowY(true);
0295 }
0296 _cTiming_DigivsLaserMon_SubdetPM.initialize(
0297 _name,
0298 "Timing_DigivsLaserMon",
0299 hcaldqm::hashfunctions::fSubdetPM,
0300 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTime_ns_250_coarse),
0301 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTime_ns_250_coarse),
0302 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN),
0303 0);
0304 _cTiming_DigivsLaserMon_SubdetPM.showOverflowX(true);
0305 _cTiming_DigivsLaserMon_SubdetPM.showOverflowY(true);
0306
0307 _xTimingRefLMSum.initialize(hcaldqm::hashfunctions::fDChannel);
0308 _xTimingRefLMSum2.initialize(hcaldqm::hashfunctions::fDChannel);
0309 _xNBadTimingRefLM.initialize(hcaldqm::hashfunctions::fFED);
0310 _xNChs.initialize(hcaldqm::hashfunctions::fFED);
0311 _xMissingLaserMon = 0;
0312
0313 std::vector<int> vFEDs = hcaldqm::utilities::getFEDList(_emap);
0314 std::vector<int> vFEDsVME = hcaldqm::utilities::getFEDVMEList(_emap);
0315 std::vector<int> vFEDsuTCA = hcaldqm::utilities::getFEDuTCAList(_emap);
0316 for (std::vector<int>::const_iterator it = vFEDsVME.begin(); it != vFEDsVME.end(); ++it)
0317 _vhashFEDs.push_back(
0318 HcalElectronicsId(constants::FIBERCH_MIN, FIBER_VME_MIN, SPIGOT_MIN, (*it) - FED_VME_MIN).rawId());
0319 for (std::vector<int>::const_iterator it = vFEDsuTCA.begin(); it != vFEDsuTCA.end(); ++it) {
0320 std::pair<uint16_t, uint16_t> cspair = utilities::fed2crate(*it);
0321 _vhashFEDs.push_back(HcalElectronicsId(cspair.first, cspair.second, FIBER_uTCA_MIN1, FIBERCH_MIN, false).rawId());
0322 }
0323
0324 _cSummaryvsLS_FED.initialize(_name,
0325 "SummaryvsLS",
0326 hcaldqm::hashfunctions::fFED,
0327 new hcaldqm::quantity::LumiSection(_maxLS),
0328 new hcaldqm::quantity::FlagQuantity(_vflags),
0329 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fState),
0330 0);
0331 _cSummaryvsLS.initialize(_name,
0332 "SummaryvsLS",
0333 new hcaldqm::quantity::LumiSection(_maxLS),
0334 new hcaldqm::quantity::FEDQuantity(vFEDs),
0335 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fState),
0336 0);
0337 }
0338
0339
0340 _cSignalMean_Subdet.book(ib, _emap, _subsystem);
0341 _cSignalRMS_Subdet.book(ib, _emap, _subsystem);
0342 _cTimingMean_Subdet.book(ib, _emap, _subsystem);
0343 _cTimingRMS_Subdet.book(ib, _emap, _subsystem);
0344
0345 _cSignalMean_depth.book(ib, _emap, _subsystem);
0346 _cSignalRMS_depth.book(ib, _emap, _subsystem);
0347 _cTimingMean_depth.book(ib, _emap, _subsystem);
0348 _cTimingRMS_depth.book(ib, _emap, _subsystem);
0349
0350 if (_ptype == fLocal) {
0351 _cTimingvsEvent_SubdetPM.book(ib, _emap, _subsystem);
0352 _cSignalvsEvent_SubdetPM.book(ib, _emap, _subsystem);
0353 } else {
0354 _cTimingvsLS_SubdetPM.book(ib, _emap, _subsystem);
0355 _cSignalvsLS_SubdetPM.book(ib, _emap, _subsystem);
0356 _cTimingvsBX_SubdetPM.book(ib, _emap, _subsystem);
0357 _cSignalvsBX_SubdetPM.book(ib, _emap, _subsystem);
0358 }
0359
0360 if (_ptype == fLocal) {
0361 _cSignalMean_FEDuTCA.book(ib, _emap, _filter_VME, _subsystem);
0362 _cSignalRMS_FEDuTCA.book(ib, _emap, _filter_VME, _subsystem);
0363 _cTimingMean_FEDuTCA.book(ib, _emap, _filter_VME, _subsystem);
0364 _cTimingRMS_FEDuTCA.book(ib, _emap, _filter_VME, _subsystem);
0365 _cMissing_FEDuTCA.book(ib, _emap, _filter_VME, _subsystem);
0366 _cShapeCut_FEDSlot.book(ib, _emap, _subsystem);
0367 }
0368 _cADC_SubdetPM.book(ib, _emap, _subsystem);
0369
0370 _cMissing_depth.book(ib, _emap, _subsystem);
0371
0372 _xSignalSum.book(_emap);
0373 _xSignalSum2.book(_emap);
0374 _xEntries.book(_emap);
0375 _xTimingSum.book(_emap);
0376 _xTimingSum2.book(_emap);
0377
0378 if (_ptype == fOnline || _ptype == fLocal) {
0379 _cADCvsTS_SubdetPM.book(ib, _emap, _subsystem);
0380 _cLaserMonSumQ.book(ib, _subsystem);
0381 _cLaserMonTiming.book(ib, _subsystem);
0382 if (_ptype == fOnline) {
0383 _cLaserMonSumQ_LS.book(ib, _subsystem);
0384 _cLaserMonTiming_LS.book(ib, _subsystem);
0385 _cTimingDiffLS_SubdetPM.book(ib, _emap, _subsystem);
0386 } else if (_ptype == fLocal) {
0387 _cLaserMonSumQ_Event.book(ib, _subsystem);
0388 _cLaserMonTiming_Event.book(ib, _subsystem);
0389 _cTimingDiffEvent_SubdetPM.book(ib, _emap, _subsystem);
0390 }
0391 _cTiming_DigivsLaserMon_SubdetPM.book(ib, _emap, _subsystem);
0392 _xTimingRefLMSum.book(_emap);
0393 _xTimingRefLMSum2.book(_emap);
0394 _xNBadTimingRefLM.book(_emap);
0395 _xNChs.book(_emap);
0396 _cSummaryvsLS_FED.book(ib, _emap, _subsystem);
0397 _cSummaryvsLS.book(ib, _subsystem);
0398 }
0399
0400 _ehashmap.initialize(_emap, electronicsmap::fD2EHashMap);
0401 }
0402
0403 void LaserTask::_resetMonitors(hcaldqm::UpdateFreq uf) { DQTask::_resetMonitors(uf); }
0404
0405 void LaserTask::_dump() {
0406 _cSignalMean_Subdet.reset();
0407 _cSignalRMS_Subdet.reset();
0408 _cTimingMean_Subdet.reset();
0409 _cTimingRMS_Subdet.reset();
0410 _cSignalMean_depth.reset();
0411 _cSignalRMS_depth.reset();
0412 _cTimingMean_depth.reset();
0413 _cTimingRMS_depth.reset();
0414
0415 if (_ptype == fLocal) {
0416 _cSignalMean_FEDuTCA.reset();
0417 _cSignalRMS_FEDuTCA.reset();
0418 _cTimingMean_FEDuTCA.reset();
0419 _cTimingRMS_FEDuTCA.reset();
0420 }
0421 if (_ptype != fOffline) {
0422 _xNChs.reset();
0423 }
0424
0425 std::vector<HcalGenericDetId> dids = _emap->allPrecisionId();
0426 for (std::vector<HcalGenericDetId>::const_iterator it = dids.begin(); it != dids.end(); ++it) {
0427 if (!it->isHcalDetId())
0428 continue;
0429 HcalDetId did = HcalDetId(it->rawId());
0430 HcalElectronicsId eid(_ehashmap.lookup(*it));
0431 int n = _xEntries.get(did);
0432
0433 if (n == 0) {
0434 _cMissing_depth.fill(did);
0435 if (_ptype == fLocal) {
0436 if (!eid.isVMEid())
0437 _cMissing_FEDuTCA.fill(eid);
0438 }
0439 continue;
0440 }
0441
0442 _xNChs.get(eid)++;
0443
0444 double msig = _xSignalSum.get(did) / n;
0445 double mtim = _xTimingSum.get(did) / n;
0446 double rsig = sqrt(_xSignalSum2.get(did) / n - msig * msig);
0447 double rtim = sqrt(_xTimingSum2.get(did) / n - mtim * mtim);
0448
0449 _cSignalMean_Subdet.fill(did, msig);
0450 _cSignalMean_depth.fill(did, msig);
0451 _cSignalRMS_Subdet.fill(did, rsig);
0452 _cSignalRMS_depth.fill(did, rsig);
0453 _cTimingMean_Subdet.fill(did, mtim);
0454 _cTimingMean_depth.fill(did, mtim);
0455 _cTimingRMS_Subdet.fill(did, rtim);
0456 _cTimingRMS_depth.fill(did, rtim);
0457 if (_ptype == fLocal) {
0458 if (!eid.isVMEid()) {
0459 _cSignalMean_FEDuTCA.fill(eid, msig);
0460 _cSignalRMS_FEDuTCA.fill(eid, rsig);
0461 _cTimingMean_FEDuTCA.fill(eid, mtim);
0462 _cTimingRMS_FEDuTCA.fill(eid, rtim);
0463 }
0464 }
0465
0466
0467
0468 double timingreflm_mean = _xTimingRefLMSum.get(did) / n;
0469
0470
0471 if ((timingreflm_mean < _thresh_timingreflm[did.subdet()].first) ||
0472 (timingreflm_mean > _thresh_timingreflm[did.subdet()].second)) {
0473 _xNBadTimingRefLM.get(eid)++;
0474 }
0475 }
0476 if (_ptype != fOffline) {
0477 for (std::vector<uint32_t>::const_iterator it = _vhashFEDs.begin(); it != _vhashFEDs.end(); ++it) {
0478 hcaldqm::flag::Flag fSum("LASER");
0479 HcalElectronicsId eid = HcalElectronicsId(*it);
0480 std::vector<uint32_t>::const_iterator jt = std::find(_vcdaqEids.begin(), _vcdaqEids.end(), (*it));
0481 if (jt == _vcdaqEids.end()) {
0482
0483 for (uint32_t iflag = 0; iflag < _vflags.size(); iflag++)
0484 _cSummaryvsLS_FED.setBinContent(eid, _currentLS, int(iflag), int(hcaldqm::flag::fNCDAQ));
0485 _cSummaryvsLS.setBinContent(eid, _currentLS, int(hcaldqm::flag::fNCDAQ));
0486 continue;
0487 }
0488
0489 if (hcaldqm::utilities::isFEDHBHE(eid) || hcaldqm::utilities::isFEDHO(eid) || hcaldqm::utilities::isFEDHF(eid)) {
0490 int min_nchs = 10;
0491 if (_xNChs.get(eid) < min_nchs) {
0492 _vflags[fBadTiming]._state = hcaldqm::flag::fNA;
0493 } else {
0494 double frbadtimingreflm = double(_xNBadTimingRefLM.get(eid)) / double(_xNChs.get(eid));
0495 if (frbadtimingreflm > _thresh_frac_timingreflm) {
0496 _vflags[fBadTiming]._state = hcaldqm::flag::fBAD;
0497 } else {
0498 _vflags[fBadTiming]._state = hcaldqm::flag::fGOOD;
0499 }
0500 }
0501 if (_xMissingLaserMon) {
0502 _vflags[fMissingLaserMon]._state = hcaldqm::flag::fBAD;
0503 } else {
0504 _vflags[fMissingLaserMon]._state = hcaldqm::flag::fGOOD;
0505 }
0506 }
0507
0508
0509 int iflag = 0;
0510 for (std::vector<hcaldqm::flag::Flag>::iterator ft = _vflags.begin(); ft != _vflags.end(); ++ft) {
0511 _cSummaryvsLS_FED.setBinContent(eid, _currentLS, iflag, int(ft->_state));
0512 fSum += (*ft);
0513 iflag++;
0514 ft->reset();
0515 }
0516 _cSummaryvsLS.setBinContent(eid, _currentLS, fSum._state);
0517 }
0518
0519
0520 _xNBadTimingRefLM.reset();
0521 _xMissingLaserMon = 0;
0522 }
0523 }
0524
0525 void LaserTask::_process(edm::Event const& e, edm::EventSetup const& es) {
0526 edm::Handle<QIE11DigiCollection> c_QIE11;
0527 edm::Handle<HODigiCollection> c_ho;
0528 edm::Handle<QIE10DigiCollection> c_QIE10;
0529
0530 if (!e.getByToken(_tokQIE11, c_QIE11))
0531 _logger.dqmthrow("Collection QIE11DigiCollection isn't available " + _tagQIE11.label() + " " +
0532 _tagQIE11.instance());
0533 if (!e.getByToken(_tokHO, c_ho))
0534 _logger.dqmthrow("Collection HODigiCollection isn't available " + _tagHO.label() + " " + _tagHO.instance());
0535 if (!e.getByToken(_tokQIE10, c_QIE10))
0536 _logger.dqmthrow("Collection QIE10DigiCollection isn't available " + _tagQIE10.label() + " " +
0537 _tagQIE10.instance());
0538
0539
0540 int bx = e.bunchCrossing();
0541
0542
0543 edm::Handle<QIE10DigiCollection> cLaserMon;
0544 if (!e.getByToken(_tokLaserMon, cLaserMon)) {
0545 _logger.dqmthrow("QIE10DigiCollection for laserMonDigis isn't available.");
0546 }
0547
0548 QIE11DataFrame laserMonDigi;
0549 processLaserMon(c_QIE11, laserMonDigi);
0550
0551 double laserMonSumQ = hcaldqm::utilities::sumQ_v10<QIE11DataFrame>(laserMonDigi, 0, 0, laserMonDigi.samples() - 1);
0552 double laserMonTiming = 0.;
0553
0554 if (laserMonSumQ > 0.) {
0555 laserMonTiming = hcaldqm::utilities::aveTS_v10<QIE11DataFrame>(laserMonDigi, 0, 0, laserMonDigi.samples() - 1);
0556 } else {
0557 ++_xMissingLaserMon;
0558 }
0559
0560 if (laserMonSumQ > _laserMonThreshold) {
0561 _cLaserMonSumQ.fill(laserMonSumQ);
0562 _cLaserMonTiming.fill(laserMonTiming);
0563 if (_ptype == fOnline) {
0564 _cLaserMonSumQ_LS.fill(_currentLS, laserMonSumQ);
0565 _cLaserMonTiming_LS.fill(_currentLS, laserMonTiming);
0566 } else if (_ptype == fLocal) {
0567 int currentEvent = e.eventAuxiliary().id().event();
0568 _cLaserMonSumQ_Event.fill(currentEvent, laserMonSumQ);
0569 _cLaserMonTiming_Event.fill(currentEvent, laserMonTiming);
0570 }
0571 }
0572
0573 for (QIE11DigiCollection::const_iterator it = c_QIE11->begin(); it != c_QIE11->end(); ++it) {
0574 const QIE11DataFrame digi = static_cast<const QIE11DataFrame>(*it);
0575 HcalDetId const& did = digi.detid();
0576 if ((did.subdet() != HcalBarrel) && (did.subdet() != HcalEndcap)) {
0577 continue;
0578 }
0579 uint32_t rawid = _ehashmap.lookup(did);
0580 HcalElectronicsId const& eid(rawid);
0581
0582 CaloSamples digi_fC = hcaldqm::utilities::loadADC2fCDB<QIE11DataFrame>(_dbService, did, digi);
0583
0584 double sumQ = hcaldqm::utilities::sumQDB<QIE11DataFrame>(_dbService, digi_fC, did, digi, 0, digi.samples() - 1);
0585 if (sumQ < _lowHBHE)
0586 continue;
0587
0588
0589 double aveTS = hcaldqm::utilities::aveTSDB<QIE11DataFrame>(_dbService, digi_fC, did, digi, 0, digi.samples() - 1);
0590 _xSignalSum.get(did) += sumQ;
0591 _xSignalSum2.get(did) += sumQ * sumQ;
0592 _xTimingSum.get(did) += aveTS;
0593 _xTimingSum2.get(did) += aveTS * aveTS;
0594 _xEntries.get(did)++;
0595
0596 for (int i = 0; i < digi.samples(); i++) {
0597 if (_ptype == fLocal) {
0598 _cShapeCut_FEDSlot.fill(
0599 eid, i, hcaldqm::utilities::adc2fCDBMinusPedestal<QIE11DataFrame>(_dbService, digi_fC, did, digi, i));
0600 }
0601 _cADC_SubdetPM.fill(did, digi[i].adc());
0602 }
0603
0604
0605 double digiTimingSOI = (aveTS - digi.presamples()) * 25.;
0606 double deltaTiming = digiTimingSOI - laserMonTiming;
0607 _cTiming_DigivsLaserMon_SubdetPM.fill(did, laserMonTiming, digiTimingSOI);
0608 _xTimingRefLMSum.get(did) += deltaTiming;
0609 _xTimingRefLMSum2.get(did) += deltaTiming * deltaTiming;
0610 if (_ptype == fLocal) {
0611 int currentEvent = e.eventAuxiliary().id().event();
0612 _cTimingvsEvent_SubdetPM.fill(did, currentEvent, aveTS);
0613 _cSignalvsEvent_SubdetPM.fill(did, currentEvent, sumQ);
0614 _cTimingDiffEvent_SubdetPM.fill(did, currentEvent, hcaldqm::utilities::getRBX(did.iphi()), deltaTiming);
0615 } else {
0616 _cTimingvsLS_SubdetPM.fill(did, _currentLS, aveTS);
0617 _cSignalvsLS_SubdetPM.fill(did, _currentLS, sumQ);
0618 _cTimingvsBX_SubdetPM.fill(did, bx, aveTS);
0619 _cSignalvsBX_SubdetPM.fill(did, bx, sumQ);
0620 _cTimingDiffLS_SubdetPM.fill(did, _currentLS, hcaldqm::utilities::getRBX(did.iphi()), deltaTiming);
0621 }
0622
0623 if (_ptype == fOnline || _ptype == fLocal) {
0624 for (int iTS = 0; iTS < digi.samples(); ++iTS) {
0625 _cADCvsTS_SubdetPM.fill(did, iTS, digi[iTS].adc());
0626 }
0627 }
0628 }
0629 for (HODigiCollection::const_iterator it = c_ho->begin(); it != c_ho->end(); ++it) {
0630 const HODataFrame digi = (const HODataFrame)(*it);
0631 double sumQ = hcaldqm::utilities::sumQ<HODataFrame>(digi, 8.5, 0, digi.size() - 1);
0632 if (sumQ < _lowHO)
0633 continue;
0634 HcalDetId did = digi.id();
0635 HcalElectronicsId eid = digi.elecId();
0636
0637 double aveTS = hcaldqm::utilities::aveTS<HODataFrame>(digi, 8.5, 0, digi.size() - 1);
0638 _xSignalSum.get(did) += sumQ;
0639 _xSignalSum2.get(did) += sumQ * sumQ;
0640 _xTimingSum.get(did) += aveTS;
0641 _xTimingSum2.get(did) += aveTS * aveTS;
0642 _xEntries.get(did)++;
0643
0644 for (int i = 0; i < digi.size(); i++) {
0645 if (_ptype == fLocal) {
0646 _cShapeCut_FEDSlot.fill(eid, i, digi.sample(i).nominal_fC() - 8.5);
0647 }
0648 _cADC_SubdetPM.fill(did, digi.sample(i).adc());
0649 }
0650
0651
0652 double digiTimingSOI = (aveTS - digi.presamples()) * 25.;
0653 double deltaTiming = digiTimingSOI - laserMonTiming;
0654 _cTiming_DigivsLaserMon_SubdetPM.fill(did, laserMonTiming, digiTimingSOI);
0655 _xTimingRefLMSum.get(did) += deltaTiming;
0656 _xTimingRefLMSum2.get(did) += deltaTiming * deltaTiming;
0657 if (_ptype == fLocal) {
0658 int currentEvent = e.eventAuxiliary().id().event();
0659 _cTimingvsEvent_SubdetPM.fill(did, currentEvent, aveTS);
0660 _cSignalvsEvent_SubdetPM.fill(did, currentEvent, sumQ);
0661 _cTimingDiffEvent_SubdetPM.fill(did, currentEvent, hcaldqm::utilities::getRBX(did.iphi()), deltaTiming);
0662 } else {
0663 _cTimingvsLS_SubdetPM.fill(did, _currentLS, aveTS);
0664 _cSignalvsLS_SubdetPM.fill(did, _currentLS, sumQ);
0665 _cTimingvsBX_SubdetPM.fill(did, bx, aveTS);
0666 _cSignalvsBX_SubdetPM.fill(did, bx, sumQ);
0667 _cTimingDiffLS_SubdetPM.fill(did, _currentLS, hcaldqm::utilities::getRBX(did.iphi()), deltaTiming);
0668 }
0669
0670 if (_ptype == fOnline || _ptype == fLocal) {
0671 for (int iTS = 0; iTS < digi.size(); ++iTS) {
0672 _cADCvsTS_SubdetPM.fill(did, iTS, digi.sample(iTS).adc());
0673 }
0674 }
0675 }
0676 for (QIE10DigiCollection::const_iterator it = c_QIE10->begin(); it != c_QIE10->end(); ++it) {
0677 const QIE10DataFrame digi = (const QIE10DataFrame)(*it);
0678 HcalDetId did = digi.detid();
0679 if (did.subdet() != HcalForward) {
0680 continue;
0681 }
0682 HcalElectronicsId eid = HcalElectronicsId(_ehashmap.lookup(did));
0683
0684 CaloSamples digi_fC = hcaldqm::utilities::loadADC2fCDB<QIE10DataFrame>(_dbService, did, digi);
0685 double sumQ = hcaldqm::utilities::sumQDB<QIE10DataFrame>(_dbService, digi_fC, did, digi, 0, digi.samples() - 1);
0686
0687 if (sumQ < _lowHF)
0688 continue;
0689
0690
0691 double aveTS = hcaldqm::utilities::aveTSDB<QIE10DataFrame>(_dbService, digi_fC, did, digi, 0, digi.samples() - 1);
0692
0693 _xSignalSum.get(did) += sumQ;
0694 _xSignalSum2.get(did) += sumQ * sumQ;
0695 _xTimingSum.get(did) += aveTS;
0696 _xTimingSum2.get(did) += aveTS * aveTS;
0697 _xEntries.get(did)++;
0698
0699 for (int i = 0; i < digi.samples(); i++) {
0700 if (_ptype == fLocal) {
0701 _cShapeCut_FEDSlot.fill(
0702 eid, (int)i, hcaldqm::utilities::adc2fCDBMinusPedestal<QIE10DataFrame>(_dbService, digi_fC, did, digi, i));
0703 }
0704 _cADC_SubdetPM.fill(did, digi[i].adc());
0705 }
0706
0707
0708 double digiTimingSOI = (aveTS - digi.presamples()) * 25.;
0709 double deltaTiming = digiTimingSOI - laserMonTiming;
0710 _cTiming_DigivsLaserMon_SubdetPM.fill(did, laserMonTiming, digiTimingSOI);
0711 _xTimingRefLMSum.get(did) += deltaTiming;
0712 _xTimingRefLMSum2.get(did) += deltaTiming * deltaTiming;
0713 if (_ptype == fLocal) {
0714 int currentEvent = e.eventAuxiliary().id().event();
0715 _cTimingvsEvent_SubdetPM.fill(did, currentEvent, aveTS);
0716 _cSignalvsEvent_SubdetPM.fill(did, currentEvent, sumQ);
0717 _cTimingDiffEvent_SubdetPM.fill(did, currentEvent, hcaldqm::utilities::getRBX(did.iphi()), deltaTiming);
0718 } else {
0719 _cTimingvsLS_SubdetPM.fill(did, _currentLS, aveTS);
0720 _cSignalvsLS_SubdetPM.fill(did, _currentLS, sumQ);
0721 _cTimingvsBX_SubdetPM.fill(did, bx, aveTS);
0722 _cSignalvsBX_SubdetPM.fill(did, bx, sumQ);
0723 _cTimingDiffLS_SubdetPM.fill(did, _currentLS, hcaldqm::utilities::getRBX(did.iphi()), deltaTiming);
0724 }
0725 if (_ptype == fOnline || _ptype == fLocal) {
0726 for (int iTS = 0; iTS < digi.samples(); ++iTS) {
0727 _cADCvsTS_SubdetPM.fill(did, iTS, digi[iTS].adc());
0728 }
0729 }
0730 }
0731 }
0732
0733 void LaserTask::processLaserMon(edm::Handle<QIE11DigiCollection>& col, QIE11DataFrame& iLaserMonDigi) {
0734 for (QIE11DigiCollection::const_iterator it = col->begin(); it != col->end(); ++it) {
0735 const QIE11DataFrame digi = static_cast<const QIE11DataFrame>(*it);
0736 HcalCalibDetId hcdid(digi.id());
0737
0738 if (hcdid.rawId() != constants::HBLasMon.rawId())
0739 continue;
0740
0741 iLaserMonDigi = digi;
0742 }
0743 }
0744
0745 void LaserTask::globalEndLuminosityBlock(edm::LuminosityBlock const& lb, edm::EventSetup const& es) {
0746 auto lumiCache = luminosityBlockCache(lb.index());
0747 _currentLS = lumiCache->currentLS;
0748
0749 if (_ptype == fLocal)
0750 return;
0751 this->_dump();
0752
0753 DQTask::globalEndLuminosityBlock(lb, es);
0754 }
0755
0756 bool LaserTask::_isApplicable(edm::Event const& e) {
0757 if (_ptype != fOnline)
0758 return true;
0759 else {
0760
0761 edm::Handle<HcalUMNioDigi> cumn;
0762 if (!e.getByToken(_tokuMN, cumn))
0763 return false;
0764
0765
0766
0767 uint32_t laserType = cumn->valueUserWord(0);
0768 if (laserType != _laserType)
0769 return false;
0770
0771
0772 int eventflag_uHTR = -1;
0773 edm::Handle<FEDRawDataCollection> craw;
0774 if (!e.getByToken(_tokFEDs, craw))
0775 _logger.dqmthrow("Collection FEDRawDataCollection isn't available " + _tagFEDs.label() + " " +
0776 _tagFEDs.instance());
0777
0778 for (int fed = FEDNumbering::MINHCALFEDID; fed <= FEDNumbering::MAXHCALuTCAFEDID && !(eventflag_uHTR >= 0); fed++) {
0779 if ((fed > FEDNumbering::MAXHCALFEDID && fed < FEDNumbering::MINHCALuTCAFEDID) ||
0780 fed > FEDNumbering::MAXHCALuTCAFEDID)
0781 continue;
0782 FEDRawData const& raw = craw->FEDData(fed);
0783 if (raw.size() < constants::RAW_EMPTY)
0784 continue;
0785
0786 hcal::AMC13Header const* hamc13 = (hcal::AMC13Header const*)raw.data();
0787 if (!hamc13)
0788 continue;
0789
0790 for (int iamc = 0; iamc < hamc13->NAMC(); iamc++) {
0791 HcalUHTRData uhtr(hamc13->AMCPayload(iamc), hamc13->AMCSize(iamc));
0792 if (static_cast<int>(uhtr.crateId()) == 25 && static_cast<int>(uhtr.slot()) == 11) {
0793 eventflag_uHTR = uhtr.getEventType();
0794 break;
0795 }
0796 }
0797 }
0798 if (eventflag_uHTR == constants::EVENTTYPE_LASER) {
0799 return true;
0800 }
0801 }
0802
0803 return false;
0804 }
0805
0806 DEFINE_FWK_MODULE(LaserTask);