CMS 3D CMS Logo

EventTimeDistribution.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: SiStripTools
4 // Class: EventTimeDistribution
5 //
13 //
14 // Original Author: Andrea Venturi
15 // Created: Tue Jul 19 11:56:00 CEST 2009
16 //
17 //
18 
19 // system include files
20 #include <memory>
21 
22 // user include files
23 #include <string>
24 #include <vector>
25 
26 #include "TH1F.h"
27 #include "TH2F.h"
28 
31 
35 
37 
39 
41 
43 
46 
48 //
49 // class decleration
50 //
51 
52 class EventTimeDistribution : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
53 public:
55  ~EventTimeDistribution() override;
56 
57 private:
58  void beginJob() override;
59  void beginRun(const edm::Run&, const edm::EventSetup&) override;
60  void endRun(const edm::Run&, const edm::EventSetup&) override;
61  void analyze(const edm::Event&, const edm::EventSetup&) override;
62  void endJob() override;
63 
64  // ----------member data ---------------------------
65 
69  const bool _wantdbxvsbxincycle;
70  const bool _wantdbxvsbx;
71  const bool _wantbxincyclevsbx;
73  unsigned int _nevents;
74  const unsigned int m_maxLS;
75  const unsigned int m_LSfrac;
76  const bool m_ewhdepthHisto;
77 
79 
80  TH1F** _dbx;
81  std::vector<TH1F**> m_dbxhistos;
82  std::vector<std::pair<unsigned int, unsigned int> > m_dbxindices;
83  TH1F** _bx;
84  TH1F** _bxincycle;
85  TH1F** _orbit;
87  TH2F** _dbxvsbx;
90  TH1F** m_ewhdepth;
91 };
92 
93 //
94 // constants, enums and typedefs
95 //
96 
97 //
98 // static data member definitions
99 //
100 
101 //
102 // constructors and destructor
103 //
105  : _historyProductToken(consumes<EventWithHistory>(iConfig.getParameter<edm::InputTag>("historyProduct"))),
106  _apvphasecollToken(consumes<APVCyclePhaseCollection>(iConfig.getParameter<edm::InputTag>("apvPhaseCollection"))),
107  _phasepart(iConfig.getUntrackedParameter<std::string>("phasePartition", "None")),
108  _wantdbxvsbxincycle(iConfig.getUntrackedParameter<bool>("wantDBXvsBXincycle", false)),
109  _wantdbxvsbx(iConfig.getUntrackedParameter<bool>("wantDBXvsBX", false)),
110  _wantbxincyclevsbx(iConfig.getUntrackedParameter<bool>("wantBXincyclevsBX", false)),
111  _wantorbitvsbxincycle(iConfig.getUntrackedParameter<bool>("wantOrbitvsBXincycle", false)),
112  _nevents(0),
113  m_maxLS(iConfig.getUntrackedParameter<unsigned int>("maxLSBeforeRebin", 100)),
114  m_LSfrac(iConfig.getUntrackedParameter<unsigned int>("startingLSFraction", 4)),
115  m_ewhdepthHisto(iConfig.getUntrackedParameter<bool>("wantEWHDepthHisto", false)),
116  _rhm(consumesCollector()),
117  _dbxvsbxincycle(nullptr),
118  _dbxvsbx(nullptr),
119  _bxincyclevsbx(nullptr),
120  _orbitvsbxincycle(nullptr),
121  m_ewhdepth(nullptr) {
122  //now do what ever initialization is needed
123 
124  std::vector<edm::ParameterSet> dbxhistoparams(iConfig.getUntrackedParameter<std::vector<edm::ParameterSet> >(
125  "dbxHistosParams", std::vector<edm::ParameterSet>()));
126 
127  for (std::vector<edm::ParameterSet>::const_iterator params = dbxhistoparams.begin(); params != dbxhistoparams.end();
128  ++params) {
129  m_dbxindices.push_back(std::pair<unsigned int, unsigned int>(params->getParameter<unsigned int>("firstEvent"),
130  params->getParameter<unsigned int>("secondEvent")));
131  char hname[300];
132  sprintf(hname,
133  "dbx_%d_%d",
134  params->getParameter<unsigned int>("firstEvent"),
135  params->getParameter<unsigned int>("secondEvent"));
136  char htitle[300];
137  sprintf(htitle,
138  "dbx(%d,%d)",
139  params->getParameter<unsigned int>("firstEvent"),
140  params->getParameter<unsigned int>("secondEvent"));
141 
142  m_dbxhistos.push_back(_rhm.makeTH1F(hname,
143  htitle,
144  params->getParameter<int>("nbins"),
145  params->getParameter<double>("min"),
146  params->getParameter<double>("max")));
147  LogDebug("DBXHistoPreBooking") << "Booked DBX histo named " << hname << " untitled " << htitle;
148  }
149 
150  _dbx = _rhm.makeTH1F("dbx", "dbx", 1000, -0.5, 999.5);
151  _bx = _rhm.makeTH1F("bx", "BX number", 3564, -0.5, 3563.5);
152  _bxincycle = _rhm.makeTH1F("bxcycle", "bxcycle", 70, -0.5, 69.5);
153  _orbit = _rhm.makeTH1F("orbit", "orbit", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
155  _dbxvsbxincycle = _rhm.makeTH2F("dbxvsbxincycle", "dbxvsbxincycle", 70, -0.5, 69.5, 1000, -0.5, 999.5);
156  if (_wantdbxvsbx)
157  _dbxvsbx = _rhm.makeTH2F("dbxvsbx", "dbxvsbx", 3564, -0.5, 3563.5, 1000, -0.5, 999.5);
158  if (_wantbxincyclevsbx)
159  _bxincyclevsbx = _rhm.makeTH2F("bxincyclevsbx", "bxincyclevsbx", 3564, -0.5, 3563.5, 70, -0.5, 69.5);
162  _rhm.makeTH2F("orbitvsbxincycle", "orbitvsbxincycle", 70, -0.5, 69.5, m_maxLS, 0, m_maxLS * 262144);
163  if (m_ewhdepthHisto)
164  m_ewhdepth = _rhm.makeTH1F("ewhdepth", "EventWithHistory Depth", 11, -0.5, 10.5);
165 
166  edm::LogInfo("UsedAPVCyclePhaseCollection")
167  << " APVCyclePhaseCollection " << iConfig.getParameter<edm::InputTag>("apvPhaseCollection") << " used";
168 }
169 
171  // do anything here that needs to be done at desctruction time
172  // (e.g. close files, deallocate resources etc.)
173 }
174 
175 //
176 // member functions
177 //
178 
179 // ------------ method called to for each event ------------
181  using namespace edm;
182 
183  _nevents++;
184 
186  iEvent.getByToken(_historyProductToken, he);
187 
188  // improve the matchin between default and actual partitions
189 
190  (*_dbx)->Fill(he->deltaBX());
191  std::vector<std::pair<unsigned int, unsigned int> >::const_iterator indices = m_dbxindices.begin();
192  for (std::vector<TH1F**>::const_iterator dbxhist = m_dbxhistos.begin(); dbxhist != m_dbxhistos.end();
193  ++dbxhist, ++indices) {
194  (*(*dbxhist))->Fill(he->deltaBX(indices->first, indices->second));
195  }
196 
197  (*_bx)->Fill(iEvent.bunchCrossing() % 3564);
198  (*_orbit)->Fill(iEvent.orbitNumber());
199  if (_dbxvsbx && *_dbxvsbx)
200  (*_dbxvsbx)->Fill(iEvent.bunchCrossing() % 3564, he->deltaBX());
201  if (m_ewhdepth && *m_ewhdepth)
202  (*m_ewhdepth)->Fill(he->depth());
203 
205  iEvent.getByToken(_apvphasecollToken, apvphase);
206 
207  long long tbx = he->absoluteBX();
208  if (apvphase.isValid() && !apvphase.failedToGet()) {
209  const int thephase = apvphase->getPhase(_phasepart);
212  tbx -= thephase;
213  (*_bxincycle)->Fill(tbx % 70);
215  (*_dbxvsbxincycle)->Fill(tbx % 70, he->deltaBX());
217  (*_bxincyclevsbx)->Fill(iEvent.bunchCrossing() % 3564, tbx % 70);
219  (*_orbitvsbxincycle)->Fill(tbx % 70, iEvent.orbitNumber());
220 
221  } else {
222  LogDebug("InvalidPhase") << "Invalid APVCyclePhase value : " << _phasepart << " " << thephase;
223  }
224  }
225 }
226 
228  _rhm.beginRun(iRun);
229 
230  if (*_dbx) {
231  (*_dbx)->GetXaxis()->SetTitle("#DeltaBX");
232  }
233 
234  LogDebug("NomberOfHistos") << m_dbxhistos.size();
235  for (std::vector<TH1F**>::const_iterator dbxhist = m_dbxhistos.begin(); dbxhist != m_dbxhistos.end(); ++dbxhist) {
236  LogDebug("HistoPointer") << *dbxhist;
237  if (*(*dbxhist)) {
238  (*(*dbxhist))->GetXaxis()->SetTitle("#DeltaBX");
239  }
240  }
241  LogDebug("LabelDone") << "all labels set";
242 
243  if (*_bx) {
244  (*_bx)->GetXaxis()->SetTitle("BX");
245  }
246 
247  if (*_bxincycle) {
248  (*_bxincycle)->GetXaxis()->SetTitle("Event BX mod(70)");
249  }
250 
251  if (*_orbit) {
252  (*_orbit)->SetCanExtend(TH1::kXaxis);
253  (*_orbit)->GetXaxis()->SetTitle("time [Orb#]");
254  }
255 
256  LogDebug("StdPlotsDone") << "all labels in std plots set";
257 
259  (*_dbxvsbxincycle)->GetXaxis()->SetTitle("Event BX mod(70)");
260  (*_dbxvsbxincycle)->GetYaxis()->SetTitle("#DeltaBX");
261  }
262 
263  if (_dbxvsbx && *_dbxvsbx) {
264  (*_dbxvsbx)->GetXaxis()->SetTitle("BX");
265  (*_dbxvsbx)->GetYaxis()->SetTitle("#DeltaBX");
266  }
267 
268  if (_bxincyclevsbx && *_bxincyclevsbx) {
269  (*_bxincyclevsbx)->GetXaxis()->SetTitle("BX");
270  (*_bxincyclevsbx)->GetYaxis()->SetTitle("Event BX mod(70)");
271  }
272 
274  (*_orbitvsbxincycle)->SetCanExtend(TH1::kYaxis);
275  (*_orbitvsbxincycle)->GetXaxis()->SetTitle("Event BX mod(70)");
276  (*_orbitvsbxincycle)->GetYaxis()->SetTitle("time [Orb#]");
277  }
278 
279  if (m_ewhdepth && *m_ewhdepth) {
280  (*m_ewhdepth)->GetXaxis()->SetTitle("Depth");
281  }
282 }
283 
285 
286 // ------------ method called once each job just before starting event loop ------------
288 
289 // ------------ method called once each job just after ending the event loop ------------
290 void EventTimeDistribution::endJob() { edm::LogInfo("EndOfJob") << _nevents << " analyzed events"; }
291 
292 //define this as a plug-in
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
void beginRun(const edm::Run &, const edm::EventSetup &) override
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
edm::EDGetTokenT< APVCyclePhaseCollection > _apvphasecollToken
T getUntrackedParameter(std::string const &, T const &) const
edm::EDGetTokenT< EventWithHistory > _historyProductToken
bool failedToGet() const
Definition: HandleBase.h:72
void endRun(const edm::Run &, const edm::EventSetup &) override
TH1F ** makeTH1F(const char *name, const char *title, const unsigned int nbinx, const double xmin, const double xmax)
int iEvent
Definition: GenABIO.cc:224
std::vector< std::pair< unsigned int, unsigned int > > m_dbxindices
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
const unsigned int m_LSfrac
RunHistogramManager _rhm
const std::string _phasepart
EventTimeDistribution(const edm::ParameterSet &)
Log< level::Info, false > LogInfo
void beginRun(const edm::Run &iRun)
bool isValid() const
Definition: HandleBase.h:70
const unsigned int m_maxLS
const int getPhase(const std::string partition) const
HLT enums.
void analyze(const edm::Event &, const edm::EventSetup &) override
TH2F ** makeTH2F(const char *name, const char *title, const unsigned int nbinx, const double xmin, const double xmax, const unsigned int nbiny, const double ymin, const double ymax)
std::vector< TH1F ** > m_dbxhistos
Definition: Run.h:45
#define LogDebug(id)