CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
BtlLocalRecoHarvester.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: Validation/MtdValidation
4 // Class: BtlLocalRecoHarvester
5 //
14 #include <string>
15 
20 
23 
25 
27 public:
28  explicit BtlLocalRecoHarvester(const edm::ParameterSet& iConfig);
29  ~BtlLocalRecoHarvester() override;
30 
31  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
32 
33 protected:
35 
36 private:
38 
39  // --- Histograms
41 };
42 
43 // ------------ constructor and destructor --------------
45  : folder_(iConfig.getParameter<std::string>("folder")) {}
46 
48 
49 // ------------ endjob tasks ----------------------------
51  // --- Get the monitoring histograms
52  MonitorElement* meBtlHitLogEnergy = igetter.get(folder_ + "BtlHitLogEnergy");
53  MonitorElement* meNevents = igetter.get(folder_ + "BtlNevents");
54 
55  if (!meBtlHitLogEnergy || !meNevents) {
56  edm::LogError("BtlLocalRecoHarvester") << "Monitoring histograms not found!" << std::endl;
57  return;
58  }
59 
60  // --- Get the number of BTL crystals and the number of processed events
61  const float NBtlCrystals = BTLDetId::kCrystalsBTL;
62  const float Nevents = meNevents->getEntries();
63  const float scale = (Nevents > 0 ? 1. / (Nevents * NBtlCrystals) : 1.);
64 
65  // --- Book the cumulative histogram
66  ibook.cd(folder_);
67  meHitOccupancy_ = ibook.book1D("BtlHitOccupancy",
68  "BTL cell occupancy vs RECO hit energy;log_{10}(E_{RECO} [MeV]); Occupancy per event",
69  meBtlHitLogEnergy->getNbinsX(),
70  meBtlHitLogEnergy->getTH1()->GetXaxis()->GetXmin(),
71  meBtlHitLogEnergy->getTH1()->GetXaxis()->GetXmax());
72 
73  // --- Calculate the cumulative histogram
74  double bin_sum = meBtlHitLogEnergy->getBinContent(meBtlHitLogEnergy->getNbinsX() + 1);
75  for (int ibin = meBtlHitLogEnergy->getNbinsX(); ibin >= 1; --ibin) {
76  bin_sum += meBtlHitLogEnergy->getBinContent(ibin);
77  meHitOccupancy_->setBinContent(ibin, scale * bin_sum);
78  }
79 }
80 
81 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
84 
85  desc.add<std::string>("folder", "MTD/BTL/LocalReco/");
86 
87  descriptions.add("btlLocalRecoPostProcessor", desc);
88 }
89 
static constexpr uint32_t kCrystalsBTL
Definition: BTLDetId.h:43
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
std::string folder_
Log< level::Error, false > LogError
BtlLocalRecoHarvester(const edm::ParameterSet &iConfig)
MonitorElement * meHitOccupancy_
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
virtual double getEntries() const
get # of entries
virtual void setBinContent(int binx, double content)
set content of bin (1-D)
void add(std::string const &label, ParameterSetDescription const &psetDescription)
virtual MonitorElement * get(std::string const &fullpath) const
Definition: DQMStore.cc:697
virtual TH1 * getTH1() const
virtual int getNbinsX() const
get # of bins in X-axis
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
void dqmEndJob(DQMStore::IBooker &, DQMStore::IGetter &) override
virtual double getBinContent(int binx) const
get content of bin (1-D)