Skip to content

Commit

Permalink
feat: support for writing hist derived profiles (#1000)
Browse files Browse the repository at this point in the history 
* Adds way to write TProfiles to file, but misses computations of fTsumwx, fTsumwx2, fTsumwy, fTsumw2.

* style: pre-commit fixes

* Fixes test name and imports.

* Changes data and values for hist and adds fSumw2 to metadata - needed by root to not throw segmentation error.

* style: pre-commit fixes

* Defines missing vars.

* Fix error throwing.

---------

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
  • Loading branch information
@ioanaif @pre-commit-ci
ioanaif and pre-commit-ci[bot] committed Jan 18, 2024
1 parent 0913f36 commit 4c2a5e5
Show file tree
Hide file tree
Showing 3 changed files with 247 additions and 87 deletions.
4 changes: 3 additions & 1 deletion src/uproot/behaviors/TProfile.py
View file
Expand Up @@ -297,7 +297,8 @@ def to_boost(self, metadata=boost_metadata, axis_metadata=boost_axis_metadata):
boost_histogram = uproot.extras.boost_histogram()

effective_counts = self.counts(flow=True)
values, errors = self._values_errors(True, self.member("fErrorMode"))
_, errors = self._values_errors(True, self.member("fErrorMode"))
values = self._bases[0]._bases[-1]
variances = numpy.square(errors)
sum_of_bin_weights = numpy.asarray(self.member("fBinEntries"))

Expand All @@ -314,6 +315,7 @@ def to_boost(self, metadata=boost_metadata, axis_metadata=boost_axis_metadata):
variances = variances[1:]
sum_of_bin_weights = sum_of_bin_weights[1:]

out.metadata = {"fSumw2": self.member("fSumw2")}
view = out.view(flow=True)

# https://github.com/root-project/root/blob/ffc7c588ac91aca30e75d356ea971129ee6a836a/hist/hist/src/TProfileHelper.h#L668-L671
Expand Down
218 changes: 132 additions & 86 deletions src/uproot/writing/identify.py
View file
Expand Up @@ -245,11 +245,7 @@ def to_writable(obj):
):
import boost_histogram

if obj.kind == "MEAN":
raise NotImplementedError(
"PlottableHistogram with kind='MEAN' (i.e. profile plots) not supported yet"
)
elif obj.kind != "COUNT":
if obj.kind != "COUNT" and obj.kind != "MEAN":
raise ValueError(
"PlottableHistogram can only be converted to ROOT TH* if kind='COUNT' or 'MEAN'"
)
Expand Down Expand Up @@ -347,91 +343,141 @@ def to_writable(obj):

# make TH1, TH2, TH3 types independently
if len(axes) == 1:
fTsumw, fTsumw2, fTsumwx, fTsumwx2 = _root_stats_1d(
obj.values(flow=False), obj.axes[0].edges
)
return to_TH1x(
fName=None,
fTitle=title,
data=data,
fEntries=fEntries,
fTsumw=fTsumw,
fTsumw2=fTsumw2,
fTsumwx=fTsumwx,
fTsumwx2=fTsumwx2,
fSumw2=fSumw2,
fXaxis=axes[0],
)
if obj.kind == "MEAN":
if hasattr(obj, "storage_type"):
if "fSumw2" in obj.metadata.keys():
fSumw2 = obj.metadata["fSumw2"]
else:
raise ValueError(f"fSumw2 has not been set for {obj}")
return to_TProfile(
fName=None,
fTitle=title,
data=obj.values(flow=True),
fEntries=obj.size + 1,
fTsumw=obj.sum()["sum_of_weights"],
fTsumw2=obj.sum()["sum_of_weights_squared"],
fTsumwx=0,
fTsumwx2=0,
fTsumwy=0,
fTsumwy2=0,
fSumw2=fSumw2,
fBinEntries=obj.counts(flow=True),
fBinSumw2=numpy.asarray([], numpy.float64),
fXaxis=axes[0],
)
else:
return to_TProfile(
fName=None,
fTitle=title,
data=obj._bases[0]._bases[-1],
fEntries=obj.member("fEntries"),
fTsumw=obj.member("fTsumw"),
fTsumw2=obj.member("fTsumw2"),
fTsumwx=obj.member("fTsumwx"),
fTsumwx2=obj.member("fTsumwx2"),
fTsumwy=obj.member("fTsumwy"),
fTsumwy2=obj.member("fTsumwy2"),
fSumw2=obj.member("fSumw2"),
fBinEntries=obj.member("fBinEntries"),
fBinSumw2=obj.member("fBinSumw2"),
fXaxis=axes[0],
)
else:
fTsumw, fTsumw2, fTsumwx, fTsumwx2 = _root_stats_1d(
obj.values(flow=False), obj.axes[0].edges
)
return to_TH1x(
fName=None,
fTitle=title,
data=data,
fEntries=fEntries,
fTsumw=fTsumw,
fTsumw2=fTsumw2,
fTsumwx=fTsumwx,
fTsumwx2=fTsumwx2,
fSumw2=fSumw2,
fXaxis=axes[0],
)

elif len(axes) == 2:
(
fTsumw,
fTsumw2,
fTsumwx,
fTsumwx2,
fTsumwy,
fTsumwy2,
fTsumwxy,
) = _root_stats_2d(
obj.values(flow=False), obj.axes[0].edges, obj.axes[1].edges
)
return to_TH2x(
fName=None,
fTitle=title,
data=data,
fEntries=fEntries,
fTsumw=fTsumw,
fTsumw2=fTsumw2,
fTsumwx=fTsumwx,
fTsumwx2=fTsumwx2,
fTsumwy=fTsumwy,
fTsumwy2=fTsumwy2,
fTsumwxy=fTsumwxy,
fSumw2=fSumw2,
fXaxis=axes[0],
fYaxis=axes[1],
)
if obj.kind == "MEAN":
raise NotImplementedError(
"2D PlottableHistogram with kind='MEAN' (i.e. 2D profile plots) not supported yet"
)
else:
(
fTsumw,
fTsumw2,
fTsumwx,
fTsumwx2,
fTsumwy,
fTsumwy2,
fTsumwxy,
) = _root_stats_2d(
obj.values(flow=False), obj.axes[0].edges, obj.axes[1].edges
)
return to_TH2x(
fName=None,
fTitle=title,
data=data,
fEntries=fEntries,
fTsumw=fTsumw,
fTsumw2=fTsumw2,
fTsumwx=fTsumwx,
fTsumwx2=fTsumwx2,
fTsumwy=fTsumwy,
fTsumwy2=fTsumwy2,
fTsumwxy=fTsumwxy,
fSumw2=fSumw2,
fXaxis=axes[0],
fYaxis=axes[1],
)

elif len(axes) == 3:
(
fTsumw,
fTsumw2,
fTsumwx,
fTsumwx2,
fTsumwy,
fTsumwy2,
fTsumwxy,
fTsumwz,
fTsumwz2,
fTsumwxz,
fTsumwyz,
) = _root_stats_3d(
obj.values(flow=False),
obj.axes[0].edges,
obj.axes[1].edges,
obj.axes[2].edges,
)
return to_TH3x(
fName=None,
fTitle=title,
data=data,
fEntries=fEntries,
fTsumw=fTsumw,
fTsumw2=fTsumw2,
fTsumwx=fTsumwx,
fTsumwx2=fTsumwx2,
fTsumwy=fTsumwy,
fTsumwy2=fTsumwy2,
fTsumwxy=fTsumwxy,
fTsumwz=fTsumwz,
fTsumwz2=fTsumwz2,
fTsumwxz=fTsumwxz,
fTsumwyz=fTsumwyz,
fSumw2=fSumw2,
fXaxis=axes[0],
fYaxis=axes[1],
fZaxis=axes[2],
)
if obj.kind == "MEAN":
raise NotImplementedError(
"3D PlottableHistogram with kind='MEAN' (i.e. 3D profile plots) not supported yet"
)
else:
(
fTsumw,
fTsumw2,
fTsumwx,
fTsumwx2,
fTsumwy,
fTsumwy2,
fTsumwxy,
fTsumwz,
fTsumwz2,
fTsumwxz,
fTsumwyz,
) = _root_stats_3d(
obj.values(flow=False),
obj.axes[0].edges,
obj.axes[1].edges,
obj.axes[2].edges,
)
return to_TH3x(
fName=None,
fTitle=title,
data=data,
fEntries=fEntries,
fTsumw=fTsumw,
fTsumw2=fTsumw2,
fTsumwx=fTsumwx,
fTsumwx2=fTsumwx2,
fTsumwy=fTsumwy,
fTsumwy2=fTsumwy2,
fTsumwxy=fTsumwxy,
fTsumwz=fTsumwz,
fTsumwz2=fTsumwz2,
fTsumwxz=fTsumwxz,
fTsumwyz=fTsumwyz,
fSumw2=fSumw2,
fXaxis=axes[0],
fYaxis=axes[1],
fZaxis=axes[2],
)

elif (
isinstance(obj, (tuple, list))
Expand Down
112 changes: 112 additions & 0 deletions tests/test_1000-write-TProfiles.py
View file
@@ -0,0 +1,112 @@
# BSD 3-Clause License; see https://github.com/scikit-hep/uproot5/blob/main/LICENSE

import pytest
import uproot
import os
import math
import numpy as np

pytest.importorskip("hist")
ROOT = pytest.importorskip("ROOT")


def test_write_TProfile(tmp_path):
newfile = os.path.join(tmp_path, "newfile.root")

h1 = ROOT.TProfile("h1", "title", 2, -3.14, 2.71)
h1.Fill(-4, 10)
h1.Fill(-3.1, 10)
h1.Fill(-3.1, 20)
h1.Fill(2.7, 20)
h1.Fill(3, 20)

hhist = uproot.from_pyroot(h1).to_hist()

uhist = uproot.writing.identify.to_TProfile(
fName="h1",
fTitle="title",
data=np.array([10, 30, 20, 20], np.float64),
fEntries=5.0,
fTsumw=3.0,
fTsumw2=3.0,
fTsumwx=-3.5,
fTsumwx2=26.51,
fTsumwy=50.0,
fTsumwy2=900.0,
fSumw2=np.array([100, 500, 400, 400], np.float64),
fBinEntries=np.array([1, 2, 1, 1], np.float64),
fBinSumw2=np.array([], np.float64),
fXaxis=uproot.writing.identify.to_TAxis(
fName="xaxis",
fTitle="",
fNbins=2,
fXmin=-3.14,
fXmax=2.71,
),
)

with uproot.recreate(newfile) as fin:
fin["hhist"] = hhist
fin["uhist"] = uhist

f = ROOT.TFile(newfile)
h2 = f.Get("hhist")
h3 = f.Get("uhist")

assert h1.GetEntries() == h2.GetEntries() == h3.GetEntries() == 5
assert h1.GetSumOfWeights() == h2.GetSumOfWeights() == h3.GetSumOfWeights() == 35
assert (
h1.GetBinLowEdge(1)
== h2.GetBinLowEdge(1)
== h3.GetBinLowEdge(1)
== pytest.approx(-3.14)
)
assert (
h1.GetBinWidth(1)
== h2.GetBinWidth(1)
== h3.GetBinWidth(1)
== pytest.approx((2.71 + 3.14) / 2)
)
assert (
h1.GetBinContent(0)
== h2.GetBinContent(0)
== h3.GetBinContent(0)
== pytest.approx(10)
)
assert (
h1.GetBinContent(1)
== h2.GetBinContent(1)
== h3.GetBinContent(1)
== pytest.approx(15)
)
assert (
h1.GetBinContent(2)
== h2.GetBinContent(2)
== h3.GetBinContent(2)
== pytest.approx(20)
)
assert (
h1.GetBinContent(3)
== h2.GetBinContent(3)
== h3.GetBinContent(3)
== pytest.approx(20)
)
assert (
h1.GetBinError(0) == h2.GetBinError(0) == h3.GetBinError(0) == pytest.approx(0)
)
assert (
h1.GetBinError(1)
== h2.GetBinError(1)
== h3.GetBinError(1)
== pytest.approx(np.sqrt(12.5))
)
assert (
h1.GetBinError(2) == h2.GetBinError(2) == h3.GetBinError(2) == pytest.approx(0)
)
assert (
h1.GetBinError(3) == h2.GetBinError(3) == h3.GetBinError(3) == pytest.approx(0)
)

assert hhist[0].variance == pytest.approx(12.5)

f.Close()

0 comments on commit 4c2a5e5

Please sign in to comment.