CMS 3D CMS Logo

SiStripBaselineValidator.cc
Go to the documentation of this file.
1 // Original Author: Ivan Amos Cali
2 // Created: Mon Jul 28 14:10:52 CEST 2008
3 //
4 //
5 
6 // system include files
7 #include <memory>
8 #include <iostream>
9 
10 // user include files
20 
21 //ROOT inclusion
22 #include "TH1F.h"
23 #include "TH2F.h"
24 #include "TString.h"
25 #include <cassert>
26 #include <fstream>
27 
28 class TFile;
29 
30 using namespace edm;
31 using namespace std;
32 
34  srcProcessedRawDigi_ = conf.getParameter<edm::InputTag>("srcProcessedRawDigi");
35  moduleRawDigiToken_ =
36  consumes<edm::DetSetVector<SiStripDigi> >(conf.getParameter<edm::InputTag>("srcProcessedRawDigi"));
37 }
38 
40 
42  const edm::Run& run,
43  const edm::EventSetup& es) {
45  ibooker.setCurrentFolder("SiStrip/BaselineValidator");
46 
47  h1NumbadAPVsRes_ = ibooker.book1D("ResAPVs", ";#ResAPVs", 100, 1.0, 10001);
48 
49  h1ADC_vs_strip_ = ibooker.book2D("ADCvsAPVs", ";ADCvsAPVs", 768, -0.5, 767.5, 1023, -0.5, 1022.5);
50 
51  return;
52 }
53 
56  e.getByToken(moduleRawDigiToken_, moduleRawDigi);
57  edm::DetSetVector<SiStripDigi>::const_iterator itRawDigis = moduleRawDigi->begin();
58 
59  int NumResAPVs = 0;
60  for (; itRawDigis != moduleRawDigi->end(); ++itRawDigis) {
61 
62  edm::DetSet<SiStripDigi>::const_iterator itRaw = itRawDigis->begin();
63  int strip = 0, totStripAPV = 0, apv = 0, prevapv = itRaw->strip() / 128;
64 
65  for (; itRaw != itRawDigis->end(); ++itRaw) {
66 
67  strip = itRaw->strip();
68  apv = strip / 128;
69  float adc = itRaw->adc();
70  h1ADC_vs_strip_->Fill(strip, adc);
71 
72  if (prevapv != apv) {
73  if (totStripAPV > 64) {
74  NumResAPVs++;
75  }
76  prevapv = apv;
77  totStripAPV = 0;
78  }
79  if (adc > 0)
80  ++totStripAPV;
81 
82  }
83  if (totStripAPV > 64) {
84  NumResAPVs++;
85  }
86 
87  }
88 
90 
91  h1NumbadAPVsRes_->Fill(NumResAPVs);
92 
93 }
94 
95 //define this as a plug-in
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX)
Definition: DQMStore.cc:239
T getParameter(std::string const &) const
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:525
void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:418
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
SiStripBaselineValidator(const edm::ParameterSet &)
constexpr int adc(sample_type sample)
get the ADC sample (12 bits)
iterator end()
Return the off-the-end iterator.
Definition: DetSetVector.h:325
void analyze(const edm::Event &, const edm::EventSetup &) override
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
HLT enums.
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Definition: DQMStore.cc:266
iterator begin()
Return an iterator to the first DetSet.
Definition: DetSetVector.h:314
collection_type::const_iterator const_iterator
Definition: DetSet.h:32
collection_type::const_iterator const_iterator
Definition: DetSetVector.h:102
Definition: Run.h:45