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 hasAbilityToProduceInLumis () const final
 
bool hasAbilityToProduceInRuns () 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 () 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
 
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)
 
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 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::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 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 35 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_.

18  :
19  currentPedestals_(nullptr),channelStatus_(nullptr){
20 
21  chStatusToExclude_= StringToEnumValue<EcalChannelStatusCode::Code>(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 183 of file ECALpedestalPCLHarvester.cc.

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

Referenced by dqmEndJob().

183  {
184 
186  dbstatusPtr = channelStatus_->getMap().find(id.rawId());
187  EcalChannelStatusCode::Code dbstatus = dbstatusPtr->getStatusCode();
188 
189  std::vector<int>::const_iterator res =
190  std::find( chStatusToExclude_.begin(), chStatusToExclude_.end(), dbstatus );
191  if ( res != chStatusToExclude_.end() ) return false;
192 
193  return true;
194 }
const self & getMap() const
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
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 208 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().

209  {
210 
211  uint32_t nAnomaliesEB =0;
212  uint32_t nAnomaliesEE =0;
213 
214  for (uint16_t i =0; i< EBDetId::kSizeForDenseIndexing; ++i) {
215 
217  const EcalPedestal& newped=* newPedestals.find(id.rawId());
218  const EcalPedestal& oldped=* oldPedestals.find(id.rawId());
219 
220  if (std::abs(newped.mean_x12 -oldped.mean_x12) > nSigma_ * oldped.rms_x12) nAnomaliesEB++;
221 
222  }
223 
224  for (uint16_t i =0; i< EEDetId::kSizeForDenseIndexing; ++i) {
225 
227  const EcalPedestal& newped=* newPedestals.find(id.rawId());
228  const EcalPedestal& oldped=* oldPedestals.find(id.rawId());
229 
230  if (std::abs(newped.mean_x12 -oldped.mean_x12) > nSigma_ * oldped.rms_x12) nAnomaliesEE++;
231 
232  }
233 
234  if (nAnomaliesEB > thresholdAnomalies_ * EBDetId::kSizeForDenseIndexing ||
235  nAnomaliesEE > thresholdAnomalies_ * EEDetId::kSizeForDenseIndexing)
236  return true;
237 
238  return false;
239 
240 
241 }
static EEDetId detIdFromDenseIndex(uint32_t din)
Definition: EEDetId.h:220
static EBDetId detIdFromDenseIndex(uint32_t di)
Definition: EBDetId.h:111
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Definition: DetId.h:18
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(), diffTreeTool::diff, dqmDir_, dqmPlots(), entriesEB_, entriesEE_, EcalCondObjectContainer< T >::find(), g6g1Pedestals_, DQMStore::IGetter::get(), MonitorElement::getEntries(), MonitorElement::getMean(), 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 
35 
36  // calculate pedestals and fill db record
37  EcalPedestals pedestals;
38  int kBarrelSize = 61200;
39  int kEndcapSize = 2*7324;
40  int skipped_channels_EB = 0;
41  int skipped_channels_EE = 0;
42 
43 
44  for (uint16_t i =0; i< EBDetId::kSizeForDenseIndexing; ++i) {
45  std::string hname = dqmDir_+"/EB/"+std::to_string(int(i/100))+"/eb_" + std::to_string(i);
46  MonitorElement* ch= igetter_.get(hname);
47  double mean = ch->getMean();
48  double rms = ch->getRMS();
49  entriesEB_[i] = ch->getEntries();
50 
52  EcalPedestal ped;
53  EcalPedestal oldped =* currentPedestals_->find(id.rawId());
54  EcalPedestal g6g1ped=* g6g1Pedestals_->find(id.rawId());
55 
56  ped.mean_x12=mean;
57  ped.rms_x12=rms;
58 
59  float diff = std::abs(mean-oldped.mean_x12);
60 
61  // if bad channel or low stat skip or the difference is too large wrt to previous record
62  if(ch->getEntries()< minEntries_ || !checkStatusCode(id) || diff>threshDiffEB_){
63 
64  ped.mean_x12=oldped.mean_x12;
65  ped.rms_x12=oldped.rms_x12;
66 
67  skipped_channels_EB++;
68 
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 
81  for (uint16_t i =0; i< EEDetId::kSizeForDenseIndexing; ++i) {
82 
83  std::string hname = dqmDir_+"/EE/"+std::to_string(int(i/100))+"/ee_" + std::to_string(i);
84 
85  MonitorElement* ch= igetter_.get(hname);
86  double mean = ch->getMean();
87  double rms = ch->getRMS();
88  entriesEE_[i] = ch->getEntries();
89 
91  EcalPedestal ped;
92  EcalPedestal oldped = *currentPedestals_->find(id.rawId());
93  EcalPedestal g6g1ped= *g6g1Pedestals_->find(id.rawId());
94 
95  ped.mean_x12=mean;
96  ped.rms_x12=rms;
97 
98  float diff = std::abs(mean-oldped.mean_x12);
99 
100  // if bad channel or low stat skip or the difference is too large wrt to previous record
101  if(ch->getEntries()< minEntries_ || !checkStatusCode(id)|| diff>threshDiffEE_){
102 
103  ped.mean_x12=oldped.mean_x12;
104  ped.rms_x12=oldped.rms_x12;
105 
106  skipped_channels_EE++;
107  }
108 
109  // copy g6 and g1 pedestals from corresponding record
110  ped.mean_x6= g6g1ped.mean_x6;
111  ped.rms_x6 = g6g1ped.rms_x6;
112  ped.mean_x1= g6g1ped.mean_x1;
113  ped.rms_x1 = g6g1ped.rms_x1;
114 
115  pedestals.setValue(id.rawId(),ped);
116  }
117 
118 
119  bool enough_stat = false;
120  if ((kBarrelSize-skipped_channels_EB)/kBarrelSize > threshChannelsAnalyzed_ &&
121  (kEndcapSize-skipped_channels_EE)/kEndcapSize > threshChannelsAnalyzed_ ){
122 
123  enough_stat=true;
124  }
125 
126 
127  if(enough_stat){
128 
129  dqmPlots(pedestals, ibooker_);
130 
131  }
132 
133  // check if there are large variations wrt exisiting pedstals
134 
135  if (checkAnomalies_){
136  if (checkVariation(*currentPedestals_, pedestals)) {
137  edm::LogError("Large Variations found wrt to old pedestals, no file created");
138  return;
139  }
140  }
141 
142  // write out pedestal record
144 
145  if ( !poolDbService.isAvailable() ) {
146  throw std::runtime_error("PoolDBService required.");
147  }
148  if (enough_stat)
149  poolDbService->writeOne( &pedestals, poolDbService->currentTime(),"EcalPedestalsRcd" );
150 
151 }
static EEDetId detIdFromDenseIndex(uint32_t din)
Definition: EEDetId.h:220
int entriesEE_[EEDetId::kSizeForDenseIndexing]
MonitorElement * get(const std::string &path)
Definition: DQMStore.cc:307
bool checkVariation(const EcalPedestalsMap &oldPedestals, const EcalPedestalsMap &newPedestals)
double getMean(int axis=1) const
get mean value of histogram along x, y or z axis (axis=1, 2, 3 respectively)
static EBDetId detIdFromDenseIndex(uint32_t di)
Definition: EBDetId.h:111
void setValue(const uint32_t id, const Item &item)
const EcalPedestals * g6g1Pedestals_
bool isAvailable() const
Definition: Service.h:46
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:18
double getEntries() const
get # of entries
bool checkStatusCode(const DetId &id)
const EcalPedestals * currentPedestals_
double getRMS(int axis=1) const
get RMS of histogram along x, y or z axis (axis=1, 2, 3 respectively)
const_iterator find(uint32_t rawId) const
int entriesEB_[EBDetId::kSizeForDenseIndexing]
void ECALpedestalPCLHarvester::dqmPlots ( const EcalPedestals newpeds,
DQMStore::IBooker ibooker 
)
private

Definition at line 246 of file ECALpedestalPCLHarvester.cc.

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

Referenced by dqmEndJob().

246  {
247 
248  ibooker.cd();
249  ibooker.setCurrentFolder(dqmDir_+"/Summary");
250 
251  MonitorElement * pmeb = ibooker.book2D("meaneb","Pedestal Means EB",360, 1., 361., 171, -85., 86.);
252  MonitorElement * preb = ibooker.book2D("rmseb","Pedestal RMS EB ",360, 1., 361., 171, -85., 86.);
253 
254  MonitorElement * pmeep = ibooker.book2D("meaneep","Pedestal Means EEP",100,1,101,100,1,101);
255  MonitorElement * preep = ibooker.book2D("rmseep","Pedestal RMS EEP",100,1,101,100,1,101);
256 
257  MonitorElement * pmeem = ibooker.book2D("meaneem","Pedestal Means EEM",100,1,101,100,1,101);
258  MonitorElement * preem = ibooker.book2D("rmseem","Pedestal RMS EEM",100,1,101,100,1,101);
259 
260  MonitorElement * pmebd = ibooker.book2D("meanebdiff","Abs Rel Pedestal Means Diff EB",360,1., 361., 171, -85., 86.);
261  MonitorElement * prebd = ibooker.book2D("rmsebdiff","Abs Rel Pedestal RMS Diff E ",360, 1., 361., 171, -85., 86.);
262 
263  MonitorElement * pmeepd = ibooker.book2D("meaneepdiff","Abs Rel Pedestal Means Diff EEP",100,1,101,100,1,101);
264  MonitorElement * preepd = ibooker.book2D("rmseepdiff","Abs Rel Pedestal RMS Diff EEP",100,1,101,100,1,101);
265 
266  MonitorElement * pmeemd = ibooker.book2D("meaneemdiff","Abs Rel Pedestal Means Diff EEM",100,1,101,100,1,101);
267  MonitorElement * preemd = ibooker.book2D("rmseemdiff","Abs RelPedestal RMS Diff EEM",100,1,101,100,1,101);
268 
269  MonitorElement * poeb = ibooker.book2D("occeb","Occupancy EB",360, 1., 361., 171, -85., 86.);
270  MonitorElement * poeep = ibooker.book2D("occeep","Occupancy EEP",100,1,101,100,1,101);
271  MonitorElement * poeem = ibooker.book2D("occeem","Occupancy EEM",100,1,101,100,1,101);
272 
273 
274  MonitorElement * hdiffeb = ibooker.book1D("diffeb","Pedestal Differences EB",100,-2.5,2.5);
275  MonitorElement * hdiffee = ibooker.book1D("diffee","Pedestal Differences EE",100,-2.5,2.5);
276 
278 
280  float mean= newpeds[di].mean_x12;
281  float rms = newpeds[di].rms_x12;
282 
283  float cmean = (*currentPedestals_)[di].mean_x12;
284  float crms = (*currentPedestals_)[di].rms_x12;
285 
286  if (!isGood(di) ) continue; // only good channels are plotted
287 
288  pmeb->Fill(di.iphi(),di.ieta(),mean);
289  preb->Fill(di.iphi(),di.ieta(),rms);
290  if (cmean) pmebd->Fill(di.iphi(),di.ieta(),std::abs(mean-cmean)/cmean);
291  if (crms) prebd->Fill(di.iphi(),di.ieta(),std::abs(rms-crms)/crms);
292  poeb->Fill(di.iphi(),di.ieta(),entriesEB_[hash]);
293  hdiffeb->Fill(mean-cmean);
294  }
295 
296 
298 
300  float mean= newpeds[di].mean_x12;
301  float rms = newpeds[di].rms_x12;
302  float cmean = (*currentPedestals_)[di].mean_x12;
303  float crms = (*currentPedestals_)[di].rms_x12;
304 
305  if (!isGood(di) ) continue; // only good channels are plotted
306 
307  if (di.zside() >0){
308  pmeep->Fill(di.ix(),di.iy(),mean);
309  preep->Fill(di.ix(),di.iy(),rms);
310  poeep->Fill(di.ix(),di.iy(),entriesEE_[hash]);
311  if (cmean) pmeepd->Fill(di.ix(),di.iy(),std::abs(mean-cmean)/cmean);
312  if (crms) preepd->Fill(di.ix(),di.iy(),std::abs(rms-crms)/crms);
313  } else{
314  pmeem->Fill(di.ix(),di.iy(),mean);
315  preem->Fill(di.ix(),di.iy(),rms);
316  if (cmean) pmeemd->Fill(di.ix(),di.iy(),std::abs(mean-cmean)/cmean);
317  poeem->Fill(di.ix(),di.iy(),entriesEE_[hash]);
318  if (crms)preemd->Fill(di.ix(),di.iy(),std::abs(rms-crms)/crms);
319  }
320  hdiffee->Fill(mean-cmean);
321 
322  }
323 
324 
325 
326 }
static EEDetId detIdFromDenseIndex(uint32_t din)
Definition: EEDetId.h:220
int ix() const
Definition: EEDetId.h:76
int entriesEE_[EEDetId::kSizeForDenseIndexing]
void Fill(long long x)
int iphi() const
get the crystal iphi
Definition: EBDetId.h:53
static EBDetId detIdFromDenseIndex(uint32_t di)
Definition: EBDetId.h:111
int zside() const
Definition: EEDetId.h:70
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:118
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int iy() const
Definition: EEDetId.h:82
int ieta() const
get the crystal ieta
Definition: EBDetId.h:51
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:279
MonitorElement * book2D(Args &&...args)
Definition: DQMStore.h:136
int entriesEB_[EBDetId::kSizeForDenseIndexing]
void ECALpedestalPCLHarvester::endRun ( edm::Run const &  run,
edm::EventSetup const &  isetup 
)
override

Definition at line 166 of file ECALpedestalPCLHarvester.cc.

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

166  {
167 
168 
170  isetup.get<EcalChannelStatusRcd>().get(chStatus);
171  channelStatus_=chStatus.product();
172 
174  isetup.get<EcalPedestalsRcd>().get(peds);
175  currentPedestals_=peds.product();
176 
178  isetup.get<EcalPedestalsRcd>().get(labelG6G1_,g6g1peds);
179  g6g1Pedestals_=g6g1peds.product();
180 
181 }
const EcalPedestals * g6g1Pedestals_
const EcalChannelStatus * channelStatus_
const EcalPedestals * currentPedestals_
T const * product() const
Definition: ESHandle.h:86
void ECALpedestalPCLHarvester::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

Definition at line 158 of file ECALpedestalPCLHarvester.cc.

References edm::ConfigurationDescriptions::addDefault(), and edm::ParameterSetDescription::setUnknown().

158  {
159 
161  desc.setUnknown();
162  descriptions.addDefault(desc);
163 }
void addDefault(ParameterSetDescription const &psetDescription)
bool ECALpedestalPCLHarvester::isGood ( const DetId id)
private

Definition at line 196 of file ECALpedestalPCLHarvester.cc.

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

Referenced by dqmPlots().

196  {
197 
199  dbstatusPtr = channelStatus_->getMap().find(id.rawId());
200  if (dbstatusPtr == channelStatus_->getMap().end())
201  edm::LogError("Invalid DetId supplied");
202  EcalChannelStatusCode::Code dbstatus = dbstatusPtr->getStatusCode();
203  if (dbstatus ==0 ) return true;
204  return false;
205 }
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 50 of file ECALpedestalPCLHarvester.h.

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

bool ECALpedestalPCLHarvester::checkAnomalies_
private

Definition at line 61 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

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

Definition at line 55 of file ECALpedestalPCLHarvester.h.

Referenced by checkStatusCode(), and ECALpedestalPCLHarvester().

const EcalPedestals* ECALpedestalPCLHarvester::currentPedestals_
private

Definition at line 48 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and endRun().

std::string ECALpedestalPCLHarvester::dqmDir_
private

Definition at line 64 of file ECALpedestalPCLHarvester.h.

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

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

Definition at line 59 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and dqmPlots().

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

Definition at line 60 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and dqmPlots().

const EcalPedestals* ECALpedestalPCLHarvester::g6g1Pedestals_
private

Definition at line 49 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and endRun().

std::string ECALpedestalPCLHarvester::labelG6G1_
private

Definition at line 65 of file ECALpedestalPCLHarvester.h.

Referenced by ECALpedestalPCLHarvester(), and endRun().

int ECALpedestalPCLHarvester::minEntries_
private

Definition at line 56 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

double ECALpedestalPCLHarvester::nSigma_
private

Definition at line 62 of file ECALpedestalPCLHarvester.h.

Referenced by checkVariation(), and ECALpedestalPCLHarvester().

float ECALpedestalPCLHarvester::threshChannelsAnalyzed_
private

Definition at line 68 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

float ECALpedestalPCLHarvester::threshDiffEB_
private

Definition at line 66 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

float ECALpedestalPCLHarvester::threshDiffEE_
private

Definition at line 67 of file ECALpedestalPCLHarvester.h.

Referenced by dqmEndJob(), and ECALpedestalPCLHarvester().

double ECALpedestalPCLHarvester::thresholdAnomalies_
private

Definition at line 63 of file ECALpedestalPCLHarvester.h.

Referenced by checkVariation(), and ECALpedestalPCLHarvester().