File indexing completed on 2023-03-17 10:54:36
0001 #include "DQM/EcalMonitorTasks/interface/PedestalTask.h"
0002
0003 #include <iomanip>
0004
0005 #include "DataFormats/DetId/interface/DetId.h"
0006 #include "DataFormats/EcalDigi/interface/EcalDataFrame.h"
0007
0008 #include "DQM/EcalCommon/interface/MESetMulti.h"
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 namespace ecaldqm {
0013 PedestalTask::PedestalTask() : DQWorkerTask(), gainToME_(), pnGainToME_() { std::fill_n(enable_, nDCC, false); }
0014
0015 void PedestalTask::setParams(edm::ParameterSet const& _params) {
0016 std::vector<int> MGPAGains(_params.getUntrackedParameter<std::vector<int> >("MGPAGains"));
0017 std::vector<int> MGPAGainsPN(_params.getUntrackedParameter<std::vector<int> >("MGPAGainsPN"));
0018
0019 MESet::PathReplacements repl;
0020
0021 MESetMulti& pedestal(static_cast<MESetMulti&>(MEs_.at("Pedestal")));
0022 unsigned nG(MGPAGains.size());
0023 for (unsigned iG(0); iG != nG; ++iG) {
0024 int gain(MGPAGains[iG]);
0025 if (gain != 1 && gain != 6 && gain != 12)
0026 throw cms::Exception("InvalidConfiguration") << "MGPA gain";
0027 repl["gain"] = std::to_string(gain);
0028 gainToME_[gain] = pedestal.getIndex(repl);
0029 }
0030
0031 repl.clear();
0032
0033 MESetMulti& pnPedestal(static_cast<MESetMulti&>(MEs_.at("PNPedestal")));
0034 unsigned nGPN(MGPAGainsPN.size());
0035 for (unsigned iG(0); iG != nGPN; ++iG) {
0036 int gain(MGPAGainsPN[iG]);
0037 if (gain != 1 && gain != 16)
0038 throw cms::Exception("InvalidConfiguration") << "PN MGPA gain";
0039 repl["pngain"] = std::to_string(gain);
0040 pnGainToME_[gain] = pnPedestal.getIndex(repl);
0041 }
0042 }
0043
0044 bool PedestalTask::filterRunType(short const* _runType) {
0045 bool enable(false);
0046
0047 for (int iFED(0); iFED < nDCC; iFED++) {
0048 if (_runType[iFED] == EcalDCCHeaderBlock::PEDESTAL_STD || _runType[iFED] == EcalDCCHeaderBlock::PEDESTAL_GAP) {
0049 enable = true;
0050 enable_[iFED] = true;
0051 } else
0052 enable_[iFED] = false;
0053 }
0054
0055 return enable;
0056 }
0057
0058 template <typename DigiCollection>
0059 void PedestalTask::runOnDigis(DigiCollection const& _digis) {
0060 MESet& mePedestal(MEs_.at("Pedestal"));
0061 MESet& meOccupancy(MEs_.at("Occupancy"));
0062
0063 unsigned iME(-1);
0064
0065 for (typename DigiCollection::const_iterator digiItr(_digis.begin()); digiItr != _digis.end(); ++digiItr) {
0066 DetId id(digiItr->id());
0067
0068 int iDCC(dccId(id, GetElectronicsMap()) - 1);
0069
0070 if (!enable_[iDCC])
0071 continue;
0072
0073
0074 EcalDataFrame dataFrame(*digiItr);
0075
0076 int gain(0);
0077 switch (dataFrame.sample(0).gainId()) {
0078 case 1:
0079 gain = 12;
0080 break;
0081 case 2:
0082 gain = 6;
0083 break;
0084 case 3:
0085 gain = 1;
0086 break;
0087 default:
0088 continue;
0089 }
0090
0091 if (gainToME_.find(gain) == gainToME_.end())
0092 continue;
0093
0094 if (iME != gainToME_[gain]) {
0095 iME = gainToME_[gain];
0096 static_cast<MESetMulti&>(mePedestal).use(iME);
0097 }
0098
0099 meOccupancy.fill(getEcalDQMSetupObjects(), id);
0100
0101 for (int iSample(0); iSample < EcalDataFrame::MAXSAMPLES; iSample++)
0102 mePedestal.fill(getEcalDQMSetupObjects(), id, double(dataFrame.sample(iSample).adc()));
0103 }
0104 }
0105
0106 void PedestalTask::runOnPnDigis(EcalPnDiodeDigiCollection const& _digis) {
0107 MESet& mePNPedestal(MEs_.at("PNPedestal"));
0108
0109 unsigned iME(-1);
0110
0111 for (EcalPnDiodeDigiCollection::const_iterator digiItr(_digis.begin()); digiItr != _digis.end(); ++digiItr) {
0112 EcalPnDiodeDetId id(digiItr->id());
0113
0114 int iDCC(dccId(id, GetElectronicsMap()) - 1);
0115
0116 if (!enable_[iDCC])
0117 continue;
0118
0119 int gain(0);
0120 switch (digiItr->sample(0).gainId()) {
0121 case 0:
0122 gain = 1;
0123 break;
0124 case 1:
0125 gain = 16;
0126 break;
0127 default:
0128 continue;
0129 }
0130
0131 if (pnGainToME_.find(gain) == pnGainToME_.end())
0132 continue;
0133
0134 if (iME != pnGainToME_[gain]) {
0135 iME = pnGainToME_[gain];
0136 static_cast<MESetMulti&>(mePNPedestal).use(iME);
0137 }
0138
0139 for (int iSample(0); iSample < 50; iSample++)
0140 mePNPedestal.fill(getEcalDQMSetupObjects(), id, double(digiItr->sample(iSample).adc()));
0141 }
0142 }
0143
0144 DEFINE_ECALDQM_WORKER(PedestalTask);
0145 }