yoda is hosted by Hepforge, IPPP Durham
YODA - Yet more Objects for Data Analysis  1.7.2
Binning1D.h
Go to the documentation of this file.
1 #ifndef YODA_Binned1D_h
2 #define YODA_Binned1D_h
3 
4 #include "YODA/AnalysisObject.h"
5 #include "YODA/Exceptions.h"
6 #include "YODA/Utils/MathUtils.h"
8 #include <limits>
9 #include <string>
10 
11 namespace YODA {
12 
13 
18  template <typename T>
19  class Binned1D {
20  public:
21 
23 
24 
27  : _locked(false)
28  { }
29 
31  Binned1D(const std::vector<double>& binedges, const std::vector<T>& binvalues=std::vector<T>())
32  : _locked(false)
33  {
34  addBins(binedges);
36  }
37 
39  Binned1D(const std::vector< std::pair<double,double> >& binedges, const std::vector<T>& binvalues=std::vector<T>())
40  : _locked(false)
41  {
42  addBins(binedges);
44  }
45 
47  Binned1D(size_t nbins, double lower, double upper, const std::vector<T>& binvalues=std::vector<T>())
48  : _locked(false)
49  {
50  addBins(linspace(nbins, lower, upper));
51  }
52 
54 
55 
57 
58 
60  size_t numBins() const {
61  return bins().size();
62  }
63 
65  double xMin() const {
66  if (numBins() == 0) throw RangeError("This axis contains no bins and so has no defined range");
67  return bins().front().xMin();
68  }
69 
71  double xMax() const {
72  if (numBins() == 0) throw RangeError("This axis contains no bins and so has no defined range");
73  return bins().back().xMax();
74  }
75 
76 
78  T& valueAtIndex(size_t index) {
79  if (index >= numBins()) throw RangeError("YODA::Histo1D: index out of range!");
80  return _values[index];
81  }
82 
84  ssize_t indexAt(double coord) const {
85  // Yes, this is robust even with an empty axis: there's always at least one outflow
86  return _indexes[_binsearcher.index(coord)];
87  }
88 
90  T& valueAt(double x) {
91  const ssize_t index = binIndexAt(x);
92  if (index == -1) throw RangeError("There is no bin at the specified x");
93  return valueAtIndex(index);
94  }
95 
97  const T& valueAt(double x) const {
98  const ssize_t index = binIndexAt(x);
99  if (index == -1) throw RangeError("There is no bin at the specified x");
100  return valueAtIndex(index);
101  }
102 
106  T valueAt(double x, const T& fallback) {
107  const ssize_t index = binIndexAt(x);
108  if (index == -1) return fallback;
109  return valueAtIndex(index);
110  }
111 
113 
114 
116 
117 
119  void reset() {
120  BOOST_FOREACH(Bin& bin, _bins) bin.reset();
121  _locked = false;
122  }
123 
124 
126  void addBin(double low, double high, const T& value=T()) {
127  Bins newBins(_bins);
128  newBins.push_back(Bin(low, high));
129  _updateAxis(newBins);
130  }
131 
133  void addBins(const std::vector<double>& binedges) {
134  Bins newBins(_bins);
135  if (binedges.size() == 0)
136  return;
137 
138  double low = binedges.front();
139 
140  for (size_t i = 1; i < binedges.size(); i++) {
141  double high = binedges[i];
142  newBins.push_back(Bin(low, high));
143  low = high;
144  }
145 
146  _updateAxis(newBins);
147  }
148 
150  void addBins(const std::vector<std::pair<double, double> >& binpairs) {
151  // Make a copy of the current binning
152  Bins newBins(_bins);
153 
154  // Iterate over given bins
155  for (size_t i = 0; i < binpairs.size(); i++) {
156  std::pair<double, double> b = binpairs[i];
157  newBins.push_back(Bin(b.first, b.second));
158  }
159  _updateAxis(newBins);
160  }
161 
162 
163  void addBins(const Bins& bins) {
164  Bins newBins(_bins);
165  BOOST_FOREACH(const Bin& b, bins) newBins.push_back(b);
166  _updateAxis(newBins);
167  }
168 
169 
171  void eraseBin(const size_t i) {
172  // Might as well erase from the internal bins, as we can guarantee
173  // consistency.
174  if (i >= numBins())
175  throw RangeError("Bin index is out of range");
176 
177  const bool wasLocked = _locked;
178  _locked = false;
179  _bins.erase(_bins.begin() + i);
180  _updateAxis(_bins);
181  _locked = wasLocked;
182  }
183 
184 
186  void eraseBins(const size_t from, const size_t to) {
187  if (from >= numBins())
188  throw RangeError("Initial index out of range");
189  if (to >= numBins())
190  throw RangeError("Final index out of range");
191  if (from > to)
192  throw RangeError("Final index is less than initial index");
193 
194  const bool wasLocked = _locked;
195  _locked = false;
196  _bins.erase(_bins.begin() + from, _bins.begin() + to + 1);
197  _updateAxis(_bins);
198  _locked = wasLocked;
199  }
200 
201 
203  // @todo What if somebody passes in negative scalefactor? (test idea)
204  void scaleX(double scalefactor) {
205  _dbn.scaleX(scalefactor);
206  _underflow.scaleX(scalefactor);
207  _overflow.scaleX(scalefactor);
208  for (size_t i = 0; i < _bins.size(); ++i)
209  _bins[i].scaleX(scalefactor);
210  _updateAxis(_bins);
211  }
212 
214 
215 
217  bool sameBinning(const Binned1D& other) const {
218  if (numBins() != other.numBins()) return false;
219  for (size_t i = 0; i < numBins(); i++)
221  if (!(fuzzyEquals(bin(i).xMin(), other.bin(i).xMin()) && fuzzyEquals(bin(i).xMax(), other.bin(i).xMax())))
222  return false;
223  return true;
224  }
225 
226 
227  private:
228 
230  //
233  void _updateAxis(Bins& bins) {
234  // Ensure that axis is not locked
235  if (_locked) throw LockError("Attempting to update a locked axis");
236 
237  // Define the new cuts and indexes
238  std::vector<double> edges; edges.reserve(bins.size()+1); // Nbins+1 edges
239  std::vector<long> indexes; edges.reserve(bins.size()+2); // Nbins + 2 outflows
240 
241  // Sort the bins
242  std::sort(bins.begin(), bins.end());
243 
244  // Keep a note of the last high edge
245  double last_high = -std::numeric_limits<double>::infinity();
246 
247  // Check for overlaps
248  for (size_t i = 0; i < bins.size(); ++i) {
249  Bin& currentBin = bins[i];
250  const double new_low = currentBin.xMin();
251  const double reldiff = (new_low - last_high) / currentBin.width();
252  if (reldiff < -1e-3) { //< @note If there is a "large" negative gap (i.e. overlap), throw an exception
253  std::stringstream ss;
254  ss << "Bin edges overlap: " << last_high << " -> " << new_low;
255  throw RangeError(ss.str());
256  } else if (reldiff > 1e-3) { //< @note If there is a "large" positive gap, create a bin gap
257  indexes.push_back(-1); // First index will be for underflow
258  edges.push_back(new_low); // Will include first edge
259  }
260 
261  // Bins check that they are not zero or negative width. It's okay for
262  // them to throw an exception here, as we haven't changed anything yet.
263  indexes.push_back(i);
264  edges.push_back(currentBin.xMax());
265 
266  last_high = currentBin.xMax();
267  }
268  indexes.push_back(-1); // Overflow
269 
270  // Everything was okay, so let's make our changes
271  _binsearcher = Utils::BinSearcher(edges);
272  _indexes = indexes;
273  _bins = bins;
274  }
275 
276 
277  private:
278 
280 
281 
283  std::vector<T> _values;
284 
285  // Binsearcher, for searching bins
286  Utils::BinSearcher _binsearcher;
287  // Mapping from binsearcher indices to bin indices (allowing gaps)
288  std::vector<long> _indexes;
289 
291  bool _locked;
292 
294 
295  };
296 
297 
298 }
299 
300 #endif
bool sameBinning(const Binned1D &other) const
Check if two Binned1D objects have the same bin edges.
Definition: Binning1D.h:217
Base class for bins in 1D and 2D histograms.
Definition: Bin.h:20
virtual void reset()=0
Reset this bin.
void scaleX(double scalefactor)
Scale the size of an axis by a factor.
Definition: Binning1D.h:204
Binned1D(size_t nbins, double lower, double upper, const std::vector< T > &binvalues=std::vector< T >())
Constructor with the number of bins and the axis limits.
Definition: Binning1D.h:47
void addBin(double low, double high, const T &value=T())
Add a bin, providing its low and high edge and an optional value.
Definition: Binning1D.h:126
void addBins(const std::vector< std::pair< double, double > > &binpairs)
Add a list of bins as pairs of lowEdge, highEdge.
Definition: Binning1D.h:150
Error for e.g. use of invalid bin ranges.
Definition: Exceptions.h:34
T & valueAt(double x)
Return the bin value at a given coordinate (non-const ref)
Definition: Binning1D.h:90
bool fuzzyEquals(double a, double b, double tolerance=1E-5)
Compare two floating point numbers for equality with a degree of fuzziness.
Definition: MathUtils.h:72
std::vector< double > linspace(size_t nbins, double start, double end, bool include_end=true)
Make a list of nbins + 1 values equally spaced between start and end inclusive.
Definition: MathUtils.h:252
Binned1D(const std::vector< double > &binedges, const std::vector< T > &binvalues=std::vector< T >())
Constructor accepting a list of bin edges and optional corresponding vector of values.
Definition: Binning1D.h:31
const T & valueAt(double x) const
Return the bin value at a given coordinate (const ref)
Definition: Binning1D.h:97
1D binned container of T objects, without outflows
Definition: Binning1D.h:19
void eraseBin(const size_t i)
Remove a bin.
Definition: Binning1D.h:171
void addBins(const Bins &bins)
Definition: Binning1D.h:163
Error for modification of a data object where filling has already begun.
Definition: Exceptions.h:41
Binned1D()
Empty constructor.
Definition: Binning1D.h:26
void addBins(const std::vector< double > &binedges)
Add a contiguous set of bins to an axis, via their list of edges.
Definition: Binning1D.h:133
T valueAt(double x, const T &fallback)
Definition: Binning1D.h:106
Binned1D(const std::vector< std::pair< double, double > > &binedges, const std::vector< T > &binvalues=std::vector< T >())
Constructor accepting a list of bin edge pairs and optional corresponding vector of values...
Definition: Binning1D.h:39
size_t numBins() const
Get the number of bins on the axis.
Definition: Binning1D.h:60
ssize_t indexAt(double coord) const
Returns the index of a bin at a given coord, -1 if no bin matches.
Definition: Binning1D.h:84
void reset()
Reset all the bin statistics on the axis.
Definition: Binning1D.h:119
T & valueAtIndex(size_t index)
Return the bin value at a given index (non-const)
Definition: Binning1D.h:78
void eraseBins(const size_t from, const size_t to)
Remove a bin range.
Definition: Binning1D.h:186
double xMin() const
Return the lowest-value bin edge on the axis.
Definition: Binning1D.h:65
double xMax() const
Return the highest-value bin edge on the axis.
Definition: Binning1D.h:71