Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-03-17 23:26:21

0001 #include "DQM/EcalMonitorClient/interface/LaserClient.h"
0002 
0003 #include "DataFormats/EcalDetId/interface/EcalPnDiodeDetId.h"
0004 
0005 #include "CondFormats/EcalObjects/interface/EcalDQMStatusHelper.h"
0006 
0007 #include "DQM/EcalCommon/interface/EcalDQMCommonUtils.h"
0008 #include "DQM/EcalCommon/interface/MESetMulti.h"
0009 
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 
0012 #include <cmath>
0013 
0014 namespace ecaldqm {
0015   LaserClient::LaserClient()
0016       : DQWorkerClient(),
0017         wlToME_(),
0018         minChannelEntries_(0),
0019         expectedAmplitude_(0),
0020         toleranceAmplitudeLo_(0.),
0021         toleranceAmplitudeHi_(0.),
0022         toleranceAmpRMSRatio_(0.),
0023         expectedTiming_(0),
0024         toleranceTiming_(0.),
0025         toleranceTimRMS_(0.),
0026         expectedPNAmplitude_(0),
0027         tolerancePNAmp_(0.),
0028         tolerancePNRMSRatio_(0.),
0029         forwardFactor_(0.) {}
0030 
0031   void LaserClient::setParams(edm::ParameterSet const& _params) {
0032     minChannelEntries_ = _params.getUntrackedParameter<int>("minChannelEntries");
0033     toleranceAmplitudeLo_ = _params.getUntrackedParameter<double>("toleranceAmplitudeLo");
0034     toleranceAmplitudeHi_ = _params.getUntrackedParameter<double>("toleranceAmplitudeHi");
0035     toleranceAmpRMSRatio_ = _params.getUntrackedParameter<double>("toleranceAmpRMSRatio");
0036     toleranceTiming_ = _params.getUntrackedParameter<double>("toleranceTiming");
0037     toleranceTimRMS_ = _params.getUntrackedParameter<double>("toleranceTimRMS");
0038     tolerancePNAmp_ = _params.getUntrackedParameter<double>("tolerancePNAmp");
0039     tolerancePNRMSRatio_ = _params.getUntrackedParameter<double>("tolerancePNRMSRatio");
0040     forwardFactor_ = _params.getUntrackedParameter<double>("forwardFactor");
0041 
0042     std::vector<int> laserWavelengths(_params.getUntrackedParameter<std::vector<int> >("laserWavelengths"));
0043 
0044     // wavelengths are not necessarily ordered
0045     // create a map wl -> MESet index
0046     // using Amplitude here but any multi-wavelength plot is fine
0047 
0048     MESet::PathReplacements repl;
0049 
0050     MESetMulti const& amplitude(static_cast<MESetMulti const&>(sources_.at("Amplitude")));
0051     unsigned nWL(laserWavelengths.size());
0052     for (unsigned iWL(0); iWL != nWL; ++iWL) {
0053       int wl(laserWavelengths[iWL]);
0054       if (wl <= 0 || wl >= 5)
0055         throw cms::Exception("InvalidConfiguration") << "Laser Wavelength";
0056       repl["wl"] = std::to_string(wl);
0057       wlToME_[wl] = amplitude.getIndex(repl);
0058     }
0059 
0060     expectedAmplitude_.resize(nWL);
0061     expectedTiming_.resize(nWL);
0062     expectedPNAmplitude_.resize(nWL);
0063 
0064     std::vector<double> inExpectedAmplitude(_params.getUntrackedParameter<std::vector<double> >("expectedAmplitude"));
0065     std::vector<double> inExpectedTiming(_params.getUntrackedParameter<std::vector<double> >("expectedTiming"));
0066     std::vector<double> inExpectedPNAmplitude(
0067         _params.getUntrackedParameter<std::vector<double> >("expectedPNAmplitude"));
0068 
0069     for (std::map<int, unsigned>::iterator wlItr(wlToME_.begin()); wlItr != wlToME_.end(); ++wlItr) {
0070       unsigned iME(wlItr->second);
0071       int iWL(wlItr->first - 1);
0072       expectedAmplitude_[iME] = inExpectedAmplitude[iWL];
0073       expectedTiming_[iME] = inExpectedTiming[iWL];
0074       expectedPNAmplitude_[iME] = inExpectedPNAmplitude[iWL];
0075     }
0076 
0077     qualitySummaries_.insert("Quality");
0078     qualitySummaries_.insert("QualitySummary");
0079     qualitySummaries_.insert("PNQualitySummary");
0080   }
0081 
0082   void LaserClient::producePlots(ProcessType) {
0083     uint32_t mask(1 << EcalDQMStatusHelper::LASER_MEAN_ERROR | 1 << EcalDQMStatusHelper::LASER_RMS_ERROR |
0084                   1 << EcalDQMStatusHelper::LASER_TIMING_MEAN_ERROR | 1 << EcalDQMStatusHelper::LASER_TIMING_RMS_ERROR);
0085 
0086     MESetMulti& meQuality(static_cast<MESetMulti&>(MEs_.at("Quality")));
0087     MESetMulti& meQualitySummary(static_cast<MESetMulti&>(MEs_.at("QualitySummary")));
0088     MESetMulti& meAmplitudeMean(static_cast<MESetMulti&>(MEs_.at("AmplitudeMean")));
0089     MESetMulti& meAmplitudeRMS(static_cast<MESetMulti&>(MEs_.at("AmplitudeRMS")));
0090     MESetMulti& meTimingMean(static_cast<MESetMulti&>(MEs_.at("TimingMean")));
0091     MESetMulti& meTimingRMSMap(static_cast<MESetMulti&>(MEs_.at("TimingRMSMap")));
0092     MESetMulti& meTimingRMS(static_cast<MESetMulti&>(MEs_.at("TimingRMS")));
0093     MESetMulti& mePNQualitySummary(static_cast<MESetMulti&>(MEs_.at("PNQualitySummary")));
0094 
0095     MESetMulti const& sAmplitude(static_cast<MESetMulti const&>(sources_.at("Amplitude")));
0096     MESetMulti const& sTiming(static_cast<MESetMulti const&>(sources_.at("Timing")));
0097     MESetMulti const& sPNAmplitude(static_cast<MESetMulti const&>(sources_.at("PNAmplitude")));
0098     MESet const& sCalibStatus(static_cast<MESet const&>(sources_.at("CalibStatus")));
0099 
0100     for (std::map<int, unsigned>::iterator wlItr(wlToME_.begin()); wlItr != wlToME_.end(); ++wlItr) {
0101       meQuality.use(wlItr->second);
0102       meQualitySummary.use(wlItr->second);
0103       meAmplitudeMean.use(wlItr->second);
0104       meAmplitudeRMS.use(wlItr->second);
0105       meTimingMean.use(wlItr->second);
0106       meTimingRMSMap.use(wlItr->second);
0107       meTimingRMS.use(wlItr->second);
0108       mePNQualitySummary.use(wlItr->second);
0109 
0110       sAmplitude.use(wlItr->second);
0111       sTiming.use(wlItr->second);
0112       sPNAmplitude.use(wlItr->second);
0113 
0114       MESet::iterator qEnd(meQuality.end(GetElectronicsMap()));
0115 
0116       MESet::const_iterator tItr(GetElectronicsMap(), sTiming);
0117       MESet::const_iterator aItr(GetElectronicsMap(), sAmplitude);
0118 
0119       int wl(wlItr->first - 1);
0120       bool enabled(wl < 0 ? false : sCalibStatus.getBinContent(getEcalDQMSetupObjects(), wl) > 0 ? true : false);
0121       for (MESet::iterator qItr(meQuality.beginChannel(GetElectronicsMap())); qItr != qEnd;
0122            qItr.toNextChannel(GetElectronicsMap())) {
0123         DetId id(qItr->getId());
0124 
0125         bool doMask(meQuality.maskMatches(id, mask, statusManager_, GetTrigTowerMap()));
0126 
0127         aItr = qItr;
0128 
0129         float aEntries(aItr->getBinEntries());
0130 
0131         if (aEntries < minChannelEntries_) {
0132           qItr->setBinContent(enabled ? (doMask ? kMUnknown : kUnknown) : kMUnknown);
0133           continue;
0134         }
0135 
0136         float aMean(aItr->getBinContent());
0137         float aRms(aItr->getBinError() * sqrt(aEntries));
0138 
0139         meAmplitudeMean.fill(getEcalDQMSetupObjects(), id, aMean);
0140         meAmplitudeRMS.setBinContent(getEcalDQMSetupObjects(), id, aRms);
0141 
0142         tItr = qItr;
0143 
0144         float tEntries(tItr->getBinEntries());
0145 
0146         if (tEntries < minChannelEntries_)
0147           continue;
0148 
0149         float tMean(tItr->getBinContent());
0150         float tRms(tItr->getBinError() * sqrt(tEntries));
0151 
0152         meTimingMean.fill(getEcalDQMSetupObjects(), id, tMean);
0153         meTimingRMS.fill(getEcalDQMSetupObjects(), id, tRms);
0154         meTimingRMSMap.setBinContent(getEcalDQMSetupObjects(), id, tRms);
0155 
0156         float intensity(aMean / expectedAmplitude_[wlItr->second]);
0157         if (isForward(id))
0158           intensity /= forwardFactor_;
0159 
0160         if (intensity < toleranceAmplitudeLo_ || intensity > toleranceAmplitudeHi_ ||
0161             aRms > aMean * toleranceAmpRMSRatio_ ||
0162             std::abs(tMean - expectedTiming_[wlItr->second]) > toleranceTiming_ /*|| tRms > toleranceTimRMS_*/)
0163           qItr->setBinContent(doMask ? kMBad : kBad);
0164         else
0165           qItr->setBinContent(doMask ? kMGood : kGood);
0166       }
0167 
0168       towerAverage_(meQualitySummary, meQuality, 0.2);
0169 
0170       for (unsigned iDCC(0); iDCC < nDCC; ++iDCC) {
0171         if (memDCCIndex(iDCC + 1) == unsigned(-1))
0172           continue;
0173         int subdet(0);
0174         if (iDCC >= kEBmLow && iDCC <= kEBpHigh)
0175           subdet = EcalBarrel;
0176         else
0177           subdet = EcalEndcap;
0178 
0179         for (unsigned iPN(0); iPN < 10; ++iPN) {
0180           EcalPnDiodeDetId id(subdet, iDCC + 1, iPN + 1);
0181 
0182           bool doMask(mePNQualitySummary.maskMatches(id, mask, statusManager_, GetTrigTowerMap()));
0183 
0184           float pEntries(sPNAmplitude.getBinEntries(getEcalDQMSetupObjects(), id));
0185 
0186           if (pEntries < minChannelEntries_) {
0187             mePNQualitySummary.setBinContent(getEcalDQMSetupObjects(), id, doMask ? kMUnknown : kUnknown);
0188             continue;
0189           }
0190 
0191           float pMean(sPNAmplitude.getBinContent(getEcalDQMSetupObjects(), id));
0192           float pRms(sPNAmplitude.getBinError(getEcalDQMSetupObjects(), id) * sqrt(pEntries));
0193           float intensity(pMean / expectedPNAmplitude_[wlItr->second]);
0194 
0195           if (intensity < tolerancePNAmp_ || pRms > pMean * tolerancePNRMSRatio_)
0196             mePNQualitySummary.setBinContent(getEcalDQMSetupObjects(), id, doMask ? kMBad : kBad);
0197           else
0198             mePNQualitySummary.setBinContent(getEcalDQMSetupObjects(), id, doMask ? kMGood : kGood);
0199         }
0200       }
0201     }
0202   }
0203 
0204   DEFINE_ECALDQM_WORKER(LaserClient);
0205 }  // namespace ecaldqm