No OneTemporary
Actions

Size

46 KB

Subscribers

None

View Options

	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 (11 h, 25 m)
Storage Engine: blob
Storage Format: Raw Data
Storage Handle: 4023032
Default Alt Text: (46 KB)

No OneTemporaryActions

View Options

File Metadata

Event Timeline

No OneTemporary
Actions