|
[Rivet-svn] r3230 - in branches/2011-07-aida2yoda/include: Rivet YODAblackhole at projects.hepforge.org blackhole at projects.hepforge.orgThu Jul 21 09:23:53 BST 2011
Author: mkawalec Date: Thu Jul 21 09:23:52 2011 New Revision: 3230 Log: Added Histo2D and all relevant files. Note that RivetYODA doesn't yet work and the files in YODA do need some tweaking (like comments), but they compile, so here is the code. Added: branches/2011-07-aida2yoda/include/YODA/ branches/2011-07-aida2yoda/include/YODA/.Bin2D.h.swn (contents, props changed) branches/2011-07-aida2yoda/include/YODA/.HistoBin2D.h.swo (contents, props changed) branches/2011-07-aida2yoda/include/YODA/Axis2D.h branches/2011-07-aida2yoda/include/YODA/Bin1D.h branches/2011-07-aida2yoda/include/YODA/Bin2D.h branches/2011-07-aida2yoda/include/YODA/Histo2D.h branches/2011-07-aida2yoda/include/YODA/HistoBin2D.h branches/2011-07-aida2yoda/include/YODA/Makefile.am branches/2011-07-aida2yoda/include/YODA/ProfileBin1D.h branches/2011-07-aida2yoda/include/YODA/Reader.h branches/2011-07-aida2yoda/include/YODA/ReaderYODA.h branches/2011-07-aida2yoda/include/YODA/Scatter3D.h branches/2011-07-aida2yoda/include/YODA/WriterMethods.icc branches/2011-07-aida2yoda/include/YODA/WriterYODA.h Modified: branches/2011-07-aida2yoda/include/Rivet/RivetYODA.hh Modified: branches/2011-07-aida2yoda/include/Rivet/RivetYODA.hh ============================================================================== --- branches/2011-07-aida2yoda/include/Rivet/RivetYODA.hh Wed Jul 20 19:22:37 2011 (r3229) +++ branches/2011-07-aida2yoda/include/Rivet/RivetYODA.hh Thu Jul 21 09:23:52 2011 (r3230) @@ -17,6 +17,36 @@ #include "YODA/Point2D.h" #include "YODA/ReaderAIDA.h" +// Yoda ones: +#include "YODA/AnalysisObject.h" +#include "YODA/Axis1D.h" +#include "YODA/Bin1D.h" +#include "YODA/Bin.h" +#include "YODA/Dbn1D.h" +#include "YODA/Exceptions.h" +#include "YODA/Histo1D.h" +#include "YODA/HistoBin1D.h" +#include "YODA/Point2D.h" +#include "YODA/Profile1D.h" +#include "YODA/ProfileBin1D.h" +#include "YODA/ReaderAIDA.h" +#include "YODA/Reader.h" +#include "YODA/ReaderYODA.h" +#include "YODA/Scatter2D.h" +#include "YODA/WriterAIDA.h" +#include "YODA/Writer.h" +#include "YODA/WriterYODA.h" + +//And second dimension goes here: +#include "YODA/Axis2D.h" +#include "YODA/Bin2D.h" +#include "YODA/Dbn2D.h" +#include "YODA/Histo2D.h" +#include "YODA/HistoBin2D.h" +#include "YODA/Point3D.h" +#include "YODA/Scatter3D.h" + + namespace Rivet { /// Function to get a map of all the bin edge vectors in a paper with the @@ -34,7 +64,7 @@ /// Return the integral over the histogram bins inline double integral(Histo1DPtr histo) { - double intg = 0.; + double intg = 0.0; for ( size_t i = 0; i < histo->numBins(); ++i ) { intg += histo->bin(i).area(); } Added: branches/2011-07-aida2yoda/include/YODA/.Bin2D.h.swn ============================================================================== Binary file. No diff available. Added: branches/2011-07-aida2yoda/include/YODA/.HistoBin2D.h.swo ============================================================================== Binary file. No diff available. Added: branches/2011-07-aida2yoda/include/YODA/Axis2D.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/Axis2D.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,300 @@ +#ifndef YODA_Axis2D_h +#define YODA_Axis2D_h + +#include "YODA/AnalysisObject.h" +#include "YODA/Exceptions.h" +#include "YODA/Bin.h" +#include "YODA/Utils/sortedvector.h" +#include "YODA/Utils/MathUtils.h" +#include "YODA/Dbn2D.h" +#include <string> +#include <cassert> +#include <cmath> +#include <algorithm> + +using namespace std; + +namespace YODA { + template <typename BIN> + class Axis2D { + public: + + typedef BIN Bin; + typedef typename Utils::sortedvector<Utils::sortedvector<BIN> > Bins; + + static inline pair<std::vector<double>, std::vector<double> > mkBinEdgesLinLin(double startX, double startY, double endX, double endY, size_t nbinsX, size_t nbinsY) { + return make_pair(linspace(startX, endX, nbinsX), linspace(startY, endY, nbinsY)); + } + + static inline pair<std::vector<double>, std::vector<double> > mkBinEdgesLinLog(double startX, double startY, double endX, double endY, size_t nbinsX, size_t nbinsY) { + return make_pair(linspace(startX, endX, nbinsX), logspace(startY, endY, nbinsY)); + } + + static inline pair<std::vector<double>, std::vector<double> > mkBinEdgesLogLog(double startX, double startY, double endX, double endY, size_t nbinsX, size_t nbinsY) { + return make_pair(logspace(startX, endX, nbinsX), logspace(startY, endY, nbinsY)); + } + + static inline pair<std::vector<double>, std::vector<double> > mkBinEdgesLogLin(double startX, double startY, double endX, double endY, size_t nbinsX, size_t nbinsY) { + return make_pair(logspace(startX, endX, nbinsX), linspace(startY, endY, nbinsY)); + } + + private: + + /* _binHash is of type pair< vector <double>, vector <double> >, + * it is used for indexing the bin edges + * _cachedBinEdges is of type pair<vector<double>,vector<double> > + */ + + void _mkBinHash() { + for (size_t i = 0; i < numBinsX(); i++) { + _binHash.first.insert(make_pair(_cachedBinEdges.first[i+1], i)); + } + for (size_t j = 0; j < numBinsY(); j++) { + _binHash.second.insert(make_pair(_cachedBinEdges.second[j+1], j)); + } + } + + void _mkAxis(const pair<vector<double>,vector<double> >& binedges) { + const size_t nbinsX = binedges.first.size() - 1; + const size_t nbinsY = binedges.second.size() -1; + + //Now we are pushing back the bins: + for(int i=0; i<nbinsY; i++) { + Utils::sortedvector<BIN> Temp; + for(int j=0; j<nbinsX; j++) { + Temp.push_back( BIN(binedges.first[j], binedges.second[i], binedges.first[j+1], binedges.second[i+1]) ); + } + Temp.sort(); + _bins.push_back(Temp); + } + _bins.sort(); + + //Hashing the bin edges: + _cachedBinEdges = binedges; + std::sort(_cachedBinEdges.first.begin(), _cachedBinEdges.first.end()); + std::sort(_cachedBinEdges.second.begin(), _cachedBinEdges.second.end()); + _mkBinHash(); + } + + public: + //Constructors: + + //Default constructor with + Axis2D(const pair<vector<double>,vector<double> >& binedges) { + assert(binedges.first.size() > 1 && binedges.second.size() > 1); + _mkAxis(binedges); + } + + //Most standard constructor, should be self-explanatory + Axis2D(size_t nbinsX, double lowerX, double upperX, size_t nbinsY, double lowerY, double upperY) { + _mkAxis( make_pair(linspace(nbinsX, lowerX, upperX), linspace(nbinsY, lowerY, upperY)) ); + } + + //Constructor allowing to specify a vector of vectors that house bins is not specified and + // not supported as would make filling much more computationally involving + + + //Some helper functions: + + //Functions returning different size statistics + unsigned int numBinsX() const { + return _bins[0].size(); + } + unsigned int numBinsY() const { + return _bins.size(); + } + unsigned int numBinsTotal() const { + return _bins.size()*_bins[0].size(); + } + + Bins& bins() { + return _bins; + } + + const Bins& bins() const { + return _bins; + } + + std::pair<pair<double,double>, pair<double,double> > binEdges(size_t binIdX, size_t binIdY) const { + assert(binIdX < numBinsX() && binIdY < numBinsY()); + return make_pair(make_pair(_cachedBinEdges.first[binIdX], _cachedBinEdges.first[binIdX+1]), make_pair(_cachedBinEdges.second[binIdY], _cachedBinEdges.second[binIdY+1])); + } + + //Edges: + double lowEdgeX() const { + return _bins.front().front().lowEdgeX(); + } + double lowEdgeY() const { + return _bins.front().front().lowEdgeY(); + } + double highEdgeX() const { + return _bins.front().back().highEdgeX(); + } + double highEdgeY() const { + return _bins.back().front().highEdgeY(); + } + + //Maybe throwing different errors for X and Y? + BIN& bin(size_t indexX, size_t indexY) { + if(indexX >= numBinsX() || indexY >= numBinsY()) { + throw RangeError("YODA::Histo: index out of range"); } + return _bins[indexY][indexX]; + } + const BIN& bin(size_t indexX, size_t indexY) const { + if(indexX >= numBinsX() || indexY >= numBinsY()) { + throw RangeError("YODA::Histo: index out of range");} + return _bins[indexY][indexX]; + } + BIN& bin(std::pair<size_t, size_t> index) { + bin(index.first, index.second); + } + const BIN& bin(std::pair<size_t, size_t> index) const { + bin(index.first, index.second); + } + + + BIN& binByCoord(double x, double y) { + return bin(findBinIndex(x, y)); + } + + const BIN& binByCoord(double x, double y) const { + return bin(findBinIndex(x, y)); + } + BIN& binByCoord(pair<double, double>& coords) { + return bin(findBinIndex(coords.first, coords.second)); + } + + const BIN& binByCoord(pair<double, double>& coords) const { + return bin(findBinIndex(coords.first, coords.second)); + } + + + + Dbn2D& totalDbn() { + return _dbn; + } + const Dbn2D& totalDbn() const{ + return _dbn; + } + + Dbn2D& overflow() { + return _overflow; + } + const Dbn2D& overflow() const { + return _overflow; + } + + std::pair<size_t, size_t> findBinIndex(double coordX, double coordY) const { + if (coordX < *_cachedBinEdges.first.begin() || + coordY < *_cachedBinEdges.second.begin() || + coordX >= *_cachedBinEdges.first.end() || + coordY >= *_cachedBinEdges.second.end()) { + throw RangeError("Coordinate is outside the valid range: you should request the underflow or overflow"); + } + + size_t x = _binHash.first.upper_bound(coordX) -> second; + size_t y = _binHash.second.lower_bound(coordY) -> second; + return make_pair(x, y); + } + + void reset() { + _dbn.reset(); + _underflow.reset(); + _overflow.reset(); + for (unsigned int i=0; i<_bins.size(); i++) { + for (unsigned int j=0; j < _bins[i].size(); j++) { + _bins[i][j].reset(); + } + } + } + + void scale(double scalefactor) { + /// @todo Implement! + throw std::runtime_error("Axis coordinate transformations not yet implemented!"); + } + + void scaleW(double scalefactor) { + _dbn.scaleW(scalefactor); + _underflow.scaleW(scalefactor); + _overflow.scaleW(scalefactor); + for (int i=0; i<_bins.size(); i++) { + for (int j=0; j < _bins[i].size(); j++ { + _bins[i][j].scaleW(scalefactor); + } + } + } + // Operators: + bool operator == (const Axis2D& other) const { + return + _cachedBinEdges.first == _cachedBinEdges.first; + _cachedBinEdges.second == _cachedBinEdges.second; + } + + bool operator != (const Axis2D& other) const { + return ! operator == (other); + } + + //We must use at() as bins() is a function + Axis2D<BIN>& operator += (const Axis2D<BIN>& toAdd) { + if (*this != toAdd) { + throw LogicError("YODA::Histo1D: Cannot add axes with different binnings."); + } + for (int i=0; i < bins().size(); i++) { + for (unsigned int j = 0; j < bins().at(i).size(); i++) { + bins().at(i).at(j) += toAdd.bins().at(i).at(j); + } + } + _dbn += toAdd._dbn; + _underflow += toAdd._underflow; + _overflow += toAdd._overflow; + return *this; + } + + Axis2D<BIN>& operator -= (const Axis2D<BIN>& toSubstract) { + if (*this != toSubstract) { + throw LogicError("YODA::Histo1D: Cannot add axes with different binnings."); + } + for (int i=0; i < bins().size(); i++) { + for (size_t j = 0; j < bins().at(i).size(); i++) { + bins().at(i).at(j) -= toSubstract.bins().at(i).at(j); + } + } + _dbn -= toSubstract._dbn; + _underflow -= toSubstract._underflow; + _overflow -= toSubstract._overflow; + return *this; + } + + private: + //Bins contained in this histogram: + Bins _bins; + Dbn2D _underflow; + Dbn2D _overflow; + Dbn2D _dbn; + + std::pair<std::vector<double>, std::vector<double> > _cachedBinEdges; + std::pair<std::map<double,size_t>, std::map<double,size_t> > _binHash; + + }; + + template <typename BIN> + Axis2D<BIN> operator + (const Axis2D<BIN>& first, const Axis2D<BIN>& second) { + Axis2D<BIN> tmp = first; + tmp += second; + return tmp; + } + + template <typename BIN> + Axis2D<BIN> operator - (const Axis2D<BIN>& first, const Axis2D<BIN>& second) { + Axis2D<BIN> tmp = first; + tmp -= second; + return tmp; + } + + inline bool operator < (const pair<vector<double>, vector<double> >& a, const pair<vector<double>, vector<double> >& b) { + if(a.first == b.first) return a.second < b.second; + return a.first < b.first; + } +} + +#endif Added: branches/2011-07-aida2yoda/include/YODA/Bin1D.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/Bin1D.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,163 @@ +// -*- C++ -*- +// +// This file is part of YODA -- Yet more Objects for Data Analysis +// Copyright (C) 2008-2011 The YODA collaboration (see AUTHORS for details) +// +#ifndef YODA_Bin1D_h +#define YODA_Bin1D_h + +#include "YODA/Bin.h" +#include "YODA/Dbn1D.h" +#include <string> +#include <utility> + +namespace YODA { + + + /// @brief Base class for bins in 1D normal and profile histograms. + /// The lower bin edge is inclusive. This base class provides no fill + /// method, since the signatures for standard and profile histos differ. + class Bin1D : public Bin { + public: + + /// @name Constructors, giving bin low and high edges. + //@{ + Bin1D(double lowedge, double highedge); + + Bin1D(std::pair<double,double> edges); + //@} + + + /// @name Miscellaneous + //@{ + + /// Reset this bin + virtual void reset(); + + //@} + + + public: + + /// @name X-axis info + //@{ + /// Lower limit of the bin (inclusive). + double lowEdge() const; + double xMin() const { return lowEdge(); } + + /// Upper limit of the bin (exclusive). + double highEdge() const; + double xMax() const { return highEdge(); } + + /// Get the {low,high} edges as an STL @c pair. + std::pair<double,double> edges() const; + + /// Separation of low and high edges, i.e. high-low. + double width() const; + + /// The mean position in the bin, or the midpoint if that is not available. + double focus() const; + + /// Geometric centre of the bin, i.e. high+low/2.0 + double midpoint() const; + //@} + + + public: + + /// @name X distribution statistics + //@{ + + /// Mean value of x-values in the bin. + double xMean() const; + + /// The variance of x-values in the bin. + double xVariance() const; + + /// The standard deviation (spread) of x-values in the bin. + double xStdDev() const; + + /// The standard error on the bin focus. + double xStdError() const; + //@} + + + public: + + /// @name Raw x distribution statistics + //@{ + + /// The number of entries + unsigned long numEntries() const; + + /// The sum of weights + double sumW() const; + + /// The sum of weights squared + double sumW2() const; + + /// The sum of x*weight + double sumWX() const; + + /// The sum of x^2 * weight + double sumWX2() const; + + //@} + + + public: + + /// @name Operators + //@{ + + /// Add two bins + Bin1D& operator += (const Bin1D&); + + /// Subtract one bin from another + Bin1D& operator -= (const Bin1D&); + //@} + + + protected: + + /// @name Named operators + //@{ + + /// Add two bins (internal, explicitly named version) + Bin1D& add(const Bin1D&); + + /// Subtract one bin from another (internal, explicitly named version) + Bin1D& subtract(const Bin1D&); + //@} + + + protected: + + /// The bin limits + std::pair<double,double> _edges; + + // Distribution of weighted x values + Dbn1D _xdbn; + + }; + + + /// Add two bins + Bin1D operator + (const Bin1D& a, const Bin1D& b); + + /// Subtract one bin from another + Bin1D operator - (const Bin1D& a, const Bin1D& b); + + + /// Bin1Ds are compared for axis sorting by lower edge position + /// @todo Check for overlaps somewhere... on Axis1D, I guess. + inline bool operator<(const Bin1D& a, const Bin1D& b) { + return b.edges().first > a.edges().first; + } + + +} + + + +#endif Added: branches/2011-07-aida2yoda/include/YODA/Bin2D.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/Bin2D.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,94 @@ +#ifndef YODA_Bin2D_h +#define YODA_Bin2D_h + +#include "YODA/Bin.h" +#include "YODA/Dbn2D.h" +#include <string> +#include <utility> + +namespace YODA { + + class Bin2D : public Bin { + public: + + Bin2D(double lowedgeX, double lowedgeY, double highedgeX, double highedgeY); + Bin2D(std::pair<std::pair<double, double>, std::pair<double, double> > edges); + + virtual void reset(); + + + double lowEdgeX() const; + double xMin() const { return lowEdgeX(); } + + double lowEdgeY() const; + double yMin() const { return lowEdgeY(); } + + double highEdgeX() const; + double xMax() const { return highEdgeX(); } + + double highEdgeY() const; + double yMax() const { return highEdgeY(); } + + std::pair<double, double> edgesX() const; + std::pair<double, double> edgesY() const; + + double widthX() const; + double widthY() const; + + std::pair<double, double> focus() const; + std::pair<double, double> midpoint() const; + + + //Now some distribution statistics: + double xMean() const; + double yMean() const; + std::pair<double, double> Mean() const; + + double xVariance() const; + double yVariance() const; + std::pair<double, double> Variance() const; + + double xStdDev() const; + double yStdDev() const; + std::pair<double, double> StdDev() const; + + //Standar error, previously named StdError! + double xStdErr() const; + double yStdErr() const; + std::pair<double, double> StdErr() const; + + + //Some "raw distribution statistics" + unsigned long numEntries() const; + double sumW() const; + double sumW2() const; + double sumWX() const; + double sumWY() const; + double sumWXY() const; + double sumX2() const; + double sumY2() const; + + Bin2D& operator += (const Bin2D&); + Bin2D& operator -= (const Bin2D&); + + protected: + Bin2D& add(const Bin2D&); + Bin2D& substract(const Bin2D&); + + std::pair<double,double> _edgesX; + std::pair<double,double> _edgesY; + Dbn2D _xdbn; + + }; + + Bin2D operator + (const Bin2D& a, const Bin2D& b); + Bin2D operator - (const Bin2D& a, const Bin2D& b); + + //Comparison functions that _may_ be used for sorting later, we will see... + inline bool operator<(const Bin2D& a, const Bin2D& b) { + if(a.edgesX().first == b.edgesX().first) return b.edgesY().first > a.edgesY().first; + return b.edgesX().first > a.edgesX().first; +} + + +#endif Added: branches/2011-07-aida2yoda/include/YODA/Histo2D.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/Histo2D.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,253 @@ +// -*- C++ -*- +// +// This file is part of YODA -- Yet more Objects for Data Analysis +// Copyright (C) 2008-2011 The YODA collaboration (see AUTHORS for details) +// +#ifndef YODA_Histo2D_h +#define YODA_Histo2D_h + +#include "YODA/AnalysisObject.h" +#include "YODA/HistoBin2D.h" +#include "YODA/Scatter3D.h" +#include "YODA/Axis2D.h" +#include "YODA/Exceptions.h" +#include <vector> +#include <string> +#include <map> + +namespace YODA { + + + /// Convenience typedef + typedef Axis2D<HistoBin2D> Histo2DAxis; + + + /// A one-dimensional histogram. + class Histo2D : public AnalysisObject { + + public: + + /// Convenience typedefs + typedef Histo2DAxis Axis; + typedef Axis::Bins Bins; + + + /// @name Constructors + //@{ + + /// Constructor giving range and number of bins. + /// @todo Remove binning enum stuff + Histo2D(size_t nbinsX, double lowerX, double upperX, + size_t nbinsY, double lowerY, double upperY, + const std::string& path="", const std::string& title="") + : AnalysisObject("Histo2D", path, title), + _axis(nbinsX, lowerX, upperX, + nbinsY, lowerY, upperY) + { } + + /// @brief Constructor giving explicit bin edges. + /// For n bins, binedges.size() == n+1, the last + /// one being the upper bound of the last bin + Histo2D(const std::pair<std::vector<double>, std::vector<double> >& binedges, + const std::string& path="", const std::string& title="") + : AnalysisObject("Histo2D", path, title), + _axis(binedges) + { } + + /// Copy constructor with optional new path + Histo2D(const Histo2D& h, const std::string& path=""); + + /// Constructor from a Scatter3D's binning, with optional new path + Histo2D(const Scatter3D& s, const std::string& path=""); + + + //@} + + + public: + + /// @name Persistency hooks + //@{ + + /// Get name of the analysis object type, for persisting + std::string _aotype() const { return "Histo2D"; } + + /// Set the state of the histo object, for unpersisting + /// @todo Need to set annotations (do that on AO), all-histo Dbns, and dbns for every bin. Delegate! + // void _setstate() = 0; + + //@} + + + /// @name Modifiers + //@{ + + /// Fill histo by value and weight + void fill(double x, double weight=1.0); + + /// @brief Reset the histogram. + /// Keep the binning but set all bin contents and related quantities to zero + virtual void reset() { + _axis.reset(); + } + + /// Rescale as if all fill weights had been different by factor @a scalefactor. + void scaleW(double scalefactor) { + _axis.scaleW(scalefactor); + } + + //@} + + + public: + + /// @name Bin accessors + //@{ + + /// Number of bins on this axis (not counting under/overflow) + size_t numBinsX() const { + return bins().at(0).size(); + } + + size_t numBinsY() const { + return bins().size(); + } + + /// Low edge of this histo's axis + double lowEdgeX() const { + return _axis.lowEdgeX(); + } + + double lowEdgeY() const { + return _axis.lowEdgeY(); + } + + /// High edge of this histo's axis + double highEdgeX() const { + return _axis.highEdgeX(); + } + + double highEdgeY() const { + return _axis.highEdgeY(); + } + + /// Access the bin vector + /// @todo Actually, it's a Histo + Utils::sortedvector<Utils::sortedvector<YODA::HistoBin2D> >& bins() { + return _axis.bins(); + } + + /// Access the bin vector (const version) + const Utils::sortedvector<Utils::sortedvector<YODA::HistoBin2D> >& bins() const { + return _axis.bins(); + } + + /// Access a bin by index (non-const version) + HistoBin2D& bin(size_t indexX, size_t indexY) { + return _axis.bins()[indexX][indexY]; + } + + /// Access a bin by index (const version) + const HistoBin2D& bin(size_t indexX, size_t indexY) const { + return _axis.bins()[indexX][indexY]; + } + + /// Access a bin by coordinate (non-const version) + HistoBin2D& binByCoord(double x, double y) { + return _axis.binByCoord(x, y); + } + + /// Access a bin by coordinate (const version) + const HistoBin2D& binByCoord(double x, double y) const { + return _axis.binByCoord(x, y); + } + + //@} + + + public: + + /// @name Whole histo data + //@{ + + /// Get the total area of the histogram + double integral(bool includeoverflows=true) const { + return sumW(includeoverflows); + } + + /// Get sum of weights in histo + double sumW(bool includeoverflows=true) const; + + /// Get sum of squared weights in histo + double sumW2(bool includeoverflows=true) const; + + /// Get the mean + double mean(bool includeoverflows=true) const; + + /// Get the variance + double variance(bool includeoverflows=true) const; + + /// Get the standard deviation + double stdDev(bool includeoverflows=true) const { + return std::sqrt(variance(includeoverflows)); + } + + //@} + + + public: + + /// @name Adding and subtracting histograms + //@{ + + /// Add another histogram to this + Histo2D& operator += (const Histo2D& toAdd) { + _axis += toAdd._axis; + return *this; + } + + /// Subtract another histogram from this + Histo2D& operator -= (const Histo2D& toSubtract) { + _axis -= toSubtract._axis; + return *this; + } + + //@} + + + private: + + /// @name Bin data + //@{ + + /// Definition of bin edges and contents + Axis2D<HistoBin2D> _axis; + + //@} + + }; + + + /// @name Combining histos: global operators + //@{ + + /// Add two histograms + inline Histo2D operator + (const Histo2D& first, const Histo2D& second) { + Histo2D tmp = first; + tmp += second; + return tmp; + } + + /// Subtract two histograms + inline Histo2D operator - (const Histo2D& first, const Histo2D& second) { + Histo2D tmp = first; + tmp -= second; + return tmp; + } + + //@} + + +} + +#endif Added: branches/2011-07-aida2yoda/include/YODA/HistoBin2D.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/HistoBin2D.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,52 @@ +#ifndef YODA_HistoBin2D_h +#define YODA_HistoBin2D_h + +#include "YODA/Bin2D.h" +#include "YODA/Exceptions.h" + +namespace YODA { + class HistoBin2D : public Bin2D { + public: + //Constructors: + HistoBin2D(double lowEdgeX, double highEdgeX, + double lowEdgeY, double highEdgeY); + + HistoBin2D(std::pair<std::pair<double,double>, std::pair<double,double> >& edges); + + //Modifiers: + void fill(std::pair<double,double>, double weight=1.0); + void fill(double coordX, double coordY, double weight=1.0); + + void fillBin(double weight=1.0); + void reset() { + Bin2D::reset(); + } + + void scaleW(double scalefactor) { + _xdbn.scaleW(scalefactor); + } + + //Bin content info: + + //Note that area is actually a volume, kept to keep variable names + //the same with *1D + double area() const; + double height() const; + double areaErr() const; + double heightErr() const; + + //Operators: + HistoBin2D& operator += (const HistoBin2D&); + HistoBin2D& operator -= (const HistoBin2D&); + + protected: + HistoBin2D& add(const HistoBin2D&); + HistoBin2D& substract(const HistoBin2D&); + }; + + HistoBin2D operator + (const HistoBin2D& a, const HistoBin2D& b); + HistoBin2D operator - (const HistoBin2D& a, const HistoBin2D& b); + +} + +#endif Added: branches/2011-07-aida2yoda/include/YODA/Makefile.am ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/Makefile.am Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,16 @@ +pkginclude_HEADERS = \ + Exceptions.h \ + AnalysisObject.h Bin.h \ + Axis1D.h Bin1D.h \ + Histo1D.h HistoBin1D.h \ + Profile1D.h ProfileBin1D.h \ + Scatter2D.h Point2D.h \ + Writer.h WriterAIDA.h WriterYODA.h + +nobase_pkginclude_HEADERS = \ + WriterMethods.icc + Utils/MathUtils.h + Utils/sortedvector.h \ + Utils/nvector.h \ + Config/YodaConfig.h \ + Config/BuildConfig.h Added: branches/2011-07-aida2yoda/include/YODA/ProfileBin1D.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/ProfileBin1D.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,117 @@ +// -*- C++ -*- +// +// This file is part of YODA -- Yet more Objects for Data Analysis +// Copyright (C) 2008-2011 The YODA collaboration (see AUTHORS for details) +// +#ifndef YODA_ProfileBin1D_h +#define YODA_ProfileBin1D_h + +#include "YODA/Bin1D.h" +#include "YODA/Exceptions.h" + +namespace YODA { + + + /// A Bin in a 1D profile histogram. + class ProfileBin1D : public Bin1D { + + /// Profile1D is a friend to add/subtract bins + friend class Profile1D; + + + public: + + /// @name Constructors + //@{ + + /// Constructor giving bin low and high edges. + ProfileBin1D(double lowedge, double highedge); + + /// Constructor giving bin low and high edges as a pair. + ProfileBin1D(std::pair<double,double> edges); + + //@} + + + /// @name Modifiers + //@{ + + /// Fill histo by value and weight. + void fill(double x, double d, double weight=1.0); + + /// Fill histo with @a weight at bin midpoint. + void fillBin(double d, double weight=1.0); + + /// Reset the bin. + void reset(); + + /// Rescale as if all fill weights had been different by factor @a scalefactor. + void scaleW(double scalefactor) { + _xdbn.scaleW(scalefactor); + _ydbn.scaleW(scalefactor); + } + + //@} + + + public: + + /// @name Bin content info + //@{ + + /// Make the "y" explicit? + + double mean() const; + + double stdDev() const; + + double variance() const; + + double stdErr() const; + + //@} + + + public: + + /// Add two bins (for use by Profile1D). + ProfileBin1D& operator += (const ProfileBin1D&); + + /// Subtract two bins + ProfileBin1D& operator -= (const ProfileBin1D&); + + + protected: + + /// Add two bins (internal, explicitly named version) + ProfileBin1D& add(const ProfileBin1D&); + + /// Subtract one bin from another (internal, explicitly named version) + ProfileBin1D& subtract(const ProfileBin1D&); + + + public: + + /// The sum of y*weight + double sumWY() const; + + /// The sum of y^2 * weight + double sumWY2() const; + + + private: + + // Distribution of weighted data values + Dbn1D _ydbn; + + + }; + + + ProfileBin1D operator + (const ProfileBin1D& a, const ProfileBin1D& b); + + ProfileBin1D operator - (const ProfileBin1D& a, const ProfileBin1D& b); + +} + +#endif Added: branches/2011-07-aida2yoda/include/YODA/Reader.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/Reader.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,79 @@ +// -*- C++ -*- +// +// This file is part of YODA -- Yet more Objects for Data Analysis +// Copyright (C) 2008-2011 The YODA collaboration (see AUTHORS for details) +// +#ifndef YODA_Reader_h +#define YODA_Reader_h + +#include "YODA/AnalysisObject.h" +#include "YODA/Histo1D.h" +#include "YODA/Profile1D.h" +#include "YODA/Scatter2D.h" +#include <string> +#include <fstream> +#include <vector> + +namespace YODA { + + + /// Pure virtual base class for various output writers. + class Reader { + public: + + /// Virtual destructor + virtual ~Reader() {} + + + /// @name Reading multiple analysis objects, + //@{ + + /// @brief Read in a collection of objects @a objs from output stream @a stream. + /// This version fills (actually, appends to) a supplied vector, avoiding copying, + /// and is hence CPU efficient. + virtual void read(std::istream& stream, std::vector<AnalysisObject*>& aos) = 0; + + /// @brief Read in a collection of objects from output stream @a stream. + /// This version returns a vector by value, involving copying, and is hence less + /// CPU efficient than the alternative version where a vector is filled by reference. + std::vector<AnalysisObject*> read(std::istream& stream) { + std::vector<AnalysisObject*> rtn; + read(stream, rtn); + return rtn; + } + + /// @brief Read in a collection of objects @a objs from file @a filename. + /// This version fills (actually, appends to) a supplied vector, avoiding copying, + /// and is hence CPU efficient. + void read(const std::string& filename, std::vector<AnalysisObject*>& aos) { + std::ifstream instream; + instream.open(filename.c_str()); + read(instream, aos); + instream.close(); + } + + /// @brief Read in a collection of objects from output stream @a stream. + /// This version returns a vector by value, involving copying, and is hence less + /// CPU efficient than the alternative version where a vector is filled by reference. + std::vector<AnalysisObject*> read(const std::string& filename) { + std::vector<AnalysisObject*> rtn; + read(filename, rtn); + return rtn; + } + + //@} + + + + protected: + // virtual void readGenericAO(std::istream& stream) = 0; + // // virtual void readHisto(std::istream& stream, const Histo1D& h) = 0; + // // virtual void readProfile(std::istream& stream, const Profile1D& p) = 0; + // virtual void readScatter(std::istream& stream, const Scatter2D& p) = 0; + + }; + + +} + +#endif Added: branches/2011-07-aida2yoda/include/YODA/ReaderYODA.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/ReaderYODA.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,45 @@ +// -*- C++ -*- +// +// This file is part of YODA -- Yet more Objects for Data Analysis +// Copyright (C) 2008-2011 The YODA collaboration (see AUTHORS for details) +// +#ifndef YODA_READERYODA_H +#define YODA_READERYODA_H + +#include "YODA/AnalysisObject.h" +#include "YODA/Reader.h" + +#include <vector> +#include <string> +#include <istream> + + +namespace YODA { + + + /// @brief Persistency reader from YODA flat text data format. + class ReaderYODA : public Reader { + public: + + static Reader& create() { + static ReaderYODA _instance; + return _instance; + } + + protected: + void writeHeader(std::ostream& stream); + void writeFooter(std::ostream& stream); + + void writeHisto(std::ostream& stream, const Histo1D& h); + void writeProfile(std::ostream& stream, const Profile1D& p); + void writeScatter2D(std::ostream& stream, const Scatter2D& s); + + private: + ReaderYODA() { } + + }; + + +} + +#endif Added: branches/2011-07-aida2yoda/include/YODA/Scatter3D.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/Scatter3D.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,230 @@ +// -*- C++ -*- +// +// This file is part of YODA -- Yet more Objects for Data Analysis +// Copyright (C) 2008-2011 The YODA collaboration (see AUTHORS for details) +// +#ifndef YODA_SCATTER3D_H +#define YODA_SCATTER3D_H + +#include "YODA/AnalysisObject.h" +#include "YODA/Point3D.h" +#include "YODA/Utils/sortedvector.h" +#include <vector> +#include <set> +#include <string> +#include <utility> + +namespace YODA { + + + class Histo3D; + + + /// A very generic data type which is just a collection of 3D data points with errors + class Scatter3D : public AnalysisObject { + public: + + /// Type of the native Point3D collection + typedef Utils::sortedvector<Point3D> Points; + + + /// @name Constructors + //@{ + + Scatter3D(const std::string& path="", const std::string& title="") + : AnalysisObject("Scatter3D", path, title) + { } + + + Scatter3D(const Points& points, + const std::string& path="", const std::string& title="") + : AnalysisObject("Scatter3D", path, title), + _points(points) + { } + + + /// Values with asymmetric errors on both x and y + Scatter3D(const std::vector<double>& x, const std::vector<double>& y, const std::vector<double>& z, + const std::vector<std::pair<double,double> >& ex, const std::vector<std::pair<double,double> >& ey, const std::vector<std::pair<double,double> >& ez, + const std::string& path="", const std::string& title="") + : AnalysisObject("Scatter3D", path, title) + { + assert(x.size() == y.size() && x.size() == ex.size() && x.size() == ey.size()); + for (size_t i = 0; i < x.size(); ++i) { + addPoint(Point3D(x[i], y[i], z[i], ex[i], ey[i], ez[i])); + } + } + + /// Values with no errors: + Scatter3D(const std::vector<double>& x, const std::vector<double>& y, const std::vector<double>& z, + const std::string& path="", const std::string& title="") { + std::vector<pair<double,double> > null; + null.resize(x.size()); + + for(unsigned int i=0; i < null.size(); i++) null[i] = make_pair(0.0, 0.0); + Scatter3D(x, y, z, null, null, null, path, title); + } + + /// Values with completely explicit asymmetric errors + Scatter3D(const std::vector<double>& x, const std::vector<double>& y, const std::vector<double> z, + const std::vector<double>& exminus, + const std::vector<double>& explus, + const std::vector<double>& eyminus, + const std::vector<double>& eyplus, + const std::vector<double>& ezminus, + const std::vector<double>& ezplus, + const std::string& path="", const std::string& title="") + : AnalysisObject("Scatter3D", path, title) + { + assert(x.size() == y.size() == z.size() && + x.size() == exminus.size() && x.size() == explus.size() && + y.size() == eyminus.size() && y.size() == eyplus.size() && + z.size() == ezminus.size() && z.size() == ezplus.size()); + for (size_t i = 0; i < x.size(); ++i) { + addPoint(Point3D(x[i], exminus[i], explus[i], y[i], eyminus[i], eyplus[i], z[i], ezminus[i], ezplus[i])); + } + } + + //@} + + + /// Clear all points + void reset() { + _points.clear(); + } + + + /// @name Persistency hooks + //@{ + + /// Set the state of the profile object, for unpersisting + /// @todo Need to set annotations (do that on AO), all-histo Dbns, and dbns for every bin. Delegate! + // void _setstate() = 0; + + //@} + + + /////////////////////////////////////////////////// + + /// @name Point accessors + //@{ + + size_t numPoints() const { + return _points.size(); + } + + + const Points& points() const { + return _points; + } + + + Point3D& point(size_t index) { + assert(index < numPoints()); + return _points.at(index); + } + + + const Point3D& point(size_t index) const { + assert(index < numPoints()); + return _points.at(index); + } + + //@} + + + /// @name Point adders + //@{ + + Scatter3D& addPoint(const Point3D& pt) { + _points.insert(pt); + return *this; + } + + Scatter3D& addPoint(double x, double y, double z) { + _points.insert(Point3D(x, y, z)); + return *this; + } + + Scatter3D& addPoint(double x, double y, double z, + std::pair<double,double> ex, std::pair<double,double> ey, std::pair<double,double> ez) { + _points.insert(Point3D(x, y, z, ex, ey, ez)); + return *this; + } + + Scatter3D& addPoint(double x, double exminus, double explus, + double y, double eyminus, double eyplus, + double z, double ezminus, double ezplus) { + _points.insert(Point3D(x, exminus, explus, y, eyminus, eyplus, z, ezminus, ezplus)); + return *this; + } + + Scatter3D& addPoints(Points pts) { + foreach (const Point3D& pt, pts) { + addPoint(pt); + } + return *this; + } + + //@} + + + /// @todo Better name? + Scatter3D& combineWith(const Scatter3D& other) { + addPoints(other.points()); + return *this; + } + + + /// @todo Better name? + /// @todo Convert to accept a Range or generic + Scatter3D& combineWith(const std::vector<Scatter3D>& others) { + foreach (const Scatter3D& s, others) { + combineWith(s); + } + return *this; + } + + + private: + + Points _points; + + std::string _myaotype; + + }; + + + + inline Scatter3D combine(const Scatter3D& a, const Scatter3D& b) { + Scatter3D rtn = a; + rtn.combineWith(b); + return rtn; + } + + + inline Scatter3D combine(const std::vector< Scatter3D >& scatters) { + Scatter3D rtn; + for (std::vector<Scatter3D>::const_iterator s = scatters.begin(); + s != scatters.end(); ++s) { + rtn.combineWith(*s); + } + return rtn; + } + + + ////////////////////////////////// + + + /// @name Conversion functions from other data types + //@{ + + /// Make a Scatter3D representation of a Histo1D + Scatter3D mkScatter(const Histo3D& h); + + //@} + + +} + +#endif Added: branches/2011-07-aida2yoda/include/YODA/WriterMethods.icc ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/WriterMethods.icc Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,87 @@ +// This file contains boilerplate code for static writer functions in all +// classes inheriting from Writer. These methods just forward to the methods on +// the Writer base class, but code duplication can't be avoided without a +// preprocessor hack like this, AFAIK. + + +/// @name Writing a single analysis object. +//@{ + +/// Write out object @a ao to output stream @a stream. +static void write(std::ostream& stream, const AnalysisObject& ao) { + create().write(stream, ao); +} + +/// Write out object @a ao to file @a filename. +static void write(const std::string& filename, const AnalysisObject& ao) { + create().write(filename, ao); +} + +//@} + + +/// @name Writing multiple analysis objects by collection. +//@{ + +/// Write out a collection of objects @a objs to output stream @a stream. +static void write(std::ostream& stream, const std::vector<AnalysisObject*>& aos) { + create().write(stream, aos.begin(), aos.end()); +} +/// Write out a collection of objects @a objs to file @a filename. +static void write(const std::string& filename, const std::vector<AnalysisObject*>& aos) { + create().write(filename, aos.begin(), aos.end()); +} + + +/// Write out a collection of objects @a objs to output stream @a stream. +static void write(std::ostream& stream, const std::list<AnalysisObject*>& aos) { + create().write(stream, aos.begin(), aos.end()); +} +/// Write out a collection of objects @a objs to file @a filename. +static void write(const std::string& filename, const std::list<AnalysisObject*>& aos) { + create().write(filename, aos.begin(), aos.end()); +} + + +/// Write out a collection of objects @a objs to output stream @a stream. +static void write(std::ostream& stream, const std::set<AnalysisObject*>& aos) { + create().write(stream, aos.begin(), aos.end()); +} +/// Write out a collection of objects @a objs to file @a filename. +static void write(const std::string& filename, const std::set<AnalysisObject*>& aos) { + create().write(filename, aos.begin(), aos.end()); +} + + +/// Write out a collection of objects @a objs to output stream @a stream. +static void write(std::ostream& stream, const std::deque<AnalysisObject*>& aos) { + create().write(stream, aos.begin(), aos.end()); +} +/// Write out a collection of objects @a objs to file @a filename. +static void write(const std::string& filename, const std::deque<AnalysisObject*>& aos) { + create().write(filename, aos.begin(), aos.end()); +} + +//@} + + +/// @name Writing multiple analysis objects by iterator range. +//@{ + +/// Write out the objects specified by start iterator @a begin and end +/// iterator @a end to output stream @a stream. +template <typename AOITER> +static void write(std::ostream& stream, + const AOITER& begin, const AOITER& end) { + create().write(stream, begin, end); +} + +/// Write out the objects specified by start iterator @a begin and end +/// iterator @a end to file @a filename. +template <typename AOITER> +static void write(const std::string& filename, + const AOITER& begin, const AOITER& end) { + create().write(filename, begin, end); +} + +//@} Added: branches/2011-07-aida2yoda/include/YODA/WriterYODA.h ============================================================================== --- /dev/null 00:00:00 1970 (empty, because file is newly added) +++ branches/2011-07-aida2yoda/include/YODA/WriterYODA.h Thu Jul 21 09:23:52 2011 (r3230) @@ -0,0 +1,55 @@ +// -*- C++ -*- +// +// This file is part of YODA -- Yet more Objects for Data Analysis +// Copyright (C) 2008-2011 The YODA collaboration (see AUTHORS for details) +// +#ifndef YODA_WRITERYODA_H +#define YODA_WRITERYODA_H + +#include "YODA/AnalysisObject.h" +#include "YODA/Writer.h" + +#include <vector> +#include <string> +#include <ostream> + + +namespace YODA { + + + /// @brief Persistency writer for YODA flat text format. + class WriterYODA : public Writer { + public: + + /// Singleton creation function + static Writer& create() { + static WriterYODA _instance; + return _instance; + } + + + // Include definitions of all write methods (all fulfilled by Writer::write(...)) + #include "YODA/WriterMethods.icc" + + + protected: + + void writeHeader(std::ostream& stream); + void writeFooter(std::ostream& stream); + + void writeHisto1D(std::ostream& stream, const Histo1D& h); + void writeProfile1D(std::ostream& stream, const Profile1D& p); + void writeScatter2D(std::ostream& stream, const Scatter2D& s); + + + private: + + /// Private since it's a singleton. + WriterYODA() { } + + }; + + +} + +#endif
More information about the Rivet-svn mailing list |