From 9ea0182b677b91cfb8d6d68374150cdec6b45695 Mon Sep 17 00:00:00 2001 From: Simon Date: Mon, 30 Sep 2019 23:20:06 +0200 Subject: [PATCH 01/15] Implementing a python module --- .gitignore | 6 +- .gitmodules | 3 + CMakeLists.txt | 4 +- MANIFEST.in | 4 ++ requirements.txt | 1 + setup.py | 91 +++++++++++++++++++++++++ src/CMakeLists.txt | 106 +----------------------------- src/python_wrapper/CMakeLists.txt | 8 +++ src/python_wrapper/_stk.cpp | 50 ++++++++++++++ src/stk/CMakeLists.txt | 105 +++++++++++++++++++++++++++++ stk/__init__.py | 1 + third_party/CMakeLists.txt | 4 ++ third_party/pybind11 | 1 + 13 files changed, 279 insertions(+), 105 deletions(-) create mode 100644 .gitmodules create mode 100644 MANIFEST.in create mode 100644 requirements.txt create mode 100644 setup.py create mode 100644 src/python_wrapper/CMakeLists.txt create mode 100644 src/python_wrapper/_stk.cpp create mode 100644 src/stk/CMakeLists.txt create mode 100644 stk/__init__.py create mode 160000 third_party/pybind11 diff --git a/.gitignore b/.gitignore index 8234407..5276635 100644 --- a/.gitignore +++ b/.gitignore @@ -1,5 +1,9 @@ build_debug/ build_vs/ build/ +dist/ +*.egg-info +__pycache__ +.env/ .vscode/ -.ycm_extra_conf.py +.ycm_extra_conf.py \ No newline at end of file diff --git a/.gitmodules b/.gitmodules new file mode 100644 index 0000000..9bb357f --- /dev/null +++ b/.gitmodules @@ -0,0 +1,3 @@ +[submodule "third_party/pybind11"] + path = third_party/pybind11 + url = https://github.com/pybind/pybind11.git diff --git a/CMakeLists.txt b/CMakeLists.txt index 2626c45..561fe5e 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -2,6 +2,7 @@ cmake_minimum_required(VERSION 3.8) option(STK_BUILD_EXAMPLES "Build examples" ON) option(STK_BUILD_TESTS "Build unit tests" ON) +option(STK_BUILD_PYTHON_WRAPPER "Build Python wrapper" OFF) option(STK_USE_CUDA "Enables CUDA support" OFF) option(STK_WARNINGS_ARE_ERRORS "Warnings are treated as errors" OFF) option(STK_BUILD_WITH_DEBUG_INFO "Includes debug info in release builds" OFF) @@ -109,10 +110,11 @@ set(CMAKE_C_FLAGS_RELEASE "${CMAKE_C_FLAGS_RELEASE} ${EXTRA_FLAGS_RELEASE}") set(CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} ${EXTRA_LINK_FLAGS}") set(CMAKE_EXE_LINKER_FLAGS_RELEASE "${CMAKE_EXE_LINKER_FLAGS_RELEASE} ${EXTRA_LINK_FLAGS_RELEASE}") +include_directories(src) + # TODO: Determine QNANHIBIT for NrrdIO set(QNANHIBIT 1) - if(NOT WIN32) find_package(ZLIB QUIET) find_package(NIFTI QUIET) diff --git a/MANIFEST.in b/MANIFEST.in new file mode 100644 index 0000000..6e9c09d --- /dev/null +++ b/MANIFEST.in @@ -0,0 +1,4 @@ +recursive-include src * +recursive-include third_party * +include cmake/* +include CMakeLists.txt diff --git a/requirements.txt b/requirements.txt new file mode 100644 index 0000000..296d654 --- /dev/null +++ b/requirements.txt @@ -0,0 +1 @@ +numpy \ No newline at end of file diff --git a/setup.py b/setup.py new file mode 100644 index 0000000..3fa7ff0 --- /dev/null +++ b/setup.py @@ -0,0 +1,91 @@ +import os +import sys +import platform +import subprocess + +from pprint import pprint + +from setuptools import setup, Extension +from setuptools.command.build_ext import build_ext + +with open("README.md", encoding="utf-8") as f: + readme = f.read() + +# Parse command line flags +flags = {k: 'OFF' for k in ['--debug', '--use-cuda', '--use-itk']} +for flag in flags.keys(): + if flag in sys.argv: + flags[flag] = 'ON' + sys.argv.remove(flag) + +# Command line flags forwarded to CMake +cmake_cmd_args = [] +for f in sys.argv: + if f.startswith('-D'): + cmake_cmd_args.append(f) + +for f in cmake_cmd_args: + sys.argv.remove(f) + + +class CMakeExtension(Extension): + def __init__(self, name, cmake_lists_dir=''): + Extension.__init__(self, name, sources=[]) + self.cmake_lists_dir = os.path.abspath(cmake_lists_dir) + + +class CMakeBuild(build_ext): + def run(self): + try: + subprocess.check_output(['cmake', '--version']) + except OSError: + raise RuntimeError('Cannot find CMake executable') + + for ext in self.extensions: + self.build_extension(ext) + + def build_extension(self, ext): + + extdir = os.path.abspath(os.path.dirname(self.get_ext_fullpath(ext.name))) + cfg = 'Debug' if flags['--debug'] == 'ON' else 'Release' + build_args = ['--config', cfg] + + cmake_args = [ + '-DSTK_BUILD_PYTHON_WRAPPER=ON', + '-DSTK_BUILD_TESTS=OFF', + '-DSTK_BUILD_WITH_DEBUG_INFO=%s' % ('ON' if cfg == 'Debug' else 'OFF'), + '-DCMAKE_BUILD_TYPE=%s' % cfg, + '-DSTK_USE_CUDA=%s' % flags['--use-cuda'], + '-DSTK_ITK_BRIDGE=%s' % flags['--use-itk'], + '-DPYTHON_EXECUTABLE=' + sys.executable, + '-DCMAKE_LIBRARY_OUTPUT_DIRECTORY=' + extdir, + ] + + if platform.system() == "Windows": + cmake_args += ['-DCMAKE_LIBRARY_OUTPUT_DIRECTORY_{}={}'.format(cfg.upper(), extdir)] + + cmake_args += cmake_cmd_args + pprint(cmake_args) + + if not os.path.exists(self.build_temp): + os.makedirs(self.build_temp) + + subprocess.check_call(['cmake', ext.cmake_lists_dir] + cmake_args, cwd=self.build_temp) + subprocess.check_call(['cmake', '--build', '.'] + build_args, cwd=self.build_temp) + + +setup( + name='python-stk', + version='0.1', + author='Simon Ekström', + author_email='', + description='', + long_description=readme, + long_description_content_type='text/markdown', + install_requires=['numpy'], + packages=['stk'], + ext_modules=[CMakeExtension('_stk', '.')], + cmdclass={'build_ext': CMakeBuild}, + zip_safe=False, +) + diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 65fdd29..d7e040f 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -1,105 +1,5 @@ -set(STK_SRCS - "stk/common/error.cpp" - "stk/common/error.h" - "stk/common/log.cpp" - "stk/common/log.h" - "stk/common/stack_trace.cpp" - "stk/common/stack_trace.h" - "stk/cuda/volume.h" - "stk/filters/decomposable_filter.h" - "stk/filters/decomposable_filter.inl" - "stk/filters/gaussian_filter.cpp" - "stk/filters/gaussian_filter.h" - "stk/filters/map.h" - "stk/filters/sobel.cpp" - "stk/filters/sobel.h" - "stk/filters/vector_calculus.h" - "stk/filters/vector_calculus.inl" - "stk/io/io.cpp" - "stk/io/io.h" - "stk/io/nifti.cpp" - "stk/io/nifti.h" - "stk/io/nrrd.cpp" - "stk/io/nrrd.h" - "stk/io/vtk.cpp" - "stk/io/vtk.h" - "stk/math/float2.h" - "stk/math/float3.h" - "stk/math/float4.h" - "stk/math/math.cpp" - "stk/math/math.h" - "stk/math/types.h" - "stk/image/dim3.h" - "stk/image/types.cpp" - "stk/image/types.h" - "stk/image/volume.cpp" - "stk/image/volume.h" - "stk/image/volume.inl" -) +add_subdirectory(stk) -if (STK_USE_CUDA) - set(STK_SRCS - ${STK_SRCS} - "stk/cuda/cuda.cpp" - "stk/cuda/cuda.h" - "stk/cuda/stream.cpp" - "stk/cuda/stream.h" - "stk/filters/gpu/normalize.cu" - "stk/image/gpu_volume.cpp" - "stk/image/gpu_volume.cu" - "stk/image/gpu_volume.h" - ) +if(STK_BUILD_PYTHON_WRAPPER) + add_subdirectory(python_wrapper) endif() - -if(STK_ITK_BRIDGE) - set(STK_SRCS - ${STK_SRCS} - "stk/itk_bridge/itk_convert.h" - "stk/itk_bridge/itk_io.cpp" - "stk/itk_bridge/itk_io.h" - ) -endif() - -add_library(stk STATIC ${STK_SRCS}) - -target_include_directories(stk PUBLIC - ${CMAKE_CURRENT_SOURCE_DIR} - ${CMAKE_CUDA_TOOLKIT_INCLUDE_DIRECTORIES} - ) - -if(NIFTI_FOUND) - target_include_directories(stk PUBLIC ${NIFTI_INCLUDE_DIR}) - target_link_libraries(stk PUBLIC ${NIFTI_LIBRARIES}) -else() - target_link_libraries(stk PUBLIC niftiio) -endif() - -if(NrrdIO_FOUND) - target_include_directories(stk PUBLIC ${NrrdIO_INCLUDE_DIR}) - target_link_libraries(stk PUBLIC ${NrrdIO_LIBRARIES}) -else() - target_link_libraries(stk PUBLIC NrrdIO) -endif() - -if(STK_ITK_BRIDGE) - add_definitions(-DSTK_ITK_BRIDGE) - find_package(ITK REQUIRED) - include(${ITK_USE_FILE}) - - target_include_directories(stk PRIVATE ${ITK_INCLUDE_DIRS}) - target_link_libraries(stk PRIVATE ${ITK_LIBRARIES}) -endif() - -if(STK_STACK_TRACE) - add_definitions(-DSTK_STACK_TRACE) - - if(UNIX) - find_package(Backward) - if(Backward_FOUND) - add_backward(stk) - target_link_libraries(stk PUBLIC ${BACKWARD_LIBRARIES}) - add_definitions(-DSTK_USE_BACKWARD) - endif() - endif() -endif() - diff --git a/src/python_wrapper/CMakeLists.txt b/src/python_wrapper/CMakeLists.txt new file mode 100644 index 0000000..39a33a7 --- /dev/null +++ b/src/python_wrapper/CMakeLists.txt @@ -0,0 +1,8 @@ +pybind11_add_module(_stk _stk.cpp) + +target_include_directories(_stk PUBLIC ${CMAKE_CURRENT_SOURCE_DIR} + ${PYTHON_INCLUDE_DIRS} ${PYBIND11_INCLUDE_DIR} ${pybind11_INCLUDE_DIR}) +target_link_libraries(_stk PUBLIC stk) + +set(STK_PYTHON_PATH "$ENV{PYTHONPATH}${SEP}${CMAKE_CURRENT_BINARY_DIR}" CACHE STRING "Python path for Sphinx") + diff --git a/src/python_wrapper/_stk.cpp b/src/python_wrapper/_stk.cpp new file mode 100644 index 0000000..9bd628c --- /dev/null +++ b/src/python_wrapper/_stk.cpp @@ -0,0 +1,50 @@ +#include +#include +#include +#include + +#include + +#include +#include +#include + +namespace py = pybind11; + +class VolumeW + +py::array test() +{ + printf("test\n"); + return py::array_t({32,32,32}); +} + +// Volume read_volume(const std::string& filename); + +// // Attempts to write the given volume to the specified file. +// // Uses the extension of the filename to identify the target file format. +// // Triggers a fatal error if write failed (e.g. invalid file extension). +// void write_volume(const std::string&, const Volume& vol); + +PYBIND11_MODULE(_stk, m) +{ + m.def("test", &test, ""); + // m.def("read_volume", &read_volume, ""); + // m.def("write_volume", &write_volume, ""); + + // Translate relevant exception types. The exceptions not handled + // here will be translated autmatically according to pybind11's rules. + py::register_exception_translator([](std::exception_ptr p) { + try { + if (p) { + std::rethrow_exception(p); + } + } + catch (const stk::FatalException &e) { + // Map stk::FatalException to Python RutimeError + // +13 to remove the "Fatal error: " from the message + PyErr_SetString(PyExc_RuntimeError, e.what() + 13); + } + }); +} + diff --git a/src/stk/CMakeLists.txt b/src/stk/CMakeLists.txt new file mode 100644 index 0000000..aa3b4cd --- /dev/null +++ b/src/stk/CMakeLists.txt @@ -0,0 +1,105 @@ +set(STK_SRCS + "common/error.cpp" + "common/error.h" + "common/log.cpp" + "common/log.h" + "common/stack_trace.cpp" + "common/stack_trace.h" + "cuda/volume.h" + "filters/decomposable_filter.h" + "filters/decomposable_filter.inl" + "filters/gaussian_filter.cpp" + "filters/gaussian_filter.h" + "filters/map.h" + "filters/sobel.cpp" + "filters/sobel.h" + "filters/vector_calculus.h" + "filters/vector_calculus.inl" + "io/io.cpp" + "io/io.h" + "io/nifti.cpp" + "io/nifti.h" + "io/nrrd.cpp" + "io/nrrd.h" + "io/vtk.cpp" + "io/vtk.h" + "math/float2.h" + "math/float3.h" + "math/float4.h" + "math/math.cpp" + "math/math.h" + "math/types.h" + "image/dim3.h" + "image/types.cpp" + "image/types.h" + "image/volume.cpp" + "image/volume.h" + "image/volume.inl" +) + +if (STK_USE_CUDA) + set(STK_SRCS + ${STK_SRCS} + "cuda/cuda.cpp" + "cuda/cuda.h" + "cuda/stream.cpp" + "cuda/stream.h" + "filters/gpu/normalize.cu" + "image/gpu_volume.cpp" + "image/gpu_volume.cu" + "image/gpu_volume.h" + ) +endif() + +if(STK_ITK_BRIDGE) + set(STK_SRCS + ${STK_SRCS} + "itk_bridge/itk_convert.h" + "itk_bridge/itk_io.cpp" + "itk_bridge/itk_io.h" + ) +endif() + +add_library(stk STATIC ${STK_SRCS}) + +target_include_directories(stk PUBLIC + ${CMAKE_CURRENT_SOURCE_DIR} + ${CMAKE_CUDA_TOOLKIT_INCLUDE_DIRECTORIES} + ) + +if(NIFTI_FOUND) + target_include_directories(stk PUBLIC ${NIFTI_INCLUDE_DIR}) + target_link_libraries(stk PUBLIC ${NIFTI_LIBRARIES}) +else() + target_link_libraries(stk PUBLIC niftiio) +endif() + +if(NrrdIO_FOUND) + target_include_directories(stk PUBLIC ${NrrdIO_INCLUDE_DIR}) + target_link_libraries(stk PUBLIC ${NrrdIO_LIBRARIES}) +else() + target_link_libraries(stk PUBLIC NrrdIO) +endif() + +if(STK_ITK_BRIDGE) + add_definitions(-DSTK_ITK_BRIDGE) + find_package(ITK REQUIRED) + include(${ITK_USE_FILE}) + + target_include_directories(stk PRIVATE ${ITK_INCLUDE_DIRS}) + target_link_libraries(stk PRIVATE ${ITK_LIBRARIES}) +endif() + +if(STK_STACK_TRACE) + add_definitions(-DSTK_STACK_TRACE) + + if(UNIX) + find_package(Backward) + if(Backward_FOUND) + add_backward(stk) + target_link_libraries(stk PUBLIC ${BACKWARD_LIBRARIES}) + add_definitions(-DSTK_USE_BACKWARD) + endif() + endif() +endif() + diff --git a/stk/__init__.py b/stk/__init__.py new file mode 100644 index 0000000..a5a299c --- /dev/null +++ b/stk/__init__.py @@ -0,0 +1 @@ +from _stk import * \ No newline at end of file diff --git a/third_party/CMakeLists.txt b/third_party/CMakeLists.txt index 18daeee..62d76b8 100644 --- a/third_party/CMakeLists.txt +++ b/third_party/CMakeLists.txt @@ -61,3 +61,7 @@ if(NOT NrrdIO_FOUND) endif() endif() +if(STK_BUILD_PYTHON_WRAPPER) + add_subdirectory(pybind11) +endif() + diff --git a/third_party/pybind11 b/third_party/pybind11 new file mode 160000 index 0000000..7ec2ddf --- /dev/null +++ b/third_party/pybind11 @@ -0,0 +1 @@ +Subproject commit 7ec2ddfc95f65d1e986d359466a6c254aa514ef3 From 130c079deccaa1747fb38bfaed3b40866c82397b Mon Sep 17 00:00:00 2001 From: Simon Date: Tue, 1 Oct 2019 12:04:04 +0200 Subject: [PATCH 02/15] WIP: Volume wrapper --- src/python_wrapper/_stk.cpp | 242 +++++++++++++++++++++++++++++++++++- 1 file changed, 239 insertions(+), 3 deletions(-) diff --git a/src/python_wrapper/_stk.cpp b/src/python_wrapper/_stk.cpp index 9bd628c..fc4d1ec 100644 --- a/src/python_wrapper/_stk.cpp +++ b/src/python_wrapper/_stk.cpp @@ -4,6 +4,7 @@ #include #include +#include #include #include @@ -11,12 +12,232 @@ namespace py = pybind11; -class VolumeW -py::array test() +/*! + * \brief Get the stk::Type associated to a numpy image. + */ +stk::Type get_stk_type(const py::array& a) { + stk::Type base_type = stk::Type_Unknown; + + if (py::isinstance>(a)) { + base_type = stk::Type_Char; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Char; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_UChar; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Short; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_UShort; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Int; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_UInt; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Float; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Double; + } + else { + throw std::invalid_argument("Unsupported type"); + } + + // NOTE: the value of ndim can be ambiguous, e.g. + // ndim == 3 may be a scalar volume or a vector 2D image... + return stk::build_type(base_type, a.ndim() == 4 ? 3 : 1); +} + +void get_numpy_type_and_shape(stk::Type type, const dim3& size) { + stk::Type base_type = stk::base_type(type); + + stk::Type base_type = stk::Type_Unknown; + + if (py::isinstance>(a)) { + base_type = stk::Type_Char; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Char; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_UChar; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Short; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_UShort; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Int; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_UInt; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Float; + } + else if (py::isinstance>(a)) { + base_type = stk::Type_Double; + } + else { + throw std::invalid_argument("Unsupported type"); + } + + // NOTE: the value of ndim can be ambiguous, e.g. + // ndim == 3 may be a scalar volume or a vector 2D image... + return stk::build_type(base_type, a.ndim() == 4 ? 3 : 1); +} + + +/*! + * \brief Convert a numpy array to a stk::Volume. + * + * The new volume creates a copy of the input + * image data. + * + * @note The numpy array must be C-contiguous, with + * [z,y,x] indexing. + * + * @param image Array representing a volume image. + * @param origin Vector of length 3, containing + * the (x, y,z) coordinates of the + * volume origin. + * @param spacing Vector of length 3, containing + * the (x, y,z) spacing of the + * volume. + * @param direction Vector of length 9, representing + * the cosine direction matrix in + * row-major order. + * @return A volume representing the same image. + */ +stk::Volume image_to_volume( + const py::array image, + const std::vector& origin, + const std::vector& spacing, + const std::vector& direction + ) +{ + if (image.flags() & py::array::f_style) { + throw std::invalid_argument("The arrays must be C-contiguous."); + } + + float3 origin_ { + float(origin[0]), + float(origin[1]), + float(origin[2]), + }; + float3 spacing_ { + float(spacing[0]), + float(spacing[1]), + float(spacing[2]), + }; + Matrix3x3f direction_ {{ + {float(direction[0]), float(direction[1]), float(direction[2])}, + {float(direction[3]), float(direction[4]), float(direction[5])}, + {float(direction[6]), float(direction[7]), float(direction[8])}, + }}; + dim3 size { + std::uint32_t(image.shape(2)), + std::uint32_t(image.shape(1)), + std::uint32_t(image.shape(0)), + }; + stk::Volume volume {size, get_stk_type(image), image.data()}; + volume.set_origin(origin_); + volume.set_spacing(spacing_); + volume.set_direction(direction_); + return volume; +} + + +class PyVolume +{ +public: + PyVolume() {} + ~PyVolume() {} + + py::tuple origin() { + return py::make_tuple( + _volume.origin().x, + _volume.origin().y, + _volume.origin().z + ); + } + + bool valid() { + return _volume.valid(); + } + + void set_origin(py::tuple origin) { + if (origin.size() != 3) { + throw py::value_error("Expected an 3-tuple"); + } + + _volume.set_origin(float3{ + origin[0].cast(), + origin[1].cast(), + origin[2].cast() + }); + } + + py::tuple spacing() { + return py::make_tuple( + _volume.spacing().x, + _volume.spacing().y, + _volume.spacing().z + ); + } + void set_spacing(py::tuple spacing) { + if (spacing.size() != 3) { + throw py::value_error("Expected an 3-tuple"); + } + + _volume.set_spacing(float3{ + spacing[0].cast(), + spacing[1].cast(), + spacing[2].cast() + }); + } + + py::array_t direction() { + return py::array_t( + {3, 3}, + &_volume.direction()._rows[0].x + ); + } + void set_direction(py::array_t dir) { + if (dir.ndim() != 2 + || dir.shape(0) != 3 + || dir.shape(1) != 3) { + throw py::value_error("Expected an 3x3 matrix"); + } + + _volume.set_direction(Matrix3x3f{ + float3{dir.at(0, 0), dir.at(0, 1), dir.at(0, 2)}, + float3{dir.at(1, 0), dir.at(1, 1), dir.at(1, 2)}, + float3{dir.at(2, 0), dir.at(2, 1), dir.at(2, 2)} + }); + } + + py::array data() { + + } + +private: + stk::Volume _volume; +}; + +PyVolume test() { printf("test\n"); - return py::array_t({32,32,32}); + return PyVolume(); } // Volume read_volume(const std::string& filename); @@ -29,6 +250,21 @@ py::array test() PYBIND11_MODULE(_stk, m) { m.def("test", &test, ""); + + py::class_(m, "Volume") + .def(py::init<>()) + .def("valid", &PyVolume::valid) + .def_property("origin", &PyVolume::origin, &PyVolume::set_origin) + .def_property("spacing", &PyVolume::spacing, &PyVolume::set_spacing) + .def_property("direction", &PyVolume::direction, &PyVolume::set_direction) + .def_property_readonly("data", &PyVolume::data); + + // .def("go_for_a_walk", &Pet::go_for_a_walk) + // .def("get_hunger", &Pet::get_hunger) + // .def("get_name", &Pet::get_name) + // .def_property_readonly("hunger", &Pet::get_hunger) + // .def_; + // m.def("read_volume", &read_volume, ""); // m.def("write_volume", &write_volume, ""); From cf5149222913d8757bccacfe74dae1787d429339 Mon Sep 17 00:00:00 2001 From: Simon Date: Tue, 1 Oct 2019 16:23:58 +0200 Subject: [PATCH 03/15] Remove obsolete(?) set_direction --- src/stk/image/gpu_volume.cpp | 5 ----- src/stk/image/gpu_volume.h | 1 - src/stk/image/volume.cpp | 5 ----- src/stk/image/volume.h | 1 - 4 files changed, 12 deletions(-) diff --git a/src/stk/image/gpu_volume.cpp b/src/stk/image/gpu_volume.cpp index 5bc1074..54a0ec0 100644 --- a/src/stk/image/gpu_volume.cpp +++ b/src/stk/image/gpu_volume.cpp @@ -449,11 +449,6 @@ void GpuVolume::set_direction(const Matrix3x3f& direction) _direction = direction; _inverse_direction = _direction.inverse(); } -void GpuVolume::set_direction(const std::initializer_list direction) -{ - _direction.set(direction); - _inverse_direction = _direction.inverse(); -} const float3& GpuVolume::origin() const { diff --git a/src/stk/image/gpu_volume.h b/src/stk/image/gpu_volume.h index 857ddab..45d63bf 100644 --- a/src/stk/image/gpu_volume.h +++ b/src/stk/image/gpu_volume.h @@ -110,7 +110,6 @@ namespace stk void set_origin(const float3& origin); void set_spacing(const float3& spacing); void set_direction(const Matrix3x3f& direction); - void set_direction(const std::initializer_list direction); const float3& origin() const; const float3& spacing() const; diff --git a/src/stk/image/volume.cpp b/src/stk/image/volume.cpp index 28085cf..a22a534 100644 --- a/src/stk/image/volume.cpp +++ b/src/stk/image/volume.cpp @@ -255,11 +255,6 @@ void Volume::set_direction(const Matrix3x3f& direction) _direction = direction; _inverse_direction = _direction.inverse(); } -void Volume::set_direction(const std::initializer_list direction) -{ - _direction.set(direction); - _inverse_direction = _direction.inverse(); -} const float3& Volume::origin() const { return _origin; diff --git a/src/stk/image/volume.h b/src/stk/image/volume.h index 9950136..a09f2ce 100644 --- a/src/stk/image/volume.h +++ b/src/stk/image/volume.h @@ -150,7 +150,6 @@ namespace stk void set_origin(const float3& origin); void set_spacing(const float3& spacing); void set_direction(const Matrix3x3f& direction); - void set_direction(const std::initializer_list direction); const float3& origin() const; const float3& spacing() const; From 0a6a8507a0a72b0e165d699d8ac933149f671f76 Mon Sep 17 00:00:00 2001 From: Simon Date: Tue, 1 Oct 2019 16:24:23 +0200 Subject: [PATCH 04/15] Functional python module for reading/writing volumes --- src/python_wrapper/_stk.cpp | 460 ++++++++++++++++++------------------ 1 file changed, 226 insertions(+), 234 deletions(-) diff --git a/src/python_wrapper/_stk.cpp b/src/python_wrapper/_stk.cpp index fc4d1ec..9be7ff8 100644 --- a/src/python_wrapper/_stk.cpp +++ b/src/python_wrapper/_stk.cpp @@ -1,3 +1,4 @@ +#include #include #include #include @@ -5,268 +6,259 @@ #include #include +#include #include #include #include +#include namespace py = pybind11; -/*! - * \brief Get the stk::Type associated to a numpy image. - */ -stk::Type get_stk_type(const py::array& a) { - stk::Type base_type = stk::Type_Unknown; +namespace pybind11 { namespace detail { + template <> struct type_caster { + public: + PYBIND11_TYPE_CASTER(float3, _("float3")); + + bool load(handle src, bool) { + try { + tuple t = py::cast(src); + if (t.size() != 3) { + return false; + } + value = float3{ + t[0].cast(), + t[1].cast(), + t[2].cast() + }; + } catch (std::exception& e) { + return false; + } + return true; + } + static handle cast(float3 src, return_value_policy /* policy */, handle /* parent */) { + return make_tuple( + src.x, + src.y, + src.z + ).release(); + } + }; + template <> struct type_caster { + public: + PYBIND11_TYPE_CASTER(dim3, _("dim3")); + + bool load(handle src, bool) { + try { + tuple t = py::cast(src); + if (t.size() != 3) { + return false; + } + value = dim3{ + t[0].cast(), + t[1].cast(), + t[2].cast() + }; + } catch (std::exception& e) { + return false; + } + return true; + } + static handle cast(dim3 src, return_value_policy /* policy */, handle /* parent */) { + return make_tuple( + src.x, + src.y, + src.z + ).release(); + } + }; + template <> struct type_caster { + public: + PYBIND11_TYPE_CASTER(int3, _("int3")); + + bool load(handle src, bool) { + try { + tuple t = py::cast(src); + if (t.size() != 3) { + return false; + } + value = int3{ + t[0].cast(), + t[1].cast(), + t[2].cast() + }; + } catch (std::exception& e) { + return false; + } + return true; + } + static handle cast(int3 src, return_value_policy /* policy */, handle /* parent */) { + return make_tuple( + src.x, + src.y, + src.z + ).release(); + } + }; + template <> struct type_caster { + public: + PYBIND11_TYPE_CASTER(Matrix3x3f, _("Matrix3x3f")); + + bool load(handle src, bool) { + try { + array_t a = py::cast>(src); + if (a.ndim() != 2 + || a.shape(0) != 3 + || a.shape(1) != 3) { + return false; + } + value = Matrix3x3f{ + float3{a.at(0, 0), a.at(0, 1), a.at(0, 2)}, + float3{a.at(1, 0), a.at(1, 1), a.at(1, 2)}, + float3{a.at(2, 0), a.at(2, 1), a.at(2, 2)} + }; + } catch (std::exception& e) { + return false; + } + return true; + } + static handle cast(Matrix3x3f src, return_value_policy /* policy */, handle /* parent */) { + return py::array_t( + {3, 3}, + &src._rows[0].x + ).release(); + } + }; - if (py::isinstance>(a)) { - base_type = stk::Type_Char; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Char; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_UChar; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Short; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_UShort; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Int; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_UInt; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Float; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Double; +}} // namespace pybind11::detail + +std::vector get_shape(const stk::Volume& vol) +{ + dim3 size = vol.size(); + size_t ncomp = stk::num_components(vol.voxel_type()); + if (ncomp == 1) { + // Flip directions since size is in (width, height, depth), + // numpy expects (depth, height, width) + return { + size.z, + size.y, + size.x + }; } else { - throw std::invalid_argument("Unsupported type"); + return { + size.z, + size.y, + size.x, + ncomp + }; } - - // NOTE: the value of ndim can be ambiguous, e.g. - // ndim == 3 may be a scalar volume or a vector 2D image... - return stk::build_type(base_type, a.ndim() == 4 ? 3 : 1); } -void get_numpy_type_and_shape(stk::Type type, const dim3& size) { - stk::Type base_type = stk::base_type(type); - - stk::Type base_type = stk::Type_Unknown; - - if (py::isinstance>(a)) { - base_type = stk::Type_Char; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Char; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_UChar; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Short; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_UShort; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Int; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_UInt; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Float; - } - else if (py::isinstance>(a)) { - base_type = stk::Type_Double; +std::vector get_strides(const stk::Volume& vol) +{ + const size_t* strides = vol.strides(); + size_t ncomp = stk::num_components(vol.voxel_type()); + + // Reverse, because of numpy ordering + if (ncomp == 1) { + return { + strides[2], + strides[1], + strides[0] + }; } else { - throw std::invalid_argument("Unsupported type"); + return { + strides[2], + strides[1], + strides[0], + strides[0]/ncomp + }; } - - // NOTE: the value of ndim can be ambiguous, e.g. - // ndim == 3 may be a scalar volume or a vector 2D image... - return stk::build_type(base_type, a.ndim() == 4 ? 3 : 1); } - -/*! - * \brief Convert a numpy array to a stk::Volume. - * - * The new volume creates a copy of the input - * image data. - * - * @note The numpy array must be C-contiguous, with - * [z,y,x] indexing. - * - * @param image Array representing a volume image. - * @param origin Vector of length 3, containing - * the (x, y,z) coordinates of the - * volume origin. - * @param spacing Vector of length 3, containing - * the (x, y,z) spacing of the - * volume. - * @param direction Vector of length 9, representing - * the cosine direction matrix in - * row-major order. - * @return A volume representing the same image. - */ -stk::Volume image_to_volume( - const py::array image, - const std::vector& origin, - const std::vector& spacing, - const std::vector& direction - ) +const char* format_descriptor(stk::Type base_type) { - if (image.flags() & py::array::f_style) { - throw std::invalid_argument("The arrays must be C-contiguous."); - } - - float3 origin_ { - float(origin[0]), - float(origin[1]), - float(origin[2]), - }; - float3 spacing_ { - float(spacing[0]), - float(spacing[1]), - float(spacing[2]), - }; - Matrix3x3f direction_ {{ - {float(direction[0]), float(direction[1]), float(direction[2])}, - {float(direction[3]), float(direction[4]), float(direction[5])}, - {float(direction[6]), float(direction[7]), float(direction[8])}, - }}; - dim3 size { - std::uint32_t(image.shape(2)), - std::uint32_t(image.shape(1)), - std::uint32_t(image.shape(0)), - }; - stk::Volume volume {size, get_stk_type(image), image.data()}; - volume.set_origin(origin_); - volume.set_spacing(spacing_); - volume.set_direction(direction_); - return volume; + if (base_type == stk::Type_Char) return "b"; + else if (base_type == stk::Type_UChar) return "B"; + else if (base_type == stk::Type_Short) return "h"; + else if (base_type == stk::Type_UShort) return "H"; + else if (base_type == stk::Type_Int) return "i"; + else if (base_type == stk::Type_UInt) return "I"; + else if (base_type == stk::Type_Float) return "f"; + else if (base_type == stk::Type_Double) return "d"; + return ""; } - -class PyVolume +py::buffer_info get_buffer_info(stk::Volume& v) { -public: - PyVolume() {} - ~PyVolume() {} - - py::tuple origin() { - return py::make_tuple( - _volume.origin().x, - _volume.origin().y, - _volume.origin().z - ); - } - - bool valid() { - return _volume.valid(); - } - - void set_origin(py::tuple origin) { - if (origin.size() != 3) { - throw py::value_error("Expected an 3-tuple"); - } - - _volume.set_origin(float3{ - origin[0].cast(), - origin[1].cast(), - origin[2].cast() - }); - } - - py::tuple spacing() { - return py::make_tuple( - _volume.spacing().x, - _volume.spacing().y, - _volume.spacing().z - ); - } - void set_spacing(py::tuple spacing) { - if (spacing.size() != 3) { - throw py::value_error("Expected an 3-tuple"); - } - - _volume.set_spacing(float3{ - spacing[0].cast(), - spacing[1].cast(), - spacing[2].cast() - }); - } - - py::array_t direction() { - return py::array_t( - {3, 3}, - &_volume.direction()._rows[0].x - ); - } - void set_direction(py::array_t dir) { - if (dir.ndim() != 2 - || dir.shape(0) != 3 - || dir.shape(1) != 3) { - throw py::value_error("Expected an 3x3 matrix"); - } - - _volume.set_direction(Matrix3x3f{ - float3{dir.at(0, 0), dir.at(0, 1), dir.at(0, 2)}, - float3{dir.at(1, 0), dir.at(1, 1), dir.at(1, 2)}, - float3{dir.at(2, 0), dir.at(2, 1), dir.at(2, 2)} - }); - } - - py::array data() { - - } - -private: - stk::Volume _volume; -}; - -PyVolume test() -{ - printf("test\n"); - return PyVolume(); + stk::Type base_type = stk::base_type(v.voxel_type()); + int ncomp = stk::num_components(v.voxel_type()); + + return py::buffer_info( + v.ptr(), /* Pointer to buffer */ + stk::type_size(base_type), /* Size of one scalar */ + format_descriptor(base_type), /* Python struct-style format descriptor */ + ncomp == 1 ? 3 : 4, /* Number of dimensions */ + + get_shape(v), /* Buffer dimensions */ + get_strides(v) /* Strides (in bytes) for each index */ + ); } -// Volume read_volume(const std::string& filename); - -// // Attempts to write the given volume to the specified file. -// // Uses the extension of the filename to identify the target file format. -// // Triggers a fatal error if write failed (e.g. invalid file extension). -// void write_volume(const std::string&, const Volume& vol); - PYBIND11_MODULE(_stk, m) { - m.def("test", &test, ""); - - py::class_(m, "Volume") + py::enum_(m, "Type") + .value("Unknown", stk::Type_Unknown) + .value("Char", stk::Type_Char) + .value("Char2", stk::Type_Char2) + .value("Char3", stk::Type_Char3) + .value("Char4", stk::Type_Char4) + .value("UChar", stk::Type_UChar) + .value("UChar2", stk::Type_UChar2) + .value("UChar3", stk::Type_UChar3) + .value("UChar4", stk::Type_UChar4) + .value("Short", stk::Type_Short) + .value("Short2", stk::Type_Short2) + .value("Short3", stk::Type_Short3) + .value("Short4", stk::Type_Short4) + .value("UShort", stk::Type_UShort) + .value("UShort2", stk::Type_UShort2) + .value("UShort3", stk::Type_UShort3) + .value("UShort4", stk::Type_UShort4) + .value("Int", stk::Type_Int) + .value("Int2", stk::Type_Int2) + .value("Int3", stk::Type_Int3) + .value("Int4", stk::Type_Int4) + .value("UInt", stk::Type_UInt) + .value("UInt2", stk::Type_UInt2) + .value("UInt3", stk::Type_UInt3) + .value("UInt4", stk::Type_UInt4) + .value("Float", stk::Type_Float) + .value("Float2", stk::Type_Float2) + .value("Float3", stk::Type_Float3) + .value("Float4", stk::Type_Float4) + .value("Double", stk::Type_Double) + .value("Double2", stk::Type_Double2) + .value("Double3", stk::Type_Double3) + .value("Double4", stk::Type_Double4) + .export_values(); + + py::class_(m, "Volume", py::buffer_protocol()) + .def_buffer(get_buffer_info) .def(py::init<>()) - .def("valid", &PyVolume::valid) - .def_property("origin", &PyVolume::origin, &PyVolume::set_origin) - .def_property("spacing", &PyVolume::spacing, &PyVolume::set_spacing) - .def_property("direction", &PyVolume::direction, &PyVolume::set_direction) - .def_property_readonly("data", &PyVolume::data); - - // .def("go_for_a_walk", &Pet::go_for_a_walk) - // .def("get_hunger", &Pet::get_hunger) - // .def("get_name", &Pet::get_name) - // .def_property_readonly("hunger", &Pet::get_hunger) - // .def_; - - // m.def("read_volume", &read_volume, ""); - // m.def("write_volume", &write_volume, ""); + .def("valid", &stk::Volume::valid) + .def_property_readonly("size", &stk::Volume::size) + .def_property_readonly("type", &stk::Volume::voxel_type) + .def_property("origin", &stk::Volume::origin, &stk::Volume::set_origin) + .def_property("spacing", &stk::Volume::spacing, &stk::Volume::set_spacing) + .def_property("direction", &stk::Volume::direction, &stk::Volume::set_direction) + ; + + m.def("read_volume", &stk::read_volume, ""); + m.def("write_volume", &stk::write_volume, ""); // Translate relevant exception types. The exceptions not handled // here will be translated autmatically according to pybind11's rules. From 0c8f11c195ea944134479a8b393cafef0f2817d4 Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 10:31:29 +0200 Subject: [PATCH 05/15] Fix --- src/stk/CMakeLists.txt | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/stk/CMakeLists.txt b/src/stk/CMakeLists.txt index aa3b4cd..35aa4e0 100644 --- a/src/stk/CMakeLists.txt +++ b/src/stk/CMakeLists.txt @@ -61,9 +61,10 @@ if(STK_ITK_BRIDGE) endif() add_library(stk STATIC ${STK_SRCS}) +target_link_libraries(stk PUBLIC pybind11) target_include_directories(stk PUBLIC - ${CMAKE_CURRENT_SOURCE_DIR} + "${CMAKE_CURRENT_SOURCE_DIR}/.." ${CMAKE_CUDA_TOOLKIT_INCLUDE_DIRECTORIES} ) From 36a412b7948b6eee8e21f4c19d45c8d3b96c3507 Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 11:48:08 +0200 Subject: [PATCH 06/15] Numpy to stk volume conversion --- src/python_wrapper/_stk.cpp | 86 +++++++++++++++++++++++++++++++++++++ src/stk/math/math.cpp | 5 +++ src/stk/math/types.h | 2 + 3 files changed, 93 insertions(+) diff --git a/src/python_wrapper/_stk.cpp b/src/python_wrapper/_stk.cpp index 9be7ff8..e2b0d97 100644 --- a/src/python_wrapper/_stk.cpp +++ b/src/python_wrapper/_stk.cpp @@ -208,6 +208,85 @@ py::buffer_info get_buffer_info(stk::Volume& v) ); } +stk::Volume numpy_to_volume(const py::array arr) +{ + // Ensure ordering + py::array buf = py::array::ensure(arr, py::array::c_style | py::array::forcecast); + + if (!(buf.ndim() == 3 || buf.ndim() == 4)) { + throw std::invalid_argument("Expected a 3D or a 4D array"); + } + + dim3 size { + std::uint32_t(buf.shape(2)), + std::uint32_t(buf.shape(1)), + std::uint32_t(buf.shape(0)), + }; + + stk::Type base_type = stk::Type_Unknown; + if (py::isinstance>(arr)) { + base_type = stk::Type_Char; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_Char; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_UChar; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_Short; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_UShort; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_Int; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_UInt; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_Float; + } + else if (py::isinstance>(arr)) { + base_type = stk::Type_Double; + } + else { + throw std::invalid_argument("Unsupported type"); + } + + int ncomp = 1; + if (buf.ndim() == 4) { + ncomp = buf.shape(3); + } + + if (ncomp < 1 || ncomp > 4) { + throw std::invalid_argument("Unsupported number of channels"); + } + + // NOTE: the value of ndim can be ambiguous, e.g. + // ndim == 3 may be a scalar volume or a vector 2D image... + + return stk::Volume( + size, + stk::build_type(base_type, ncomp), + buf.data() + ); +} + +stk::Volume make_volume( + py::array arr, + const float3& origin, + const float3& spacing, + const Matrix3x3f& direction +) { + stk::Volume vol = numpy_to_volume(arr); + vol.set_origin(origin); + vol.set_spacing(spacing); + vol.set_direction(direction); + return vol; +} + PYBIND11_MODULE(_stk, m) { py::enum_(m, "Type") @@ -249,7 +328,14 @@ PYBIND11_MODULE(_stk, m) py::class_(m, "Volume", py::buffer_protocol()) .def_buffer(get_buffer_info) .def(py::init<>()) + .def(py::init(&make_volume), + py::arg("array"), + py::arg("origin") = float3{0.0, 0.0, 0.0}, + py::arg("spacing") = float3{1.0, 1.0, 1.0}, + py::arg("direction") = Matrix3x3f::Identity + ) .def("valid", &stk::Volume::valid) + .def("copy_meta_from", &stk::Volume::copy_meta_from) .def_property_readonly("size", &stk::Volume::size) .def_property_readonly("type", &stk::Volume::voxel_type) .def_property("origin", &stk::Volume::origin, &stk::Volume::set_origin) diff --git a/src/stk/math/math.cpp b/src/stk/math/math.cpp index 66b41cd..e75cd22 100644 --- a/src/stk/math/math.cpp +++ b/src/stk/math/math.cpp @@ -1,5 +1,10 @@ #include "math.h" +Matrix3x3f Matrix3x3f::Identity { + float3{1, 0, 0}, + float3{0, 1, 0}, + float3{0, 0, 1} +}; namespace stk { diff --git a/src/stk/math/types.h b/src/stk/math/types.h index 0d91dbe..8a12fff 100644 --- a/src/stk/math/types.h +++ b/src/stk/math/types.h @@ -136,6 +136,8 @@ struct STK_ALIGN(16) double4 struct Matrix3x3f { + static Matrix3x3f Identity; + enum : unsigned int { rows = 3, cols = 3 From ec0cb50b6f6ff92d1e78295bfe84a866f56498b1 Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 11:48:13 +0200 Subject: [PATCH 07/15] Tests --- test/test_volume.py | 85 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 85 insertions(+) create mode 100644 test/test_volume.py diff --git a/test/test_volume.py b/test/test_volume.py new file mode 100644 index 0000000..3f8bade --- /dev/null +++ b/test/test_volume.py @@ -0,0 +1,85 @@ +import os +import random +import unittest +import shutil, tempfile + +import numpy as np + +import stk + +class Test_Volume(unittest.TestCase): + def setUp(self): + self.temp_dir = tempfile.mkdtemp() + + def tearDown(self): + shutil.rmtree(self.temp_dir) + + def test_constructor(self): + arr = np.zeros((4, 5, 6)).astype(np.int32) + vol = stk.Volume(arr) + + np.testing.assert_equal(vol.size, (6, 5, 4)) + np.testing.assert_equal(vol.origin, (0, 0, 0)) + np.testing.assert_equal(vol.spacing, (1, 1, 1)) + np.testing.assert_equal(vol.direction, ((1, 0, 0), (0, 1, 0), (0, 0, 1))) + np.testing.assert_equal(np.array(vol), arr) + + d = ((10, 0, 0), (0, 11, 0), (0, 0, 12)) + vol = stk.Volume(arr, origin=(1,2,3), spacing=(4,5,6), direction=d) + + np.testing.assert_equal(vol.size, (6, 5, 4)) + np.testing.assert_equal(vol.origin, (1, 2, 3)) + np.testing.assert_equal(vol.spacing, (4, 5, 6)) + np.testing.assert_equal(vol.direction, d) + np.testing.assert_equal(np.array(vol), arr) + + def test_meta(self): + vol = stk.Volume() + vol.origin = (2, 3, 4) + np.testing.assert_equal(vol.origin, (2, 3, 4)) + + vol = stk.Volume() + vol.spacing = (2, 3, 4) + np.testing.assert_equal(vol.spacing, (2, 3, 4)) + + d = ((10, 0, 0), (0, 11, 0), (0, 0, 12)) + vol = stk.Volume() + vol.direction = d + np.testing.assert_equal(vol.direction, d) + + def test_copy_meta_from(self): + vol1 = stk.Volume() + vol2 = stk.Volume() + + vol1.origin = (2,3,4) + vol1.spacing = (5,6,7) + d = ((10, 0, 0), (0, 11, 0), (0, 0, 12)) + vol1.direction = d + + vol2.copy_meta_from(vol1) + + np.testing.assert_equal(vol2.origin, (2, 3, 4)) + np.testing.assert_equal(vol2.spacing, (5, 6, 7)) + np.testing.assert_equal(vol2.direction, d) + + def test_type(self): + data = np.zeros((5,5,5)).astype(np.float32) + vol = stk.Volume(data) + self.assertEqual(vol.type, stk.Type.Float) + + data = np.zeros((5,5,5,3)).astype(np.float32) + vol = stk.Volume(data) + self.assertEqual(vol.type, stk.Type.Float3) + + def test_io(self): + data = np.zeros((5,5,5)).astype(np.float32) + vol1 = stk.Volume(data, origin=(2,3,4), spacing=(5,6,7)) + + f = os.path.join(self.temp_dir, 'test.nrrd') + stk.write_volume(f, vol1) + + vol2 = stk.read_volume(f) + np.testing.assert_equal(vol2.origin, vol1.origin) + np.testing.assert_equal(vol2.spacing, vol1.spacing) + np.testing.assert_equal(vol2.direction, vol1.direction) + np.testing.assert_equal(np.array(vol2), np.array(vol1)) From 18fe50355729d9bf2bc1e4bdb058c58bb2609129 Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 12:29:51 +0200 Subject: [PATCH 08/15] Updated CI --- .appveyor.yml | 29 +++++++++++++-------- .travis.yml | 57 ++++++++++++++++++++++++++++-------------- src/stk/CMakeLists.txt | 1 - 3 files changed, 57 insertions(+), 30 deletions(-) diff --git a/.appveyor.yml b/.appveyor.yml index 3743e67..a1503a6 100644 --- a/.appveyor.yml +++ b/.appveyor.yml @@ -17,18 +17,27 @@ build: matrix: fast_finish: true -before_build: - - cmd: if "%PLATFORM%"=="x86" call "C:\Program Files (x86)\Microsoft Visual Studio\2017\Community\VC\Auxiliary\Build\vcvars32.bat" - - cmd: if "%PLATFORM%"=="x64" call "C:\Program Files (x86)\Microsoft Visual Studio\2017\Community\VC\Auxiliary\Build\vcvars64.bat" - - mkdir build - - cd build +environment: + matrix: + - PYTHON: "C:\\Python36-x64" + +install: +- git submodule update --init --recursive +- set PATH=%PYTHON%;%PYTHON%\Scripts;%PATH% +- pip install -r requirements.txt - - cmd: if "%PLATFORM%"=="x86" set GENERATOR="Visual Studio 15 2017" - - cmd: if "%PLATFORM%"=="x64" set GENERATOR="Visual Studio 15 2017 Win64" - - cmd: cmake .. -DSTK_BUILD_TESTS=1 -DSTK_WARNINGS_ARE_ERRORS=1 -DCMAKE_BUILD_TYPE=%CONFIGURATION% -G%GENERATOR% +before_build: +- mkdir build +- cd build +- cmake .. -DCMAKE_BUILD_TYPE=%CONFIGURATION% -DCMAKE_GENERATOR_PLATFORM=%PLATFORM% -DSTK_BUILD_TESTS=1 -DSTK_WARNINGS_ARE_ERRORS=1 build_script: - - cmd: cmake --build . --config %CONFIGURATION% +- cmake --build . --config %CONFIGURATION% +- cd .. +- python setup.py install -DCMAKE_GENERATOR_PLATFORM=%PLATFORM% test_script: - - ctest --output-on-failure \ No newline at end of file +- cd build +- ctest --output-on-failure +- cd .. +- python -m unittest discover test diff --git a/.travis.yml b/.travis.yml index 6b82c07..cc7bebb 100644 --- a/.travis.yml +++ b/.travis.yml @@ -1,37 +1,56 @@ -language: cpp +language: python +python: +- 3.6 matrix: - fast_finish: true include: - - os: linux + - name: trusty + os: linux dist: trusty - - os: linux + env: + - PYTHON_EXECUTABLE=~/virtualenv/python3.6/bin/python + - name: xenial + os: linux dist: xenial - - os: linux + env: + - PYTHON_EXECUTABLE=~/virtualenv/python3.6/bin/python + - name: g++8 + os: linux addons: apt: sources: - - ubuntu-toolchain-r-test + - ubuntu-toolchain-r-test packages: - - g++-6 - env: - - MATRIX_EVAL="CC=gcc-6 && CXX=g++-6" - - os: linux - addons: - apt: - sources: - - ubuntu-toolchain-r-test - packages: - - g++-7 - env: - - MATRIX_EVAL="CC=gcc-7 && CXX=g++-7" + - g++-8 + before_install: + - sudo update-alternatives --install /usr/bin/g++ g++ /usr/bin/g++-8 90 + - sudo update-alternatives --install /usr/bin/gcc gcc /usr/bin/gcc-8 90 + env: + - PYTHON_EXECUTABLE=~/virtualenv/python3.6/bin/python - name: osx os: osx language: generic + before_install: + - brew install llvm libomp + - python3 -m pip install virtualenv + - virtualenv venv -p python3 + - source venv/bin/activate + env: + - PYTHON_EXECUTABLE=../venv/bin/python + - HOMEBREW_NO_AUTO_UPDATE=1 + - CC=/usr/local/opt/llvm/bin/clang + - CXX=/usr/local/opt/llvm/bin/clang++ + +install: +- pip install -r requirements.txt script: - mkdir build - cd build - cmake .. -DCMAKE_BUILD_TYPE=Release -DSTK_BUILD_TESTS=1 -DSTK_WARNINGS_ARE_ERRORS=1 + -DSTK_BUILD_PYTHON_WRAPPER=ON -DPYTHON_EXECUTABLE=$PYTHON_EXECUTABLE - cmake --build . --config Release -- ctest --output-on-failure \ No newline at end of file +- ctest --output-on-failure +- cd .. +- python setup.py install +- python -m unittest discover test diff --git a/src/stk/CMakeLists.txt b/src/stk/CMakeLists.txt index 35aa4e0..7e0bac3 100644 --- a/src/stk/CMakeLists.txt +++ b/src/stk/CMakeLists.txt @@ -61,7 +61,6 @@ if(STK_ITK_BRIDGE) endif() add_library(stk STATIC ${STK_SRCS}) -target_link_libraries(stk PUBLIC pybind11) target_include_directories(stk PUBLIC "${CMAKE_CURRENT_SOURCE_DIR}/.." From 614d229a9cb1f4be1acdb2774ac08a74488e89da Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 12:47:16 +0200 Subject: [PATCH 09/15] GCC warning thanks to pybind11 --- CMakeLists.txt | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 561fe5e..c1093d2 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -57,7 +57,9 @@ if (MSVC) set(EXTRA_FLAGS "${extra_flags} /fp:fast") endif() elseif(CMAKE_CXX_COMPILER_ID MATCHES "GNU") - set(EXTRA_FLAGS "-Wall -Wextra -pedantic -Wno-missing-field-initializers") + set(EXTRA_FLAGS "-Wall -Wextra -pedantic -Wno-missing-field-initializers -Wno-attributes") + + # -Wno-attributes thanks to pybind11 if (STK_WARNINGS_ARE_ERRORS) set(EXTRA_FLAGS "${extra_flags} -Werror") From 751d2208ef292e8c19a9eb42a0a61c92b6101ee7 Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 13:04:49 +0200 Subject: [PATCH 10/15] Warnings are not errors --- .travis.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.travis.yml b/.travis.yml index cc7bebb..dbf23e1 100644 --- a/.travis.yml +++ b/.travis.yml @@ -47,7 +47,7 @@ install: script: - mkdir build - cd build -- cmake .. -DCMAKE_BUILD_TYPE=Release -DSTK_BUILD_TESTS=1 -DSTK_WARNINGS_ARE_ERRORS=1 +- cmake .. -DCMAKE_BUILD_TYPE=Release -DSTK_BUILD_TESTS=1 -DSTK_BUILD_PYTHON_WRAPPER=ON -DPYTHON_EXECUTABLE=$PYTHON_EXECUTABLE - cmake --build . --config Release - ctest --output-on-failure From fc042d1cfa76486c709043c59cb8882e7a3de2ab Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 13:27:43 +0200 Subject: [PATCH 11/15] Some vector calculus functions --- src/python_wrapper/_stk.cpp | 125 ++++++++++++++++++++++++++++++++++++ 1 file changed, 125 insertions(+) diff --git a/src/python_wrapper/_stk.cpp b/src/python_wrapper/_stk.cpp index e2b0d97..8df55cc 100644 --- a/src/python_wrapper/_stk.cpp +++ b/src/python_wrapper/_stk.cpp @@ -5,6 +5,7 @@ #include #include +#include #include #include @@ -287,6 +288,115 @@ stk::Volume make_volume( return vol; } + +std::string divergence_docstring = +R"(Compute the divergence of a displacement field. + +The divergence of a vector field + +.. math:: + f(\boldsymbol{x}) = + (f_1(\boldsymbol{x}), + f_2(\boldsymbol{x}), + f_3(\boldsymbol{x})) + +with :math:`\boldsymbol{x} = (x_1, x_2, x_3)` +is defined as + +.. math:: + \nabla \cdot f (\boldsymbol{x}) = + \sum_{i=1}^3 + \frac{\partial f_i}{\partial x_i} (\boldsymbol{x}) + +.. note:: + All the arrays must be C-contiguous. + +Parameters +---------- +displacement: np.ndarray + Displacement field used to resample the image. + +origin: np.ndarray + Origin of the displacement field. + +spacing: np.ndarray + Spacing of the displacement field. + +direction: Tuple[Int] + Cosine direction matrix of the displacement field. + +Returns +------- +np.ndarray + Scalar volume image containing the divergence of + the input displacement. +)"; + + + +std::string rotor_docstring = +R"(Compute the rotor of a displacement field. + +The rotor of a 3D 3-vector field + +.. math:: + f(\boldsymbol{x}) = + (f_1(\boldsymbol{x}), + f_2(\boldsymbol{x}), + f_3(\boldsymbol{x})) + +with :math:`\boldsymbol{x} = (x_1, x_2, x_3)` +is defined as + +.. math:: + \nabla \times f(\boldsymbol{x}) = + \left( + \frac{\partial f_3}{\partial x_2} - + \frac{\partial f_2}{\partial x_3}, + \frac{\partial f_1}{\partial x_3} - + \frac{\partial f_3}{\partial x_1}, + \frac{\partial f_2}{\partial x_1} - + \frac{\partial f_1}{\partial x_2} + \right) + +.. note:: + All the arrays must be C-contiguous. + +Parameters +---------- +displacement: stk.Volume + Displacement field used to resample the image. + +Returns +------- +stk.Volume + Vector volume image containing the rotor of + the input displacement. +)"; + + + +std::string circulation_density_docstring = +R"(Compute the circulation density of a displacement field. + +The circulation density for a 3D 3-vector field is defined as the +norm of the rotor. + +.. note:: + All the arrays must be C-contiguous. + +Parameters +---------- +displacement: stk.Volume + Displacement field used to resample the image. + +Returns +------- +stk.Volume + Vector volume image containing the circulation + density of the input displacement. +)"; + PYBIND11_MODULE(_stk, m) { py::enum_(m, "Type") @@ -346,6 +456,21 @@ PYBIND11_MODULE(_stk, m) m.def("read_volume", &stk::read_volume, ""); m.def("write_volume", &stk::write_volume, ""); + m.def("divergence", + &stk::divergence, + divergence_docstring.c_str() + ); + + m.def("rotor", + &stk::rotor, + rotor_docstring.c_str() + ); + + m.def("circulation_density", + &stk::circulation_density, + circulation_density_docstring.c_str() + ); + // Translate relevant exception types. The exceptions not handled // here will be translated autmatically according to pybind11's rules. py::register_exception_translator([](std::exception_ptr p) { From 6a1a383d87ae6d41537a06207f6bafc951d01d7b Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 13:50:24 +0200 Subject: [PATCH 12/15] Implict cast from numpy array to stk::Volume --- src/python_wrapper/_stk.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/python_wrapper/_stk.cpp b/src/python_wrapper/_stk.cpp index 8df55cc..529b18f 100644 --- a/src/python_wrapper/_stk.cpp +++ b/src/python_wrapper/_stk.cpp @@ -453,6 +453,8 @@ PYBIND11_MODULE(_stk, m) .def_property("direction", &stk::Volume::direction, &stk::Volume::set_direction) ; + py::implicitly_convertible(); + m.def("read_volume", &stk::read_volume, ""); m.def("write_volume", &stk::write_volume, ""); From 1d87a5b978315373e5914f662086eab38de53177 Mon Sep 17 00:00:00 2001 From: Simon Date: Wed, 2 Oct 2019 14:43:45 +0200 Subject: [PATCH 13/15] CUDA: Missing C++11 option --- CMakeLists.txt | 1 + 1 file changed, 1 insertion(+) diff --git a/CMakeLists.txt b/CMakeLists.txt index c1093d2..9987a62 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -91,6 +91,7 @@ endif() if (STK_USE_CUDA) add_definitions(-DSTK_USE_CUDA) enable_language(CUDA) + set(CMAKE_CUDA_STANDARD 11) if (STK_ENABLE_FAST_MATH) set(CMAKE_CUDA_FLAGS "${CMAKE_CUDA_FLAGS} --use_fast_math") From cff8c1270e7d277e5e0f85adabb9fd896357beb9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Simon=20Ekstr=C3=B6m?= Date: Wed, 2 Oct 2019 15:07:35 +0200 Subject: [PATCH 14/15] Update README.md --- README.md | 34 ++++++++++++++++++++++++++++++++++ 1 file changed, 34 insertions(+) diff --git a/README.md b/README.md index 835e66c..f61b17d 100644 --- a/README.md +++ b/README.md @@ -23,3 +23,37 @@ CMake was not compatible with the one required by CUDA, it is possible to specify a different executable with `-DCMAKE_CUDA_FLAGS="-ccbin gcc-XX"`, where `gcc-XX` is a version of `gcc` compatible with your CUDA version. +# Python API + +A minimalistic Python API is also provided. + +## Install + +```bash +python setup.py install +``` + +## Example usage + +```python +import stk +import numpy as np + +# Create volume directly from numpy +vol = stk.Volume(np.zeros((5,5,5)).astype(np.float32), spacing=(2,2,2)) +# or read volume from file +vol = stk.read_volume('test.nrrd') + +# Modify data (numpy array points to the volume data) +data = np.array(vol, copy=False) +data[0:10] = 0.0 + +# Access meta data +vol.origin = (2, 2, 2) +vol.spacing = (3, 3, 3) + +# Write volume +stk.write_volume('test-out.nrrd', vol) +``` + + From 865b0a8df895ce59840253c61a0efab1f3a896f9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Simon=20Ekstr=C3=B6m?= Date: Wed, 2 Oct 2019 15:35:55 +0200 Subject: [PATCH 15/15] Bump package version to 0.4 --- setup.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/setup.py b/setup.py index 3fa7ff0..31dd519 100644 --- a/setup.py +++ b/setup.py @@ -76,7 +76,7 @@ def build_extension(self, ext): setup( name='python-stk', - version='0.1', + version='0.4', author='Simon Ekström', author_email='', description='',