File indexing completed on 2024-04-06 12:06:35
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 <string>
0024 #include <vector>
0025
0026 #include "TH1F.h"
0027 #include "TH2F.h"
0028
0029 #include "FWCore/Framework/interface/Frameworkfwd.h"
0030 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0031
0032 #include "FWCore/Framework/interface/Event.h"
0033 #include "FWCore/Framework/interface/MakerMacros.h"
0034 #include "FWCore/Framework/interface/Run.h"
0035
0036 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0037
0038 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0039
0040 #include "FWCore/ServiceRegistry/interface/Service.h"
0041
0042 #include "FWCore/Utilities/interface/InputTag.h"
0043
0044 #include "DPGAnalysis/SiStripTools/interface/EventWithHistory.h"
0045 #include "DPGAnalysis/SiStripTools/interface/APVCyclePhaseCollection.h"
0046
0047 #include "DPGAnalysis/SiStripTools/interface/RunHistogramManager.h"
0048
0049
0050
0051
0052 class EventTimeDistribution : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
0053 public:
0054 explicit EventTimeDistribution(const edm::ParameterSet&);
0055 ~EventTimeDistribution() override;
0056
0057 private:
0058 void beginJob() override;
0059 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0060 void endRun(const edm::Run&, const edm::EventSetup&) override;
0061 void analyze(const edm::Event&, const edm::EventSetup&) override;
0062 void endJob() override;
0063
0064
0065
0066 edm::EDGetTokenT<EventWithHistory> _historyProductToken;
0067 edm::EDGetTokenT<APVCyclePhaseCollection> _apvphasecollToken;
0068 const std::string _phasepart;
0069 const bool _wantdbxvsbxincycle;
0070 const bool _wantdbxvsbx;
0071 const bool _wantbxincyclevsbx;
0072 const bool _wantorbitvsbxincycle;
0073 unsigned int _nevents;
0074 const unsigned int m_maxLS;
0075 const unsigned int m_LSfrac;
0076 const bool m_ewhdepthHisto;
0077
0078 RunHistogramManager _rhm;
0079
0080 TH1F** _dbx;
0081 std::vector<TH1F**> m_dbxhistos;
0082 std::vector<std::pair<unsigned int, unsigned int> > m_dbxindices;
0083 TH1F** _bx;
0084 TH1F** _bxincycle;
0085 TH1F** _orbit;
0086 TH2F** _dbxvsbxincycle;
0087 TH2F** _dbxvsbx;
0088 TH2F** _bxincyclevsbx;
0089 TH2F** _orbitvsbxincycle;
0090 TH1F** m_ewhdepth;
0091 };
0092
0093
0094
0095
0096
0097
0098
0099
0100
0101
0102
0103
0104 EventTimeDistribution::EventTimeDistribution(const edm::ParameterSet& iConfig)
0105 : _historyProductToken(consumes<EventWithHistory>(iConfig.getParameter<edm::InputTag>("historyProduct"))),
0106 _apvphasecollToken(consumes<APVCyclePhaseCollection>(iConfig.getParameter<edm::InputTag>("apvPhaseCollection"))),
0107 _phasepart(iConfig.getUntrackedParameter<std::string>("phasePartition", "None")),
0108 _wantdbxvsbxincycle(iConfig.getUntrackedParameter<bool>("wantDBXvsBXincycle", false)),
0109 _wantdbxvsbx(iConfig.getUntrackedParameter<bool>("wantDBXvsBX", false)),
0110 _wantbxincyclevsbx(iConfig.getUntrackedParameter<bool>("wantBXincyclevsBX", false)),
0111 _wantorbitvsbxincycle(iConfig.getUntrackedParameter<bool>("wantOrbitvsBXincycle", false)),
0112 _nevents(0),
0113 m_maxLS(iConfig.getUntrackedParameter<unsigned int>("maxLSBeforeRebin", 100)),
0114 m_LSfrac(iConfig.getUntrackedParameter<unsigned int>("startingLSFraction", 4)),
0115 m_ewhdepthHisto(iConfig.getUntrackedParameter<bool>("wantEWHDepthHisto", false)),
0116 _rhm(consumesCollector()),
0117 _dbxvsbxincycle(nullptr),
0118 _dbxvsbx(nullptr),
0119 _bxincyclevsbx(nullptr),
0120 _orbitvsbxincycle(nullptr),
0121 m_ewhdepth(nullptr) {
0122
0123
0124 std::vector<edm::ParameterSet> dbxhistoparams(iConfig.getUntrackedParameter<std::vector<edm::ParameterSet> >(
0125 "dbxHistosParams", std::vector<edm::ParameterSet>()));
0126
0127 for (std::vector<edm::ParameterSet>::const_iterator params = dbxhistoparams.begin(); params != dbxhistoparams.end();
0128 ++params) {
0129 m_dbxindices.push_back(std::pair<unsigned int, unsigned int>(params->getParameter<unsigned int>("firstEvent"),
0130 params->getParameter<unsigned int>("secondEvent")));
0131 char hname[300];
0132 sprintf(hname,
0133 "dbx_%d_%d",
0134 params->getParameter<unsigned int>("firstEvent"),
0135 params->getParameter<unsigned int>("secondEvent"));
0136 char htitle[300];
0137 sprintf(htitle,
0138 "dbx(%d,%d)",
0139 params->getParameter<unsigned int>("firstEvent"),
0140 params->getParameter<unsigned int>("secondEvent"));
0141
0142 m_dbxhistos.push_back(_rhm.makeTH1F(hname,
0143 htitle,
0144 params->getParameter<int>("nbins"),
0145 params->getParameter<double>("min"),
0146 params->getParameter<double>("max")));
0147 LogDebug("DBXHistoPreBooking") << "Booked DBX histo named " << hname << " untitled " << htitle;
0148 }
0149
0150 _dbx = _rhm.makeTH1F("dbx", "dbx", 1000, -0.5, 999.5);
0151 _bx = _rhm.makeTH1F("bx", "BX number", 3564, -0.5, 3563.5);
0152 _bxincycle = _rhm.makeTH1F("bxcycle", "bxcycle", 70, -0.5, 69.5);
0153 _orbit = _rhm.makeTH1F("orbit", "orbit", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
0154 if (_wantdbxvsbxincycle)
0155 _dbxvsbxincycle = _rhm.makeTH2F("dbxvsbxincycle", "dbxvsbxincycle", 70, -0.5, 69.5, 1000, -0.5, 999.5);
0156 if (_wantdbxvsbx)
0157 _dbxvsbx = _rhm.makeTH2F("dbxvsbx", "dbxvsbx", 3564, -0.5, 3563.5, 1000, -0.5, 999.5);
0158 if (_wantbxincyclevsbx)
0159 _bxincyclevsbx = _rhm.makeTH2F("bxincyclevsbx", "bxincyclevsbx", 3564, -0.5, 3563.5, 70, -0.5, 69.5);
0160 if (_wantorbitvsbxincycle)
0161 _orbitvsbxincycle =
0162 _rhm.makeTH2F("orbitvsbxincycle", "orbitvsbxincycle", 70, -0.5, 69.5, m_maxLS, 0, m_maxLS * 262144);
0163 if (m_ewhdepthHisto)
0164 m_ewhdepth = _rhm.makeTH1F("ewhdepth", "EventWithHistory Depth", 11, -0.5, 10.5);
0165
0166 edm::LogInfo("UsedAPVCyclePhaseCollection")
0167 << " APVCyclePhaseCollection " << iConfig.getParameter<edm::InputTag>("apvPhaseCollection") << " used";
0168 }
0169
0170 EventTimeDistribution::~EventTimeDistribution() {
0171
0172
0173 }
0174
0175
0176
0177
0178
0179
0180 void EventTimeDistribution::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0181 using namespace edm;
0182
0183 _nevents++;
0184
0185 edm::Handle<EventWithHistory> he;
0186 iEvent.getByToken(_historyProductToken, he);
0187
0188
0189
0190 (*_dbx)->Fill(he->deltaBX());
0191 std::vector<std::pair<unsigned int, unsigned int> >::const_iterator indices = m_dbxindices.begin();
0192 for (std::vector<TH1F**>::const_iterator dbxhist = m_dbxhistos.begin(); dbxhist != m_dbxhistos.end();
0193 ++dbxhist, ++indices) {
0194 (*(*dbxhist))->Fill(he->deltaBX(indices->first, indices->second));
0195 }
0196
0197 (*_bx)->Fill(iEvent.bunchCrossing() % 3564);
0198 (*_orbit)->Fill(iEvent.orbitNumber());
0199 if (_dbxvsbx && *_dbxvsbx)
0200 (*_dbxvsbx)->Fill(iEvent.bunchCrossing() % 3564, he->deltaBX());
0201 if (m_ewhdepth && *m_ewhdepth)
0202 (*m_ewhdepth)->Fill(he->depth());
0203
0204 edm::Handle<APVCyclePhaseCollection> apvphase;
0205 iEvent.getByToken(_apvphasecollToken, apvphase);
0206
0207 long long tbx = he->absoluteBX();
0208 if (apvphase.isValid() && !apvphase.failedToGet()) {
0209 const int thephase = apvphase->getPhase(_phasepart);
0210 if (thephase != APVCyclePhaseCollection::invalid && thephase != APVCyclePhaseCollection::multiphase &&
0211 thephase != APVCyclePhaseCollection::nopartition) {
0212 tbx -= thephase;
0213 (*_bxincycle)->Fill(tbx % 70);
0214 if (_dbxvsbxincycle && *_dbxvsbxincycle)
0215 (*_dbxvsbxincycle)->Fill(tbx % 70, he->deltaBX());
0216 if (_bxincyclevsbx && *_bxincyclevsbx)
0217 (*_bxincyclevsbx)->Fill(iEvent.bunchCrossing() % 3564, tbx % 70);
0218 if (_orbitvsbxincycle && *_orbitvsbxincycle)
0219 (*_orbitvsbxincycle)->Fill(tbx % 70, iEvent.orbitNumber());
0220
0221 } else {
0222 LogDebug("InvalidPhase") << "Invalid APVCyclePhase value : " << _phasepart << " " << thephase;
0223 }
0224 }
0225 }
0226
0227 void EventTimeDistribution::beginRun(const edm::Run& iRun, const edm::EventSetup&) {
0228 _rhm.beginRun(iRun);
0229
0230 if (*_dbx) {
0231 (*_dbx)->GetXaxis()->SetTitle("#DeltaBX");
0232 }
0233
0234 LogDebug("NomberOfHistos") << m_dbxhistos.size();
0235 for (std::vector<TH1F**>::const_iterator dbxhist = m_dbxhistos.begin(); dbxhist != m_dbxhistos.end(); ++dbxhist) {
0236 LogDebug("HistoPointer") << *dbxhist;
0237 if (*(*dbxhist)) {
0238 (*(*dbxhist))->GetXaxis()->SetTitle("#DeltaBX");
0239 }
0240 }
0241 LogDebug("LabelDone") << "all labels set";
0242
0243 if (*_bx) {
0244 (*_bx)->GetXaxis()->SetTitle("BX");
0245 }
0246
0247 if (*_bxincycle) {
0248 (*_bxincycle)->GetXaxis()->SetTitle("Event BX mod(70)");
0249 }
0250
0251 if (*_orbit) {
0252 (*_orbit)->SetCanExtend(TH1::kXaxis);
0253 (*_orbit)->GetXaxis()->SetTitle("time [Orb#]");
0254 }
0255
0256 LogDebug("StdPlotsDone") << "all labels in std plots set";
0257
0258 if (_dbxvsbxincycle && *_dbxvsbxincycle) {
0259 (*_dbxvsbxincycle)->GetXaxis()->SetTitle("Event BX mod(70)");
0260 (*_dbxvsbxincycle)->GetYaxis()->SetTitle("#DeltaBX");
0261 }
0262
0263 if (_dbxvsbx && *_dbxvsbx) {
0264 (*_dbxvsbx)->GetXaxis()->SetTitle("BX");
0265 (*_dbxvsbx)->GetYaxis()->SetTitle("#DeltaBX");
0266 }
0267
0268 if (_bxincyclevsbx && *_bxincyclevsbx) {
0269 (*_bxincyclevsbx)->GetXaxis()->SetTitle("BX");
0270 (*_bxincyclevsbx)->GetYaxis()->SetTitle("Event BX mod(70)");
0271 }
0272
0273 if (_orbitvsbxincycle && *_orbitvsbxincycle) {
0274 (*_orbitvsbxincycle)->SetCanExtend(TH1::kYaxis);
0275 (*_orbitvsbxincycle)->GetXaxis()->SetTitle("Event BX mod(70)");
0276 (*_orbitvsbxincycle)->GetYaxis()->SetTitle("time [Orb#]");
0277 }
0278
0279 if (m_ewhdepth && *m_ewhdepth) {
0280 (*m_ewhdepth)->GetXaxis()->SetTitle("Depth");
0281 }
0282 }
0283
0284 void EventTimeDistribution::endRun(const edm::Run& iRun, const edm::EventSetup&) {}
0285
0286
0287 void EventTimeDistribution::beginJob() {}
0288
0289
0290 void EventTimeDistribution::endJob() { edm::LogInfo("EndOfJob") << _nevents << " analyzed events"; }
0291
0292
0293 DEFINE_FWK_MODULE(EventTimeDistribution);