File indexing completed on 2024-04-06 12:08:47
0001
0002
0003
0004
0005
0006
0007 #include <cassert>
0008 #include <fstream>
0009 #include <iostream>
0010 #include <memory>
0011
0012
0013 #include "DQM/SiStripMonitorDigi/interface/SiStripBaselineValidator.h"
0014 #include "DQMServices/Core/interface/DQMStore.h"
0015 #include "DataFormats/Common/interface/DetSet.h"
0016 #include "FWCore/Framework/interface/ESHandle.h"
0017 #include "FWCore/Framework/interface/Frameworkfwd.h"
0018 #include "FWCore/Framework/interface/MakerMacros.h"
0019 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0020 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0021 #include "FWCore/ServiceRegistry/interface/Service.h"
0022
0023
0024 #include "TH1F.h"
0025 #include "TH2F.h"
0026 #include "TString.h"
0027
0028 class TFile;
0029
0030 using namespace edm;
0031 using namespace std;
0032
0033 SiStripBaselineValidator::SiStripBaselineValidator(const edm::ParameterSet& conf) {
0034 srcProcessedRawDigi_ = conf.getParameter<edm::InputTag>("srcProcessedRawDigi");
0035 moduleRawDigiToken_ =
0036 consumes<edm::DetSetVector<SiStripDigi> >(conf.getParameter<edm::InputTag>("srcProcessedRawDigi"));
0037 }
0038
0039 void SiStripBaselineValidator::bookHistograms(DQMStore::IBooker& ibooker,
0040 const edm::Run& run,
0041 const edm::EventSetup& es) {
0042
0043 ibooker.setCurrentFolder("SiStrip/BaselineValidator");
0044
0045 h1NumbadAPVsRes_ = ibooker.book1D("ResAPVs", ";#ResAPVs", 100, 1.0, 10001);
0046
0047 h1ADC_vs_strip_ = ibooker.book2D("ADCvsAPVs", ";ADCvsAPVs", 768, -0.5, 767.5, 1023, -0.5, 1022.5);
0048
0049 return;
0050 }
0051
0052 void SiStripBaselineValidator::analyze(const edm::Event& e, const edm::EventSetup& es) {
0053 edm::Handle<edm::DetSetVector<SiStripDigi> > moduleRawDigi;
0054 e.getByToken(moduleRawDigiToken_, moduleRawDigi);
0055 edm::DetSetVector<SiStripDigi>::const_iterator itRawDigis = moduleRawDigi->begin();
0056
0057 int NumResAPVs = 0;
0058 for (; itRawDigis != moduleRawDigi->end(); ++itRawDigis) {
0059
0060 edm::DetSet<SiStripDigi>::const_iterator itRaw = itRawDigis->begin();
0061 int strip = 0, totStripAPV = 0, apv = 0, prevapv = itRaw->strip() / 128;
0062
0063 for (; itRaw != itRawDigis->end(); ++itRaw) {
0064
0065 strip = itRaw->strip();
0066 apv = strip / 128;
0067 float adc = itRaw->adc();
0068 h1ADC_vs_strip_->Fill(strip, adc);
0069
0070 if (prevapv != apv) {
0071 if (totStripAPV > 64) {
0072 NumResAPVs++;
0073 }
0074 prevapv = apv;
0075 totStripAPV = 0;
0076 }
0077 if (adc > 0)
0078 ++totStripAPV;
0079
0080 }
0081 if (totStripAPV > 64) {
0082 NumResAPVs++;
0083 }
0084
0085 }
0086
0087
0088
0089 h1NumbadAPVsRes_->Fill(NumResAPVs);
0090
0091 }
0092
0093
0094 DEFINE_FWK_MODULE(SiStripBaselineValidator);