EcalMonitorClient/src/TimingClient.cc

0001 #include "DQM/EcalMonitorClient/interface/TimingClient.h"
0002
0003 #include "DQM/EcalCommon/interface/EcalDQMCommonUtils.h"
0004
0005 #include "CondFormats/EcalObjects/interface/EcalDQMStatusHelper.h"
0006 #include "CondFormats/EcalObjects/interface/EcalChannelStatusCode.h"
0007
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009
0010 #include <cmath>
0011
0012 namespace ecaldqm {
0013   TimingClient::TimingClient()
0014       : DQWorkerClient(),
0015         ebtoleranceMean_(0.),
0016         eetoleranceMean_(0.),
0017         toleranceMeanFwd_(0.),
0018         toleranceRMS_(0.),
0019         toleranceRMSFwd_(0.),
0020         minChannelEntries_(0),
0021         minChannelEntriesFwd_(0),
0022         minTowerEntries_(0),
0023         minTowerEntriesFwd_(0),
0024         tailPopulThreshold_(0.) {
0025     qualitySummaries_.insert("Quality");
0026     qualitySummaries_.insert("QualitySummary");
0027   }
0028
0029   void TimingClient::setParams(edm::ParameterSet const& _params) {
0030     ebtoleranceMean_ = _params.getUntrackedParameter<double>("ebtoleranceMean");
0031     eetoleranceMean_ = _params.getUntrackedParameter<double>("eetoleranceMean");
0032     toleranceMeanFwd_ = _params.getUntrackedParameter<double>("toleranceMeanFwd");
0033     toleranceRMS_ = _params.getUntrackedParameter<double>("toleranceRMS");
0034     toleranceRMSFwd_ = _params.getUntrackedParameter<double>("toleranceRMSFwd");
0035     minChannelEntries_ = _params.getUntrackedParameter<int>("minChannelEntries");
0036     minChannelEntriesFwd_ = _params.getUntrackedParameter<int>("minChannelEntriesFwd");
0037     minTowerEntries_ = _params.getUntrackedParameter<int>("minTowerEntries");
0038     minTowerEntriesFwd_ = _params.getUntrackedParameter<int>("minChannelEntriesFwd");
0039     tailPopulThreshold_ = _params.getUntrackedParameter<double>("tailPopulThreshold");
0040   }
0041
0042   void TimingClient::producePlots(ProcessType) {
0043     MESet& meQuality(MEs_.at("Quality"));
0044     MESet& meMeanSM(MEs_.at("MeanSM"));
0045     MESet& meMeanAll(MEs_.at("MeanAll"));
0046     MESet& meFwdBkwdDiff(MEs_.at("FwdBkwdDiff"));
0047     MESet& meFwdvBkwd(MEs_.at("FwdvBkwd"));
0048     MESet& meRMSMap(MEs_.at("RMSMap"));
0049     MESet& meRMSAll(MEs_.at("RMSAll"));
0050     MESet& meProjEta(MEs_.at("ProjEta"));
0051     MESet& meProjPhi(MEs_.at("ProjPhi"));
0052     MESet& meQualitySummary(MEs_.at("QualitySummary"));
0053
0054     MESet const& sTimeAllMap(sources_.at("TimeAllMap"));
0055     MESet const& sTimeMap(sources_.at("TimeMap"));
0056     MESet const& sTimeMapByLS(sources_.at("TimeMapByLS"));
0057     MESet const& sChStatus(sources_.at("ChStatus"));
0058
0059     uint32_t mask(1 << EcalDQMStatusHelper::PHYSICS_BAD_CHANNEL_WARNING);
0060
0061     MESet::iterator qEnd(meQuality.end(GetElectronicsMap()));
0062
0063     MESet::iterator rItr(GetElectronicsMap(), meRMSMap);
0064     MESet::const_iterator tItr(GetElectronicsMap(), sTimeMap);
0065     MESet::const_iterator tLSItr(GetElectronicsMap(), sTimeMapByLS);
0066
0067     float EBentries(0.), EEentries(0.);
0068     float EBmean(0.), EEmean(0.);
0069     float EBrms(0.), EErms(0.);
0070     for (MESet::iterator qItr(meQuality.beginChannel(GetElectronicsMap())); qItr != qEnd;
0071          qItr.toNextChannel(GetElectronicsMap())) {
0072       tItr = qItr;
0073       rItr = qItr;
0074
0075       DetId id(qItr->getId());
0076
0077       int minChannelEntries(minChannelEntries_);
0078       float meanThresh;
0079       float rmsThresh(toleranceRMS_);
0080
0081       if (id.subdetId() == EcalBarrel)
0082         meanThresh = ebtoleranceMean_;
0083       else
0084         meanThresh = eetoleranceMean_;
0085
0086       if (isForward(id)) {
0087         minChannelEntries = minChannelEntriesFwd_;
0088         meanThresh = toleranceMeanFwd_;
0089         rmsThresh = toleranceRMSFwd_;
0090       }
0091
0092       bool doMask(meQuality.maskMatches(id, mask, statusManager_, GetTrigTowerMap()));
0093
0094       float entries(tItr->getBinEntries());
0095
0096       if (entries < minChannelEntries) {
0097         qItr->setBinContent(doMask ? kMUnknown : kUnknown);
0098         rItr->setBinContent(-1.);
0099         continue;
0100       }
0101
0102       float mean(tItr->getBinContent());
0103       float rms(tItr->getBinError() * sqrt(entries));
0104
0105       meMeanSM.fill(getEcalDQMSetupObjects(), id, mean);
0106       meMeanAll.fill(getEcalDQMSetupObjects(), id, mean);
0107       meProjEta.fill(getEcalDQMSetupObjects(), id, mean);
0108       meProjPhi.fill(getEcalDQMSetupObjects(), id, mean);
0109       meRMSAll.fill(getEcalDQMSetupObjects(), id, rms);
0110       rItr->setBinContent(rms);
0111
0112       bool negative(false);
0113       float posTime(0.);
0114
0115       if (id.subdetId() == EcalBarrel) {
0116         EBDetId ebid(id);
0117         if (ebid.zside() < 0) {
0118           negative = true;
0119           EBDetId posId(EBDetId::switchZSide(ebid));
0120           posTime = sTimeMap.getBinContent(getEcalDQMSetupObjects(), posId);
0121         }
0122       } else {
0123         EEDetId eeid(id);
0124         if (eeid.zside() < 0) {
0125           negative = true;
0126           EEDetId posId(EEDetId::switchZSide(eeid));
0127           posTime = sTimeMap.getBinContent(getEcalDQMSetupObjects(), posId);
0128         }
0129       }
0130       if (negative) {
0131         meFwdBkwdDiff.fill(getEcalDQMSetupObjects(), id, posTime - mean);
0132         meFwdvBkwd.fill(getEcalDQMSetupObjects(), id, mean, posTime);
0133       }
0134
0135       if (std::abs(mean) > meanThresh || rms > rmsThresh)
0136         qItr->setBinContent(doMask ? kMBad : kBad);
0137       else
0138         qItr->setBinContent(doMask ? kMGood : kGood);
0139
0140       // For Trend plots:
0141       tLSItr = qItr;
0142       float entriesLS(tLSItr->getBinEntries());
0143       float meanLS(tLSItr->getBinContent());
0144       float rmsLS(tLSItr->getBinError() * sqrt(entriesLS));
0145       int chStatus = static_cast<int>(sChStatus.getBinContent(getEcalDQMSetupObjects(), id));
0146
0147       if (entriesLS < minChannelEntries)
0148         continue;
0149       if (chStatus != EcalChannelStatusCode::kOk)
0150         continue;  // exclude problematic channels
0151
0152       // Keep running count of timing mean, rms, and N_hits
0153       if (id.subdetId() == EcalBarrel) {
0154         EBmean += meanLS;
0155         EBrms += rmsLS;
0156         EBentries += entriesLS;
0157       } else {
0158         EEmean += meanLS;
0159         EErms += rmsLS;
0160         EEentries += entriesLS;
0161       }
0162
0163     }  // channel loop
0164
0165     // Fill Timing Trend plots at each LS
0166     MESet& meTrendMean(MEs_.at("TrendMean"));
0167     MESet& meTrendRMS(MEs_.at("TrendRMS"));
0168     if (EBentries > 0.) {
0169       if (std::abs(EBmean) > 0.)
0170         meTrendMean.fill(getEcalDQMSetupObjects(), EcalBarrel, double(timestamp_.iLumi), EBmean / EBentries);
0171       if (std::abs(EBrms) > 0.)
0172         meTrendRMS.fill(getEcalDQMSetupObjects(), EcalBarrel, double(timestamp_.iLumi), EBrms / EBentries);
0173     }
0174     if (EEentries > 0.) {
0175       if (std::abs(EEmean) > 0.)
0176         meTrendMean.fill(getEcalDQMSetupObjects(), EcalEndcap, double(timestamp_.iLumi), EEmean / EEentries);
0177       if (std::abs(EErms) > 0.)
0178         meTrendRMS.fill(getEcalDQMSetupObjects(), EcalEndcap, double(timestamp_.iLumi), EErms / EEentries);
0179     }
0180
0181     MESet::iterator qsEnd(meQualitySummary.end(GetElectronicsMap()));
0182
0183     for (MESet::iterator qsItr(meQualitySummary.beginChannel(GetElectronicsMap())); qsItr != qsEnd;
0184          qsItr.toNextChannel(GetElectronicsMap())) {
0185       DetId tId(qsItr->getId());
0186
0187       std::vector<DetId> ids;
0188
0189       if (tId.subdetId() == EcalTriggerTower)
0190         ids = GetTrigTowerMap()->constituentsOf(EcalTrigTowerDetId(tId));
0191       else
0192         ids = scConstituents(EcalScDetId(tId));
0193
0194       int minTowerEntries(minTowerEntries_);
0195       float meanThresh;
0196       float rmsThresh(toleranceRMS_);
0197
0198       if (tId.subdetId() == EcalBarrel)
0199         meanThresh = ebtoleranceMean_;
0200       else
0201         meanThresh = eetoleranceMean_;
0202
0203       if (isForward(tId)) {
0204         minTowerEntries = minTowerEntriesFwd_;
0205         meanThresh = toleranceMeanFwd_;
0206         rmsThresh = toleranceRMSFwd_;
0207       }
0208
0209       // tower entries != sum(channel entries) because of the difference in timing cut at the source
0210       float summaryEntries(sTimeAllMap.getBinEntries(getEcalDQMSetupObjects(), tId));
0211
0212       float towerEntries(0.);
0213       float towerMean(0.);
0214       float towerMean2(0.);
0215
0216       bool doMask(false);
0217
0218       for (std::vector<DetId>::iterator idItr(ids.begin()); idItr != ids.end(); ++idItr) {
0219         DetId& id(*idItr);
0220
0221         doMask |= meQuality.maskMatches(id, mask, statusManager_, GetTrigTowerMap());
0222
0223         MESet::const_iterator tmItr(GetElectronicsMap(), sTimeMap, id);
0224
0225         float entries(tmItr->getBinEntries());
0226         if (entries < 0.)
0227           continue;
0228         towerEntries += entries;
0229         float mean(tmItr->getBinContent());
0230         towerMean += mean * entries;
0231         float rms(tmItr->getBinError() * sqrt(entries));
0232         towerMean2 += (rms * rms + mean * mean) * entries;
0233       }
0234
0235       double quality(doMask ? kMUnknown : kUnknown);
0236       if (towerEntries / ids.size() > minTowerEntries / 25.) {
0237         if (summaryEntries < towerEntries * (1. - tailPopulThreshold_))  // large timing deviation
0238           quality = doMask ? kMBad : kBad;
0239         else {
0240           towerMean /= towerEntries;
0241           towerMean2 /= towerEntries;
0242
0243           float towerRMS(sqrt(towerMean2 - towerMean * towerMean));
0244
0245           if (std::abs(towerMean) > meanThresh || towerRMS > rmsThresh)
0246             quality = doMask ? kMBad : kBad;
0247           else
0248             quality = doMask ? kMGood : kGood;
0249         }
0250       }
0251       qsItr->setBinContent(quality);
0252     }
0253   }
0254
0255   DEFINE_ECALDQM_WORKER(TimingClient);
0256 }  // namespace ecaldqm