From 2c881413933ed7bdbc9444ae433e00e5d6a2c2cd Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Franz=20P=C3=B6schel?= Date: Tue, 1 Aug 2023 15:25:25 +0200 Subject: [PATCH] Add simple Python bindings and an example --- CMakeLists.txt | 2 ++ examples/14_custom_hierarchy.py | 48 ++++++++++++++++++++++++++ src/binding/python/Container.cpp | 6 ++++ src/binding/python/CustomHierarchy.cpp | 36 +++++++++++++++++++ src/binding/python/Iteration.cpp | 22 ++++-------- src/binding/python/openPMD.cpp | 2 ++ 6 files changed, 101 insertions(+), 15 deletions(-) create mode 100755 examples/14_custom_hierarchy.py create mode 100644 src/binding/python/CustomHierarchy.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 672340fc2f..9a13e8135e 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -602,6 +602,7 @@ if(openPMD_HAVE_PYTHON) src/binding/python/BaseRecordComponent.cpp src/binding/python/ChunkInfo.cpp src/binding/python/Container.cpp + src/binding/python/CustomHierarchy.cpp src/binding/python/Dataset.cpp src/binding/python/Datatype.cpp src/binding/python/Error.cpp @@ -773,6 +774,7 @@ set(openPMD_PYTHON_EXAMPLE_NAMES 11_particle_dataframe 12_span_write 13_write_dynamic_configuration + 14_custom_hierarchy ) if(openPMD_USE_INVASIVE_TESTS) diff --git a/examples/14_custom_hierarchy.py b/examples/14_custom_hierarchy.py new file mode 100755 index 0000000000..b3eff208a9 --- /dev/null +++ b/examples/14_custom_hierarchy.py @@ -0,0 +1,48 @@ +import numpy as np +import openpmd_api as io + + +def main(): + if "bp" in io.file_extensions: + filename = "../samples/custom_hierarchy.bp" + else: + filename = "../samples/custom_hierarchy.json" + s = io.Series(filename, io.Access.create) + it = s.write_iterations()[100] + + # write openPMD part + temp = it.meshes["temperature"] + temp.axis_labels = ["x", "y"] + temp.unit_dimension = {io.Unit_Dimension.T: 1} + temp.position = [0.5, 0.5] + temp.grid_spacing = [1, 1] + temp.grid_global_offset = [0, 0] + temp.reset_dataset(io.Dataset(np.dtype("double"), [5, 5])) + temp[()] = np.zeros((5, 5)) + + # write NeXus part + nxentry = it["Scan"] + nxentry.set_attribute("NX_class", "NXentry") + nxentry.set_attribute("default", "data") + + data = nxentry["data"] + data.set_attribute("NX_class", "NXdata") + data.set_attribute("signal", "counts") + data.set_attribute("axes", ["two_theta"]) + data.set_attribute("two_theta_indices", [0]) + + counts = data.as_container_of_datasets()["counts"] + counts.set_attribute("units", "counts") + counts.set_attribute("long_name", "photodiode counts") + counts.reset_dataset(io.Dataset(np.dtype("int"), [15])) + counts[()] = np.zeros(15, dtype=np.dtype("int")) + + two_theta = data.as_container_of_datasets()["two_theta"] + two_theta.set_attribute("units", "degrees") + two_theta.set_attribute("long_name", "two_theta (degrees)") + two_theta.reset_dataset(io.Dataset(np.dtype("double"), [15])) + two_theta[()] = np.zeros(15) + + +if __name__ == "__main__": + main() diff --git a/src/binding/python/Container.cpp b/src/binding/python/Container.cpp index fbdb96815f..6da0b1eb3d 100644 --- a/src/binding/python/Container.cpp +++ b/src/binding/python/Container.cpp @@ -26,6 +26,7 @@ #include +#include "openPMD/CustomHierarchy.hpp" #include "openPMD/Iteration.hpp" #include "openPMD/Mesh.hpp" #include "openPMD/ParticlePatches.hpp" @@ -52,6 +53,7 @@ using PyPatchRecordContainer = Container; using PyRecordComponentContainer = Container; using PyMeshRecordComponentContainer = Container; using PyPatchRecordComponentContainer = Container; +using PyCustomHierarchyContainer = Container; PYBIND11_MAKE_OPAQUE(PyIterationContainer) PYBIND11_MAKE_OPAQUE(PyMeshContainer) PYBIND11_MAKE_OPAQUE(PyPartContainer) @@ -61,6 +63,7 @@ PYBIND11_MAKE_OPAQUE(PyPatchRecordContainer) PYBIND11_MAKE_OPAQUE(PyRecordComponentContainer) PYBIND11_MAKE_OPAQUE(PyMeshRecordComponentContainer) PYBIND11_MAKE_OPAQUE(PyPatchRecordComponentContainer) +PYBIND11_MAKE_OPAQUE(PyCustomHierarchyContainer) void init_Container(py::module &m) { @@ -85,4 +88,7 @@ void init_Container(py::module &m) ::detail::create_and_bind_container< PyPatchRecordComponentContainer, Attributable>(m, "Patch_Record_Component_Container"); + ::detail:: + create_and_bind_container( + m, "Custom_Hierarchy_Container"); } diff --git a/src/binding/python/CustomHierarchy.cpp b/src/binding/python/CustomHierarchy.cpp new file mode 100644 index 0000000000..cf6d37f2ef --- /dev/null +++ b/src/binding/python/CustomHierarchy.cpp @@ -0,0 +1,36 @@ + + +#include "openPMD/CustomHierarchy.hpp" +#include "openPMD/ParticleSpecies.hpp" +#include "openPMD/RecordComponent.hpp" +#include "openPMD/backend/Attributable.hpp" +#include + +namespace py = pybind11; +using namespace openPMD; + +void init_CustomHierarchy(py::module &m) +{ + py::class_, Attributable>( + m, "CustomHierarchy") + .def( + "as_container_of_datasets", + &CustomHierarchy::asContainerOf) + .def("as_container_of_meshes", &CustomHierarchy::asContainerOf) + .def( + "as_container_of_particles", + &CustomHierarchy::asContainerOf) + + .def_readwrite( + "meshes", + &CustomHierarchy::meshes, + py::return_value_policy::copy, + // garbage collection: return value must be freed before Iteration + py::keep_alive<1, 0>()) + .def_readwrite( + "particles", + &CustomHierarchy::particles, + py::return_value_policy::copy, + // garbage collection: return value must be freed before Iteration + py::keep_alive<1, 0>()); +} diff --git a/src/binding/python/Iteration.cpp b/src/binding/python/Iteration.cpp index b54c81c626..7087f672b7 100644 --- a/src/binding/python/Iteration.cpp +++ b/src/binding/python/Iteration.cpp @@ -21,6 +21,7 @@ #include #include +#include "openPMD/CustomHierarchy.hpp" #include "openPMD/Iteration.hpp" #include @@ -32,7 +33,11 @@ using namespace openPMD; void init_Iteration(py::module &m) { - py::class_(m, "Iteration") + py::class_< + Iteration, + CustomHierarchy, + Container, + Attributable>(m, "Iteration") .def(py::init()) .def( @@ -83,18 +88,5 @@ void init_Iteration(py::module &m) // TODO remove in future versions (deprecated) .def("set_time", &Iteration::setTime) .def("set_dt", &Iteration::setDt) - .def("set_time_unit_SI", &Iteration::setTimeUnitSI) - - .def_readwrite( - "meshes", - &Iteration::meshes, - py::return_value_policy::copy, - // garbage collection: return value must be freed before Iteration - py::keep_alive<1, 0>()) - .def_readwrite( - "particles", - &Iteration::particles, - py::return_value_policy::copy, - // garbage collection: return value must be freed before Iteration - py::keep_alive<1, 0>()); + .def("set_time_unit_SI", &Iteration::setTimeUnitSI); } diff --git a/src/binding/python/openPMD.cpp b/src/binding/python/openPMD.cpp index e784d26f3a..e71a55ccee 100644 --- a/src/binding/python/openPMD.cpp +++ b/src/binding/python/openPMD.cpp @@ -41,6 +41,7 @@ void init_Dataset(py::module &); void init_Datatype(py::module &); void init_Error(py::module &); void init_Helper(py::module &); +void init_CustomHierarchy(py::module &); void init_Iteration(py::module &); void init_IterationEncoding(py::module &); void init_Mesh(py::module &); @@ -94,6 +95,7 @@ PYBIND11_MODULE(openpmd_api_cxx, m) init_Dataset(m); init_Datatype(m); init_Helper(m); + init_CustomHierarchy(m); init_Iteration(m); init_IterationEncoding(m); init_BaseRecordComponent(m);