diff --git a/examples/example_files.yml b/examples/example_files.yml index 8f0bf0f..67ea2aa 100644 --- a/examples/example_files.yml +++ b/examples/example_files.yml @@ -1,10 +1,11 @@ 'data.root': type: data legend: 'Data' + pretty-name: 'Data' 'MC_sample1.root': type: mc - sample-name: "TTJets_TuneCUETP8M1_amcatnloFXFX_25ns_15c5fa1_c106444" + pretty-name: 'MC sample 1' cross-section: 245.8 generated-events: 2167 fill-color: "#D95B43" @@ -14,6 +15,7 @@ 'MC_sample2.root': type: mc + pretty-name: 'MC sample 2' cross-section: 666.3 generated-events: 2404 fill-color: '#53777A' diff --git a/external/build-external.sh b/external/build-external.sh index 1c7b1b2..eb307fa 100755 --- a/external/build-external.sh +++ b/external/build-external.sh @@ -2,10 +2,10 @@ # YAML -curl -O https://yaml-cpp.googlecode.com/files/yaml-cpp-0.5.1.tar.gz -tar xf yaml-cpp-0.5.1.tar.gz +curl -O -J -L https://github.com/jbeder/yaml-cpp/archive/release-0.5.3.tar.gz +tar xf yaml-cpp-release-0.5.3.tar.gz -cd yaml-cpp-0.5.1 +cd yaml-cpp-release-0.5.3 mkdir build cd build @@ -15,7 +15,7 @@ make -j4 make install cd ../.. -rm yaml-cpp-0.5.1.tar.gz +rm yaml-cpp-release-0.5.3.tar.gz # TCLAP curl -L "http://downloads.sourceforge.net/project/tclap/{tclap-1.2.1.tar.gz}?r=http%3A%2F%2Fsourceforge.net%2Fprojects%2Ftclap%2Ffiles%2F&ts=1431017326&use_mirror=freefr" -o "#1" diff --git a/include/TH1Plotter.h b/include/TH1Plotter.h index d971dad..e2d2b8e 100644 --- a/include/TH1Plotter.h +++ b/include/TH1Plotter.h @@ -9,7 +9,7 @@ namespace plotIt { plotter(plotIt) { } - virtual bool plot(TCanvas& c, Plot& plot); + virtual boost::optional plot(TCanvas& c, Plot& plot); virtual bool supports(TObject& object); private: diff --git a/include/colors.h b/include/colors.h new file mode 100644 index 0000000..6dea6f4 --- /dev/null +++ b/include/colors.h @@ -0,0 +1,20 @@ +#pragma once + +#include + +namespace Color { + enum Code { + RESET = 0, + FG_RED = 31, + FG_GREEN = 32, + FG_YELLOW = 33, + FG_BLUE = 34, + FG_MAGENTA = 35, + FG_CYAN = 36, + FG_WHITE = 37 + }; + + std::ostream& operator<<(std::ostream& os, Code code) { + return os << "\033[" << static_cast(code) << "m"; + } +} diff --git a/include/plotter.h b/include/plotter.h index 9e4db0f..3a402fe 100644 --- a/include/plotter.h +++ b/include/plotter.h @@ -1,6 +1,9 @@ #pragma once #include +#include + +#include class TCanvas; class TObject; @@ -15,7 +18,7 @@ namespace plotIt { } - virtual bool plot(TCanvas& c, Plot& plot) = 0; + virtual boost::optional plot(TCanvas& c, Plot& plot) = 0; virtual bool supports(TObject& object) = 0; protected: diff --git a/include/plotters.h b/include/plotters.h index d6ca853..3efc649 100644 --- a/include/plotters.h +++ b/include/plotters.h @@ -1,6 +1,9 @@ #pragma once #include +#include + +#include namespace plotIt { static std::vector> s_plotters; @@ -8,12 +11,12 @@ namespace plotIt { s_plotters.push_back(std::make_shared(plotIt)); } - bool plot(const File& file, TCanvas& c, Plot& plot) { + boost::optional plot(const File& file, TCanvas& c, Plot& plot) { for (auto& plotter: s_plotters) { if (plotter->supports(*file.object)) return plotter->plot(c, plot); } - return false; + return boost::none; } } diff --git a/include/summary.h b/include/summary.h new file mode 100644 index 0000000..c19fb73 --- /dev/null +++ b/include/summary.h @@ -0,0 +1,46 @@ +#pragma once + +#include + +#include +#include +#include + +namespace plotIt { + + struct SummaryItem { + std::string name; + + float events = 0; + float events_uncertainty = 0; + + float efficiency = 0; + float efficiency_uncertainty = 0; + }; + + class Summary { + public: + void add(const Type& type, const SummaryItem& item); + void addSystematics(const Type& type, const SummaryItem& item); + + std::vector get(const Type& type) const; + std::vector getSystematics(const Type& type) const; + + private: + std::map> m_items; + std::map> m_systematics_items; + }; + + class SummaryPrinter { + public: + virtual void print(const Summary& summary) const = 0; + }; + + class ConsoleSummaryPrinter: public SummaryPrinter { + public: + virtual void print(const Summary& summary) const override; + + private: + void printItems(const Type& type, const Summary& summary) const; + }; +} diff --git a/include/types.h b/include/types.h index 637677d..1fe5c53 100644 --- a/include/types.h +++ b/include/types.h @@ -9,6 +9,8 @@ #include #include +#include + #include #include #include @@ -21,6 +23,21 @@ namespace plotIt { DATA }; + inline std::string type_to_string(const Type& type) { + switch (type) { + case MC: + return "Simulation"; + + case SIGNAL: + return "Signal"; + + case DATA: + return "Data"; + } + + return "Unknown"; + } + enum ErrorsType { Normal = 0, Poisson = 1, @@ -59,25 +76,11 @@ namespace plotIt { VERTICAL }; - struct Summary { - float n_events = 0; - float n_events_error = 0; - - float efficiency = 0; - float efficiency_error = 0; - - void clear() { - n_events = n_events_error = efficiency = efficiency_error = 0; - } - }; - struct Systematic { std::string path; TObject* object = nullptr; std::map objects; std::shared_ptr handle; - - Summary summary; }; struct LineStyle { @@ -111,6 +114,7 @@ namespace plotIt { struct File { std::string path; + std::string pretty_name; // For MC and Signal float cross_section = 1.; @@ -130,7 +134,6 @@ namespace plotIt { std::vector systematics; int16_t order = std::numeric_limits::min(); - Summary summary; std::shared_ptr chain; std::shared_ptr handle; diff --git a/include/utilities.h b/include/utilities.h index 262d6ef..30ed01e 100644 --- a/include/utilities.h +++ b/include/utilities.h @@ -164,4 +164,14 @@ namespace plotIt { globfree(&glob_result); return ret; } + + inline std::string truncate(const std::string& str, size_t max_len) { + if (str.length() > max_len - 1) { + std::string ret = str; + ret.resize(max_len - 1); + return ret + u8"…"; + } else { + return str; + } + } } diff --git a/src/TH1Plotter.cc b/src/TH1Plotter.cc index a930e04..caf4072 100644 --- a/src/TH1Plotter.cc +++ b/src/TH1Plotter.cc @@ -1,6 +1,7 @@ #include #include +#include #include #include #include @@ -61,9 +62,11 @@ namespace plotIt { return object.InheritsFrom("TH1"); } - bool TH1Plotter::plot(TCanvas& c, Plot& plot) { + boost::optional TH1Plotter::plot(TCanvas& c, Plot& plot) { c.cd(); + Summary global_summary; + // Rescale and style histograms for (File& file: m_plotIt.getFiles()) { setHistogramStyle(file); @@ -77,20 +80,35 @@ namespace plotIt { factor *= m_plotIt.getConfiguration().scale * file.scale; } - float n_entries = h->Integral(); - file.summary.efficiency = n_entries / file.generated_events; - file.summary.efficiency_error = sqrt( (file.summary.efficiency * (1 - file.summary.efficiency)) / file.generated_events ); + SummaryItem summary; + summary.name = file.pretty_name; - file.summary.n_events = n_entries * factor; - file.summary.n_events_error = m_plotIt.getConfiguration().luminosity * file.cross_section * file.branching_ratio * file.summary.efficiency_error; - if (! m_plotIt.getConfiguration().ignore_scales) { - file.summary.n_events_error *= m_plotIt.getConfiguration().scale * file.scale; - } + double integral_error = 0; + double integral = h->IntegralAndError(h->GetXaxis()->GetFirst(), h->GetXaxis()->GetLast(), integral_error); + + summary.events = integral * factor; + summary.events_uncertainty = integral_error * factor; + + // FIXME: Probably invalid in case of weights... + + // Bayesian efficiency + // Taken from https://root.cern.ch/doc/master/TEfficiency_8cxx_source.html#l02428 + + // Use a flat prior (equivalent to Beta(1, 1)) + float alpha = 1.; + float beta = 1.; + + summary.efficiency = TEfficiency::BetaMean(integral + alpha, file.generated_events - integral + beta); + summary.efficiency_uncertainty = TEfficiency::Bayesian(file.generated_events, integral, 0.68, alpha, beta, true) - summary.efficiency; + + global_summary.add(file.type, summary); h->Scale(factor); } else { - file.summary.n_events = h->Integral(); - file.summary.n_events_error = 0; + SummaryItem summary; + summary.name = file.pretty_name; + summary.events = h->Integral(); + global_summary.add(file.type, summary); } h->Rebin(plot.rebin); @@ -225,8 +243,8 @@ namespace plotIt { for (Systematic& syst: file.systematics) { // This histogram should contains syst errors // in percent - float total_syst_error = 0; TH1* h = dynamic_cast(syst.object); + float total_syst_error = 0; for (uint32_t i = 1; i <= (uint32_t) mc_histo_syst_only->GetNbinsX(); i++) { float total_error = mc_histo_syst_only->GetBinError(i); float syst_error_percent = h->GetBinError(i); @@ -237,9 +255,13 @@ namespace plotIt { mc_histo_syst_only->SetBinError(i, std::sqrt(total_error * total_error + syst_error * syst_error)); } - syst.summary.n_events = file.summary.n_events; - syst.summary.n_events_error = std::sqrt(total_syst_error); + SummaryItem summary; + summary.name = fs::path(syst.path).stem().native(); + summary.events_uncertainty = std::sqrt(total_syst_error); + + global_summary.addSystematics(file.type, summary); } + } // Propagate syst errors to the stat + syst histogram @@ -265,7 +287,7 @@ namespace plotIt { if (!toDraw.size()) { std::cerr << "Error: nothing to draw." << std::endl; - return false; + return boost::none; }; // Sort object by minimum @@ -610,7 +632,7 @@ namespace plotIt { if (hi_pad.get()) hi_pad->cd(); - return true; + return global_summary; } void TH1Plotter::setHistogramStyle(const File& file) { diff --git a/src/plotIt.cc b/src/plotIt.cc index 3848aee..e9467f6 100644 --- a/src/plotIt.cc +++ b/src/plotIt.cc @@ -35,6 +35,8 @@ #include #include +#include + namespace fs = boost::filesystem; using std::setw; @@ -262,6 +264,12 @@ namespace plotIt { YAML::Node node = it->second; + if (node["pretty-name"]) { + file.pretty_name = node["pretty-name"].as(); + } else { + file.pretty_name = path.stem().native(); + } + if (node["type"]) { std::string type = node["type"].as(); if (type == "signal") @@ -746,79 +754,14 @@ namespace plotIt { // Create canvas TCanvas c("canvas", "canvas", m_config.width, m_config.height); - bool success = ::plotIt::plot(m_files[0], c, plot); + boost::optional summary = ::plotIt::plot(m_files[0], c, plot); - auto printSummary = [&](Type type) { - float sum_n_events = 0; - float sum_n_events_error = 0; - - auto truncate = [](const std::string& str, size_t max_len) -> std::string { - if (str.length() > max_len - 1) { - std::string ret = str; - ret.resize(max_len - 1); - return ret + u8"…"; - } else { - return str; - } - }; - - printf("%50s%18s ± %11s%4s%10s ± %10s\n", " ", "N", u8"ΔN", " ", u8"ε", u8"Δε"); - for (File& file: m_files) { - if (file.type == type) { - fs::path path(file.path); - printf("%50s%18.2f ± %10.2f%5s%3.5f%% ± %3.5f%%\n", truncate(path.stem().native(), 50).c_str(), file.summary.n_events, file.summary.n_events_error, " ", file.summary.efficiency * 100, file.summary.efficiency_error * 100); - - sum_n_events += file.summary.n_events; - sum_n_events_error += file.summary.n_events_error * file.summary.n_events_error; - } - } - - if (sum_n_events) { - float systematics = 0; - if (type == MC && m_config.luminosity_error_percent > 0) { - std::cout << "------------------------------------------" << std::endl; - std::cout << "Systematic uncertainties" << std::endl; - printf("%50s%18s ± %10.2f\n", "Luminosity", " ", sum_n_events * m_config.luminosity_error_percent); - systematics = sum_n_events * m_config.luminosity_error_percent; - } - for (File& file: m_files) { - if (file.type == type) { - for (Systematic& s: file.systematics) { - fs::path path(s.path); - printf("%50s%18s ± %10.2f\n", truncate(path.stem().native(), 50).c_str(), " ", s.summary.n_events_error); - - sum_n_events_error += s.summary.n_events_error * s.summary.n_events_error; - } - } - } - std::cout << "------------------------------------------" << std::endl; - printf("%50s%18.2f ± %10.2f\n", " ", sum_n_events, sqrt(sum_n_events_error + systematics * systematics)); - } - }; - - if (! success) + if (! summary) return false; if (m_config.verbose) { - std::cout << "Summary: " << std::endl; - - if (hasData) { - std::cout << "Data" << std::endl; - printSummary(DATA); - std::cout << std::endl; - } - - if (hasMC) { - std::cout << "MC: " << std::endl; - printSummary(MC); - std::cout << std::endl; - } - - if (hasSignal) { - std::cout << std::endl << "Signal: " << std::endl; - printSummary(SIGNAL); - std::cout << std::endl; - } + ConsoleSummaryPrinter printer; + printer.print(*summary); } if (plot.log_y) @@ -926,11 +869,6 @@ namespace plotIt { group.second.added = false; } - // Clear summary - for (auto& file: m_files) { - file.summary.clear(); - } - return true; } diff --git a/src/summary.cc b/src/summary.cc new file mode 100644 index 0000000..c7dab84 --- /dev/null +++ b/src/summary.cc @@ -0,0 +1,85 @@ +#include +#include +#include + +#include + +namespace plotIt { + void Summary::add(const Type& type, const SummaryItem& item) { + m_items[type].push_back(item); + } + + void Summary::addSystematics(const Type& type, const SummaryItem& item) { + m_systematics_items[type].push_back(item); + } + + std::vector Summary::get(const Type& type) const { + auto it = m_items.find(type); + if (it == m_items.end()) + return std::vector(); + + return it->second; + } + + std::vector Summary::getSystematics(const Type& type) const { + auto it = m_systematics_items.find(type); + if (it == m_systematics_items.end()) + return std::vector(); + + return it->second; + } + + void ConsoleSummaryPrinter::print(const Summary& summary) const { + printItems(DATA, summary); + printItems(MC, summary); + printItems(SIGNAL, summary); + } + + void ConsoleSummaryPrinter::printItems(const Type& type, const Summary& summary) const { + using namespace boost; + + std::vector nominal = summary.get(type); + std::vector systematics = summary.getSystematics(type); + + float nominal_events = 0; + float nominal_events_uncertainty = 0; + + if (nominal.size() == 0) + return; + + std::cout << Color::FG_MAGENTA << type_to_string(type) << Color::RESET << std::endl; + + // Print header + std::cout << format(u8"%|1$50| %|1$9|N ± %|1$6|ΔN") % u8" "; + if (type != DATA) { + std::cout << format(" %|1$7|ε ± %|1$6|Δε") % " "; + } + std::cout << std::endl; + + for (const auto& n: nominal) { + std::cout << Color::FG_YELLOW << format("%|50|") % truncate(n.name, 50) << Color::RESET << " " << format("%|10.2f| ± %|8.2f|") % n.events % n.events_uncertainty; + if (type != DATA) { + std::cout << " " << format("%|8.4f| ± %|8.4f|") % (n.efficiency * 100) % (n.efficiency_uncertainty * 100); + } + std::cout << std::endl; + + nominal_events += n.events; + nominal_events_uncertainty += n.events_uncertainty * n.events_uncertainty; + } + + if (systematics.size()) { + // Print systematics + std::cout << format("%|1$50| ---------------------") % " " << std::endl; + for (const auto& n: systematics) { + std::cout << Color::FG_YELLOW << format("%|50|") % truncate(n.name, 50) << Color::RESET << " " << format(" ± %|8.2f|") % n.events_uncertainty; + std::cout << std::endl; + + nominal_events_uncertainty += n.events_uncertainty * n.events_uncertainty; + } + } + + // Print total sum + std::cout << format("%|1$50| ---------------------") % " " << std::endl; + std::cout << format("%|50t| %|10.2f| ± %|8.2f|") % nominal_events % std::sqrt(nominal_events_uncertainty) << std::endl; + } +} diff --git a/test/unit_tests/configuration.py b/test/unit_tests/configuration.py index df08610..f6f27d8 100644 --- a/test/unit_tests/configuration.py +++ b/test/unit_tests/configuration.py @@ -53,7 +53,7 @@ def conf_get_files(): f = { 'data.root': {'type': 'data', 'legend': 'Data'}, 'MC_sample1.root': {'type': 'mc', 'legend': 'MC 1', 'cross-section': 245.8, 'generated-events': 2167, 'fill-color': '#D95B43', 'order': 1, 'group': 'mygroup'}, - 'MC_sample2.root': {'type': 'mc', 'legend': 'MC 2', 'cross-section': 666.3, 'generated-events': 2403, 'fill-color': '#53777A', 'order': 0, 'group': 'mygroup'}, + 'MC_sample2.root': {'type': 'mc', 'legend': 'MC 2', 'cross-section': 666.3, 'generated-events': 2404, 'fill-color': '#53777A', 'order': 0, 'group': 'mygroup'}, } return {'files': f}