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/EDProducer.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/MessageLogger/interface/MessageLogger.h"
0031
0032 #include "FWCore/Utilities/interface/InputTag.h"
0033
0034 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0035 #include "FWCore/Utilities/interface/Exception.h"
0036
0037 #include <map>
0038 #include <vector>
0039 #include <string>
0040
0041 #include "TH1F.h"
0042 #include "TProfile.h"
0043
0044 #include "DataFormats/Scalers/interface/L1AcceptBunchCrossing.h"
0045 #include "DataFormats/TCDS/interface/TCDSRecord.h"
0046 #include "DPGAnalysis/SiStripTools/interface/APVCyclePhaseCollection.h"
0047
0048 #include "DPGAnalysis/SiStripTools/interface/RunHistogramManager.h"
0049
0050
0051
0052
0053
0054 class APVCyclePhaseProducerFromL1ABC : public edm::one::EDProducer<edm::one::WatchRuns> {
0055 public:
0056 explicit APVCyclePhaseProducerFromL1ABC(const edm::ParameterSet&);
0057 ~APVCyclePhaseProducerFromL1ABC() override;
0058
0059 private:
0060 void beginJob() override;
0061 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0062 void endRun(const edm::Run&, const edm::EventSetup&) override;
0063 void produce(edm::Event&, const edm::EventSetup&) override;
0064 void endJob() override;
0065
0066
0067
0068 edm::EDGetTokenT<L1AcceptBunchCrossingCollection> _l1abccollectionToken;
0069 edm::EDGetTokenT<TCDSRecord> _tcdsRecordToken;
0070 const std::vector<std::string> _defpartnames;
0071 const std::vector<int> _defphases;
0072 const int _orbitoffsetSOR;
0073 const bool _wantHistos;
0074 const bool _forceSCAL;
0075 RunHistogramManager m_rhm;
0076
0077 TH1F** _hbx;
0078 TH1F** _hdbx;
0079 TH1F** _hdorbit;
0080 const unsigned int _firstgoodrun;
0081 std::map<edm::EventNumber_t, long long> _offsets;
0082 long long _curroffset;
0083 edm::EventNumber_t _curroffevent;
0084 };
0085
0086
0087
0088
0089
0090
0091
0092
0093
0094
0095
0096
0097 APVCyclePhaseProducerFromL1ABC::APVCyclePhaseProducerFromL1ABC(const edm::ParameterSet& iConfig)
0098 : _l1abccollectionToken(
0099 mayConsume<L1AcceptBunchCrossingCollection>(iConfig.getParameter<edm::InputTag>("l1ABCCollection"))),
0100 _tcdsRecordToken(mayConsume<TCDSRecord>(iConfig.getParameter<edm::InputTag>("tcdsRecordLabel"))),
0101 _defpartnames(iConfig.getParameter<std::vector<std::string> >("defaultPartitionNames")),
0102 _defphases(iConfig.getParameter<std::vector<int> >("defaultPhases")),
0103 _orbitoffsetSOR(iConfig.getParameter<int>("StartOfRunOrbitOffset")),
0104 _wantHistos(iConfig.getUntrackedParameter<bool>("wantHistos", false)),
0105 _forceSCAL(iConfig.getParameter<bool>("forceSCAL")),
0106 m_rhm(consumesCollector()),
0107 _hbx(nullptr),
0108 _hdbx(nullptr),
0109 _hdorbit(nullptr),
0110 _firstgoodrun(110878),
0111 _offsets(),
0112 _curroffset(0),
0113 _curroffevent(0) {
0114 produces<APVCyclePhaseCollection, edm::InEvent>();
0115
0116
0117
0118 if (_wantHistos) {
0119 _hbx = m_rhm.makeTH1F("l1abcbx", "BX number from TCDS (or L1ABC as fallback)", 4096, -0.5, 4095.5);
0120 _hdbx = m_rhm.makeTH1F("dbx", "BX number difference", 4096 * 2 - 1, -4095.5, 4095.5);
0121 _hdorbit = m_rhm.makeTH1F("dorbit", "Orbit Number difference", 9999, -4999.5, 4999.5);
0122 }
0123 }
0124
0125 APVCyclePhaseProducerFromL1ABC::~APVCyclePhaseProducerFromL1ABC() {
0126
0127
0128 }
0129
0130
0131
0132
0133
0134
0135 void APVCyclePhaseProducerFromL1ABC::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup)
0136
0137 {
0138
0139
0140 _offsets.clear();
0141 edm::LogInfo("AbsoluteBXOffsetReset") << "Absolute BX offset map reset";
0142
0143 if (_wantHistos) {
0144 m_rhm.beginRun(iRun);
0145
0146 if (_hbx && *_hbx) {
0147 (*_hbx)->GetXaxis()->SetTitle("BX");
0148 (*_hbx)->GetYaxis()->SetTitle("Events");
0149 }
0150
0151 if (_hdbx && *_hdbx) {
0152 (*_hdbx)->GetXaxis()->SetTitle("#DeltaBX");
0153 (*_hdbx)->GetYaxis()->SetTitle("Events");
0154 }
0155
0156 if (_hdorbit && *_hdorbit) {
0157 (*_hdorbit)->GetXaxis()->SetTitle("#Deltaorbit");
0158 (*_hdorbit)->GetYaxis()->SetTitle("Events");
0159 }
0160 }
0161
0162 if (iRun.run() < _firstgoodrun) {
0163 edm::LogInfo("UnreliableMissingL1AcceptBunchCrossingCollection")
0164 << "In this run L1AcceptBunchCrossingCollection is missing or unreliable: default phases will be used";
0165 }
0166 }
0167
0168 void APVCyclePhaseProducerFromL1ABC::endRun(const edm::Run&, const edm::EventSetup&) {
0169
0170
0171 edm::LogInfo("L1AcceptBunchCrossingAbsoluteBXOffsetSummary") << "Absolute BX offset summary:";
0172 for (std::map<edm::EventNumber_t, long long>::const_iterator offset = _offsets.begin(); offset != _offsets.end();
0173 ++offset) {
0174 edm::LogVerbatim("L1AcceptBunchCrossingAbsoluteBXOffsetSummary") << offset->first << " " << offset->second;
0175 }
0176 }
0177
0178 void APVCyclePhaseProducerFromL1ABC::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0179 using namespace edm;
0180
0181 std::unique_ptr<APVCyclePhaseCollection> apvphases(new APVCyclePhaseCollection());
0182
0183 const std::vector<int>& phases = _defphases;
0184 const std::vector<std::string>& partnames = _defpartnames;
0185
0186
0187
0188 int phasechange = 0;
0189
0190 if (iEvent.run() >= _firstgoodrun) {
0191 Handle<L1AcceptBunchCrossingCollection> pIn;
0192 iEvent.getByToken(_l1abccollectionToken, pIn);
0193 Handle<TCDSRecord> tcds_pIn;
0194 iEvent.getByToken(_tcdsRecordToken, tcds_pIn);
0195 bool useTCDS(tcds_pIn.isValid() && !_forceSCAL);
0196 const auto* tcdsRecord = useTCDS ? tcds_pIn.product() : nullptr;
0197
0198
0199 long long orbitoffset = _orbitoffsetSOR;
0200 int bxoffset = 0;
0201
0202 if (useTCDS) {
0203 if (tcdsRecord->getL1aHistoryEntry(0).getIndex() == 0) {
0204 if (tcdsRecord->getEventType() != 0) {
0205 orbitoffset = (long long)tcdsRecord->getOrbitNr() - (long long)iEvent.orbitNumber();
0206
0207 bxoffset = iEvent.bunchCrossing() - tcdsRecord->getBXID();
0208
0209 if (_wantHistos) {
0210 if (_hbx && *_hbx)
0211 (*_hbx)->Fill(tcdsRecord->getBXID());
0212 if (_hdbx && *_hdbx)
0213 (*_hdbx)->Fill(bxoffset);
0214 if (_hdorbit && *_hdorbit)
0215 (*_hdorbit)->Fill(orbitoffset);
0216 }
0217 }
0218 } else {
0219 edm::LogWarning("L1AcceptBunchCrossingNoType") << "L1AcceptBunchCrossing with no type found: ";
0220 edm::LogPrint("L1AcceptBunchCrossingNoType") << *tcdsRecord;
0221 }
0222 } else {
0223 for (L1AcceptBunchCrossingCollection::const_iterator l1abc = pIn->begin(); l1abc != pIn->end(); ++l1abc) {
0224 if (l1abc->l1AcceptOffset() == 0) {
0225 if (l1abc->eventType() != 0) {
0226 orbitoffset = (long long)iEvent.orbitNumber() - (long long)l1abc->orbitNumber();
0227 bxoffset =
0228 iEvent.bunchCrossing() -
0229 l1abc->bunchCrossing();
0230
0231 if (_wantHistos) {
0232 if (_hbx && *_hbx)
0233 (*_hbx)->Fill(l1abc->bunchCrossing());
0234 if (_hdbx && *_hdbx)
0235 (*_hdbx)->Fill(bxoffset);
0236 if (_hdorbit && *_hdorbit)
0237 (*_hdorbit)->Fill(orbitoffset);
0238 }
0239 } else {
0240 edm::LogWarning("L1AcceptBunchCrossingNoType") << "L1AcceptBunchCrossing with no type found: ";
0241 for (L1AcceptBunchCrossingCollection::const_iterator debu = pIn->begin(); debu != pIn->end(); ++debu) {
0242 edm::LogPrint("L1AcceptBunchCrossingNoType") << *debu;
0243 }
0244 }
0245 }
0246 }
0247 }
0248
0249 long long absbxoffset = orbitoffset * 3564 + bxoffset;
0250
0251 if (orbitoffset != _orbitoffsetSOR)
0252 phasechange = (orbitoffset * 3564) % 70;
0253
0254 if (_offsets.empty()) {
0255 _curroffset = absbxoffset;
0256 _curroffevent = iEvent.id().event();
0257 _offsets[iEvent.id().event()] = absbxoffset;
0258 } else {
0259 if (_curroffset != absbxoffset || iEvent.id().event() < _curroffevent) {
0260 if (_curroffset != absbxoffset) {
0261 edm::LogInfo("L1AcceptBunchCrossingAbsoluteBXOffsetChanged")
0262 << "Absolute BX offset changed from " << _curroffset << " to " << absbxoffset << " at orbit "
0263 << iEvent.orbitNumber() << " and BX " << iEvent.bunchCrossing();
0264 if (useTCDS) {
0265 edm::LogVerbatim("AbsoluteBXOffsetChanged") << *tcdsRecord;
0266 } else {
0267 for (L1AcceptBunchCrossingCollection::const_iterator l1abc = pIn->begin(); l1abc != pIn->end(); ++l1abc) {
0268 edm::LogVerbatim("L1AcceptBunchCrossingAbsoluteBXOffsetChanged") << *l1abc;
0269 }
0270 }
0271 }
0272
0273 _curroffset = absbxoffset;
0274 _curroffevent = iEvent.id().event();
0275 _offsets[iEvent.id().event()] = absbxoffset;
0276 }
0277 }
0278 }
0279
0280 if (phases.size() < partnames.size()) {
0281
0282 throw cms::Exception("InvalidAPVCyclePhases")
0283 << " Inconsistent phases/partitions vector sizes: " << phases.size() << " " << partnames.size();
0284 }
0285
0286 for (unsigned int ipart = 0; ipart < partnames.size(); ++ipart) {
0287 if (phases[ipart] >= 0) {
0288 apvphases->get()[partnames[ipart]] = (phases[ipart] + phasechange) % 70;
0289 }
0290 }
0291
0292 iEvent.put(std::move(apvphases));
0293 }
0294
0295
0296 void APVCyclePhaseProducerFromL1ABC::beginJob() {}
0297
0298
0299 void APVCyclePhaseProducerFromL1ABC::endJob() {}
0300
0301
0302 DEFINE_FWK_MODULE(APVCyclePhaseProducerFromL1ABC);