Skip to content

Commit

Permalink
Manually add a Model for TMatrixTSym<double>. (#484)
Browse files Browse the repository at this point in the history
* Manually add a Model for TMatrixTSym<double>.

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* Add test and fix Cursor index, for use in collections.

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
  • Loading branch information
jpivarski and pre-commit-ci[bot] authored Jun 21, 2022
1 parent c70975e commit e60db1c
Show file tree
Hide file tree
Showing 4 changed files with 109 additions and 0 deletions.
2 changes: 2 additions & 0 deletions src/uproot/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -144,6 +144,8 @@
import uproot.models.RNTuple
import uproot.models.TH
import uproot.models.TGraph
import uproot.models.TMatrixT

from uproot.models.TTree import num_entries

from uproot.containers import STLVector
Expand Down
1 change: 1 addition & 0 deletions src/uproot/model.py
Original file line number Diff line number Diff line change
Expand Up @@ -102,6 +102,7 @@ def reset_classes():
reload(uproot.models.RNTuple)
reload(uproot.models.TH)
reload(uproot.models.TGraph)
reload(uproot.models.TMatrixT)


_classname_regularize = re.compile(r"\s*(<|>|::)\s*")
Expand Down
88 changes: 88 additions & 0 deletions src/uproot/models/TMatrixT.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,88 @@
# BSD 3-Clause License; see https://github.com/scikit-hep/uproot4/blob/main/LICENSE

"""
This module defines versioned models for ``TLeaf`` and its subclasses.
"""


import struct

import numpy

import uproot
import uproot._util
import uproot.model


class Model_TMatrixTSym_3c_double_3e__v5(uproot.model.VersionedModel):
"""
A :doc:`uproot.model.VersionedModel` for ``TMatrixTSym<double>`` version 2,
which shows up as version 5 because it's reading the ``TMatrixTBase<double>``
header.
"""

def read_members(self, chunk, cursor, context, file):
if self.is_memberwise:
raise NotImplementedError(
"memberwise serialization of {}\nin file {}".format(
type(self).__name__, self.file.file_path
)
)
self._bases.append(
file.class_named("TObject", 1).read(
chunk,
cursor,
context,
file,
self._file,
self._parent,
concrete=self.concrete,
)
)
(
self._members["fNrows"],
self._members["fNcols"],
self._members["fRowLwb"],
self._members["fColLwb"],
self._members["fNelems"],
self._members["fNrowIndex"],
self._members["fTol"],
) = cursor.fields(chunk, self._format0, context)

num_elements = self.member("fNrows") * (self.member("fNcols") + 1) // 2
self._members["fElements"] = cursor.array(
chunk, num_elements, self._dtype0, context
)
self._num_bytes += self._members["fElements"].nbytes

_format0 = struct.Struct(">iiiiiid")
_format_memberwise0 = struct.Struct(">i")
_format_memberwise1 = struct.Struct(">i")
_format_memberwise2 = struct.Struct(">i")
_format_memberwise3 = struct.Struct(">i")
_format_memberwise4 = struct.Struct(">i")
_format_memberwise5 = struct.Struct(">i")
_format_memberwise6 = struct.Struct(">d")
_dtype0 = numpy.dtype(">f8")
base_names_versions = [("TObject", 1)]
member_names = [
"fNrows",
"fNcols",
"fRowLwb",
"fColLwb",
"fNelems",
"fNrowIndex",
"fTol",
]
class_flags = {}


class Model_TMatrixTSym_3c_double_3e_(uproot.model.DispatchByVersion):
"""
A :doc:`uproot.model.DispatchByVersion` for ``TMatrixTSym<double>``.
"""

known_versions = {5: Model_TMatrixTSym_3c_double_3e__v5}


uproot.classes["TMatrixTSym<double>"] = Model_TMatrixTSym_3c_double_3e_
18 changes: 18 additions & 0 deletions tests/test_0484-manually-add-model-for-TMatrixTSym_double_.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
# BSD 3-Clause License; see https://github.com/scikit-hep/uproot4/blob/main/LICENSE


import io
import os

import numpy as np
import pytest
import skhep_testdata

import uproot


def test():
with uproot.open(skhep_testdata.data_path("uproot-issue-359.root")) as file:
matrix = file["covmatrixOCratio"]
num_elements = matrix.member("fNrows") * (matrix.member("fNcols") + 1) // 2
assert len(matrix.member("fElements")) == num_elements

0 comments on commit e60db1c

Please sign in to comment.