Page MenuHomeHEPForge

No OneTemporary

diff --git a/include/YODA/Point1D.h b/include/YODA/Point1D.h
--- a/include/YODA/Point1D.h
+++ b/include/YODA/Point1D.h
@@ -1,366 +1,369 @@
// -*- C++ -*-
//
// This file is part of YODA -- Yet more Objects for Data Analysis
// Copyright (C) 2008-2018 The YODA collaboration (see AUTHORS for details)
//
#ifndef YODA_POINT1D_H
#define YODA_POINT1D_H
#include "YODA/Point.h"
#include "YODA/Exceptions.h"
#include "YODA/Utils/MathUtils.h"
#include <utility>
namespace YODA {
/// A 1D data point to be contained in a Scatter1D
class Point1D : public Point {
public:
/// @name Constructors
//@{
// Default constructor
Point1D() { }
/// Constructor from values with optional symmetric errors
Point1D(double x, double ex=0.0, std::string source="")
: _x(x)
{
_ex[source] = std::make_pair(ex, ex);
}
/// Constructor from values with explicit asymmetric errors
Point1D(double x, double exminus, double explus, std::string source="")
: _x(x)
{
_ex[source] = std::make_pair(exminus, explus);
}
/// Constructor from values with asymmetric errors
Point1D(double x, const std::pair<double,double>& ex, std::string source="")
: _x(x)
{
_ex[source] = ex;
}
/// Copy constructor
Point1D(const Point1D& p)
: _x(p._x), _ex(p._ex)
{
this->setParentAO( p.getParentAO());
}
/// Copy assignment
Point1D& operator = (const Point1D& p) {
_x = p._x;
_ex = p._ex;
this->setParentAO( p.getParentAO());
return *this;
}
//@}
public:
/// Space dimension of the point
size_t dim() { return 1; }
/// @name Value accessors
//@{
/// Get x value
double x() const { return _x; }
/// Set x value
void setX(double x) { _x = x; }
/// @todo Uniform "coords" accessor across all Scatters: returning fixed-size tuple?
//@}
/// @name x error accessors
//@{
/// Get x-error values
const std::pair<double,double>& xErrs( std::string source="") const {
+ if (source!="") getVariationsFromParent();
if (!_ex.count(source)) throw RangeError("xErrs has no such key: "+source);
return _ex.at(source);
}
/// Get negative x-error value
double xErrMinus( std::string source="") const {
+ if (source!="") getVariationsFromParent();
if (!_ex.count(source)) throw RangeError("xErrs has no such key: "+source);
return _ex.at(source).first;
}
/// Get positive x-error value
double xErrPlus( std::string source="") const {
+ if (source!="") getVariationsFromParent();
if (!_ex.count(source)) throw RangeError("xErrs has no such key: "+source);
return _ex.at(source).second;
}
/// Get average x-error value
double xErrAvg( std::string source="") const {
+ if (source!="") getVariationsFromParent();
if (!_ex.count(source)) throw RangeError("xErrs has no such key: "+source);
return (_ex.at(source).first + _ex.at(source).second)/2.0;
}
/// Set negative x error
void setXErrMinus(double exminus, std::string source="") {
- if (source!="") getVariationsFromParent();
if (!_ex.count(source)) _ex[source] = std::make_pair(0.,0.);
_ex.at(source).first = exminus;
}
/// Set positive x error
void setXErrPlus(double explus, std::string source="") {
- if (source!="") getVariationsFromParent();
if (!_ex.count(source)) _ex[source] = std::make_pair(0.,0.);
_ex.at(source).second = explus;
}
/// Set symmetric x error
void setXErr(double ex, std::string source="") {
setXErrMinus(ex, source);
setXErrPlus(ex, source);
}
/// Set symmetric x error (alias)
void setXErrs(double ex, std::string source="") {
setXErr(ex, source);
}
/// Set asymmetric x error
void setXErrs(double exminus, double explus, std::string source="") {
setXErrMinus(exminus, source);
setXErrPlus(explus, source);
}
/// Set asymmetric x error
void setXErrs(const std::pair<double,double>& ex, std::string source="") {
- if (source!="") getVariationsFromParent();
_ex[source] = ex;
}
/// Get value minus negative x-error
double xMin(std::string source="") const {
+ if (source!="") getVariationsFromParent();
if (!_ex.count(source)) throw RangeError("xErrs has no such key: "+source);
return _x - _ex.at(source).first;
}
/// Get value plus positive x-error
double xMax(std::string source="") const {
+ if (source!="") getVariationsFromParent();
if (!_ex.count(source)) throw RangeError("xErrs has no such key: "+source);
return _x + _ex.at(source).second;
}
//@}
/// @name Combined x value and error setters
//@{
/// Set x value and symmetric error
void setX(double x, double ex, std::string source="") {
setX(x);
setXErr(ex, source);
}
/// Set x value and asymmetric error
void setX(double x, double exminus, double explus, std::string source="") {
setX(x);
setXErrs(exminus, explus, source);
}
/// Set x value and asymmetric error
void setX(double x, std::pair<double,double>& ex, std::string source="") {
setX(x);
setXErrs(ex, source);
}
//@}
// @name Manipulations
//@{
/// Scaling of x axis
void scaleX(double scalex) {
setX(x()*scalex);
for (const auto &source : _ex){
setXErrs(xErrMinus()*scalex, xErrPlus()*scalex, source.first);
}
}
//@}
/// @name Integer axis accessor equivalents
//@{
/// Get the point value for direction @a i
double val(size_t i) const {
if (i == 0 || i > 1) throw RangeError("Invalid axis int, must be in range 1..dim");
return x();
}
/// Set the point value for direction @a i
void setVal(size_t i, double val) {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setX(val);
}
/// Get error map for direction @a i
const std::map< std::string, std::pair<double,double>> & errMap() const {
getVariationsFromParent();
return _ex;
}
// Parse the variations from the parent AO if it exists
void getVariationsFromParent() const;
/// Get error values for direction @a i
const std::pair<double,double>& errs(size_t i, std::string source="") const {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
return xErrs(source);
}
/// Get negative error value for direction @a i
double errMinus(size_t i, std::string source="") const {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
return xErrMinus(source);
}
/// Get positive error value for direction @a i
double errPlus(size_t i, std::string source="") const {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
return xErrPlus(source);
}
/// Get average error value for direction @a i
double errAvg(size_t i, std::string source="") const {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
return xErrAvg(source);
}
/// Set negative error for direction @a i
void setErrMinus(size_t i, double eminus, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setXErrMinus(eminus, source);
}
/// Set positive error for direction @a i
void setErrPlus(size_t i, double eplus, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setXErrPlus(eplus, source);
}
/// Set symmetric error for direction @a i
void setErr(size_t i, double e, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setXErr(e, source);
}
/// Set asymmetric error for direction @a i
void setErrs(size_t i, double eminus, double eplus, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setXErrs(eminus, eplus, source);
}
/// Set asymmetric error for direction @a i
void setErrs(size_t i, std::pair<double,double>& e, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setXErrs(e, source);
}
/// Set value and symmetric error for direction @a i
void set(size_t i, double val, double e, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setX(val, e, source);
}
/// Set value and asymmetric error for direction @a i
void set(size_t i, double val, double eminus, double eplus, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setX(val, eminus, eplus, source);
}
/// Set value and asymmetric error for direction @a i
void set(size_t i, double val, std::pair<double,double>& e, std::string source="") {
if (i != 1) throw RangeError("Invalid axis int, must be in range 1..dim");
setX(val, e, source);
}
//@}
protected:
/// @name Value and error variables
//@{
double _x;
// a map of the errors for each source. Nominal stored under ""
// to ensure backward compatibility
std::map< std::string, std::pair<double,double> > _ex;
//@}
};
/// @name Comparison operators
//@{
/// Equality test of x characteristics only
/// @todo Base on a named fuzzyEquals(a,b,tol=1e-3) unbound function
inline bool operator==(const YODA::Point1D& a, const YODA::Point1D& b) {
if (!YODA::fuzzyEquals(a.x(), b.x()) ||
!YODA::fuzzyEquals(a.xErrMinus(), b.xErrMinus()) ||
!YODA::fuzzyEquals(a.xErrPlus(), b.xErrPlus()) ) return false;
return true;
}
/// Equality test of x characteristics only
inline bool operator != (const YODA::Point1D& a, const YODA::Point1D& b) {
return !(a == b);
}
/// Less-than operator used to sort bins by x-ordering
inline bool operator < (const YODA::Point1D& a, const YODA::Point1D& b) {
if (!YODA::fuzzyEquals(a.x(), b.x())) {
return a.x() < b.x();
}
if (!YODA::fuzzyEquals(a.xErrMinus(), b.xErrMinus())) {
return a.xErrMinus() < b.xErrMinus();
}
if (!YODA::fuzzyEquals(a.xErrPlus(), b.xErrPlus())) {
return a.xErrPlus() < b.xErrPlus();
}
return false;
}
/// Less-than-or-equals operator used to sort bins by x-ordering
inline bool operator <= (const YODA::Point1D& a, const YODA::Point1D& b) {
if (a == b) return true;
return a < b;
}
/// Greater-than operator used to sort bins by x-ordering
inline bool operator > (const YODA::Point1D& a, const YODA::Point1D& b) {
return !(a <= b);
}
/// Greater-than-or-equals operator used to sort bins by x-ordering
inline bool operator >= (const YODA::Point1D& a, const YODA::Point1D& b) {
return !(a < b);
}
//@}
}
#endif
diff --git a/include/YODA/Point3D.h b/include/YODA/Point3D.h
--- a/include/YODA/Point3D.h
+++ b/include/YODA/Point3D.h
@@ -1,695 +1,692 @@
// -*- C++ -*-
//
// This file is part of YODA -- Yet more Objects for Data Analysis
// Copyright (C) 2008-2018 The YODA collaboration (see AUTHORS for details)
//
#ifndef YODA_POINT3D_H
#define YODA_POINT3D_H
#include "YODA/Point.h"
#include "YODA/Exceptions.h"
#include "YODA/Utils/MathUtils.h"
#include <utility>
namespace YODA {
/// A 3D data point to be contained in a Scatter3D
class Point3D : public Point {
public:
/// @name Constructors
//@{
// Default constructor
Point3D() { }
/// Constructor from values with optional symmetric errors
Point3D(double x, double y, double z, double ex=0.0, double ey=0.0, double ez=0.0, std::string source="")
: _x(x), _y(y), _z(z)
{
_ex = std::make_pair(ex, ex);
_ey = std::make_pair(ey, ey);
_ez[source] = std::make_pair(ez, ez);
}
/// Constructor from values with explicit asymmetric errors
Point3D(double x, double y, double z,
double exminus,
double explus,
double eyminus,
double eyplus,
double ezminus,
double ezplus, std::string source="")
: _x(x), _y(y), _z(z)
{
_ex = std::make_pair(exminus, explus);
_ey = std::make_pair(eyminus, eyplus);
_ez[source] = std::make_pair(ezminus, ezplus);
}
/// Constructor from asymmetric errors given as vectors
Point3D(double x, double y, double z,
const std::pair<double,double>& ex,
const std::pair<double,double>& ey,
const std::pair<double,double>& ez, std::string source="")
: _x(x), _y(y), _z(z),
_ex(ex), _ey(ey)
{
_ez[source] = ez;
}
/// Copy constructor
Point3D(const Point3D& p)
: _x(p._x), _y(p._y), _z(p._z),
_ex(p._ex), _ey(p._ey), _ez(p._ez)
{
this->setParentAO( p.getParentAO());
}
/// Copy assignment
Point3D& operator = (const Point3D& p) {
_x = p._x;
_y = p._y;
_z = p._z;
_ex = p._ex;
_ey = p._ey;
_ez = p._ez;
this->setParentAO( p.getParentAO());
return *this;
}
//@}
public:
/// Space dimension of the point
size_t dim() { return 3; }
/// @name Value and error accessors
//@{
/// Get x value
double x() const { return _x; }
/// Set x value
void setX(double x) { _x = x; }
/// Get y value
double y() const { return _y; }
/// Set y value
void setY(double y) { _y = y; }
/// Get z value
double z() const { return _z;}
/// Set z value
void setZ(double z) { _z = z;}
/// @todo Uniform "coords" accessor across all Scatters: returning fixed-size tuple?
// /// Get x,y,z value tuple
// triple<double,double,double> xyz() const { return std::triple(_x, _y, _z); }
/// Set x, y and z values
void setXYZ(double x, double y, double z) { setX(x); setY(y); setZ(z); }
// /// Set x and y values
// void setXY(triple<double,double,double> xyz) { setX(xy.first); setY(xy.second); setZ(xy.third); }
//@}
/// @name x error accessors
//@{
/// Get x-error values
const std::pair<double,double>& xErrs() const {
return _ex;
}
/// Get negative x-error value
double xErrMinus() const {
return _ex.first;
}
/// Get positive x-error value
double xErrPlus() const {
return _ex.second;
}
/// Get average x-error value
double xErrAvg() const {
return (_ex.first + _ex.second)/2.0;
}
/// Set negative x error
void setXErrMinus(double exminus) {
_ex.first = exminus;
}
/// Set positive x error
void setXErrPlus(double explus) {
_ex.second = explus;
}
/// Set symmetric x error
void setXErr(double ex) {
setXErrMinus(ex);
setXErrPlus(ex);
}
/// Set symmetric x error (alias)
void setXErrs(double ex) {
setXErr(ex);
}
/// Set asymmetric x error
void setXErrs(double exminus, double explus) {
setXErrMinus(exminus);
setXErrPlus(explus);
}
/// Set asymmetric x error
void setXErrs(const std::pair<double,double>& ex) {
_ex = ex;
}
/// Get value minus negative x-error
double xMin() const {
return _x - _ex.first;
}
/// Get value plus positive x-error
double xMax() const {
return _x + _ex.second;
}
//@}
/// @name y error accessors
//@{
/// Get y-error values
const std::pair<double,double>& yErrs() const {
return _ey;
}
/// Get negative y-error value
double yErrMinus() const {
return _ey.first;
}
/// Get positive y-error value
double yErrPlus() const {
return _ey.second;
}
/// Get average y-error value
double yErrAvg() const {
return (_ey.first + _ey.second)/2.0;
}
/// Set negative y error
void setYErrMinus(double eyminus) {
_ey.first = eyminus;
}
/// Set positive y error
void setYErrPlus(double eyplus) {
_ey.second = eyplus;
}
/// Set symmetric y error
void setYErr(double ey) {
setYErrMinus(ey);
setYErrPlus(ey);
}
/// Set symmetric y error (alias)
void setYErrs(double ey) {
setYErr(ey);
}
/// Set asymmetric y error
void setYErrs(double eyminus, double eyplus) {
setYErrMinus(eyminus);
setYErrPlus(eyplus);
}
/// Set asymmetric y error
void setYErrs(const std::pair<double,double>& ey) {
_ey = ey;
}
/// Get value minus negative y-error
double yMin() const {
return _y - _ey.first;
}
/// Get value plus positive y-error
double yMax() const {
return _y + _ey.second;
}
//@}
/// @name z error accessors
//@{
/// Get z-error values
const std::pair<double,double>& zErrs( std::string source="") const {
if (source!="") getVariationsFromParent();
if (!_ez.count(source)) throw RangeError("zErrs has no such key: "+source);
return _ez.at(source);
}
/// Get negative z-error value
double zErrMinus( std::string source="") const {
if (source!="") getVariationsFromParent();
if (!_ez.count(source)) throw RangeError("zErrs has no such key: "+source);
return _ez.at(source).first;
}
/// Get positive z-error value
double zErrPlus( std::string source="") const {
if (source!="") getVariationsFromParent();
if (!_ez.count(source)) throw RangeError("zErrs has no such key: "+source);
return _ez.at(source).second;
}
/// Get average z-error value
double zErrAvg( std::string source="") const {
if (source!="") getVariationsFromParent();
if (!_ez.count(source)) throw RangeError("zErrs has no such key: "+source);
return (_ez.at(source).first + _ez.at(source).second)/2.0;
}
/// Set negative z error
void setZErrMinus(double ezminus, std::string source="") {
- if (source!="") getVariationsFromParent();
if (!_ez.count(source)) _ez[source] = std::make_pair(0.,0.);
_ez.at(source).first = ezminus;
}
/// Set positive z error
void setZErrPlus(double ezplus, std::string source="") {
- if (source!="") getVariationsFromParent();
if (!_ez.count(source)) _ez[source] = std::make_pair(0.,0.);
_ez.at(source).second = ezplus;
}
/// Set symmetric z error
void setZErr(double ez, std::string source="") {
setZErrMinus(ez, source);
setZErrPlus(ez, source);
}
/// Set symmetric z error (alias)
void setZErrs(double ez, std::string source="") {
setZErr(ez, source);
}
/// Set asymmetric z error
void setZErrs(double ezminus, double ezplus, std::string source="") {
setZErrMinus(ezminus, source);
setZErrPlus(ezplus, source);
}
/// Set asymmetric z error
void setZErrs(const std::pair<double,double>& ez, std::string source="") {
- if (source!="") getVariationsFromParent();
_ez[source] = ez;
}
/// Get value minus negative z-error
double zMin( std::string source="") const {
if (source!="") getVariationsFromParent();
if (!_ez.count(source)) throw RangeError("zErrs has no such key: "+source);
return _z - _ez.at(source).first;
}
/// Get value plus positive z-error
double zMax( std::string source="") const {
if (source!="") getVariationsFromParent();
if (!_ez.count(source)) throw RangeError("zErrs has no such key: "+source);
return _z + _ez.at(source).second;
}
//@}
/// @name Combined x/y value and error setters
//@{
/// Set x value and symmetric error
void setX(double x, double ex) {
setX(x);
setXErrs(ex);
}
/// Set x value and asymmetric error
void setX(double x, double exminus, double explus) {
setX(x);
setXErrs(exminus, explus);
}
/// Set x value and asymmetric error
void setX(double x, std::pair<double,double>& ex) {
setX(x);
setXErrs(ex);
}
/// Set y value and symmetric error
void setY(double y, double ey) {
setY(y);
setYErrs(ey);
}
/// Set y value and asymmetric error
void setY(double y, double eyminus, double eyplus) {
setY(y);
setYErrs(eyminus, eyplus);
}
/// Set y value and asymmetric error
void setY(double y, std::pair<double,double>& ey) {
setY(y);
setYErrs(ey);
}
/// Set z value and symmetric error
void setZ(double z, double ez, std::string source="") {
setZ(z);
setZErrs(ez, source);
}
/// Set z value and asymmetric error
void setZ(double z, double ezminus, double ezplus, std::string source="") {
setZ(z);
setZErrs(ezminus, ezplus, source);
}
/// Set z value and asymmetric error
void setZ(double z, std::pair<double,double>& ez, std::string source="") {
setZ(z);
setZErrs(ez, source);
}
//@}
// @name Manipulations
//@{
/// Scaling of x axis
void scaleX(double scalex) {
setX(x()*scalex);
setXErrs(xErrMinus()*scalex, xErrPlus()*scalex);
}
/// Scaling of y axis
void scaleY(double scaley) {
setY(y()*scaley);
setYErrs(yErrMinus()*scaley, yErrPlus()*scaley);
}
/// Scaling of z axis
void scaleZ(double scalez) {
setZ(z()*scalez);
for (const auto &source : _ez){
setZErrs(zErrMinus()*scalez, zErrPlus()*scalez, source.first);
}
}
/// Scaling of all three axes
void scaleXYZ(double scalex, double scaley, double scalez) {
scaleX(scalex);
scaleY(scaley);
scaleZ(scalez);
}
/// Scaling of both axes
/// @deprecated Use scaleXYZ
void scale(double scalex, double scaley, double scalez) {
scaleXYZ(scalex, scaley, scalez);
}
//@}
/// @name Integer axis accessor equivalents
//@{
/// Get the point value for direction @a i
double val(size_t i) const {
switch (i) {
case 1: return x();
case 2: return y();
case 3: return z();
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set the point value for direction @a i
void setVal(size_t i, double val) {
switch (i) {
case 1: setX(val); break;
case 2: setY(val); break;
case 3: setZ(val); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Get error map for direction @a i
const std::map< std::string, std::pair<double,double>> & errMap() const {
getVariationsFromParent();
return _ez;
}
// Parse the variations from the parent AO if it exists
void getVariationsFromParent() const;
/// Get error values for direction @a i
const std::pair<double,double>& errs(size_t i, std::string source="") const {
switch (i) {
case 1: return xErrs();
case 2: return yErrs();
case 3: return zErrs(source);
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Get negative error value for direction @a i
double errMinus(size_t i, std::string source="") const {
switch (i) {
case 1: return xErrMinus();
case 2: return yErrMinus();
case 3: return zErrMinus(source);
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Get positive error value for direction @a i
double errPlus(size_t i, std::string source="") const {
switch (i) {
case 1: return xErrPlus();
case 2: return yErrPlus();
case 3: return zErrPlus(source);
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Get average error value for direction @a i
double errAvg(size_t i, std::string source="") const {
switch (i) {
case 1: return xErrAvg();
case 2: return yErrAvg();
case 3: return zErrAvg(source);
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set negative error for direction @a i
void setErrMinus(size_t i, double eminus, std::string source="") {
switch (i) {
case 1: setXErrMinus(eminus); break;
case 2: setYErrMinus(eminus); break;
case 3: setZErrMinus(eminus, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set positive error for direction @a i
void setErrPlus(size_t i, double eplus, std::string source="") {
switch (i) {
case 1: setXErrPlus(eplus); break;
case 2: setYErrPlus(eplus); break;
case 3: setZErrPlus(eplus, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set symmetric error for direction @a i
void setErr(size_t i, double e, std::string source="") {
switch (i) {
case 1: setXErrs(e); break;
case 2: setYErrs(e); break;
case 3: setZErrs(e, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set asymmetric error for direction @a i
void setErrs(size_t i, double eminus, double eplus, std::string source="") {
switch (i) {
case 1: setXErrs(eminus, eplus); break;
case 2: setYErrs(eminus, eplus); break;
case 3: setZErrs(eminus, eplus, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set asymmetric error for direction @a i
void setErrs(size_t i, std::pair<double,double>& e, std::string source="") {
switch (i) {
case 1: setXErrs(e); break;
case 2: setYErrs(e); break;
case 3: setZErrs(e, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set value and symmetric error for direction @a i
void set(size_t i, double val, double e, std::string source="") {
switch (i) {
case 1: setX(val, e); break;
case 2: setY(val, e); break;
case 3: setZ(val, e, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set value and asymmetric error for direction @a i
void set(size_t i, double val, double eminus, double eplus, std::string source="") {
switch (i) {
case 1: setX(val, eminus, eplus); break;
case 2: setY(val, eminus, eplus); break;
case 3: setZ(val, eminus, eplus, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
/// Set value and asymmetric error for direction @a i
void set(size_t i, double val, std::pair<double,double>& e, std::string source="") {
switch (i) {
case 1: setX(val, e); break;
case 2: setY(val, e); break;
case 3: setZ(val, e, source); break;
default: throw RangeError("Invalid axis int, must be in range 1..dim");
}
}
//@}
protected:
/// @name Value and error variables
//@{
double _x;
double _y;
double _z;
std::pair<double,double> _ex;
std::pair<double,double> _ey;
// a map of the errors for each source. Nominal stored under ""
// to ensure backward compatibility
std::map< std::string, std::pair<double,double> >_ez;
//@}
};
/// @name Comparison operators
//@{
/// Equality test of x, y & z characteristics only
/// @todo Base on a named fuzzyEquals(a,b,tol=1e-3) unbound function
inline bool operator==(const Point3D& a, const YODA::Point3D& b) {
if (!YODA::fuzzyEquals(a.x(), b.x()) ||
!YODA::fuzzyEquals(a.xErrMinus(), b.xErrMinus()) ||
!YODA::fuzzyEquals(a.xErrPlus(), b.xErrPlus()) ) return false;
if (!YODA::fuzzyEquals(a.y(), b.y()) ||
!YODA::fuzzyEquals(a.yErrMinus(), b.yErrMinus()) ||
!YODA::fuzzyEquals(a.yErrPlus(), b.yErrPlus()) ) return false;
if (!YODA::fuzzyEquals(a.z(), b.z()) ||
!YODA::fuzzyEquals(a.zErrMinus(), b.zErrMinus()) ||
!YODA::fuzzyEquals(a.zErrPlus(), b.zErrPlus()) ) return false;
return true;
const bool same_val = fuzzyEquals(a.x(), b.x()) && fuzzyEquals(a.y(), b.y());
const bool same_eminus = fuzzyEquals(a.xErrMinus(), b.xErrMinus()) &&
fuzzyEquals(a.yErrMinus(), b.yErrMinus());
const bool same_eplus = fuzzyEquals(a.xErrPlus(), b.xErrPlus()) &&
fuzzyEquals(a.yErrPlus(), b.yErrPlus());
return same_val && same_eminus && same_eplus;
}
/// Inequality operator
inline bool operator != (const Point3D& a, const YODA::Point3D& b) {
return !(a == b);
}
/// Less-than operator used to sort bins by x-first ordering
inline bool operator < (const Point3D& a, const YODA::Point3D& b) {
if (! fuzzyEquals(a.x(), b.x())) {
return a.x() < b.x();
}
if (!fuzzyEquals(a.y(), b.y())) {
return a.y() < b.y();
}
if (! fuzzyEquals(a.xErrMinus(), b.xErrMinus())) {
return a.xErrMinus() < b.xErrMinus();
}
if (!fuzzyEquals(a.yErrMinus(), b.yErrMinus())) {
return a.yErrMinus() < b.yErrMinus();
}
if (! fuzzyEquals(a.xErrPlus(), b.xErrPlus())) {
return a.xErrPlus() < b.xErrPlus();
}
if (!fuzzyEquals(a.yErrPlus(), b.yErrPlus())) {
return a.yErrPlus() < b.yErrPlus();
}
return false;
}
/// Less-than-or-equals operator
inline bool operator <= (const Point3D& a, const YODA::Point3D& b) {
if (a == b) return true;
return a < b;
}
/// Greater-than operator
inline bool operator > (const Point3D& a, const YODA::Point3D& b) {
return !(a <= b);
}
/// Greater-than-or-equals operator
inline bool operator >= (const Point3D& a, const YODA::Point3D& b) {
return !(a < b);
}
//@}
}
#endif
diff --git a/src/WriterYODA.cc b/src/WriterYODA.cc
--- a/src/WriterYODA.cc
+++ b/src/WriterYODA.cc
@@ -1,385 +1,328 @@
// -*- C++ -*-
//
// This file is part of YODA -- Yet more Objects for Data Analysis
// Copyright (C) 2008-2018 The YODA collaboration (see AUTHORS for details)
//
#include "YODA/WriterYODA.h"
#include "yaml-cpp/yaml.h"
#ifdef YAML_NAMESPACE
#define YAML YAML_NAMESPACE
#endif
#include <iostream>
#include <iomanip>
using namespace std;
namespace YODA {
/// Singleton creation function
Writer& WriterYODA::create() {
static WriterYODA _instance;
_instance.setPrecision(6);
return _instance;
}
// Format version:
// - V1/empty = make-plots annotations style
// - V2 = YAML annotations
static const int YODA_FORMAT_VERSION = 2;
// Version-formatting helper function
inline string _iotypestr(const string& baseiotype) {
ostringstream os;
os << "YODA_" << Utils::toUpper(baseiotype) << "_V" << YODA_FORMAT_VERSION;
return os.str();
}
void WriterYODA::_writeAnnotations(std::ostream& os, const AnalysisObject& ao) {
os << scientific << setprecision(_precision);
for (const string& a : ao.annotations()) {
if (a.empty()) continue;
/// @todo Write out floating point annotations as scientific notation
string ann = ao.annotation(a);
// remove stpurious line returns at the end of a string so that we don't
// end up with two line returns.
ann.erase(std::remove(ann.begin(), ann.end(), '\n'), ann.end());
os << a << ": " << ann << "\n";
}
os << "---\n";
}
void WriterYODA::writeCounter(std::ostream& os, const Counter& c) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
os << "BEGIN " << _iotypestr("COUNTER") << " " << c.path() << "\n";
_writeAnnotations(os, c);
os << "# sumW\t sumW2\t numEntries\n";
os << c.sumW() << "\t" << c.sumW2() << "\t" << c.numEntries() << "\n";
os << "END " << _iotypestr("COUNTER") << "\n\n";
os.flags(oldflags);
}
void WriterYODA::writeHisto1D(std::ostream& os, const Histo1D& h) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
os << "BEGIN " << _iotypestr("HISTO1D") << " " << h.path() << "\n";
_writeAnnotations(os, h);
try {
//if ( h.totalDbn().effNumEntries() > 0 ) {
os << "# Mean: " << h.xMean() << "\n";
os << "# Area: " << h.integral() << "\n";
} catch (LowStatsError& e) {
//
}
os << "# ID\t ID\t sumw\t sumw2\t sumwx\t sumwx2\t numEntries\n";
os << "Total \tTotal \t";
os << h.totalDbn().sumW() << "\t" << h.totalDbn().sumW2() << "\t";
os << h.totalDbn().sumWX() << "\t" << h.totalDbn().sumWX2() << "\t";
os << h.totalDbn().numEntries() << "\n";
os << "Underflow\tUnderflow\t";
os << h.underflow().sumW() << "\t" << h.underflow().sumW2() << "\t";
os << h.underflow().sumWX() << "\t" << h.underflow().sumWX2() << "\t";
os << h.underflow().numEntries() << "\n";
os << "Overflow\tOverflow\t";
os << h.overflow().sumW() << "\t" << h.overflow().sumW2() << "\t";
os << h.overflow().sumWX() << "\t" << h.overflow().sumWX2() << "\t";
os << h.overflow().numEntries() << "\n";
os << "# xlow\t xhigh\t sumw\t sumw2\t sumwx\t sumwx2\t numEntries\n";
for (const HistoBin1D& b : h.bins()) {
os << b.xMin() << "\t" << b.xMax() << "\t";
os << b.sumW() << "\t" << b.sumW2() << "\t";
os << b.sumWX() << "\t" << b.sumWX2() << "\t";
os << b.numEntries() << "\n";
}
os << "END " << _iotypestr("HISTO1D") << "\n\n";
os.flags(oldflags);
}
void WriterYODA::writeHisto2D(std::ostream& os, const Histo2D& h) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
os << "BEGIN " << _iotypestr("HISTO2D") << " " << h.path() << "\n";
_writeAnnotations(os, h);
try {
//if ( h.totalDbn().numEntries() > 0 )
os << "# Mean: (" << h.xMean() << ", " << h.yMean() << ")\n";
os << "# Volume: " << h.integral() << "\n";
} catch (LowStatsError& e) {
//
}
os << "# ID\t ID\t sumw\t sumw2\t sumwx\t sumwx2\t sumwy\t sumwy2\t sumwxy\t numEntries\n";
// Total distribution
const Dbn2D& td = h.totalDbn();
os << "Total \tTotal \t";
os << td.sumW() << "\t" << td.sumW2() << "\t";
os << td.sumWX() << "\t" << td.sumWX2() << "\t";
os << td.sumWY() << "\t" << td.sumWY2() << "\t";
os << td.sumWXY() << "\t";
os << td.numEntries() << "\n";
// Outflows
/// @todo Disabled for now, reinstate with a *full* set of outflow info to allow marginalisation
os << "# 2D outflow persistency not currently supported until API is stable\n";
// for (int ix = -1; ix <= 1; ++ix) {
// for (int iy = -1; iy <= 1; ++iy) {
// if (ix == 0 && iy == 0) continue;
// os << "Outflow\t" << ix << ":" << iy << "\t";
// const Dbn2D& d = h.outflow(ix, iy);
// os << d.sumW() << "\t" << d.sumW2() << "\t";
// os << d.sumWX() << "\t" << d.sumWX2() << "\t";
// os << d.sumWY() << "\t" << d.sumWY2() << "\t";
// os << d.sumWXY() << "\t";
// os << d.numEntries() << "\n";
// }
// }
// Bins
os << "# xlow\t xhigh\t ylow\t yhigh\t sumw\t sumw2\t sumwx\t sumwx2\t sumwy\t sumwy2\t sumwxy\t numEntries\n";
for (const HistoBin2D& b : h.bins()) {
os << b.xMin() << "\t" << b.xMax() << "\t";
os << b.yMin() << "\t" << b.yMax() << "\t";
os << b.sumW() << "\t" << b.sumW2() << "\t";
os << b.sumWX() << "\t" << b.sumWX2() << "\t";
os << b.sumWY() << "\t" << b.sumWY2() << "\t";
os << b.sumWXY() << "\t";
os << b.numEntries() << "\n";
}
os << "END " << _iotypestr("HISTO2D") << "\n\n";
os.flags(oldflags);
}
void WriterYODA::writeProfile1D(std::ostream& os, const Profile1D& p) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
os << "BEGIN " << _iotypestr("PROFILE1D") << " " << p.path() << "\n";
_writeAnnotations(os, p);
os << "# ID\t ID\t sumw\t sumw2\t sumwx\t sumwx2\t sumwy\t sumwy2\t numEntries\n";
os << "Total \tTotal \t";
os << p.totalDbn().sumW() << "\t" << p.totalDbn().sumW2() << "\t";
os << p.totalDbn().sumWX() << "\t" << p.totalDbn().sumWX2() << "\t";
os << p.totalDbn().sumWY() << "\t" << p.totalDbn().sumWY2() << "\t";
os << p.totalDbn().numEntries() << "\n";
os << "Underflow\tUnderflow\t";
os << p.underflow().sumW() << "\t" << p.underflow().sumW2() << "\t";
os << p.underflow().sumWX() << "\t" << p.underflow().sumWX2() << "\t";
os << p.underflow().sumWY() << "\t" << p.underflow().sumWY2() << "\t";
os << p.underflow().numEntries() << "\n";
os << "Overflow\tOverflow\t";
os << p.overflow().sumW() << "\t" << p.overflow().sumW2() << "\t";
os << p.overflow().sumWX() << "\t" << p.overflow().sumWX2() << "\t";
os << p.overflow().sumWY() << "\t" << p.overflow().sumWY2() << "\t";
os << p.overflow().numEntries() << "\n";
os << "# xlow\t xhigh\t sumw\t sumw2\t sumwx\t sumwx2\t sumwy\t sumwy2\t numEntries\n";
for (const ProfileBin1D& b : p.bins()) {
os << b.xMin() << "\t" << b.xMax() << "\t";
os << b.sumW() << "\t" << b.sumW2() << "\t";
os << b.sumWX() << "\t" << b.sumWX2() << "\t";
os << b.sumWY() << "\t" << b.sumWY2() << "\t";
os << b.numEntries() << "\n";
}
os << "END " << _iotypestr("PROFILE1D") << "\n\n";
os.flags(oldflags);
}
void WriterYODA::writeProfile2D(std::ostream& os, const Profile2D& p) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
os << "BEGIN " << _iotypestr("PROFILE2D") << " " << p.path() << "\n";
_writeAnnotations(os, p);
os << "# sumw\t sumw2\t sumwx\t sumwx2\t sumwy\t sumwy2\t sumwz\t sumwz2\t sumwxy\t numEntries\n";
// Total distribution
const Dbn3D& td = p.totalDbn();
os << "Total \tTotal \t";
os << td.sumW() << "\t" << td.sumW2() << "\t";
os << td.sumWX() << "\t" << td.sumWX2() << "\t";
os << td.sumWY() << "\t" << td.sumWY2() << "\t";
os << td.sumWZ() << "\t" << td.sumWZ2() << "\t";
os << td.sumWXY() << "\t"; // << td.sumWXZ() << "\t" << td.sumWYZ() << "\t";
os << td.numEntries() << "\n";
// Outflows
/// @todo Disabled for now, reinstate with a *full* set of outflow info to allow marginalisation
os << "# 2D outflow persistency not currently supported until API is stable\n";
// for (int ix = -1; ix <= 1; ++ix) {
// for (int iy = -1; iy <= 1; ++iy) {
// if (ix == 0 && iy == 0) continue;
// os << "Outflow\t" << ix << ":" << iy << "\t";
// const Dbn3D& d = p.outflow(ix, iy);
// os << d.sumW() << "\t" << d.sumW2() << "\t";
// os << d.sumWX() << "\t" << d.sumWX2() << "\t";
// os << d.sumWY() << "\t" << d.sumWY2() << "\t";
// os << d.sumWZ() << "\t" << d.sumWZ2() << "\t";
// os << d.sumWXY() << "\t"; // << d.sumWXZ() << "\t" << d.sumWYZ() << "\t";
// os << d.numEntries() << "\n";
// }
// }
// Bins
os << "# xlow\t xhigh\t ylow\t yhigh\t sumw\t sumw2\t sumwx\t sumwx2\t sumwy\t sumwy2\t sumwz\t sumwz2\t sumwxy\t numEntries\n";
for (const ProfileBin2D& b : p.bins()) {
os << b.xMin() << "\t" << b.xMax() << "\t";
os << b.yMin() << "\t" << b.yMax() << "\t";
os << b.sumW() << "\t" << b.sumW2() << "\t";
os << b.sumWX() << "\t" << b.sumWX2() << "\t";
os << b.sumWY() << "\t" << b.sumWY2() << "\t";
os << b.sumWZ() << "\t" << b.sumWZ2() << "\t";
os << b.sumWXY() << "\t"; // << b.sumWXZ() << "\t" << b.sumWYZ() << "\t";
os << b.numEntries() << "\n";
}
os << "END " << _iotypestr("PROFILE2D") << "\n\n";
os.flags(oldflags);
}
void WriterYODA::writeScatter1D(std::ostream& os, const Scatter1D& s) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
os << "BEGIN " << _iotypestr("SCATTER1D") << " " << s.path() << "\n";
- //first write the Variations, a dummy annotation which
- //contains the additional columns which will be written out
- //for sytematic variations
- YAML::Emitter out;
- out << YAML::Flow ;
- out << s.variations();
- //os << "Variations" << ": " << out.c_str() << "\n";
- // then write the regular annotations
_writeAnnotations(os, s);
- std::vector<std::string> variations= s.variations();
-
//write headers
- std::string headers="# xval\t ";
- for (const auto &source : variations){
- headers+=" xerr-"+source+"\t xerr+"+source+"\t";
- }
+ std::string headers="# xval\t xerr-\t xerr+\t";
os << headers << "\n";
//write points
for (const Point1D& pt : s.points()) {
// fill central value
- os << pt.x();
- // fill errors for variations. The first should always be "" which is nominal.
- // Assumes here that all points in the Scatter have the same
- // variations... if not a range error will get thrown from
- // the point when the user tries to access a variation it
- // doesn't have... @todo maybe better way to do this?
- for (const auto &source : variations){
- os << "\t" << pt.xErrMinus(source) << "\t" << pt.xErrPlus(source) ;
- }
+ os << pt.x() << "\t" << pt.xErrMinus() << "\t" << pt.xErrPlus() ;
os << "\n";
}
os << "END " << _iotypestr("SCATTER1D") << "\n\n";
os << flush;
os.flags(oldflags);
}
void WriterYODA::writeScatter2D(std::ostream& os, const Scatter2D& s) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
os << "BEGIN " << _iotypestr("SCATTER2D") << " " << s.path() << "\n";
- //first write the Variations, a dummy annotation which
- //contains the additional columns which will be written out
- //for sytematic variations
- YAML::Emitter out;
- out << YAML::Flow << YAML::BeginMap;
- int counter=0;
- std::vector<std::string> variations= s.variations();
// write annotations
_writeAnnotations(os, s);
//write headers
/// @todo Change ordering to {vals} {errs} {errs} ...
std::string headers="# xval\t xerr-\t xerr+\t yval\t yerr-\t yerr+\t";
- //for (const auto &source : variations){
- // headers+=" yerr-"+source+"\t yerr+"+source+"\t";
- //}
os << headers << "\n";
//write points
for (const Point2D& pt : s.points()) {
/// @todo Change ordering to {vals} {errs} {errs} ...
// fill central value
os << pt.x() << "\t" << pt.xErrMinus() << "\t" << pt.xErrPlus() << "\t";
- os << pt.y();
- // fill errors for variations. The first should always be "" which is nominal.
- // Assumes here that all points in the Scatter have the same
- // variations... if not a range error will get thrown from
- // the point when the user tries to access a variation it
- // doesn't have... @todo maybe better way to do this?
- //for (const auto &source : variations){
- os << "\t" << pt.yErrMinus() << "\t" << pt.yErrPlus() ;
- // }
+ os << pt.y() << "\t" << pt.yErrMinus() << "\t" << pt.yErrPlus() ;
os << "\n";
}
os << "END " << _iotypestr("SCATTER2D") << "\n\n";
os << flush;
os.flags(oldflags);
}
void WriterYODA::writeScatter3D(std::ostream& os, const Scatter3D& s) {
ios_base::fmtflags oldflags = os.flags();
os << scientific << showpoint << setprecision(_precision);
-
os << "BEGIN " << _iotypestr("SCATTER3D") << " " << s.path() << "\n";
- //first write the Variations, a dummy annotation which
- //contains the additional columns which will be written out
- //for sytematic variations
- YAML::Emitter out;
- out << YAML::Flow ;
- out << s.variations();
- // then write the regular annotations
+ // write the annotations
_writeAnnotations(os, s);
- std::vector<std::string> variations= s.variations();
+ // std::vector<std::string> variations= s.variations();
//write headers
/// @todo Change ordering to {vals} {errs} {errs} ...
std::string headers="# xval\t xerr-\t xerr+\t yval\t yerr-\t yerr+\t zval\t ";
- for (const auto &source : variations){
- headers+=" zerr-"+source+"\t zerr+"+source+"\t";
- }
os << headers << "\n";
//write points
for (const Point3D& pt : s.points()) {
/// @todo Change ordering to {vals} {errs} {errs} ...
// fill central value
os << pt.x() << "\t" << pt.xErrMinus() << "\t" << pt.xErrPlus() << "\t";
os << pt.y() << "\t" << pt.yErrMinus() << "\t" << pt.yErrPlus() << "\t";
- os << pt.z();
- // fill errors for variations. The first should always be "" which is nominal.
- // Assumes here that all points in the Scatter have the same
- // variations... if not a range error will get thrown from
- // the point when the user tries to access a variation it
- // doesn't have... @todo maybe better way to do this?
- for (const auto &source : variations){
- os << "\t" << pt.zErrMinus(source) << "\t" << pt.zErrPlus(source) ;
- }
+ os << pt.z() << "\t" << pt.zErrMinus() << "\t" << pt.zErrPlus() ;
os << "\n";
}
os << "END " << _iotypestr("SCATTER3D") << "\n\n";
os << flush;
os.flags(oldflags);
}
}

File Metadata

Mime Type
text/x-diff
Expires
Sat, Dec 21, 1:53 PM (15 h, 26 m)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4023032
Default Alt Text
(46 KB)

Event Timeline