File indexing completed on 2024-09-07 04:37:24
0001 #ifndef ConfigurableAnalysis_ConfigurableHisto_H
0002 #define ConfigurableAnalysis_ConfigurableHisto_H
0003
0004 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0005 #include "CommonTools/Utils/interface/TFileDirectory.h"
0006
0007 #include "TH1.h"
0008 #include "TH1F.h"
0009 #include "TH2F.h"
0010 #include "TProfile.h"
0011 #include "THStack.h"
0012
0013 #include "PhysicsTools/UtilAlgos/interface/VariableHelper.h"
0014 #include "PhysicsTools/UtilAlgos/interface/CachingVariable.h"
0015
0016 class ConfigurableAxis {
0017 public:
0018 ConfigurableAxis() {}
0019 ConfigurableAxis(const edm::ParameterSet& par) : nBin_(0), Min_(0), Max_(0), Label_("") {
0020 Label_ = par.getParameter<std::string>("Label");
0021
0022 if (par.exists("nBins")) {
0023 nBin_ = par.getParameter<unsigned int>("nBins");
0024 Min_ = par.getParameter<double>("Min");
0025 Max_ = par.getParameter<double>("Max");
0026 } else {
0027 if (par.exists("vBins"))
0028 vBins_ = par.getParameter<std::vector<double> >("vBins");
0029 else {
0030 Min_ = par.getParameter<double>("Min");
0031 Max_ = par.getParameter<double>("Max");
0032 }
0033 }
0034 }
0035 bool variableSize() { return !vBins_.empty(); }
0036 unsigned int nBin() {
0037 if (variableSize())
0038 return vBins_.size() - 1;
0039 else
0040 return nBin_;
0041 }
0042 double Min() {
0043 if (variableSize())
0044 return vBins_.front();
0045 else
0046 return Min_;
0047 }
0048 double Max() {
0049 if (variableSize())
0050 return vBins_.back();
0051 else
0052 return Max_;
0053 }
0054 const std::string& Label() { return Label_; }
0055 const double* xBins() {
0056 if (!vBins_.empty())
0057 return &(vBins_.front());
0058 else
0059 return nullptr;
0060 }
0061
0062 private:
0063 std::vector<double> vBins_;
0064 unsigned int nBin_;
0065 double Min_;
0066 double Max_;
0067 std::string Label_;
0068 };
0069
0070 class ConfigurableHisto {
0071 public:
0072 enum HType { h1, h2, prof };
0073 ConfigurableHisto(HType t, std::string name, edm::ParameterSet& iConfig)
0074 : type_(t), h_(nullptr), name_(name), conf_(iConfig), x_(nullptr), y_(nullptr), z_(nullptr), w_(nullptr) {}
0075
0076 virtual ~ConfigurableHisto() {}
0077
0078 virtual ConfigurableHisto* clone() const { return new ConfigurableHisto(*this); }
0079
0080 virtual void book(TFileDirectory* dir) {
0081 std::string title = conf_.getParameter<std::string>("title");
0082 edm::ParameterSet xAxisPSet = conf_.getParameter<edm::ParameterSet>("xAxis");
0083 ConfigurableAxis xAxis(xAxisPSet);
0084 x_ = edm::Service<VariableHelperService>()->get().variable(xAxisPSet.getParameter<std::string>("var"));
0085
0086 std::string yLabel = "";
0087 bool yVBin = false;
0088 ConfigurableAxis yAxis;
0089 if (conf_.exists("yAxis")) {
0090 edm::ParameterSet yAxisPSet = conf_.getParameter<edm::ParameterSet>("yAxis");
0091 yAxis = ConfigurableAxis(yAxisPSet);
0092 yLabel = yAxis.Label();
0093
0094 if (yAxisPSet.exists("var"))
0095 y_ = edm::Service<VariableHelperService>()->get().variable(yAxisPSet.getParameter<std::string>("var"));
0096 yVBin = yAxis.variableSize();
0097 }
0098
0099 if (conf_.exists("zAxis")) {
0100 throw;
0101 }
0102
0103 bool xVBin = xAxis.variableSize();
0104
0105 if (type() == h1) {
0106
0107 if (xVBin)
0108 h_ = dir->make<TH1F>(name_.c_str(), title.c_str(), xAxis.nBin(), xAxis.xBins());
0109 else
0110 h_ = dir->make<TH1F>(name_.c_str(), title.c_str(), xAxis.nBin(), xAxis.Min(), xAxis.Max());
0111 } else if (type() == h2) {
0112
0113 if (xVBin) {
0114 if (yVBin)
0115 h_ = dir->make<TH2F>(name_.c_str(), title.c_str(), xAxis.nBin(), xAxis.xBins(), yAxis.nBin(), yAxis.xBins());
0116 else
0117 h_ = dir->make<TH2F>(
0118 name_.c_str(), title.c_str(), xAxis.nBin(), xAxis.xBins(), yAxis.nBin(), yAxis.Min(), yAxis.Max());
0119 } else {
0120 if (yVBin)
0121 h_ = dir->make<TH2F>(
0122 name_.c_str(), title.c_str(), xAxis.nBin(), xAxis.Min(), xAxis.Max(), yAxis.nBin(), yAxis.xBins());
0123 else
0124 h_ = dir->make<TH2F>(name_.c_str(),
0125 title.c_str(),
0126 xAxis.nBin(),
0127 xAxis.Min(),
0128 xAxis.Max(),
0129 yAxis.nBin(),
0130 yAxis.Min(),
0131 yAxis.Max());
0132 }
0133 } else if (type() == prof) {
0134
0135 std::string pFopt = "";
0136 if (conf_.exists("Option"))
0137 pFopt = conf_.getParameter<std::string>("Option");
0138
0139 if (xVBin)
0140 h_ = dir->make<TProfile>(
0141 name_.c_str(), title.c_str(), xAxis.nBin(), xAxis.xBins(), yAxis.Min(), yAxis.Max(), pFopt.c_str());
0142 else
0143 h_ = dir->make<TProfile>(name_.c_str(),
0144 title.c_str(),
0145 xAxis.nBin(),
0146 xAxis.Min(),
0147 xAxis.Max(),
0148 yAxis.Min(),
0149 yAxis.Max(),
0150 pFopt.c_str());
0151
0152 } else {
0153 edm::LogError("ConfigurableHisto") << "cannot book: " << name_ << "\n" << conf_.dump();
0154 throw;
0155 }
0156
0157
0158 h_->GetXaxis()->SetTitle(xAxis.Label().c_str());
0159 h_->SetYTitle(yLabel.c_str());
0160
0161 if (conf_.exists("weight")) {
0162 w_ = edm::Service<VariableHelperService>()->get().variable(conf_.getParameter<std::string>("weight"));
0163 }
0164 }
0165
0166 virtual void fill(const edm::Event& iEvent) {
0167 if (!h_)
0168 return;
0169
0170 double weight = 1.0;
0171 if (w_) {
0172 if (w_->compute(iEvent))
0173 weight = (*w_)(iEvent);
0174 else {
0175 edm::LogInfo("ConfigurableHisto") << "could not compute the weight for: " << name_ << " with config:\n"
0176 << conf_.dump() << " default to 1.0";
0177 }
0178 }
0179
0180 TProfile* pcast(nullptr);
0181 TH2* h2cast(nullptr);
0182 switch (type_) {
0183 case h1:
0184 if (!h_)
0185 throw;
0186 if (x_->compute(iEvent))
0187 h_->Fill((*x_)(iEvent), weight);
0188 else {
0189 edm::LogInfo("ConfigurableHisto") << "could not fill: " << name_ << " with config:\n" << conf_.dump();
0190 }
0191 break;
0192 case prof:
0193 pcast = dynamic_cast<TProfile*>(h_);
0194 if (!pcast)
0195 throw;
0196 if (x_->compute(iEvent) && y_->compute(iEvent))
0197 pcast->Fill((*x_)(iEvent), (*y_)(iEvent), weight);
0198 else {
0199 edm::LogInfo("ConfigurableHisto") << "could not fill: " << name_ << " with config:\n" << conf_.dump();
0200 }
0201 break;
0202 case h2:
0203 h2cast = dynamic_cast<TH2*>(h_);
0204 if (!h2cast)
0205 throw;
0206 if (x_->compute(iEvent) && y_->compute(iEvent))
0207 h2cast->Fill((*x_)(iEvent), (*y_)(iEvent), weight);
0208 else {
0209 edm::LogInfo("ConfigurableHisto") << "could not fill: " << name_ << " with config:\n" << conf_.dump();
0210 }
0211 break;
0212 }
0213 }
0214 const HType& type() { return type_; }
0215
0216 void complete() {}
0217 TH1* h() { return h_; }
0218
0219 protected:
0220 ConfigurableHisto(const ConfigurableHisto& master) {
0221 type_ = master.type_;
0222 h_ = nullptr;
0223 name_ = master.name_;
0224 conf_ = master.conf_;
0225 x_ = master.x_;
0226 y_ = master.y_;
0227 z_ = master.z_;
0228 w_ = master.w_;
0229 }
0230 HType type_;
0231 TH1* h_;
0232 std::string name_;
0233 edm::ParameterSet conf_;
0234
0235 const CachingVariable* x_;
0236 const CachingVariable* y_;
0237 const CachingVariable* z_;
0238 const CachingVariable* w_;
0239 };
0240
0241 class SplittingConfigurableHisto : public ConfigurableHisto {
0242 public:
0243 SplittingConfigurableHisto(HType t, std::string name, edm::ParameterSet& pset)
0244 : ConfigurableHisto(t, name, pset), splitter_(nullptr) {
0245 std::string title = pset.getParameter<std::string>("title");
0246
0247
0248 if (pset.exists("splitters")) {
0249
0250 std::vector<std::string> splitters = pset.getParameter<std::vector<std::string> >("splitters");
0251 for (unsigned int s = 0; s != splitters.size(); ++s) {
0252 const CachingVariable* v = edm::Service<VariableHelperService>()->get().variable(splitters[s]);
0253 const Splitter* splitter = dynamic_cast<const Splitter*>(v);
0254 if (!splitter) {
0255 edm::LogError("SplittingConfigurableHisto")
0256 << "for: " << name_ << " the splitting variable: " << splitters[s] << " is not a Splitter";
0257 continue;
0258 }
0259
0260
0261 std::vector<ConfigurableHisto*>& insertedHisto = subHistoMap_[splitter];
0262
0263 unsigned int mSlots = splitter->maxSlots();
0264 for (unsigned int i = 0; i != mSlots; ++i) {
0265
0266 const std::string& slabel = splitter->shortLabel(i);
0267 const std::string& label = splitter->label(i);
0268 edm::ParameterSet mPset = pset;
0269 edm::Entry e("string", title + " for " + label, true);
0270 mPset.insert(true, "title", e);
0271 insertedHisto.push_back(new ConfigurableHisto(t, name + slabel, mPset));
0272 }
0273 }
0274
0275 }
0276 else {
0277
0278
0279 const CachingVariable* v =
0280 edm::Service<VariableHelperService>()->get().variable(pset.getParameter<std::string>("splitter"));
0281 splitter_ = dynamic_cast<const Splitter*>(v);
0282 if (!splitter_) {
0283 edm::LogError("SplittingConfigurableHisto")
0284 << "for: " << name_ << " the splitting variable: " << v->name() << " is not a Splitter";
0285 } else {
0286
0287 unsigned int mSlots = splitter_->maxSlots();
0288 for (unsigned int i = 0; i != mSlots; i++) {
0289 const std::string& slabel = splitter_->shortLabel(i);
0290 const std::string& label = splitter_->label(i);
0291 edm::ParameterSet mPset = pset;
0292 edm::Entry e("string", title + " for " + label, true);
0293 mPset.insert(true, "title", e);
0294 subHistos_.push_back(new ConfigurableHisto(t, name + slabel, mPset));
0295 }
0296 }
0297 }
0298 }
0299
0300 void book(TFileDirectory* dir) override {
0301
0302 ConfigurableHisto::book(dir);
0303
0304 if (!subHistoMap_.empty()) {
0305 SubHistoMap::iterator i = subHistoMap_.begin();
0306 SubHistoMap::iterator i_end = subHistoMap_.end();
0307 for (; i != i_end; ++i) {
0308 for (unsigned int h = 0; h != i->second.size(); ++h) {
0309 i->second[h]->book(dir);
0310 }
0311
0312 std::string sName = name_ + "_" + i->first->name();
0313 std::string sTitle = "Stack histogram of " + name_ + " for splitter " + i->first->name();
0314 subHistoStacks_[i->first] = dir->make<THStack>(sName.c_str(), sTitle.c_str());
0315 }
0316 } else {
0317 for (unsigned int h = 0; h != subHistos_.size(); h++) {
0318 subHistos_[h]->book(dir);
0319 }
0320
0321 std::string sName = name_ + "_" + splitter_->name();
0322 std::string sTitle = "Stack histogram of " + name_ + " for splitter " + splitter_->name();
0323 stack_ = dir->make<THStack>(sName.c_str(), sTitle.c_str());
0324 }
0325 }
0326
0327 ConfigurableHisto* clone() const override { return new SplittingConfigurableHisto(*this); }
0328
0329 void fill(const edm::Event& e) override {
0330
0331 ConfigurableHisto::fill(e);
0332
0333 if (!subHistoMap_.empty()) {
0334 SubHistoMap::iterator i = subHistoMap_.begin();
0335 SubHistoMap::iterator i_end = subHistoMap_.end();
0336 for (; i != i_end; ++i) {
0337 const Splitter* splitter = i->first;
0338 if (!splitter)
0339 continue;
0340 if (!splitter->compute(e))
0341 continue;
0342 unsigned int slot = (unsigned int)(*splitter)(e);
0343 if (slot >= i->second.size()) {
0344 edm::LogError("SplittingConfigurableHisto")
0345 << "slot index: " << slot << " is bigger than slots size: " << i->second.size()
0346 << " from variable value: " << (*splitter)(e);
0347 continue;
0348 }
0349
0350 i->second[slot]->fill(e);
0351 }
0352 } else {
0353
0354 if (!splitter_)
0355 return;
0356 if (!splitter_->compute(e)) {
0357 return;
0358 }
0359 unsigned int slot = (unsigned int)(*splitter_)(e);
0360 if (slot >= subHistos_.size()) {
0361 edm::LogError("SplittingConfigurableHisto")
0362 << "slot index: " << slot << " is bigger than slots size: " << subHistos_.size()
0363 << " from variable value: " << (*splitter_)(e);
0364 return;
0365 }
0366 subHistos_[slot]->fill(e);
0367 }
0368 }
0369
0370 void complete() {
0371 if (!subHistoMap_.empty()) {
0372
0373 SubHistoMap::iterator i = subHistoMap_.begin();
0374 SubHistoMap::iterator i_end = subHistoMap_.end();
0375 for (; i != i_end; ++i) {
0376 for (unsigned int h = 0; h != i->second.size(); h++) {
0377
0378 subHistoStacks_[i->first]->Add(i->second[h]->h(), i->first->label(h).c_str());
0379 }
0380 }
0381
0382 } else {
0383
0384 for (unsigned int i = 0; i != subHistos_.size(); i++) {
0385 stack_->Add(subHistos_[i]->h(), splitter_->label(i).c_str());
0386 }
0387 }
0388 }
0389
0390 private:
0391 SplittingConfigurableHisto(const SplittingConfigurableHisto& master) : ConfigurableHisto(master) {
0392 splitter_ = master.splitter_;
0393 if (!master.subHistoMap_.empty()) {
0394 SubHistoMap::const_iterator i = master.subHistoMap_.begin();
0395 SubHistoMap::const_iterator i_end = master.subHistoMap_.end();
0396 for (; i != i_end; ++i) {
0397 const std::vector<ConfigurableHisto*>& masterHistos = i->second;
0398 std::vector<ConfigurableHisto*>& clonedHistos = subHistoMap_[i->first];
0399 for (unsigned int i = 0; i != masterHistos.size(); i++) {
0400 clonedHistos.push_back(masterHistos[i]->clone());
0401 }
0402 }
0403 } else {
0404 for (unsigned int i = 0; i != master.subHistos_.size(); i++) {
0405 subHistos_.push_back(master.subHistos_[i]->clone());
0406 }
0407 }
0408 }
0409
0410 typedef std::map<const Splitter*, std::vector<ConfigurableHisto*> > SubHistoMap;
0411 typedef std::map<const Splitter*, THStack*> SubHistoStacks;
0412 SubHistoStacks subHistoStacks_;
0413 SubHistoMap subHistoMap_;
0414
0415 const Splitter* splitter_;
0416 std::vector<ConfigurableHisto*> subHistos_;
0417
0418 THStack* stack_;
0419 };
0420
0421 #endif