yoda is hosted by Hepforge, IPPP Durham
YODA - Yet more Objects for Data Analysis  1.7.2
ROOTCnv.h
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // This file is part of YODA -- Yet more Objects for Data Analysis
4 // Copyright (C) 2008-2017 The YODA collaboration (see AUTHORS for details)
5 //
6 #ifndef YODA_ROOTCnv_h
7 #define YODA_ROOTCnv_h
8 
9 #include "YODA/Histo1D.h"
10 #include "YODA/Histo2D.h"
11 #include "YODA/Profile1D.h"
12 
13 #include "TH1.h"
14 #include "TH2.h"
15 #include "TProfile.h"
16 #include "TGraphAsymmErrors.h"
17 #include "TVectorF.h"
18 
19 namespace YODA {
20 
21 
23 
24 
26 
28 
30 
31 
32  // /// @brief Convert a ROOT 1D histogram to a YODA Histo1D
33  // ///
34  // /// Note that ROOT's histograms do not contain enough information to properly rebuild
35  // /// @a x distributions within bins, in underflow and overflow bins, or across the whole histogram.
36  // inline Histo1D toHisto1D(const TH1& th1) {
37  // std::vector<HistoBin1D> bins;
38  // TArrayD sumw2s = *th1.GetSumw2();
39  // Dbn1D dbn_uflow, dbn_oflow;
40  // double sumWtot(0), sumW2tot(0)
41  // for (int i = 0; i =< th1.GetNbinsX()+1; ++i) {
42  // Dbn1D dbn(static_cast<unsigned long>(th1.GetBinContent(i)), th1.GetBinContent(i), sumw2s[i], 0, 0);
43  // // th1.GetBinContent(i)*th1.GetBinCenter(i), th1.GetBinContent(i)*sqr(th1.GetBinCenter(i)));
44  // if (i == 0) dbn_uflow = dbn;
45  // else if (i == th1.GetNbinsX()+1) dbn_oflow = dbn;
46  // else bins.push_back(HistoBin1D(std::make_pair(th1.GetBinLowEdge(i), th1.GetBinLowEdge(i+1)), dbn));
47  // sumWtot += th1.GetBinContent(i);
48  // sumW2tot += sumw2s[i];
49  // }
50  // Dbn1D dbn_tot(static_cast<unsigned long>(th1.GetEntries()), sumWtot, sumW2tot, 0, 0);
51 
52  // Histo1D rtn(bins, dbn_tot, dbn_uflow, const Dbn1D& dbn_oflow, th1.GetName(), th1.GetTitle());
53  // rtn.addAnnotation("XLabel", th1.GetXaxis->GetTitle());
54  // rtn.addAnnotation("YLabel", th1.GetYaxis->GetTitle());
55  // return rtn;
56  // }
57 
58 
59  // /// @brief Convert a ROOT 1D histogram to a YODA Histo1D
60  // ///
61  // /// Note that ROOT's histograms do not contain enough information to properly rebuild
62  // /// @a x distributions within bins, in underflow and overflow bins, or across the whole histogram.
63  // inline Histo1D toHisto1D(const TH1* th1) {
64  // return toHisto1D(*th1);
65  // }
66 
67 
69 
70 
74  inline Scatter2D toScatter2D(const TH1& th1, bool scalebywidth=true) {
75  Scatter2D rtn;
76  for (int ix = 1; ix <= th1.GetNbinsX(); ++ix) {
77  const TAxis& axx = *th1.GetXaxis();
78  const double x = axx.GetBinCenter(ix);
79  const double exminus = x - axx.GetBinLowEdge(ix);
80  const double explus = axx.GetBinUpEdge(ix) - x;
81  const double xwidth = axx.GetBinWidth(ix);
82  //
83  const double val = th1.GetBinContent(ix) / (scalebywidth ? xwidth : 1);
84  const double evalminus = th1.GetBinErrorLow(ix) / (scalebywidth ? xwidth : 1);
85  const double evalplus = th1.GetBinErrorUp(ix) / (scalebywidth ? xwidth : 1);
86  rtn.addPoint(x, val,
87  exminus, explus,
88  evalminus, evalplus);
89  }
90  rtn.setPath(th1.GetName());
91  rtn.setTitle(th1.GetTitle());
92  rtn.addAnnotation("XLabel", th1.GetXaxis()->GetTitle());
93  rtn.addAnnotation("YLabel", th1.GetYaxis()->GetTitle());
94  return rtn;
95  }
96 
97 
99  inline Scatter2D toScatter2D(const TH1* th1, bool scalebywidth=true) {
100  if (th1 == NULL) throw UserError("Null TH1 pointer passed as argument");
101  return toScatter2D(*th1, scalebywidth);
102  }
103 
104 
106  inline Scatter2D* toNewScatter2D(const TH1& th1, bool scalebywidth=true) {
107  return toScatter2D(th1, scalebywidth).newclone();
108  }
110  inline Scatter2D* toNewScatter2D(const TH1* th1, bool scalebywidth=true) {
111  return toScatter2D(th1, scalebywidth).newclone();
112  }
113 
114 
116 
117 
119  inline Scatter2D toScatter2D(const TProfile& tp1) {
120  return toScatter2D((TH1&)tp1, false);
121  }
122 
124  inline Scatter2D toScatter2D(const TProfile* tp1) {
125  if (tp1 == NULL) throw UserError("Null TProfile pointer passed as argument");
126  return toScatter2D(*tp1); //< TProfile inherits from TH1, so this just uses the function above
127  }
128 
129 
131  inline Scatter2D* toNewScatter2D(const TProfile& tp1) {
132  return toScatter2D(tp1).newclone();
133  }
135  inline Scatter2D* toNewScatter2D(const TProfile* tp1) {
136  return toScatter2D(tp1).newclone();
137  }
138 
139 
141 
142 
146  inline Scatter3D toScatter3D(const TH2& th2, bool scalebyarea=true) {
147  Scatter3D rtn;
148  for (int ix = 1; ix <= th2.GetNbinsX(); ++ix) {
149  for (int iy = 1; iy <= th2.GetNbinsY(); ++iy) {
150  const TAxis& axx = *th2.GetXaxis();
151  const double x = axx.GetBinCenter(ix);
152  const double exminus = x - axx.GetBinLowEdge(ix);
153  const double explus = axx.GetBinUpEdge(ix) - x;
154  const double xwidth = axx.GetBinWidth(ix);
155  //
156  const TAxis& axy = *th2.GetYaxis();
157  const double y = axy.GetBinCenter(iy);
158  const double eyminus = y - axy.GetBinLowEdge(iy);
159  const double eyplus = axy.GetBinUpEdge(iy) - y;
160  const double ywidth = axy.GetBinWidth(iy);
161  //
162  const int ixy = th2.GetBin(ix, iy);
163  const double area = xwidth*ywidth;
164  const double val = th2.GetBinContent(ixy) / (scalebyarea ? area : 1);
165  const double evalminus = th2.GetBinErrorLow(ixy) / (scalebyarea ? area : 1);
166  const double evalplus = th2.GetBinErrorUp(ixy) / (scalebyarea ? area : 1);
167  rtn.addPoint(x, y, val,
168  exminus, explus,
169  eyminus, eyplus,
170  evalminus, evalplus);
171  }
172  }
173  rtn.setPath(th2.GetName());
174  rtn.setTitle(th2.GetTitle());
175  rtn.addAnnotation("XLabel", th2.GetXaxis()->GetTitle());
176  rtn.addAnnotation("YLabel", th2.GetYaxis()->GetTitle());
177  rtn.addAnnotation("ZLabel", th2.GetZaxis()->GetTitle());
178  return rtn;
179  }
180 
181 
183  inline Scatter3D toScatter3D(const TH2* th2, bool scalebyarea=true) {
184  if (th2 == NULL) throw UserError("Null TH2 pointer passed as argument");
185  return toScatter3D(*th2, scalebyarea);
186  }
187 
188 
190  inline Scatter3D* toNewScatter3D(const TH2& th2, bool scalebyarea=true) {
191  return toScatter3D(th2, scalebyarea).newclone();
192  }
194  inline Scatter3D* toNewScatter3D(const TH2* th2, bool scalebyarea=true) {
195  return toScatter3D(th2, scalebyarea).newclone();
196  }
197 
198 
200 
201 
203 
204 
206 
207 
208 
210 
211 
212 
214 
215 
219  inline TH1D toTH1D(const Histo1D& h) {
220  // Work out bin edges first
221  std::vector<double> edges;
222  edges.reserve(h.numBins()+1);
223  edges.push_back(h.bin(0).xMin());
224  for (size_t i = 0; i < h.numBins(); ++i) {
225  const HistoBin1D& b = h.bin(i);
226  if (!fuzzyEquals(edges.back(), b.xMin())) edges.push_back(b.xMin());
227  if (!fuzzyEquals(edges.back(), b.xMax())) edges.push_back(b.xMax());
228  }
229  // Book ROOT histogram
230  TH1D rtn(h.path().c_str(), h.title().c_str(), edges.size()-1, &edges[0]);
231  rtn.Sumw2();
232  TArrayD& sumw2s = *rtn.GetSumw2();
233  for (int i = 1; i <= rtn.GetNbinsX(); ++i) {
234  try {
235  const HistoBin1D& b = h.binAt(rtn.GetBinCenter(i)); // throws if in a gap
236  rtn.SetBinContent(i, b.sumW());
237  sumw2s[i] = b.sumW2();
238  } catch (const Exception& e) { }
239  }
240  // Overflows
241  rtn.SetBinContent(0, h.underflow().sumW());
242  rtn.SetBinContent(rtn.GetNbinsX()+1, h.overflow().sumW());
243  sumw2s[0] = h.underflow().sumW2();
244  sumw2s[rtn.GetNbinsX()+1] = h.overflow().sumW2();
245  // Labels
246  if (h.hasAnnotation("XLabel")) rtn.SetXTitle(h.annotation("XLabel").c_str());
247  if (h.hasAnnotation("YLabel")) rtn.SetYTitle(h.annotation("YLabel").c_str());
248  return rtn;
249  }
250 
252  inline TH1D* toNewTH1D(const Histo1D& h, const std::string& rootname) {
253  return (TH1D*) toTH1D(h).Clone(rootname.c_str());
254  }
255 
256 
257 
261  inline TH2D toTH2D(const Histo2D& h) {
262  // Work out bin edges first
263  std::vector<double> xedges, yedges;
264  xedges.reserve(h.numBins()+1);
265  yedges.reserve(h.numBins()+1);
266  xedges.push_back(h.bin(0).xMin());
267  yedges.push_back(h.bin(0).yMin());
268  for (size_t i = 0; i < h.numBins(); ++i) {
269  const HistoBin2D& b = h.bin(i);
270  xedges.push_back(b.xMin());
271  xedges.push_back(b.xMax());
272  yedges.push_back(b.yMin());
273  yedges.push_back(b.yMax());
274  }
275  // Sort and remove (fuzzy) duplicate edges
276  std::sort(xedges.begin(), xedges.end());
277  std::sort(yedges.begin(), yedges.end());
278  const std::vector<double>::iterator xlast = std::unique(xedges.begin(), xedges.end());
279  const std::vector<double>::iterator ylast = std::unique(yedges.begin(), yedges.end());
280  xedges.erase(xlast, xedges.end());
281  yedges.erase(ylast, yedges.end());
282 
283  // Book ROOT histogram
284  TH2D rtn(h.path().c_str(), h.title().c_str(), xedges.size()-1, &xedges[0], yedges.size()-1, &yedges[0]);
285  rtn.Sumw2();
286  TArrayD& sumw2s = *rtn.GetSumw2();
287  for (int ix = 1; ix <= rtn.GetNbinsX(); ++ix) {
288  for (int iy = 1; iy <= rtn.GetNbinsY(); ++iy) {
289  const int i = rtn.GetBin(ix, iy);
290  try {
291  //const HistoBin2D& b = h.binAt(rtn.GetBinCenter(ix), rtn.GetBinCenter(iy)); // throws if in a gap
292  const HistoBin2D& b = h.binAt(rtn.GetXaxis()->GetBinCenter(ix), rtn.GetYaxis()->GetBinCenter(iy)); // throws if in a gap
293  rtn.SetBinContent(i, b.sumW());
294  sumw2s[i] = b.sumW2();
295  } catch (const Exception& e) { }
296  }
297  }
298  // Overflows
300  // rtn.SetBinContent(0, h.underflow().sumW());
301  // rtn.SetBinContent(rtn.GetNbinsX()+1, h.overflow().sumW());
302  // sumw2s[0] = h.underflow().sumW2();
303  // sumw2s[rtn.GetNbinsX()+1] = h.overflow().sumW2();
304  // Labels
305  if (h.hasAnnotation("XLabel")) rtn.SetXTitle(h.annotation("XLabel").c_str());
306  if (h.hasAnnotation("YLabel")) rtn.SetYTitle(h.annotation("YLabel").c_str());
307  if (h.hasAnnotation("ZLabel")) rtn.SetZTitle(h.annotation("ZLabel").c_str());
308  return rtn;
309  }
310 
312  inline TH2D* toNewTH2D(const Histo2D& h, const std::string& rootname) {
313  return (TH2D*) toTH2D(h).Clone(rootname.c_str());
314  }
315 
316 
317 
319 
320 
321 
327  inline TProfile toTProfile(const Profile1D& p) {
328  // Work out bin edges first
329  std::vector<double> edges;
330  edges.reserve(p.numBins()+1);
331  edges.push_back(p.bin(0).xMin());
332  for (size_t i = 0; i < p.numBins(); ++i) {
333  const ProfileBin1D& b = p.bin(i);
334  if (!fuzzyEquals(edges.back(), b.xMin())) edges.push_back(b.xMin());
335  if (!fuzzyEquals(edges.back(), b.xMax())) edges.push_back(b.xMax());
336  }
337  // Book ROOT histogram
338  TProfile rtn(p.path().c_str(), p.title().c_str(), edges.size()-1, &edges[0]);
339  rtn.Sumw2();
340  Double_t* sumwys = rtn.GetArray(); //< YUCK!!!
341  TArrayD& sumwy2s = *rtn.GetSumw2(); //< YUCK!!!
342  for (int i = 1; i <= rtn.GetNbinsX(); ++i) {
343  try {
344  const ProfileBin1D& b = p.binAt(rtn.GetBinCenter(i)); // throws if in a gap
345  // rtn.SetBinContent(i, b.mean());
346  // rtn.SetBinError(i, b.stdErr());
348  // - for sum(y*y): TProfile::GetSumw2()
349  // - for sum(y): TProfile::GetArray()
350  // - for sum(w): TProfile::SetBinEntries(bin, w)
351  // Clearly, the names of accessors/methods are confusing...
352  rtn.SetBinEntries(i, b.sumW());
353  sumwys[i] = b.sumWY();
354  sumwy2s[i] = b.sumWY2();
355  } catch (const Exception& e) { }
356  }
357  // Overflows
358  rtn.SetBinEntries(0, p.underflow().sumW());
359  rtn.SetBinEntries(rtn.GetNbinsX()+1, p.overflow().sumW());
360  sumwys[0] = p.underflow().sumWY();
361  sumwys[0] = p.underflow().sumWY();
362  sumwy2s[rtn.GetNbinsX()+1] = p.overflow().sumWY2();
363  sumwy2s[rtn.GetNbinsX()+1] = p.overflow().sumWY2();
364  // Labels
365  if (p.hasAnnotation("XLabel")) rtn.SetXTitle(p.annotation("XLabel").c_str());
366  if (p.hasAnnotation("YLabel")) rtn.SetYTitle(p.annotation("YLabel").c_str());
367  return rtn;
368  }
369 
370 
371 
373  inline TH1D toTH1D(const Profile1D& p) {
374  // Work out bin edges first
375  std::vector<double> edges;
376  edges.reserve(p.numBins()+1);
377  edges.push_back(p.bin(0).xMin());
378  for (size_t i = 0; i < p.numBins(); ++i) {
379  const ProfileBin1D& b = p.bin(i);
380  if (!fuzzyEquals(edges.back(), b.xMin())) edges.push_back(b.xMin());
381  if (!fuzzyEquals(edges.back(), b.xMax())) edges.push_back(b.xMax());
382  }
383  // Book ROOT histogram
384  TH1D rtn(p.path().c_str(), p.title().c_str(), edges.size()-1, &edges[0]);
385  for (int i = 1; i <= rtn.GetNbinsX(); ++i) {
386  try {
387  const ProfileBin1D& b = p.binAt(rtn.GetBinCenter(i)); // throws if in a gap
388  rtn.SetBinContent(i, b.mean());
389  rtn.SetBinError(i, b.stdErr());
390  } catch (const Exception& e) { }
391  }
392  // Overflows
393  rtn.SetBinContent(0, p.underflow().yMean());
394  rtn.SetBinContent(rtn.GetNbinsX()+1, p.overflow().yMean());
395  rtn.SetBinError(0, p.underflow().yStdErr());
396  rtn.SetBinError(rtn.GetNbinsX()+1, p.overflow().yStdErr());
397  // Labels
398  if (p.hasAnnotation("XLabel")) rtn.SetXTitle(p.annotation("XLabel").c_str());
399  if (p.hasAnnotation("YLabel")) rtn.SetYTitle(p.annotation("YLabel").c_str());
400  return rtn;
401  }
402 
404  inline TProfile* toNewTProfile(const Profile1D& p, const std::string& rootname) {
405  return (TProfile*) toTProfile(p).Clone(rootname.c_str());
406  }
407 
408 
409 
411 
412 
413 
417  inline TGraphAsymmErrors toTGraph(const Scatter2D& s) {
418  TVectorF xs(s.numPoints()), ys(s.numPoints());
419  TVectorF exls(s.numPoints()), exhs(s.numPoints());
420  TVectorF eyls(s.numPoints()), eyhs(s.numPoints());
421  for (size_t i = 0; i < s.numPoints(); ++i) {
422  const Point2D& p = s.point(i);
423  xs[i] = p.x();
424  ys[i] = p.y();
425  exls[i] = p.xErrMinus();
426  exhs[i] = p.xErrPlus();
427  eyls[i] = p.yErrMinus();
428  eyhs[i] = p.yErrPlus();
429  }
430  // Make the ROOT object... mm, the constructors don't take name+title, unlike all this histos!
431  TGraphAsymmErrors rtn(xs, ys, exls, exhs, eyls, eyhs);
432  rtn.SetName(s.path().c_str());
433  rtn.SetTitle(s.title().c_str());
434  // Labels
435  if (s.hasAnnotation("XLabel")) rtn.GetXaxis()->SetTitle(s.annotation("XLabel").c_str());
436  if (s.hasAnnotation("YLabel")) rtn.GetYaxis()->SetTitle(s.annotation("YLabel").c_str());
437  return rtn;
438  }
439 
441  inline TGraphAsymmErrors* toNewTGraph(const Scatter2D& s, const std::string& rootname) {
442  return (TGraphAsymmErrors*) toTGraph(s).Clone(rootname.c_str());
443  }
444 
445 
446 
448  inline TGraphAsymmErrors toTGraph(const Histo1D& h) {
449  return toTGraph(mkScatter(h));
450  }
451 
453  inline TGraphAsymmErrors* toNewTGraph(const Histo1D& h, const std::string& rootname) {
454  return (TGraphAsymmErrors*) toTGraph(h).Clone(rootname.c_str());
455  }
456 
457 
458 
460  inline TGraphAsymmErrors toTGraph(const Profile1D& p) {
461  return toTGraph(mkScatter(p));
462  }
463 
465  inline TGraphAsymmErrors* toNewTGraph(const Profile1D& p, const std::string& rootname) {
466  return (TGraphAsymmErrors*) toTGraph(p).Clone(rootname.c_str());
467  }
468 
469 
470 
472 
473 
474 
476 
477 
478 }
479 
480 #endif
TProfile toTProfile(const Profile1D &p)
Convert a YODA Scatter2D to a ROOT TH1D.
Definition: ROOTCnv.h:327
double sumW2() const
The sum of weights squared.
Definition: Bin2D.h:313
TH1D toTH1D(const Histo1D &h)
Convert a YODA Histo1D to a ROOT 1D histogram.
Definition: ROOTCnv.h:219
double x() const
Get x value.
Definition: Point2D.h:106
double yMax() const
Upper y limit of the bin (exclusive).
Definition: Bin2D.h:163
TProfile * toNewTProfile(const Profile1D &p, const std::string &rootname)
Convert a YODA Profile1D to a ROOT TProfile as new&#39;d pointer.
Definition: ROOTCnv.h:404
Scatter2D * newclone() const
Make a copy on the heap, via &#39;new&#39;.
Definition: Scatter2D.h:127
double xMin() const
Lower x limit of the bin (inclusive).
Definition: Bin2D.h:142
Scatter2D toScatter2D(const TH1 &th1, bool scalebywidth=true)
Definition: ROOTCnv.h:74
double sumWY2() const
The sum of y^2*weight.
Definition: Dbn2D.h:201
const std::string path() const
Get the AO path.
void addPoint(const Point2D &pt)
Insert a new point.
Definition: Scatter2D.h:216
double sumWY() const
The sum of y*weight.
Definition: Dbn2D.h:196
A very generic data type which is just a collection of 2D data points with errors.
Definition: Scatter2D.h:24
const ProfileBin1D & binAt(double x) const
Access a bin by x-coordinate (const version)
Definition: Profile1D.h:250
double y() const
Get y value.
Definition: Point2D.h:112
double yMin() const
Lower y limit of the bin (inclusive).
Definition: Bin2D.h:158
size_t numBins() const
Number of bins on this axis (not counting under/overflow)
Definition: Profile1D.h:219
Generic unspecialised YODA runtime error.
Definition: Exceptions.h:20
void setPath(const std::string &path)
double sumW() const
The sum of weights.
Definition: Dbn1D.h:152
Dbn1D & underflow()
Access underflow (non-const version)
Definition: Histo1D.h:241
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
Dbn2D & overflow()
Access overflow (non-const version)
Definition: Profile1D.h:270
double xMax() const
Upper x limit of the bin (exclusive).
Definition: Bin2D.h:147
Scatter3D * newclone() const
Make a copy on the heap, via &#39;new&#39;.
Definition: Scatter3D.h:139
double sumWY() const
The sum of y*weight.
Definition: ProfileBin1D.h:147
Dbn1D & overflow()
Access overflow (non-const version)
Definition: Histo1D.h:249
A one-dimensional histogram.
Definition: Histo1D.h:26
double yErrPlus(std::string source="") const
Get positive y-error value.
Definition: Point2D.h:223
double yMean() const
Weighted mean, , of distribution.
Definition: Dbn2D.h:133
const HistoBin2D & binAt(double x, double y) const
Access a bin by coordinate (const version)
Definition: Histo2D.h:277
A Bin1D specialised for handling histogram-type information.
Definition: HistoBin1D.h:21
double yErrMinus(std::string source="") const
Get negative y-error value.
Definition: Point2D.h:217
double sumW2() const
The sum of weights squared.
Definition: Bin1D.h:223
Error for problems introduced outside YODA, to put it nicely.
Definition: Exceptions.h:100
Point2D & point(size_t index)
Get a reference to the point with index index (non-const)
Definition: Scatter2D.h:197
size_t numPoints() const
Number of points in the scatter.
Definition: Scatter2D.h:179
TGraphAsymmErrors * toNewTGraph(const Scatter2D &s, const std::string &rootname)
Convert a YODA Scatter2D to a ROOT TGraphAsymmErrors as new&#39;d pointer.
Definition: ROOTCnv.h:441
double sumW() const
The sum of weights.
Definition: Dbn2D.h:176
TH2D * toNewTH2D(const Histo2D &h, const std::string &rootname)
Convert a YODA Histo2D to a ROOT 2D histogram as new&#39;d pointer.
Definition: ROOTCnv.h:312
double yStdErr() const
Weighted standard error on the mean, , of distribution.
Definition: Dbn2D.h:151
void addAnnotation(const std::string &name, const T &value)
Add or set an annotation by name.
const std::string & annotation(const std::string &name) const
Get an annotation by name (as a string)
double stdErr() const
The standard error on the mean.
Definition: ProfileBin1D.h:126
TGraphAsymmErrors toTGraph(const Scatter2D &s)
Convert a YODA Scatter2D to a ROOT TGraphAsymmErrors.
Definition: ROOTCnv.h:417
Scatter1D mkScatter(const Counter &c)
Make a Scatter1D representation of a Histo1D.
Definition: Scatter1D.cc:9
const std::string title() const
Get the AO title.
A Bin1D specialised for handling profile-type information.
Definition: ProfileBin1D.h:23
A very generic data type which is just a collection of 3D data points with errors.
Definition: Scatter3D.h:24
A Bin2D specialised for handling histogram-type information.
Definition: HistoBin2D.h:21
HistoBin1D & bin(size_t index)
Access a bin by index (non-const version)
Definition: Histo1D.h:217
size_t numBins() const
Number of bins on this axis (not counting under/overflow)
Definition: Histo1D.h:195
double sumW() const
The sum of weights.
Definition: Bin2D.h:308
bool hasAnnotation(const std::string &name) const
Check if an annotation is defined.
void setTitle(const std::string &title)
Set the AO title.
double xMin() const
Lower limit of the bin (inclusive).
Definition: Bin1D.h:126
HistoBin2D & bin(size_t index)
Access a bin by index (non-const version)
Definition: Histo2D.h:266
TH2D toTH2D(const Histo2D &h)
Convert a YODA Histo2D to a ROOT 2D histogram.
Definition: ROOTCnv.h:261
double sumW() const
The sum of weights.
Definition: Bin1D.h:218
double sumW2() const
The sum of weights squared.
Definition: Dbn1D.h:157
double mean() const
The mean of the y distribution.
Definition: ProfileBin1D.h:111
const HistoBin1D & binAt(double x) const
Access a bin by coordinate (const version)
Definition: Histo1D.h:229
Scatter2D * toNewScatter2D(const TH1 &th1, bool scalebywidth=true)
Convert a ROOT 1D histogram (excluding TProfile) to a new&#39;d YODA Scatter2D.
Definition: ROOTCnv.h:106
A two-dimensional histogram.
Definition: Histo2D.h:31
Scatter3D * toNewScatter3D(const TH2 &th2, bool scalebyarea=true)
Convert a ROOT 2D histogram to a new&#39;d YODA Scatter3D.
Definition: ROOTCnv.h:190
A one-dimensional profile histogram.
Definition: Profile1D.h:33
double xMax() const
Upper limit of the bin (exclusive).
Definition: Bin1D.h:131
size_t numBins() const
Number of bins.
Definition: Histo2D.h:283
TH1D * toNewTH1D(const Histo1D &h, const std::string &rootname)
Convert a YODA Histo1D to a ROOT 1D histogram as new&#39;d pointer.
Definition: ROOTCnv.h:252
void addPoint(const Point3D &pt)
Insert a new point.
Definition: Scatter3D.h:233
double xErrMinus() const
Get negative x-error value.
Definition: Point2D.h:141
double xErrPlus() const
Get positive x-error value.
Definition: Point2D.h:146
double sumWY2() const
The sum of y^2 * weight.
Definition: ProfileBin1D.h:152
ProfileBin1D & bin(size_t index)
Access a bin by index (non-const version)
Definition: Profile1D.h:242
Dbn2D & underflow()
Access underflow (non-const version)
Definition: Profile1D.h:262
A 2D data point to be contained in a Scatter2D.
Definition: Point2D.h:18
Scatter3D toScatter3D(const TH2 &th2, bool scalebyarea=true)
Definition: ROOTCnv.h:146