|
[Rivet-svn] r4330 - in trunk: . include/Rivet src/Analyses src/Coreblackhole at projects.hepforge.org blackhole at projects.hepforge.orgThu May 30 13:00:39 BST 2013
Author: buckley Date: Thu May 30 13:00:39 2013 New Revision: 4330 Log: Reverting bookScatter2D behaviour to never look at ref data, and updating a few affected analyses. This should fix some bugs with doubled datapoints introduced by the previous behaviour+addPoint. Modified: trunk/ChangeLog trunk/include/Rivet/Analysis.hh trunk/src/Analyses/MC_JetAnalysis.cc trunk/src/Analyses/MC_JetSplittings.cc trunk/src/Analyses/MC_XS.cc trunk/src/Core/Analysis.cc Modified: trunk/ChangeLog ============================================================================== --- trunk/ChangeLog Thu May 30 12:59:06 2013 (r4329) +++ trunk/ChangeLog Thu May 30 13:00:39 2013 (r4330) @@ -1,3 +1,11 @@ +2013-05-30 Andy Buckley <andy.buckley at cern.ch> + + * Reverting bookScatter2D behaviour to never look at ref data, and + updating a few affected analyses. This should fix some bugs with + doubled datapoints introduced by the previous behaviour+addPoint. + + * Adding a couple of minor Utils.hh and MathUtils.hh features. + 2013-05-29 Andy Buckley <andy.buckley at cern.ch> * Removing Constraints.hh header. Modified: trunk/include/Rivet/Analysis.hh ============================================================================== --- trunk/include/Rivet/Analysis.hh Thu May 30 12:59:06 2013 (r4329) +++ trunk/include/Rivet/Analysis.hh Thu May 30 13:00:39 2013 (r4330) @@ -539,28 +539,43 @@ /// @name 2D scatter booking //@{ - /// Book a 2-dimensional data point set, using the binnings in the reference data histogram. + /// @brief Book a 2-dimensional data point set with the given name. + /// + /// @note Unlike histogram booking, scatter booking makes no attempt to use + /// reference data to pre-fill the data object. If you want this (why!?) + /// then use the @a refData() function to retrieve the reference scatter and + /// copy its contents. Scatter2DPtr bookScatter2D(const std::string& name, const std::string& title="", const std::string& xtitle="", const std::string& ytitle=""); - /// Book a 2-dimensional data point set, using the binnings in the reference data histogram. + /// @brief Book a 2-dimensional data point set, using the binnings in the reference data histogram. /// /// The paper, dataset and x/y-axis IDs will be used to build the histo name in the HepData standard way. + /// + /// + /// @note Unlike histogram booking, scatter booking makes no attempt to use + /// reference data to pre-fill the data object. If you want this (why!?) + /// then use the @a refData() function to retrieve the reference scatter and + /// copy its contents. Scatter2DPtr bookScatter2D(size_t datasetId, size_t xAxisId, size_t yAxisId, const std::string& title="", const std::string& xtitle="", const std::string& ytitle=""); - /// Book a 2-dimensional data point set with equally spaced points in a range. + /// @brief Book a 2-dimensional data point set with equally spaced x-points in a range. + /// + /// The y values and errors will be set to 0. Scatter2DPtr bookScatter2D(const std::string& name, size_t npts, double lower, double upper, const std::string& title="", const std::string& xtitle="", const std::string& ytitle=""); - /// Book a 2-dimensional data point set based on provided contiguous bin edges. + /// @brief Book a 2-dimensional data point set based on provided contiguous "bin edges". + /// + /// The y values and errors will be set to 0. Scatter2DPtr bookScatter2D(const std::string& hname, const std::vector<double>& binedges, const std::string& title, Modified: trunk/src/Analyses/MC_JetAnalysis.cc ============================================================================== --- trunk/src/Analyses/MC_JetAnalysis.cc Thu May 30 12:59:06 2013 (r4329) +++ trunk/src/Analyses/MC_JetAnalysis.cc Thu May 30 13:00:39 2013 (r4330) @@ -60,12 +60,11 @@ _h_jet_multi_exclusive = bookHisto1D("jet_multi_exclusive", m_njet+3, -0.5, m_njet+3-0.5); _h_jet_multi_inclusive = bookHisto1D("jet_multi_inclusive", m_njet+3, -0.5, m_njet+3-0.5); - _h_jet_multi_ratio = bookScatter2D("jet_multi_ratio", m_njet+2, 0.5, m_njet+3-0.5); + _h_jet_multi_ratio = bookScatter2D("jet_multi_ratio"); _h_jet_HT = bookHisto1D("jet_HT", logspace(50, m_jetptcut, sqrtS()/GeV/2.0)); } - // Do the analysis void MC_JetAnalysis::analyze(const Event & e) { const double weight = e.weight(); @@ -143,36 +142,25 @@ scale(_h_rap_jet[i], crossSection()/sumOfWeights()); // Create eta/rapidity ratio plots - /// @todo This should be neater: there should be a way to book an empty scatter by name. - Scatter2DPtr jet_eta_pmratio( new Scatter2D(histoPath("jet_eta_pmratio_" + to_str(i+1))) ); addPlot(jet_eta_pmratio); - Scatter2DPtr jet_rap_pmratio( new Scatter2D(histoPath("jet_y_pmratio_" + to_str(i+1))) ); addPlot(jet_rap_pmratio); - divide(*_h_eta_jet_plus[i], *_h_eta_jet_minus[i], jet_eta_pmratio); - divide(*_h_rap_jet_plus[i], *_h_rap_jet_minus[i], jet_rap_pmratio); + divide(*_h_eta_jet_plus[i], *_h_eta_jet_minus[i], bookScatter2D("jet_eta_pmratio_" + to_str(i+1))); + divide(*_h_rap_jet_plus[i], *_h_rap_jet_minus[i], bookScatter2D("jet_y_pmratio_" + to_str(i+1))); } - // Scale the d{eta,R} histograms - map<pair<size_t, size_t>, Histo1DPtr>::iterator it; - for (it = _h_deta_jets.begin(); it != _h_deta_jets.end(); ++it) { - scale(it->second, crossSection()/sumOfWeights()); - } - for (it = _h_dphi_jets.begin(); it != _h_dphi_jets.end(); ++it) { - scale(it->second, crossSection()/sumOfWeights()); - } - for (it = _h_dR_jets.begin(); it != _h_dR_jets.end(); ++it) { - scale(it->second, crossSection()/sumOfWeights()); - } + // Scale the d{eta,phi,R} histograms + typedef map<pair<size_t, size_t>, Histo1DPtr> HistMap; + foreach (HistMap::value_type& it, _h_deta_jets) scale(it.second, crossSection()/sumOfWeights()); + foreach (HistMap::value_type& it, _h_dphi_jets) scale(it.second, crossSection()/sumOfWeights()); + foreach (HistMap::value_type& it, _h_dR_jets) scale(it.second, crossSection()/sumOfWeights()); // Fill inclusive jet multi ratio int Nbins = _h_jet_multi_inclusive->numBins(); for (int i = 0; i < Nbins-1; ++i) { - if (_h_jet_multi_inclusive->bin(i).area() > 0.0 && _h_jet_multi_inclusive->bin(i+1).area() > 0.0) { - double ratio = _h_jet_multi_inclusive->bin(i+1).area()/_h_jet_multi_inclusive->bin(i).area(); - double relerr_i = _h_jet_multi_inclusive->bin(i).areaErr()/_h_jet_multi_inclusive->bin(i).area(); - double relerr_j = _h_jet_multi_inclusive->bin(i+1).areaErr()/_h_jet_multi_inclusive->bin(i+1).area(); - double err = ratio * (relerr_i + relerr_j); - Point2D & pt = _h_jet_multi_ratio->point(i); - pt.setY(ratio); - pt.setYErr(err); + if (_h_jet_multi_inclusive->bin(i).sumW() > 0.0) { //< Don't add a point if the division is invalid + const double ratio = _h_jet_multi_inclusive->bin(i+1).sumW()/_h_jet_multi_inclusive->bin(i).sumW(); + const double relerr_i = _h_jet_multi_inclusive->bin(i).relErr(); + const double relerr_j = _h_jet_multi_inclusive->bin(i+1).relErr(); + const double err = ratio * (relerr_i + relerr_j); + _h_jet_multi_ratio->addPoint(i, ratio, 0.0, err); } } Modified: trunk/src/Analyses/MC_JetSplittings.cc ============================================================================== --- trunk/src/Analyses/MC_JetSplittings.cc Thu May 30 12:59:06 2013 (r4329) +++ trunk/src/Analyses/MC_JetSplittings.cc Thu May 30 13:00:39 2013 (r4330) @@ -18,20 +18,14 @@ // Book histograms void MC_JetSplittings::init() { - - for (size_t i=0; i < m_njet; ++i) { - stringstream dname; - dname << "log10_d_" << i << i+1; - - _h_log10_d[i] = bookHisto1D(dname.str(), 100, 0.2, log10(0.5*sqrtS())); - - stringstream Rname; - Rname << "log10_R_" << i; - _h_log10_R[i] = bookScatter2D(Rname.str(), 50, 0.2, log10(0.5*sqrtS())); + for (size_t i = 0; i < m_njet; ++i) { + string dname = "log10_d_" + to_str(i) + to_str(i+1); + _h_log10_d[i] = bookHisto1D(dname, 100, 0.2, log10(0.5*sqrtS())); + string Rname = "log10_R_" + to_str(i); + _h_log10_R[i] = bookScatter2D(Rname, 50, 0.2, log10(0.5*sqrtS())); } - stringstream Rname; - Rname << "log10_R_" << m_njet; - _h_log10_R[m_njet] = bookScatter2D(Rname.str(), 50, 0.2, log10(0.5*sqrtS())); + string Rname = "log10_R_" + to_str(m_njet); + _h_log10_R[m_njet] = bookScatter2D(Rname, 50, 0.2, log10(0.5*sqrtS())); } @@ -56,6 +50,7 @@ // Fill integrated jet resolution for (size_t ibin = 0; ibin < _h_log10_R[i]->numPoints(); ++ibin) { Point2D & dp = _h_log10_R[i]->point(ibin); + /// @todo Could inRange be used here (avoiding the temporary variable)? double dcut = dp.x(); if (d_ij < dcut && previous_dij > dcut) { dp.setY(dp.y() + weight); Modified: trunk/src/Analyses/MC_XS.cc ============================================================================== --- trunk/src/Analyses/MC_XS.cc Thu May 30 12:59:06 2013 (r4329) +++ trunk/src/Analyses/MC_XS.cc Thu May 30 13:00:39 2013 (r4330) @@ -27,7 +27,8 @@ /// Book histograms and initialise projections before the run void init() { - _h_XS = bookScatter2D("XS", 1, 0.0, 1.0); + /// @todo Convert to Scatter1D or Counter + _h_XS = bookScatter2D("XS"); _h_N = bookHisto1D("N", 1, 0.0, 1.0); _h_pmXS = bookHisto1D("pmXS", 2, -1.0, 1.0); _h_pmN = bookHisto1D("pmN", 2, -1.0, 1.0); @@ -38,8 +39,8 @@ /// Perform the per-event analysis void analyze(const Event& event) { _h_N->fill(0.5,1.); - _h_pmXS->fill(0.5*(event.weight()> 0 ? 1. : -1), abs(event.weight())); - _h_pmN->fill(0.5*(event.weight() > 0 ? 1. : -1), 1.); + _h_pmXS->fill(0.5*(event.weight() > 0 ? 1. : -1), abs(event.weight())); + _h_pmN ->fill(0.5*(event.weight() > 0 ? 1. : -1), 1.); #ifdef HEPMC_HAS_CROSS_SECTION _mc_xs = event.genEvent()->cross_section()->cross_section(); _mc_error = event.genEvent()->cross_section()->cross_section_error(); Modified: trunk/src/Core/Analysis.cc ============================================================================== --- trunk/src/Core/Analysis.cc Thu May 30 12:59:06 2013 (r4329) +++ trunk/src/Core/Analysis.cc Thu May 30 13:00:39 2013 (r4330) @@ -396,17 +396,17 @@ const string& title, const string& xtitle, const string& ytitle) { - const Scatter2D& refdata = refData(hname); const string path = histoPath(hname); - Scatter2DPtr s( new Scatter2D(refdata, path) ); - foreach (Point2D& p, s->points()) p.setY(0, 0); + // const Scatter2D& refdata = refData(hname); + // Scatter2DPtr s( new Scatter2D(refdata, path) ); + // foreach (Point2D& p, s->points()) p.setY(0, 0); + Scatter2DPtr s( new Scatter2D(path) ); addPlot(s); MSG_TRACE("Made scatter " << hname << " for " << name()); s->setTitle(title); s->setAnnotation("XLabel", xtitle); s->setAnnotation("YLabel", ytitle); return s; - }
More information about the Rivet-svn mailing list |