File indexing completed on 2024-04-06 12:06:34
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021
0022
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0025
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/Run.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029
0030 #include "FWCore/Framework/interface/EventSetup.h"
0031
0032 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0033
0034 #include "FWCore/Utilities/interface/InputTag.h"
0035
0036 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0037 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0038 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0039
0040 #include "FWCore/Utilities/interface/Exception.h"
0041
0042 #include <map>
0043 #include <vector>
0044 #include <utility>
0045 #include <string>
0046 #include <iostream>
0047
0048 #include "TH1F.h"
0049 #include "TProfile.h"
0050
0051 #include "DataFormats/Scalers/interface/Level1TriggerScalers.h"
0052
0053 #include "DPGAnalysis/SiStripTools/interface/RunHistogramManager.h"
0054
0055
0056
0057
0058
0059 class APVCyclePhaseDebuggerFromL1TS : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
0060 public:
0061 explicit APVCyclePhaseDebuggerFromL1TS(const edm::ParameterSet&);
0062 ~APVCyclePhaseDebuggerFromL1TS() override;
0063
0064 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0065
0066 private:
0067 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0068 void analyze(const edm::Event&, const edm::EventSetup&) override;
0069 void endRun(const edm::Run&, const edm::EventSetup&) override {}
0070
0071
0072
0073 edm::EDGetTokenT<Level1TriggerScalersCollection> _l1tscollectionToken;
0074 const unsigned int m_maxLS;
0075 const unsigned int m_LSfrac;
0076
0077 RunHistogramManager m_rhm;
0078
0079 TH1F** _hsize;
0080 TH1F** _hlresync;
0081 TH1F** _hlOC0;
0082 TH1F** _hlTE;
0083 TH1F** _hlstart;
0084 TH1F** _hlEC0;
0085 TH1F** _hlHR;
0086
0087 TH1F** _hdlec0lresync;
0088 TH1F** _hdlresynclHR;
0089
0090 long long _lastResync;
0091 long long _lastHardReset;
0092 long long _lastStart;
0093 long long _lastEventCounter0;
0094 long long _lastOrbitCounter0;
0095 long long _lastTestEnable;
0096 };
0097
0098
0099
0100
0101
0102
0103
0104
0105
0106
0107
0108
0109 APVCyclePhaseDebuggerFromL1TS::APVCyclePhaseDebuggerFromL1TS(const edm::ParameterSet& iConfig)
0110 : _l1tscollectionToken(
0111 consumes<Level1TriggerScalersCollection>(iConfig.getParameter<edm::InputTag>("l1TSCollection"))),
0112 m_maxLS(iConfig.getUntrackedParameter<unsigned int>("maxLSBeforeRebin", 250)),
0113 m_LSfrac(iConfig.getUntrackedParameter<unsigned int>("startingLSFraction", 16)),
0114 m_rhm(consumesCollector()),
0115 _hsize(nullptr),
0116 _hlresync(nullptr),
0117 _hlOC0(nullptr),
0118 _hlTE(nullptr),
0119 _hlstart(nullptr),
0120 _hlEC0(nullptr),
0121 _hlHR(nullptr),
0122 _hdlec0lresync(nullptr),
0123 _hdlresynclHR(nullptr),
0124 _lastResync(-1),
0125 _lastHardReset(-1),
0126 _lastStart(-1),
0127 _lastEventCounter0(-1),
0128 _lastOrbitCounter0(-1),
0129 _lastTestEnable(-1) {
0130
0131
0132 _hsize = m_rhm.makeTH1F("size", "Level1TriggerScalers Collection size", 20, -0.5, 19.5);
0133
0134 _hlresync = m_rhm.makeTH1F("lresync", "Orbit of last resync", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
0135 _hlOC0 = m_rhm.makeTH1F("lOC0", "Orbit of last OC0", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
0136 _hlTE = m_rhm.makeTH1F("lTE", "Orbit of last TestEnable", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
0137 _hlstart = m_rhm.makeTH1F("lstart", "Orbit of last Start", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
0138 _hlEC0 = m_rhm.makeTH1F("lEC0", "Orbit of last EC0", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
0139 _hlHR = m_rhm.makeTH1F("lHR", "Orbit of last HardReset", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
0140 _hdlec0lresync = m_rhm.makeTH1F("dlec0lresync", "Orbit difference EC0-Resync", 4000, -1999.5, 2000.5);
0141 _hdlresynclHR = m_rhm.makeTH1F("dlresynclHR", "Orbit difference Resync-HR", 4000, -1999.5, 2000.5);
0142 }
0143
0144 APVCyclePhaseDebuggerFromL1TS::~APVCyclePhaseDebuggerFromL1TS() {
0145
0146
0147 }
0148
0149
0150
0151
0152
0153
0154 void APVCyclePhaseDebuggerFromL1TS::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup)
0155
0156 {
0157
0158
0159 m_rhm.beginRun(iRun);
0160
0161 if (_hlresync && *_hlresync) {
0162 (*_hlresync)->GetXaxis()->SetTitle("Orbit");
0163 (*_hlresync)->GetYaxis()->SetTitle("Events");
0164 (*_hlresync)->SetCanExtend(TH1::kXaxis);
0165 }
0166
0167 if (_hlOC0 && *_hlOC0) {
0168 (*_hlOC0)->GetXaxis()->SetTitle("Orbit");
0169 (*_hlOC0)->GetYaxis()->SetTitle("Events");
0170 (*_hlOC0)->SetCanExtend(TH1::kXaxis);
0171 }
0172
0173 if (_hlTE && *_hlTE) {
0174 (*_hlTE)->GetXaxis()->SetTitle("Orbit");
0175 (*_hlTE)->GetYaxis()->SetTitle("Events");
0176 (*_hlTE)->SetCanExtend(TH1::kXaxis);
0177 }
0178
0179 if (_hlstart && *_hlstart) {
0180 (*_hlstart)->GetXaxis()->SetTitle("Orbit");
0181 (*_hlstart)->GetYaxis()->SetTitle("Events");
0182 (*_hlstart)->SetCanExtend(TH1::kXaxis);
0183 }
0184
0185 if (_hlEC0 && *_hlEC0) {
0186 (*_hlEC0)->GetXaxis()->SetTitle("Orbit");
0187 (*_hlEC0)->GetYaxis()->SetTitle("Events");
0188 (*_hlEC0)->SetCanExtend(TH1::kXaxis);
0189 }
0190
0191 if (_hlHR && *_hlHR) {
0192 (*_hlHR)->GetXaxis()->SetTitle("Orbit");
0193 (*_hlHR)->GetYaxis()->SetTitle("Events");
0194 (*_hlHR)->SetCanExtend(TH1::kXaxis);
0195 }
0196
0197 if (_hdlec0lresync && *_hdlec0lresync) {
0198 (*_hdlec0lresync)->GetXaxis()->SetTitle("lastEC0-lastResync");
0199 }
0200
0201 if (_hdlresynclHR && *_hdlresynclHR) {
0202 (*_hdlresynclHR)->GetXaxis()->SetTitle("lastEC0-lastResync");
0203 }
0204 }
0205
0206 void APVCyclePhaseDebuggerFromL1TS::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0207 using namespace edm;
0208
0209 Handle<Level1TriggerScalersCollection> l1ts;
0210 iEvent.getByToken(_l1tscollectionToken, l1ts);
0211
0212 if (_hsize && *_hsize)
0213 (*_hsize)->Fill(l1ts->size());
0214
0215
0216
0217 if (!l1ts->empty()) {
0218 if (_hlresync && *_hlresync)
0219 (*_hlresync)->Fill((*l1ts)[0].lastResync());
0220 if (_hlOC0 && *_hlOC0)
0221 (*_hlOC0)->Fill((*l1ts)[0].lastOrbitCounter0());
0222 if (_hlTE && *_hlTE)
0223 (*_hlTE)->Fill((*l1ts)[0].lastTestEnable());
0224 if (_hlstart && *_hlstart)
0225 (*_hlstart)->Fill((*l1ts)[0].lastStart());
0226 if (_hlEC0 && *_hlEC0)
0227 (*_hlEC0)->Fill((*l1ts)[0].lastEventCounter0());
0228 if (_hlHR && *_hlHR)
0229 (*_hlHR)->Fill((*l1ts)[0].lastHardReset());
0230
0231 if (_lastResync != (*l1ts)[0].lastResync()) {
0232 _lastResync = (*l1ts)[0].lastResync();
0233 if (_hdlec0lresync && *_hdlec0lresync)
0234 (*_hdlec0lresync)->Fill((*l1ts)[0].lastEventCounter0() - (*l1ts)[0].lastResync());
0235 LogDebug("TTCSignalReceived") << "New Resync at orbit " << _lastResync;
0236 }
0237 if (_lastHardReset != (*l1ts)[0].lastHardReset()) {
0238 _lastHardReset = (*l1ts)[0].lastHardReset();
0239 if (_hdlresynclHR && *_hdlresynclHR)
0240 (*_hdlresynclHR)->Fill((*l1ts)[0].lastResync() - (*l1ts)[0].lastHardReset());
0241 LogDebug("TTCSignalReceived") << "New HardReset at orbit " << _lastHardReset;
0242 }
0243 if (_lastTestEnable != (*l1ts)[0].lastTestEnable()) {
0244 _lastTestEnable = (*l1ts)[0].lastTestEnable();
0245
0246 }
0247 if (_lastOrbitCounter0 != (*l1ts)[0].lastOrbitCounter0()) {
0248 _lastOrbitCounter0 = (*l1ts)[0].lastOrbitCounter0();
0249 LogDebug("TTCSignalReceived") << "New OrbitCounter0 at orbit " << _lastOrbitCounter0;
0250 }
0251 if (_lastEventCounter0 != (*l1ts)[0].lastEventCounter0()) {
0252 _lastEventCounter0 = (*l1ts)[0].lastEventCounter0();
0253 LogDebug("TTCSignalReceived") << "New EventCounter0 at orbit " << _lastEventCounter0;
0254 }
0255 if (_lastStart != (*l1ts)[0].lastStart()) {
0256 _lastStart = (*l1ts)[0].lastStart();
0257 LogDebug("TTCSignalReceived") << "New Start at orbit " << _lastStart;
0258 }
0259 }
0260 }
0261
0262 void APVCyclePhaseDebuggerFromL1TS::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0263 edm::ParameterSetDescription desc;
0264 desc.add<edm::InputTag>("l1TSCollection", edm::InputTag("scalersRawToDigi"));
0265 descriptions.add("l1TSDebugger", desc);
0266 }
0267
0268
0269 DEFINE_FWK_MODULE(APVCyclePhaseDebuggerFromL1TS);