CMS 3D CMS Logo

List of all members | Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes
ECALpedestalPCLHarvester Class Reference

#include <ECALpedestalPCLHarvester.h>

Inheritance diagram for ECALpedestalPCLHarvester:
DQMEDHarvester edm::one::EDProducer< edm::Accumulator, edm::EndLuminosityBlockProducer, edm::EndRunProducer, edm::one::WatchLuminosityBlocks, edm::one::WatchRuns, edm::one::SharedResources > edm::one::EDProducerBase edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Public Member Functions

 ECALpedestalPCLHarvester (const edm::ParameterSet &ps)
 
void endRun (edm::Run const &run, edm::EventSetup const &isetup) override
 
- Public Member Functions inherited from DQMEDHarvester
void accumulate (edm::Event const &ev, edm::EventSetup const &es) final
 
void beginLuminosityBlock (edm::LuminosityBlock const &, edm::EventSetup const &) final
 
void beginRun (edm::Run const &, edm::EventSetup const &) override
 
 DQMEDHarvester ()
 
virtual void dqmEndLuminosityBlock (DQMStore::IBooker &, DQMStore::IGetter &, edm::LuminosityBlock const &, edm::EventSetup const &)
 
void endJob () final
 
void endLuminosityBlock (edm::LuminosityBlock const &, edm::EventSetup const &) final
 
void endLuminosityBlockProduce (edm::LuminosityBlock &, edm::EventSetup const &) final
 
void endRun (edm::Run const &, edm::EventSetup const &) override
 
void endRunProduce (edm::Run &run, edm::EventSetup const &setup) override
 
 ~DQMEDHarvester () override=default
 
- Public Member Functions inherited from edm::one::EDProducer< edm::Accumulator, edm::EndLuminosityBlockProducer, edm::EndRunProducer, edm::one::WatchLuminosityBlocks, edm::one::WatchRuns, edm::one::SharedResources >
 EDProducer ()=default
 
SerialTaskQueueglobalLuminosityBlocksQueue () final
 
SerialTaskQueueglobalRunsQueue () final
 
bool hasAbilityToProduceInBeginLumis () const final
 
bool hasAbilityToProduceInBeginRuns () const final
 
bool hasAbilityToProduceInEndLumis () const final
 
bool hasAbilityToProduceInEndRuns () const final
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
- Public Member Functions inherited from edm::one::EDProducerBase
 EDProducerBase ()
 
ModuleDescription const & moduleDescription () const
 
bool wantsStreamLuminosityBlocks () const
 
bool wantsStreamRuns () const
 
 ~EDProducerBase () override
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
std::vector< edm::ProductResolverIndex > const & indiciesForPutProducts (BranchType iBranchType) const
 
 ProducerBase ()
 
std::vector< edm::ProductResolverIndex > const & putTokenIndexToProductResolverIndex () const
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription const &)> registrationCallback () const
 used by the fwk to register list of products More...
 
void resolvePutIndicies (BranchType iBranchType, ModuleToResolverIndicies const &iIndicies, std::string const &moduleLabel)
 
 ~ProducerBase () noexcept(false) override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ESProxyIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
void updateLookup (eventsetup::ESRecordsToProxyIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Static Public Member Functions

static void fillDescriptions (edm::ConfigurationDescriptions &descriptions)
 
- Static Public Member Functions inherited from edm::one::EDProducerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 

Private Member Functions

bool checkStatusCode (const DetId &id)
 
bool checkVariation (const EcalPedestalsMap &oldPedestals, const EcalPedestalsMap &newPedestals)
 
void dqmEndJob (DQMStore::IBooker &ibooker_, DQMStore::IGetter &igetter_) override
 
void dqmPlots (const EcalPedestals &newpeds, DQMStore::IBooker &ibooker)
 
bool isGood (const DetId &id)
 

Private Attributes

const EcalChannelStatuschannelStatus_
 
bool checkAnomalies_
 
std::vector< int > chStatusToExclude_
 
const EcalPedestalscurrentPedestals_
 
std::string dqmDir_
 
int entriesEB_ [EBDetId::kSizeForDenseIndexing]
 
int entriesEE_ [EEDetId::kSizeForDenseIndexing]
 
const EcalPedestalsg6g1Pedestals_
 
std::string labelG6G1_
 
int minEntries_
 
double nSigma_
 
float threshChannelsAnalyzed_
 
float threshDiffEB_
 
float threshDiffEE_
 
double thresholdAnomalies_
 

Additional Inherited Members

- Public Types inherited from DQMEDHarvester
typedef dqm::harvesting::DQMStore DQMStore
 
typedef dqm::harvesting::MonitorElement MonitorElement
 
- Public Types inherited from edm::one::EDProducerBase
typedef EDProducerBase ModuleType
 
- Public Types inherited from edm::ProducerBase
using ModuleToResolverIndicies = std::unordered_multimap< std::string, std::tuple< edm::TypeID const *, const char *, edm::ProductResolverIndex >>
 
typedef ProductRegistryHelper::TypeLabelList TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Protected Member Functions inherited from edm::ProducerBase
ProducesCollector producesCollector ()
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes ()
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes (ESInputTag const &tag)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
- Protected Attributes inherited from DQMEDHarvester
edm::EDPutTokenT< DQMTokenlumiToken_
 
edm::EDPutTokenT< DQMTokenrunToken_
 

Detailed Description

Description: Fill DQM histograms with pedestals. Intended to be used on laser data from the TestEnablesEcalHcal dataset

Definition at line 33 of file ECALpedestalPCLHarvester.h.

Constructor & Destructor Documentation

ECALpedestalPCLHarvester::ECALpedestalPCLHarvester ( const edm::ParameterSet ps)
explicit

Definition at line 18 of file ECALpedestalPCLHarvester.cc.

References checkAnomalies_, chStatusToExclude_, dqmDir_, edm::ParameterSet::getParameter(), labelG6G1_, minEntries_, nSigma_, AlCaHLTBitMon_QueryRunRegistry::string, threshChannelsAnalyzed_, threshDiffEB_, threshDiffEE_, and thresholdAnomalies_.

19  : currentPedestals_(nullptr), channelStatus_(nullptr) {
20  chStatusToExclude_ = StringToEnumValue<EcalChannelStatusCode::Code>(
21  ps.getParameter<std::vector<std::string> >("ChannelStatusToExclude"));
22  minEntries_ = ps.getParameter<int>("MinEntries");
23  checkAnomalies_ = ps.getParameter<bool>("checkAnomalies");
24  nSigma_ = ps.getParameter<double>("nSigma");
25  thresholdAnomalies_ = ps.getParameter<double>("thresholdAnomalies");
26  dqmDir_ = ps.getParameter<std::string>("dqmDir");
27  labelG6G1_ = ps.getParameter<std::string>("labelG6G1");
28  threshDiffEB_ = ps.getParameter<double>("threshDiffEB");
29  threshDiffEE_ = ps.getParameter<double>("threshDiffEE");
30  threshChannelsAnalyzed_ = ps.getParameter<double>("threshChannelsAnalyzed");
31 }
T getParameter(std::string const &) const
std::vector< int > chStatusToExclude_
const EcalChannelStatus * channelStatus_
const EcalPedestals * currentPedestals_

Member Function Documentation

bool ECALpedestalPCLHarvester::checkStatusCode ( const DetId id)
private

Definition at line 171 of file ECALpedestalPCLHarvester.cc.

References channelStatus_, chStatusToExclude_, spr::find(), EcalCondObjectContainer< T >::find(), and EcalCondObjectContainer< T >::getMap().

Referenced by dqmEndJob().

171  {
173  dbstatusPtr = channelStatus_->getMap().find(id.rawId());
174  EcalChannelStatusCode::Code dbstatus = dbstatusPtr->getStatusCode();
175 
176  std::vector<int>::const_iterator res = std::find(chStatusToExclude_.begin(), chStatusToExclude_.end(), dbstatus);
177  if (res != chStatusToExclude_.end())
178  return false;
179 
180  return true;
181 }
const self & getMap() const
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
Definition: Electron.h:6
std::vector< int > chStatusToExclude_
const EcalChannelStatus * channelStatus_
std::vector< Item >::const_iterator const_iterator
const_iterator find(uint32_t rawId) const
bool ECALpedestalPCLHarvester::checkVariation ( const EcalPedestalsMap oldPedestals,
const EcalPedestalsMap newPedestals 
)
private

Definition at line 194 of file ECALpedestalPCLHarvester.cc.

References funct::abs(), EBDetId::detIdFromDenseIndex(), EEDetId::detIdFromDenseIndex(), EcalCondObjectContainer< T >::find(), mps_fire::i, EBDetId::kSizeForDenseIndexing, EEDetId::kSizeForDenseIndexing, EcalPedestal::mean_x12, nSigma_, EcalPedestal::rms_x12, and thresholdAnomalies_.

Referenced by dqmEndJob().

195  {
196  uint32_t nAnomaliesEB = 0;
197  uint32_t nAnomaliesEE = 0;
198 
199  for (uint16_t i = 0; i < EBDetId::kSizeForDenseIndexing; ++i) {
201  const EcalPedestal& newped = *newPedestals.find(id.rawId());
202  const EcalPedestal& oldped = *oldPedestals.find(id.rawId());
203 
204  if (std::abs(newped.mean_x12 - oldped.mean_x12) > nSigma_ * oldped.rms_x12)
205  nAnomaliesEB++;
206  }
207 
208  for (uint16_t i = 0; i < EEDetId::kSizeForDenseIndexing; ++i) {
210  const EcalPedestal& newped = *newPedestals.find(id.rawId());
211  const EcalPedestal& oldped = *oldPedestals.find(id.rawId());
212 
213  if (std::abs(newped.mean_x12 - oldped.mean_x12) > nSigma_ * oldped.rms_x12)
214  nAnomaliesEE++;
215  }
216 
217  if (nAnomaliesEB > thresholdAnomalies_ * EBDetId::kSizeForDenseIndexing ||
218  nAnomaliesEE > thresholdAnomalies_ * EEDetId::kSizeForDenseIndexing)
219  return true;
220 
221  return false;
222 }
static EEDetId detIdFromDenseIndex(uint32_t din)
Definition: EEDetId.h:220
static EBDetId detIdFromDenseIndex(uint32_t di)
Definition: EBDetId.h:107
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Definition: DetId.h:17
const_iterator find(uint32_t rawId) const
void ECALpedestalPCLHarvester::dqmEndJob ( DQMStore::IBooker ibooker_,
DQMStore::IGetter igetter_ 
)
overrideprivatevirtual

Implements DQMEDHarvester.

Definition at line 33 of file ECALpedestalPCLHarvester.cc.

References funct::abs(), checkAnomalies_, checkStatusCode(), checkVariation(), currentPedestals_, cond::service::PoolDBOutputService::currentTime(), EBDetId::detIdFromDenseIndex(), EEDetId::detIdFromDenseIndex(), change_name::diff, dqmDir_, dqmPlots(), entriesEB_, entriesEE_, EcalCondObjectContainer< T >::find(), g6g1Pedestals_, dqm::dqmstoreimpl::DQMStore::IGetter::get(), dqm::impl::MonitorElement::getEntries(), dqm::impl::MonitorElement::getMean(), dqm::impl::MonitorElement::getRMS(), mps_fire::i, edm::Service< T >::isAvailable(), EBDetId::kSizeForDenseIndexing, EEDetId::kSizeForDenseIndexing, SiStripPI::mean, EcalPedestal::mean_x1, EcalPedestal::mean_x12, EcalPedestal::mean_x6, minEntries_, SiStripPI::rms, EcalPedestal::rms_x1, EcalPedestal::rms_x12, EcalPedestal::rms_x6, EcalCondObjectContainer< T >::setValue(), AlCaHLTBitMon_QueryRunRegistry::string, threshChannelsAnalyzed_, threshDiffEB_, threshDiffEE_, and cond::service::PoolDBOutputService::writeOne().

33  {
34  // calculate pedestals and fill db record
35  EcalPedestals pedestals;
36  float kBarrelSize = 61200;
37  float kEndcapSize = 2 * 7324;
38  float skipped_channels_EB = 0;
39  float skipped_channels_EE = 0;
40 
41  for (uint16_t i = 0; i < EBDetId::kSizeForDenseIndexing; ++i) {
42  std::string hname = dqmDir_ + "/EB/" + std::to_string(int(i / 100)) + "/eb_" + std::to_string(i);
43  MonitorElement* ch = igetter_.get(hname);
44  if (ch == nullptr) {
45  edm::LogWarning("MissingMonitorElement") << "failed to find MonitorElement " << hname;
46  entriesEB_[i] = 0;
47  continue;
48  }
49  double mean = ch->getMean();
50  double rms = ch->getRMS();
51  entriesEB_[i] = ch->getEntries();
52 
54  EcalPedestal ped;
55  EcalPedestal oldped = *currentPedestals_->find(id.rawId());
56  EcalPedestal g6g1ped = *g6g1Pedestals_->find(id.rawId());
57 
58  ped.mean_x12 = mean;
59  ped.rms_x12 = rms;
60 
61  float diff = std::abs(mean - oldped.mean_x12);
62 
63  // if bad channel or low stat skip or the difference is too large wrt to previous record
64  if (ch->getEntries() < minEntries_ || !checkStatusCode(id) || diff > threshDiffEB_) {
65  ped.mean_x12 = oldped.mean_x12;
66  ped.rms_x12 = oldped.rms_x12;
67 
68  skipped_channels_EB++;
69  }
70 
71  // copy g6 and g1 from the corressponding record
72  ped.mean_x6 = g6g1ped.mean_x6;
73  ped.rms_x6 = g6g1ped.rms_x6;
74  ped.mean_x1 = g6g1ped.mean_x1;
75  ped.rms_x1 = g6g1ped.rms_x1;
76 
77  pedestals.setValue(id.rawId(), ped);
78  }
79 
80  for (uint16_t i = 0; i < EEDetId::kSizeForDenseIndexing; ++i) {
81  std::string hname = dqmDir_ + "/EE/" + std::to_string(int(i / 100)) + "/ee_" + std::to_string(i);
82 
83  MonitorElement* ch = igetter_.get(hname);
84  if (ch == nullptr) {
85  edm::LogWarning("MissingMonitorElement") << "failed to find MonitorElement " << hname;
86  entriesEE_[i] = 0;
87  continue;
88  }
89  double mean = ch->getMean();
90  double rms = ch->getRMS();
91  entriesEE_[i] = ch->getEntries();
92 
94  EcalPedestal ped;
95  EcalPedestal oldped = *currentPedestals_->find(id.rawId());
96  EcalPedestal g6g1ped = *g6g1Pedestals_->find(id.rawId());
97 
98  ped.mean_x12 = mean;
99  ped.rms_x12 = rms;
100 
101  float diff = std::abs(mean - oldped.mean_x12);
102 
103  // if bad channel or low stat skip or the difference is too large wrt to previous record
104  if (ch->getEntries() < minEntries_ || !checkStatusCode(id) || diff > threshDiffEE_) {
105  ped.mean_x12 = oldped.mean_x12;
106  ped.rms_x12 = oldped.rms_x12;
107 
108  skipped_channels_EE++;
109  }
110 
111  // copy g6 and g1 pedestals from corresponding record
112  ped.mean_x6 = g6g1ped.mean_x6;
113  ped.rms_x6 = g6g1ped.rms_x6;
114  ped.mean_x1 = g6g1ped.mean_x1;
115  ped.rms_x1 = g6g1ped.rms_x1;
116 
117  pedestals.setValue(id.rawId(), ped);
118  }
119 
120  bool enough_stat = false;
121  if ((kBarrelSize - skipped_channels_EB) / kBarrelSize > threshChannelsAnalyzed_ &&
122  (kEndcapSize - skipped_channels_EE) / kEndcapSize > threshChannelsAnalyzed_) {
123  enough_stat = true;
124  }
125 
126  dqmPlots(pedestals, ibooker_);
127 
128  // check if there are large variations wrt exisiting pedstals
129 
130  if (checkAnomalies_) {
131  if (checkVariation(*currentPedestals_, pedestals)) {
132  edm::LogError("Large Variations found wrt to old pedestals, no file created");
133  return;
134  }
135  }
136 
137  if (!enough_stat)
138  return;
139 
140  // write out pedestal record
142 
143  if (!poolDbService.isAvailable()) {
144  throw std::runtime_error("PoolDBService required.");
145  }
146 
147  poolDbService->writeOne(&pedestals, poolDbService->currentTime(), "EcalPedestalsRcd");
148 }
static EEDetId detIdFromDenseIndex(uint32_t din)
Definition: EEDetId.h:220
int entriesEE_[EEDetId::kSizeForDenseIndexing]
bool checkVariation(const EcalPedestalsMap &oldPedestals, const EcalPedestalsMap &newPedestals)
static EBDetId detIdFromDenseIndex(uint32_t di)
Definition: EBDetId.h:107
void setValue(const uint32_t id, const Item &item)
const EcalPedestals * g6g1Pedestals_
virtual double getEntries() const
get # of entries
bool isAvailable() const
Definition: Service.h:40
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
void writeOne(T *payload, Time_t time, const std::string &recordName, bool withlogging=false)
void dqmPlots(const EcalPedestals &newpeds, DQMStore::IBooker &ibooker)
Definition: DetId.h:17
virtual double getMean(int axis=1) const
get mean value of histogram along x, y or z axis (axis=1, 2, 3 respectively)
virtual double getRMS(int axis=1) const
get RMS of histogram along x, y or z axis (axis=1, 2, 3 respectively)
bool checkStatusCode(const DetId &id)
const EcalPedestals * currentPedestals_
const_iterator find(uint32_t rawId) const
int entriesEB_[EBDetId::kSizeForDenseIndexing]
void ECALpedestalPCLHarvester::dqmPlots ( const EcalPedestals newpeds,
DQMStore::IBooker ibooker 
)
private

Definition at line 224 of file ECALpedestalPCLHarvester.cc.

References funct::abs(), dqm::dqmstoreimpl::DQMStore::IBooker::book1D(), dqm::dqmstoreimpl::DQMStore::IBooker::book2D(), dqm::dqmstoreimpl::DQMStore::IBooker::cd(), EBDetId::detIdFromDenseIndex(), EEDetId::detIdFromDenseIndex(), dqmDir_, entriesEB_, entriesEE_, dqm::impl::MonitorElement::Fill(), cond::hash, EBDetId::ieta(), EBDetId::iphi(), isGood(), EEDetId::ix(), EEDetId::iy(), EBDetId::kSizeForDenseIndexing, EEDetId::kSizeForDenseIndexing, SiStripPI::mean, SiStripPI::rms, dqm::dqmstoreimpl::DQMStore::IBooker::setCurrentFolder(), and EEDetId::zside().

Referenced by dqmEndJob().

224  {
225  ibooker.cd();
226  ibooker.setCurrentFolder(dqmDir_ + "/Summary");
227 
228  MonitorElement* pmeb = ibooker.book2D("meaneb", "Pedestal Means EB", 360, 1., 361., 171, -85., 86.);
229  MonitorElement* preb = ibooker.book2D("rmseb", "Pedestal RMS EB ", 360, 1., 361., 171, -85., 86.);
230 
231  MonitorElement* pmeep = ibooker.book2D("meaneep", "Pedestal Means EEP", 100, 1, 101, 100, 1, 101);
232  MonitorElement* preep = ibooker.book2D("rmseep", "Pedestal RMS EEP", 100, 1, 101, 100, 1, 101);
233 
234  MonitorElement* pmeem = ibooker.book2D("meaneem", "Pedestal Means EEM", 100, 1, 101, 100, 1, 101);
235  MonitorElement* preem = ibooker.book2D("rmseem", "Pedestal RMS EEM", 100, 1, 101, 100, 1, 101);
236 
237  MonitorElement* pmebd = ibooker.book2D("meanebdiff", "Abs Rel Pedestal Means Diff EB", 360, 1., 361., 171, -85., 86.);
238  MonitorElement* prebd = ibooker.book2D("rmsebdiff", "Abs Rel Pedestal RMS Diff E ", 360, 1., 361., 171, -85., 86.);
239 
240  MonitorElement* pmeepd = ibooker.book2D("meaneepdiff", "Abs Rel Pedestal Means Diff EEP", 100, 1, 101, 100, 1, 101);
241  MonitorElement* preepd = ibooker.book2D("rmseepdiff", "Abs Rel Pedestal RMS Diff EEP", 100, 1, 101, 100, 1, 101);
242 
243  MonitorElement* pmeemd = ibooker.book2D("meaneemdiff", "Abs Rel Pedestal Means Diff EEM", 100, 1, 101, 100, 1, 101);
244  MonitorElement* preemd = ibooker.book2D("rmseemdiff", "Abs RelPedestal RMS Diff EEM", 100, 1, 101, 100, 1, 101);
245 
246  MonitorElement* poeb = ibooker.book2D("occeb", "Occupancy EB", 360, 1., 361., 171, -85., 86.);
247  MonitorElement* poeep = ibooker.book2D("occeep", "Occupancy EEP", 100, 1, 101, 100, 1, 101);
248  MonitorElement* poeem = ibooker.book2D("occeem", "Occupancy EEM", 100, 1, 101, 100, 1, 101);
249 
250  MonitorElement* hdiffeb = ibooker.book1D("diffeb", "Pedestal Differences EB", 100, -2.5, 2.5);
251  MonitorElement* hdiffee = ibooker.book1D("diffee", "Pedestal Differences EE", 100, -2.5, 2.5);
252 
253  for (int hash = 0; hash < EBDetId::kSizeForDenseIndexing; ++hash) {
255  float mean = newpeds[di].mean_x12;
256  float rms = newpeds[di].rms_x12;
257 
258  float cmean = (*currentPedestals_)[di].mean_x12;
259  float crms = (*currentPedestals_)[di].rms_x12;
260 
261  if (!isGood(di))
262  continue; // only good channels are plotted
263 
264  pmeb->Fill(di.iphi(), di.ieta(), mean);
265  preb->Fill(di.iphi(), di.ieta(), rms);
266  if (cmean)
267  pmebd->Fill(di.iphi(), di.ieta(), std::abs(mean - cmean) / cmean);
268  if (crms)
269  prebd->Fill(di.iphi(), di.ieta(), std::abs(rms - crms) / crms);
270  poeb->Fill(di.iphi(), di.ieta(), entriesEB_[hash]);
271  hdiffeb->Fill(mean - cmean);
272  }
273 
274  for (int hash = 0; hash < EEDetId::kSizeForDenseIndexing; ++hash) {
276  float mean = newpeds[di].mean_x12;
277  float rms = newpeds[di].rms_x12;
278  float cmean = (*currentPedestals_)[di].mean_x12;
279  float crms = (*currentPedestals_)[di].rms_x12;
280 
281  if (!isGood(di))
282  continue; // only good channels are plotted
283 
284  if (di.zside() > 0) {
285  pmeep->Fill(di.ix(), di.iy(), mean);
286  preep->Fill(di.ix(), di.iy(), rms);
287  poeep->Fill(di.ix(), di.iy(), entriesEE_[hash]);
288  if (cmean)
289  pmeepd->Fill(di.ix(), di.iy(), std::abs(mean - cmean) / cmean);
290  if (crms)
291  preepd->Fill(di.ix(), di.iy(), std::abs(rms - crms) / crms);
292  } else {
293  pmeem->Fill(di.ix(), di.iy(), mean);
294  preem->Fill(di.ix(), di.iy(), rms);
295  if (cmean)
296  pmeemd->Fill(di.ix(), di.iy(), std::abs(mean - cmean) / cmean);
297  poeem->Fill(di.ix(), di.iy(), entriesEE_[hash]);
298  if (crms)
299  preemd->Fill(di.ix(), di.iy(), std::abs(rms - crms) / crms);
300  }
301  hdiffee->Fill(mean - cmean);
302  }
303 }
static EEDetId detIdFromDenseIndex(uint32_t din)
Definition: EEDetId.h:220
int ix() const
Definition: EEDetId.h:77
int entriesEE_[EEDetId::kSizeForDenseIndexing]
int iphi() const
get the crystal iphi
Definition: EBDetId.h:51
static EBDetId detIdFromDenseIndex(uint32_t di)
Definition: EBDetId.h:107
void Fill(long long x)
int zside() const
Definition: EEDetId.h:71
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int iy() const
Definition: EEDetId.h:83
int ieta() const
get the crystal ieta
Definition: EBDetId.h:49
int entriesEB_[EBDetId::kSizeForDenseIndexing]
void ECALpedestalPCLHarvester::endRun ( edm::Run const &  run,
edm::EventSetup const &  isetup 
)
override

Definition at line 157 of file ECALpedestalPCLHarvester.cc.

References channelStatus_, currentPedestals_, g6g1Pedestals_, edm::EventSetup::get(), labelG6G1_, and edm::ESHandle< T >::product().

157  {
159  isetup.get<EcalChannelStatusRcd>().get(chStatus);
160  channelStatus_ = chStatus.product();
161 
163  isetup.get<EcalPedestalsRcd>().get(peds);
164  currentPedestals_ = peds.product();
165 
167  isetup.get<EcalPedestalsRcd>().get(labelG6G1_, g6g1peds);
168  g6g1Pedestals_ = g6g1peds.product();
169 }
const EcalPedestals * g6g1Pedestals_
const EcalChannelStatus * channelStatus_
const EcalPedestals * currentPedestals_
T const * product() const
Definition: ESHandle.h:86
void ECALpedestalPCLHarvester::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static
bool ECALpedestalPCLHarvester::isGood ( const DetId id)
private

Definition at line 183 of file ECALpedestalPCLHarvester.cc.

References channelStatus_, EcalCondObjectContainer< T >::end(), EcalCondObjectContainer< T >::find(), and EcalCondObjectContainer< T >::getMap().

Referenced by dqmPlots().

183  {
185  dbstatusPtr = channelStatus_->getMap().find(id.rawId());
186  if (dbstatusPtr == channelStatus_->getMap().end())
187  edm::LogError("Invalid DetId supplied");
188  EcalChannelStatusCode::Code dbstatus = dbstatusPtr->getStatusCode();
189  if (dbstatus == 0)
190  return true;
191  return false;
192 }
const self & getMap() const
const EcalChannelStatus * channelStatus_
std::vector< Item >::const_iterator const_iterator
const_iterator find(uint32_t rawId) const
const_iterator end() const

Member Data Documentation

const EcalChannelStatus* ECALpedestalPCLHarvester::channelStatus_
private

Definition at line 46 of file ECALpedestalPCLHarvester.h.

Referenced by checkStatusCode(), endRun(), and isGood().

bool ECALpedestalPCLHarvester::checkAnomalies_
private

Definition at line 56 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

std::vector<int> ECALpedestalPCLHarvester::chStatusToExclude_
private

Definition at line 51 of file ECALpedestalPCLHarvester.h.

Referenced by checkStatusCode(), and ECALpedestalPCLHarvester().

const EcalPedestals* ECALpedestalPCLHarvester::currentPedestals_
private

Definition at line 44 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and endRun().

std::string ECALpedestalPCLHarvester::dqmDir_
private

Definition at line 59 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), dqmPlots(), and ECALpedestalPCLHarvester().

int ECALpedestalPCLHarvester::entriesEB_[EBDetId::kSizeForDenseIndexing]
private

Definition at line 54 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and dqmPlots().

int ECALpedestalPCLHarvester::entriesEE_[EEDetId::kSizeForDenseIndexing]
private

Definition at line 55 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and dqmPlots().

const EcalPedestals* ECALpedestalPCLHarvester::g6g1Pedestals_
private

Definition at line 45 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and endRun().

std::string ECALpedestalPCLHarvester::labelG6G1_
private

Definition at line 60 of file ECALpedestalPCLHarvester.h.

Referenced by ECALpedestalPCLHarvester(), and endRun().

int ECALpedestalPCLHarvester::minEntries_
private

Definition at line 52 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

double ECALpedestalPCLHarvester::nSigma_
private

Definition at line 57 of file ECALpedestalPCLHarvester.h.

Referenced by checkVariation(), and ECALpedestalPCLHarvester().

float ECALpedestalPCLHarvester::threshChannelsAnalyzed_
private

Definition at line 63 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

float ECALpedestalPCLHarvester::threshDiffEB_
private

Definition at line 61 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

float ECALpedestalPCLHarvester::threshDiffEE_
private

Definition at line 62 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

double ECALpedestalPCLHarvester::thresholdAnomalies_
private

Definition at line 58 of file ECALpedestalPCLHarvester.h.

Referenced by checkVariation(), and ECALpedestalPCLHarvester().