diff --git a/.flake8 b/.flake8 index c443919c..fb450003 100644 --- a/.flake8 +++ b/.flake8 @@ -1,15 +1,17 @@ [flake8] # W503: line break before binary operator is actually considered best-practice # E203: spaces around complex variables in slices are pep-right -# F401: unused imports in __init__.py-s +# F401: unused imports in __init__.py-s and compat.py # I251: allow absolute imports in upper files # B028: !r is not supported for python<3.8 # W604: backticks in str-s are ok -ignore = W503,E203,B028,W604 +# S101: asserts are ok for now +# S102: exec in build scripts is ok +ignore = W503,E203,B028,W604,S101 per-file-ignores = - __init__.py:F401 - tests/*:I251 - benchmarks/*:I251 + setup.py,_build_utils.py:S102 + __init__.py,compat.py:F401 + tests/*,benchmarks/*:I251 max-line-length = 120 banned-modules = imops.* = Use relative imports diff --git a/.github/workflows/release.yml b/.github/workflows/release.yml index 2004088d..acd64184 100644 --- a/.github/workflows/release.yml +++ b/.github/workflows/release.yml @@ -35,7 +35,7 @@ jobs: needs: [ check_version ] strategy: matrix: - os: [ ubuntu-22.04, windows-2019, macOS-11 ] + os: [ ubuntu-22.04, windows-2019, macOS-12 ] name: Build wheels on ${{ matrix.os }} runs-on: ${{ matrix.os }} @@ -48,7 +48,7 @@ jobs: - name: Install cibuildwheel run: python -m pip install cibuildwheel==2.17.0 - name: Install gcc for mac - if: matrix.os == 'macOS-11' + if: matrix.os == 'macOS-12' run: | brew install llvm libomp echo $PATH @@ -81,7 +81,7 @@ jobs: CIBW_ENVIRONMENT_MACOS: > PATH="/usr/local/opt/llvm/bin:$PATH" LDFLAGS="-L/usr/local/opt/llvm/lib" CPPFLAGS="-I/usr/local/opt/llvm/include" CIBW_BUILD: cp37-* cp38-* cp39-* cp310-* cp311-* cp312-* - CIBW_BEFORE_BUILD_LINUX: 'if [ $(python -c "import sys; print(sys.version_info[1])") -ge 9 ]; then python -m pip install "numpy<2.0.0" --config-settings=setup-args="-Dallow-noblas=true"; fi' + CIBW_BEFORE_BUILD_LINUX: 'if [ $(python -c "import sys; print(sys.version_info[1])") -ge 9 ]; then python -m pip install "numpy<3.0.0" --config-settings=setup-args="-Dallow-noblas=true"; fi' - uses: actions/upload-artifact@v3 with: path: ./wheelhouse/*.whl diff --git a/.github/workflows/test_build.yml b/.github/workflows/test_build.yml index a6860358..9b444b19 100644 --- a/.github/workflows/test_build.yml +++ b/.github/workflows/test_build.yml @@ -9,7 +9,7 @@ jobs: build_wheels: strategy: matrix: - os: [ubuntu-22.04, windows-2019, macOS-11 ] + os: [ubuntu-22.04, windows-2019, macOS-12 ] name: Build wheels on ${{ matrix.os }} runs-on: ${{ matrix.os }} @@ -22,7 +22,7 @@ jobs: - name: Install cibuildwheel run: python -m pip install cibuildwheel==2.17.0 - name: Install gcc for mac - if: matrix.os == 'macOS-11' + if: matrix.os == 'macOS-12' run: | brew install llvm libomp echo $PATH @@ -56,4 +56,4 @@ jobs: PATH="/usr/local/opt/llvm/bin:$PATH" LDFLAGS="-L/usr/local/opt/llvm/lib" CPPFLAGS="-I/usr/local/opt/llvm/include" CIBW_BUILD: cp39-* cp312-* CIBW_SKIP: "*musllinux* *manylinux_x86_64" - CIBW_BEFORE_BUILD_LINUX: 'if [ $(python -c "import sys; print(sys.version_info[1])") -ge 9 ]; then python -m pip install "numpy<2.0.0" --config-settings=setup-args="-Dallow-noblas=true"; fi' + CIBW_BEFORE_BUILD_LINUX: 'if [ $(python -c "import sys; print(sys.version_info[1])") -ge 9 ]; then python -m pip install "numpy<3.0.0" --config-settings=setup-args="-Dallow-noblas=true"; fi' diff --git a/.gitignore b/.gitignore index b410ab7e..87e1733b 100644 --- a/.gitignore +++ b/.gitignore @@ -1,6 +1,7 @@ __pycache__/ .pytest_cache/ build/ +imops/src/*.cpp imops/src/*.c imops/src/*.so imops/src/_fast*.pyx diff --git a/_build_utils.py b/_build_utils.py index 17b51925..70e24d7b 100644 --- a/_build_utils.py +++ b/_build_utils.py @@ -80,7 +80,16 @@ def get_ext_modules(): extra_compile_args=args + cpp_args, extra_link_args=args + cpp_args, language='c++', - ) + ), + Extension( + f'{name}.src._utils', + [f'{name}/src/_utils.pyx'], + language='c++', + include_dirs=[LazyImport('numpy')], + extra_compile_args=args + cpp_args, + extra_link_args=args + cpp_args, + define_macros=[('NPY_NO_DEPRECATED_API', 'NPY_1_7_API_VERSION')], + ), ] for module in modules: libraries = [] @@ -93,6 +102,7 @@ def get_ext_modules(): if not on_windows: libraries.append('m') + # TODO: mb throw `ffast-math` away? # FIXME: import of `ffast-math` compiled modules changes global FPU state, so now `fast=True` will just # fallback to standard `-O2` compiled versions until https://github.com/neuro-ml/imops/issues/37 is resolved # for prefix, additional_args in zip(['', 'fast_'], [[], ['-ffast-math']]) @@ -101,6 +111,7 @@ def get_ext_modules(): Extension( f'{name}.src._{prefix}{module}', [f'{name}/src/_{prefix}{module}.pyx'], + language='c', include_dirs=include_dirs, library_dirs=library_dirs, libraries=libraries, diff --git a/asv.conf.json b/asv.conf.json index 85c28d0b..8dcee201 100644 --- a/asv.conf.json +++ b/asv.conf.json @@ -3,16 +3,16 @@ "project": "imops", "project_url": "https://github.com/neuro-ml/imops", "repo": ".", - "branches": ["master", "dev"], + "branches": ["dev"], "dvcs": "git", "environment_type": "conda", "show_commit_url": "https://github.com/neuro-ml/imops/commit", - "pythons": ["3.10"], + "pythons": ["3.12.3"], "conda_channels": ["defaults", "conda-forge"], "matrix": { "req": { "numpy": [], - "Cython": ["0.29.36"], + "Cython": ["3.0.10"], "scipy": [], "scikit-image": [], "numba": [], diff --git a/docs/index.md b/docs/index.md index 4e3d33b2..90ca24f5 100644 --- a/docs/index.md +++ b/docs/index.md @@ -63,3 +63,5 @@ pip install imops[numba] # additionally install Numba backend ::: imops.radon.radon ::: imops.radon.inverse_radon + +::: imops.utils.isin diff --git a/imops/__version__.py b/imops/__version__.py index a02276f7..e4e49b3b 100644 --- a/imops/__version__.py +++ b/imops/__version__.py @@ -1 +1 @@ -__version__ = '0.8.8' +__version__ = '0.9.0' diff --git a/imops/compat.py b/imops/compat.py new file mode 100644 index 00000000..07a8a1df --- /dev/null +++ b/imops/compat.py @@ -0,0 +1,14 @@ +try: + from numpy.lib.array_utils import normalize_axis_tuple +except ModuleNotFoundError: + from numpy.core.numeric import normalize_axis_tuple + +try: + from numpy.exceptions import VisibleDeprecationWarning +except ModuleNotFoundError: + from numpy import VisibleDeprecationWarning + +try: + from scipy.ndimage._morphology import _ni_support +except ImportError: + from scipy.ndimage.morphology import _ni_support diff --git a/imops/measure.py b/imops/measure.py index a19166d2..15c3b04d 100644 --- a/imops/measure.py +++ b/imops/measure.py @@ -29,6 +29,7 @@ } +# TODO: Make it work and test on immutable arrays as soon as `cc3d` package is fixed def label( label_image: np.ndarray, background: int = None, @@ -209,8 +210,8 @@ def center_of_mass( src_center_of_mass = _fast_labeled_center_of_mass if backend.fast else _labeled_center_of_mass - if array.dtype != 'float64': - array = array.astype(float) + if array.dtype not in ('float32', 'float64'): + array = array.astype(np.float32) n_dummy = 3 - ndim if n_dummy: diff --git a/imops/morphology.py b/imops/morphology.py index 9085976e..e2256ab9 100644 --- a/imops/morphology.py +++ b/imops/morphology.py @@ -5,13 +5,6 @@ from edt import edt from scipy.ndimage import distance_transform_edt as scipy_distance_transform_edt, generate_binary_structure from scipy.ndimage._nd_image import euclidean_feature_transform - - -try: - from scipy.ndimage._morphology import _ni_support -except ImportError: - from scipy.ndimage.morphology import _ni_support - from skimage.morphology import ( binary_closing as scipy_binary_closing, binary_dilation as scipy_binary_dilation, @@ -21,6 +14,7 @@ from .backend import BackendLike, Cython, Scipy, resolve_backend from .box import add_margin, box_to_shape, mask_to_box, shape_to_box +from .compat import _ni_support from .crop import crop_to_box from .pad import restore_crop from .src._fast_morphology import ( diff --git a/imops/radon.py b/imops/radon.py index 81b0e3e7..32d5f767 100644 --- a/imops/radon.py +++ b/imops/radon.py @@ -4,6 +4,7 @@ from scipy.fftpack import fft, ifft from .backend import BackendLike, resolve_backend +from .compat import normalize_axis_tuple from .numeric import copy from .src._backprojection import backprojection3d from .src._fast_backprojection import backprojection3d as fast_backprojection3d @@ -204,7 +205,7 @@ def normalize_axes(x: np.ndarray, axes): raise ValueError('For arrays of higher dimensionality the `axis` arguments is required') axes = [0, 1] - axes = np.core.numeric.normalize_axis_tuple(axes, x.ndim, 'axes') + axes = normalize_axis_tuple(axes, x.ndim, 'axes') x = np.moveaxis(x, axes, (-2, -1)) extra = x.shape[:-2] x = x.reshape(-1, *x.shape[-2:]) diff --git a/imops/src/_backprojection.pyx b/imops/src/_backprojection.pyx index 2024c552..a3e86154 100644 --- a/imops/src/_backprojection.pyx +++ b/imops/src/_backprojection.pyx @@ -19,7 +19,7 @@ ctypedef np.uint8_t uint8 @cython.boundscheck(False) @cython.wraparound(False) -cdef inline FLOAT interpolate(FLOAT x, FLOAT* ys, FLOAT radius, FLOAT right_limit) nogil: +cdef inline FLOAT interpolate(FLOAT x, const FLOAT* ys, FLOAT radius, FLOAT right_limit) noexcept nogil: cdef Py_ssize_t idx cdef FLOAT val, value @@ -40,8 +40,8 @@ cdef inline FLOAT interpolate(FLOAT x, FLOAT* ys, FLOAT radius, FLOAT right_limi @cython.boundscheck(False) @cython.wraparound(False) -cdef FLOAT accumulate(FLOAT x, FLOAT y, FLOAT* sinuses, FLOAT* cosinuses, FLOAT* ys, - Py_ssize_t size, Py_ssize_t image_size, FLOAT radius, FLOAT right_limit) nogil: +cdef FLOAT accumulate(FLOAT x, FLOAT y, const FLOAT* sinuses, const FLOAT* cosinuses, const FLOAT* ys, + Py_ssize_t size, Py_ssize_t image_size, FLOAT radius, FLOAT right_limit) noexcept nogil: cdef FLOAT accumulator = 0 cdef Py_ssize_t k @@ -52,9 +52,9 @@ cdef FLOAT accumulate(FLOAT x, FLOAT y, FLOAT* sinuses, FLOAT* cosinuses, FLOAT* @cython.boundscheck(False) @cython.wraparound(False) -cpdef FLOAT[:, :, :] backprojection3d(FLOAT[:, :, :] sinogram, FLOAT[:] theta, FLOAT[:] xs, - uint8[:, :] inside_circle, FLOAT fill_value, int image_size, int output_size, - Py_ssize_t num_threads): +cpdef FLOAT[:, :, :] backprojection3d(const FLOAT[:, :, :] sinogram, const FLOAT[:] theta, const FLOAT[:] xs, + const uint8[:, :] inside_circle, FLOAT fill_value, int image_size, + int output_size, Py_ssize_t num_threads): cdef FLOAT[:, :, :] result = np.zeros_like(sinogram, shape=(len(sinogram), output_size, output_size)) cdef Py_ssize_t slc, i, j, n_angles = len(theta), n_slices = len(sinogram) cdef FLOAT min_val = image_size // 2, right_lim = image_size - 1 diff --git a/imops/src/_measure.pyx b/imops/src/_measure.pyx index 16f4898a..2be35ede 100644 --- a/imops/src/_measure.pyx +++ b/imops/src/_measure.pyx @@ -7,11 +7,14 @@ import numpy as np +cimport cython cimport numpy as np from cython.parallel import prange +ctypedef cython.floating FLOAT + ctypedef fused LABEL: signed char short @@ -23,7 +26,7 @@ ctypedef fused LABEL: unsigned long long -cdef inline Py_ssize_t _find(LABEL num, LABEL[:] nums) nogil: +cdef inline Py_ssize_t _find(LABEL num, const LABEL[:] nums) noexcept nogil: cdef Py_ssize_t i for i in range(len(nums)): @@ -33,15 +36,16 @@ cdef inline Py_ssize_t _find(LABEL num, LABEL[:] nums) nogil: return -1 -def _labeled_center_of_mass(double[:, :, :] nums, LABEL[:, :, :] labels, LABEL[:] index) -> np.ndarray: - cdef double[:, :, ::1] contiguous_nums = np.ascontiguousarray(nums) - cdef LABEL[:, :, ::1] contiguous_labels = np.ascontiguousarray(labels) - cdef LABEL[:] contiguous_index = np.ascontiguousarray(index) +def _labeled_center_of_mass(const FLOAT[:, :, :] nums, const LABEL[:, :, :] labels, + const LABEL[:] index) -> np.ndarray: + cdef const FLOAT[:, :, ::1] contiguous_nums = np.ascontiguousarray(nums) + cdef const LABEL[:, :, ::1] contiguous_labels = np.ascontiguousarray(labels) + cdef const LABEL[:] contiguous_index = np.ascontiguousarray(index) cdef Py_ssize_t index_len = len(index) - cdef double[:, ::1] output = np.zeros((index_len, 3)) - cdef double[:] normalizers = np.zeros(index_len) + cdef FLOAT[:, ::1] output = np.zeros_like(nums, shape=(index_len, 3)) + cdef FLOAT[:] normalizers = np.zeros_like(nums, shape=(index_len,)) cdef Py_ssize_t rows = nums.shape[0], cols = nums.shape[1], dims = nums.shape[2] cdef Py_ssize_t i, j, k, pos @@ -66,11 +70,11 @@ def _labeled_center_of_mass(double[:, :, :] nums, LABEL[:, :, :] labels, LABEL[: return np.asarray(output) -def _center_of_mass(double[:, :, :] nums, Py_ssize_t num_threads) -> np.ndarray: - cdef double[:, :, ::1] contiguous_nums = np.ascontiguousarray(nums) +def _center_of_mass(const FLOAT[:, :, :] nums, Py_ssize_t num_threads) -> np.ndarray: + cdef const FLOAT[:, :, ::1] contiguous_nums = np.ascontiguousarray(nums) - cdef double output_x = 0, output_y = 0, output_z = 0 - cdef double normalizer = 0 + cdef FLOAT output_x = 0, output_y = 0, output_z = 0 + cdef FLOAT normalizer = 0 cdef Py_ssize_t rows = nums.shape[0], cols = nums.shape[1], dims = nums.shape[2] cdef Py_ssize_t i, j, k diff --git a/imops/src/_morphology.pyx b/imops/src/_morphology.pyx index c47e04bb..6a49a2bc 100644 --- a/imops/src/_morphology.pyx +++ b/imops/src/_morphology.pyx @@ -22,7 +22,7 @@ cdef struct array_iterator: int backstrides[3] -cdef int init_array_iterator(array_iterator *arr_iter, int *shape, int *strides) nogil: +cdef int init_array_iterator(array_iterator *arr_iter, const int *shape, const int *strides) noexcept nogil: cdef int i arr_iter.rank_m1 = 2 @@ -45,11 +45,11 @@ cdef struct filter_iterator: cdef int init_filter_iterator( filter_iterator *filter_iter, - int *a_shape, - int *f_shape, + const int *a_shape, + const int *f_shape, int filter_size, np.uint8_t is_dilation -) nogil: +) noexcept nogil: cdef int i, rank = 3, step, orgn filter_iter[0].strides[rank - 1] = filter_size @@ -71,12 +71,12 @@ cdef int init_filter_iterator( cdef int init_filter_offsets( - int *a_shape, int *a_strides, - np.uint8_t *footprint, int *f_shape, + const int *a_shape, const int *a_strides, + const np.uint8_t *footprint, const int *f_shape, int footprint_size, int **offsets, int *border_flag_value, np.uint8_t is_dilation -) nogil: +) noexcept nogil: cdef int i, j, k cdef int filter_size = 1, offsets_size = 1 cdef int max_size = 0, max_stride = 0 @@ -165,7 +165,7 @@ cdef int init_filter_offsets( return 0 -cdef inline int filter_iterator_offset(filter_iterator *filter_iter, int *coordinates) nogil: +cdef inline int filter_iterator_offset(const filter_iterator *filter_iter, const int *coordinates) noexcept nogil: cdef int i, position, offset = 0 for i in range(3): @@ -184,14 +184,14 @@ cdef inline int filter_iterator_offset(filter_iterator *filter_iter, int *coordi cdef inline int worker( - np.uint8_t *input, np.uint8_t *footprint, np.uint8_t *output, + const np.uint8_t *input, const np.uint8_t *footprint, np.uint8_t *output, int *a_shape, int *a_strides, int *f_shape, int footprint_size, int *offsets, int border_flag_value, int start, int end, np.uint8_t border_value, np.uint8_t is_dilation, -) nogil: +) noexcept nogil: cdef np.uint8_t _true = not is_dilation cdef np.uint8_t _false = not _true @@ -252,15 +252,15 @@ cdef inline int worker( def _binary_operation( - np.uint8_t[:, :, :] input, - np.uint8_t[:, :, :] footprint, + const np.uint8_t[:, :, :] input, + const np.uint8_t[:, :, :] footprint, np.uint8_t[:, :, ::1] out, Py_ssize_t num_threads, np.uint8_t border_value, np.uint8_t is_dilation, ) -> np.ndarray: - cdef np.uint8_t[:, :, ::1] c_input = np.ascontiguousarray(input) - cdef np.uint8_t[:, :, ::1] c_footprint = np.ascontiguousarray(footprint) + cdef const np.uint8_t[:, :, ::1] c_input = np.ascontiguousarray(input) + cdef const np.uint8_t[:, :, ::1] c_footprint = np.ascontiguousarray(footprint) cdef int f_shape[3] cdef int a_shape[3] @@ -311,8 +311,8 @@ def _binary_operation( def _binary_erosion( - np.uint8_t[:, :, :] input, - np.uint8_t[:, :, :] footprint, + const np.uint8_t[:, :, :] input, + const np.uint8_t[:, :, :] footprint, np.uint8_t[:, :, ::1] out, Py_ssize_t num_threads, ) -> np.ndarray: @@ -320,8 +320,8 @@ def _binary_erosion( def _binary_dilation( - np.uint8_t[:, :, :] input, - np.uint8_t[:, :, :] footprint, + const np.uint8_t[:, :, :] input, + const np.uint8_t[:, :, :] footprint, np.uint8_t[:, :, ::1] out, Py_ssize_t num_threads, ) -> np.ndarray: diff --git a/imops/src/_radon.pyx b/imops/src/_radon.pyx index f9f91e8f..65c6115c 100644 --- a/imops/src/_radon.pyx +++ b/imops/src/_radon.pyx @@ -18,14 +18,14 @@ ctypedef cython.integral INT # most of this code is taken from skimage and optimized for direct Radon transform -cdef inline FLOAT get_pixel2d(FLOAT* image, Py_ssize_t size, long r, long c) nogil: +cdef inline FLOAT get_pixel2d(const FLOAT* image, Py_ssize_t size, long r, long c) noexcept nogil: if (r < 0) or (r >= size) or (c < 0) or (c >= size): return 0 else: return image[r * size + c] -cdef inline FLOAT interpolate2d(FLOAT* image, Py_ssize_t size, FLOAT r, FLOAT c) nogil: +cdef inline FLOAT interpolate2d(const FLOAT* image, Py_ssize_t size, FLOAT r, FLOAT c) noexcept nogil: cdef FLOAT dr, dc cdef long minr, minc, maxr, maxc @@ -49,8 +49,8 @@ cdef inline FLOAT interpolate2d(FLOAT* image, Py_ssize_t size, FLOAT r, FLOAT c) return (1 - dr) * top + dr * bottom -cdef inline FLOAT accumulate(FLOAT* image, Py_ssize_t* size, FLOAT* sin, FLOAT* cos, - FLOAT* r_shift, FLOAT* c_shift, Py_ssize_t j, INT* limit) nogil: +cdef inline FLOAT accumulate(const FLOAT* image, const Py_ssize_t* size, const FLOAT* sin, const FLOAT* cos, + const FLOAT* r_shift, const FLOAT* c_shift, Py_ssize_t j, const INT* limit) noexcept nogil: cdef FLOAT result = 0 cdef Py_ssize_t i @@ -64,9 +64,9 @@ cdef inline FLOAT accumulate(FLOAT* image, Py_ssize_t* size, FLOAT* sin, FLOAT* return result -def radon3d(FLOAT[:, :, :] image, FLOAT[:] theta, INT[:] limits, Py_ssize_t num_threads): +def radon3d(const FLOAT[:, :, :] image, const FLOAT[:] theta, const INT[:] limits, Py_ssize_t num_threads): cdef Py_ssize_t size = image.shape[1], n_slices = image.shape[0], n_angles = len(theta) - cdef FLOAT[:, :, ::1] img = np.ascontiguousarray(image) + cdef const FLOAT[:, :, ::1] img = np.ascontiguousarray(image) cdef FLOAT[:, :, :] out = np.zeros_like(img, shape=(n_slices, size, n_angles)) sins = np.sin(theta) diff --git a/imops/src/_utils.pyx b/imops/src/_utils.pyx new file mode 100644 index 00000000..687b367a --- /dev/null +++ b/imops/src/_utils.pyx @@ -0,0 +1,46 @@ +# cython: binding=False +# cython: boundscheck=False +# cython: wraparound=False +# cython: initializedcheck=False +# cython: nonecheck=False +# cython: overflowcheck=True +# cython: overflowcheck.fold=False +# cython: embedsignature=False +# cython: embedsignature.format=c +# cython: c_api_binop_methods=True +# cython: profile=False +# cython: linetrace=False +# cython: infer_types=False +# cython: language_level=3 +# cython: cpp_locals=False + + +from libcpp.unordered_set cimport unordered_set + +import numpy as np + +cimport numpy as np + +from cython.parallel import prange + + +ctypedef fused NUM: + short + int + long long + + +cpdef void _isin(const NUM[:] elements, const NUM[:] test_elements, np.uint8_t[:] res, Py_ssize_t num_threads): + cdef unordered_set[NUM] test_elements_set + cdef Py_ssize_t i + cdef Py_ssize_t elements_len = elements.shape[0] + cdef Py_ssize_t test_elements_len = test_elements.shape[0] + + test_elements_set.reserve(test_elements_len) + + with nogil: + for i in range(test_elements_len): + test_elements_set.insert(test_elements[i]) + + for i in prange(elements_len, num_threads=num_threads): + res[i] = test_elements_set.count(elements[i]) diff --git a/imops/src/_zoom.pyx b/imops/src/_zoom.pyx index 9aabbc81..c5107ed5 100644 --- a/imops/src/_zoom.pyx +++ b/imops/src/_zoom.pyx @@ -27,12 +27,12 @@ ctypedef fused NUM: double -def _interp1d(FLOAT[:, :, :] input, +def _interp1d(const FLOAT[:, :, :] input, double[:] old_locations, double[:] new_locations, np.uint8_t bounds_error, double fill_value, np.uint8_t extrapolate, np.uint8_t assume_sorted, Py_ssize_t num_threads) -> np.ndarray: cdef Py_ssize_t rows = input.shape[0], cols = input.shape[1], dims = len(new_locations) - cdef FLOAT[:, :, ::1] contiguous_input = np.ascontiguousarray(input) + cdef const FLOAT[:, :, ::1] contiguous_input = np.ascontiguousarray(input) cdef FLOAT[:, :, ::1] interpolated = np.empty_like(contiguous_input, shape=(rows, cols, dims)) cdef double[:] dd = np.zeros(dims) @@ -120,49 +120,49 @@ def _interp1d(FLOAT[:, :, :] input, return np.asarray(interpolated) -cdef inline double get_slope(double x1, double y1, double x2, double y2) nogil: +cdef inline double get_slope(double x1, double y1, double x2, double y2) noexcept nogil: return (y2 - y1) / (x2 - x1) -cdef inline NUM get_pixel3d(NUM* input, +cdef inline NUM get_pixel3d(const NUM* input, Py_ssize_t rows, Py_ssize_t cols, Py_ssize_t dims, Py_ssize_t r, Py_ssize_t c, Py_ssize_t d, - NUM cval) nogil: + NUM cval) noexcept nogil: if (r < 0) or (r >= rows) or (c < 0) or (c >= cols) or (d < 0) or (d >= dims): return cval return input[r * cols * dims + c * dims + d] -cdef inline NUM get_pixel4d(NUM* input, +cdef inline NUM get_pixel4d(const NUM* input, Py_ssize_t stride1, Py_ssize_t stride2, Py_ssize_t stride3, Py_ssize_t stride4, Py_ssize_t dim1, Py_ssize_t dim2, Py_ssize_t dim3, Py_ssize_t dim4, Py_ssize_t c1, Py_ssize_t c2, Py_ssize_t c3, Py_ssize_t c4, - NUM cval) nogil: + NUM cval) noexcept nogil: if (c1 < 0) or (c1 >= dim1) or (c2 < 0) or (c2 >= dim2) or (c3 < 0) or (c3 >= dim3) or (c4 < 0) or (c4 >= dim4): return cval return input[c1 * stride1 + c2 * stride2 + c3 * stride3 + c4 * stride4] -cdef inline double adjusted_coef(Py_ssize_t old_n, Py_ssize_t new_n) nogil: +cdef inline double adjusted_coef(Py_ssize_t old_n, Py_ssize_t new_n) noexcept nogil: if new_n == 1: return old_n return (old_n - 1) / (new_n - 1) cdef inline double distance3d(double x1, double y1, double z1, - double x2, double y2, double z2) nogil: + double x2, double y2, double z2) noexcept nogil: return sqrt((x1 - x2) ** 2 + (y1 - y2) ** 2 + (z1 - z2) ** 2) cdef inline double distance4d(double x1, double y1, double z1, double d1, - double x2, double y2, double z2, double d2) nogil: + double x2, double y2, double z2, double d2) noexcept nogil: return sqrt((x1 - x2) ** 2 + (y1 - y2) ** 2 + (z1 - z2) ** 2 + (d1 - d2) ** 2) -cdef inline FLOAT interpolate3d_linear(FLOAT* input, +cdef inline FLOAT interpolate3d_linear(const FLOAT* input, Py_ssize_t rows, Py_ssize_t cols, Py_ssize_t dims, double r, double c, double d, - FLOAT cval) nogil: + FLOAT cval) noexcept nogil: cdef double dr, dc, dd cdef Py_ssize_t minr, minc, mind, maxr, maxc, maxd @@ -200,10 +200,10 @@ cdef inline FLOAT interpolate3d_linear(FLOAT* input, # FIXME: figure out how to avoid duplication between `linear` and `learest` without slowing down -cdef inline NUM interpolate3d_nearest(NUM* input, +cdef inline NUM interpolate3d_nearest(const NUM* input, Py_ssize_t rows, Py_ssize_t cols, Py_ssize_t dims, double r, double c, double d, - NUM cval) nogil: + NUM cval) noexcept nogil: cdef double minr, minc, mind, maxr, maxc, maxd, curr, curc, curd, distance, min_distance = 3.0 cdef short i, j, k, i_nearest = -1, j_nearest = -1, k_nearest = -1 @@ -248,11 +248,11 @@ cdef inline NUM interpolate3d_nearest(NUM* input, ) -cdef inline FLOAT interpolate4d_linear(FLOAT* input, +cdef inline FLOAT interpolate4d_linear(const FLOAT* input, Py_ssize_t stride1, Py_ssize_t stride2, Py_ssize_t stride3, Py_ssize_t stride4, Py_ssize_t dim1, Py_ssize_t dim2, Py_ssize_t dim3, Py_ssize_t dim4, double c1, double c2, double c3, double c4, - FLOAT cval) nogil: + FLOAT cval) noexcept nogil: cdef double dc1, dc2, dc3, dc4 cdef Py_ssize_t minc1, minc2, minc3, minc4, maxc1, maxc2, maxc3, maxc4 @@ -310,11 +310,11 @@ cdef inline FLOAT interpolate4d_linear(FLOAT* input, return c0_ * (1 - dc4) + c1_ * dc4 -cdef inline NUM interpolate4d_nearest(NUM* input, +cdef inline NUM interpolate4d_nearest(const NUM* input, Py_ssize_t stride1, Py_ssize_t stride2, Py_ssize_t stride3, Py_ssize_t stride4, Py_ssize_t dim1, Py_ssize_t dim2, Py_ssize_t dim3, Py_ssize_t dim4, double c1, double c2, double c3, double c4, - NUM cval) nogil: + NUM cval) noexcept nogil: cdef double minc1, minc2, minc3, minc4, maxc1, maxc2, maxc3, maxc4, curc1, curc2, curc3, curc4 cdef double distance, min_distance = 3.0 cdef short i1, i2, i3, i4, i1_nearest = -1, i2_nearest = -1, i3_nearest = -1, i4_nearest = -1 @@ -369,8 +369,8 @@ cdef inline NUM interpolate4d_nearest(NUM* input, ) -def _zoom3d_linear(FLOAT[:, :, :] input, double[:] zoom, FLOAT cval, Py_ssize_t num_threads): - cdef FLOAT[:, :, ::1] contiguous_input = np.ascontiguousarray(input) +def _zoom3d_linear(const FLOAT[:, :, :] input, double[:] zoom, FLOAT cval, Py_ssize_t num_threads): + cdef const FLOAT[:, :, ::1] contiguous_input = np.ascontiguousarray(input) cdef Py_ssize_t old_rows = input.shape[0], old_cols = input.shape[1], old_dims = input.shape[2] cdef double row_coef = zoom[0], col_coef = zoom[1], dim_coef = zoom[2] @@ -399,8 +399,8 @@ def _zoom3d_linear(FLOAT[:, :, :] input, double[:] zoom, FLOAT cval, Py_ssize_t return np.asarray(zoomed) -def _zoom3d_nearest(NUM[:, :, :] input, double[:] zoom, NUM cval, Py_ssize_t num_threads): - cdef NUM[:, :, ::1] contiguous_input = np.ascontiguousarray(input) +def _zoom3d_nearest(const NUM[:, :, :] input, double[:] zoom, NUM cval, Py_ssize_t num_threads): + cdef const NUM[:, :, ::1] contiguous_input = np.ascontiguousarray(input) cdef Py_ssize_t old_rows = input.shape[0], old_cols = input.shape[1], old_dims = input.shape[2] cdef double row_coef = zoom[0], col_coef = zoom[1], dim_coef = zoom[2] @@ -429,8 +429,8 @@ def _zoom3d_nearest(NUM[:, :, :] input, double[:] zoom, NUM cval, Py_ssize_t num return np.asarray(zoomed) -def _zoom4d_linear(FLOAT[:, :, :, :] input, double[:] zoom, FLOAT cval, Py_ssize_t num_threads): - cdef FLOAT[:, :, :, ::1] contiguous_input = np.ascontiguousarray(input) +def _zoom4d_linear(const FLOAT[:, :, :, :] input, double[:] zoom, FLOAT cval, Py_ssize_t num_threads): + cdef const FLOAT[:, :, :, ::1] contiguous_input = np.ascontiguousarray(input) cdef Py_ssize_t old_dim1 = input.shape[0] cdef Py_ssize_t old_dim2 = input.shape[1] @@ -478,8 +478,8 @@ def _zoom4d_linear(FLOAT[:, :, :, :] input, double[:] zoom, FLOAT cval, Py_ssize return np.asarray(zoomed) -def _zoom4d_nearest(NUM[:, :, :, :] input, double[:] zoom, NUM cval, Py_ssize_t num_threads): - cdef NUM[:, :, :, ::1] contiguous_input = np.ascontiguousarray(input) +def _zoom4d_nearest(const NUM[:, :, :, :] input, double[:] zoom, NUM cval, Py_ssize_t num_threads): + cdef const NUM[:, :, :, ::1] contiguous_input = np.ascontiguousarray(input) cdef Py_ssize_t old_dim1 = input.shape[0] cdef Py_ssize_t old_dim2 = input.shape[1] diff --git a/imops/testing.py b/imops/testing.py index 2da0fb4d..c9536d82 100644 --- a/imops/testing.py +++ b/imops/testing.py @@ -3,12 +3,14 @@ import numpy as np from skimage.transform import iradon as iradon_, radon as radon_ +from .compat import VisibleDeprecationWarning + def sk_iradon(xs): with warnings.catch_warnings(): warnings.filterwarnings('ignore', module='numpy') warnings.simplefilter('ignore', DeprecationWarning) - warnings.simplefilter('ignore', np.VisibleDeprecationWarning) + warnings.simplefilter('ignore', VisibleDeprecationWarning) return np.stack([iradon_(x) for x in xs]) @@ -16,7 +18,7 @@ def sk_radon(xs): with warnings.catch_warnings(): warnings.filterwarnings('ignore', module='numpy') warnings.simplefilter('ignore', DeprecationWarning) - warnings.simplefilter('ignore', np.VisibleDeprecationWarning) + warnings.simplefilter('ignore', VisibleDeprecationWarning) warnings.filterwarnings('error', '.*image must be zero.*', module='skimage') return np.stack([radon_(x) for x in xs]) diff --git a/imops/utils.py b/imops/utils.py index 92d8d5c9..5c3e633e 100644 --- a/imops/utils.py +++ b/imops/utils.py @@ -6,7 +6,9 @@ import numpy as np -from .backend import BACKEND_NAME2ENV_NUM_THREADS_VAR_NAME, SINGLE_THREADED_BACKENDS, Backend +from .backend import BACKEND_NAME2ENV_NUM_THREADS_VAR_NAME, SINGLE_THREADED_BACKENDS, Backend, Cython +from .compat import normalize_axis_tuple +from .src._utils import _isin as cython_isin AxesLike = Union[int, Sequence[int]] @@ -108,7 +110,7 @@ def axis_from_dim(axis: Union[AxesLike, None], dim: int) -> tuple: if axis is None: return tuple(range(dim)) - return np.core.numeric.normalize_axis_tuple(axis, dim, 'axis') + return normalize_axis_tuple(axis, dim, 'axis') def broadcast_axis(axis: Union[AxesLike, None], dim: int, *values: Union[AxesLike, AxesParams]): @@ -201,3 +203,52 @@ def check_len(*args) -> None: def assert_subdtype(dtype, ref_dtype, name): if not np.issubdtype(dtype, ref_dtype): raise ValueError(f'`{name}` must be of {ref_dtype.__name__} dtype, got {dtype}') + + +def isin(element: np.ndarray, test_elements: np.ndarray, num_threads: int = 1) -> np.ndarray: + """ + Calculates `element in test_elements`, broadcasting over `element` only. + Returns a boolean array of the same shape as `element` that is True where + an element of `element` is in `test_elements` and False otherwise. + + Parameters + ---------- + element: np.ndarray + n-dimensional array + test_elements: np.ndarray + 1-d array of the values against which to test each value of element + num_threads: int + the number of threads to use for computation. Default = 1. If negative value passed + cpu count + num_threads + 1 threads will be used + + Returns + ------- + isin: np.ndarray, bool + has the same shape as `element`. The values `element[isin]` are in `test_elements` + + Examples + -------- + element = 2*np.arange(4).reshape((2, 2)) + test_elements = [1, 2, 4, 8] + mask = isin(element, test_elements) + """ + if element.dtype not in ('int16', 'int32', 'int64'): + raise ValueError(f'Supported dtypes: int16, int32, int64, got {element.dtype}') + + num_threads = normalize_num_threads(num_threads, Cython(), warn_stacklevel=2) + + contiguos_element = np.ascontiguousarray(element) + test_elements = np.asarray(test_elements, dtype=element.dtype) + out = np.zeros_like(contiguos_element, dtype=bool) + + cython_isin(contiguos_element.ravel(), test_elements, out.ravel(), num_threads) + + return out + + +def make_immutable(array: np.ndarray) -> None: + array.flags.writeable = False + + +def make_mutable(array: np.ndarray) -> None: + array.flags.writeable = True diff --git a/pyproject.toml b/pyproject.toml index 7838f311..eed6e582 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,5 +1,5 @@ [build-system] -requires = ['setuptools<69.0.0', 'numpy<2.0.0', 'Cython<3.0.0', 'pybind11'] +requires = ['setuptools<69.0.0', 'numpy<3.0.0', 'Cython>=3.0.0,<4.0.0', 'pybind11'] build-backend = 'setuptools.build_meta' [project] diff --git a/requirements.txt b/requirements.txt index 387a8573..0c574477 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,5 +1,5 @@ numpy>=1.7,<2.0.0; python_version < '3.8' -numpy>=1.22,<2.0.0; python_version >= '3.8' +numpy>=1.22,<3.0.0; python_version >= '3.8' scipy>=1.0,<2.0.0 dataclasses; python_version < '3.7' scikit-image diff --git a/setup.py b/setup.py index 6f086bba..e6e3c50e 100644 --- a/setup.py +++ b/setup.py @@ -50,9 +50,10 @@ install_requires=requirements, extras_require={'numba': ['numba'], 'all': ['numba']}, setup_requires=[ - # Setuptools 18.0 properly handles Cython extensions. - 'setuptools>=18.0', - 'Cython<3.0.0', + 'setuptools<69.0.0', + 'Cython>=3.0.0,<4.0.0', + 'numpy<3.0.0', + 'pybind11', ], ext_modules=ext_modules, python_requires='>=3.7', diff --git a/tests/test_interp1d.py b/tests/test_interp1d.py index e531c194..30438a13 100644 --- a/tests/test_interp1d.py +++ b/tests/test_interp1d.py @@ -8,6 +8,7 @@ from imops._configs import interp1d_configs from imops.backend import Backend, Scipy from imops.interp1d import interp1d +from imops.utils import make_immutable np.random.seed(1337) @@ -175,6 +176,8 @@ def test_stress(backend): for i in range(2 * n_samples): shape = np.random.randint(32, 64, size=np.random.randint(1, 5)) inp = np.random.randn(*shape) + if np.random.binomial(1, 0.5): + make_immutable(inp) axis = np.random.choice(np.arange(inp.ndim)) old_locations = np.random.randn(shape[axis]) diff --git a/tests/test_interp2d.py b/tests/test_interp2d.py index 880d92e3..9179b982 100644 --- a/tests/test_interp2d.py +++ b/tests/test_interp2d.py @@ -7,6 +7,7 @@ from scipy.interpolate import griddata from imops.interp2d import Linear2DInterpolator +from imops.utils import make_immutable tests_root = Path(__file__).parent @@ -33,6 +34,14 @@ def test_test_data(num_threads): int_points = np.transpose((np.isnan(x)).nonzero()) x_points: np.ndarray = np.transpose((~np.isnan(x)).nonzero()) x_values = x[~np.isnan(x)] + + if np.random.binomial(1, 0.5): + make_immutable(x_points) + if np.random.binomial(1, 0.5): + make_immutable(x_values) + if np.random.binomial(1, 0.5): + make_immutable(int_points) + imops_values = Linear2DInterpolator(x_points, x_values, num_threads=num_threads)(int_points, fill_value=0.0) triangles = scipy.spatial.Delaunay(x_points).simplices delaunay_values = Linear2DInterpolator(x_points, x_values, num_threads=num_threads, triangles=triangles)( diff --git a/tests/test_measure.py b/tests/test_measure.py index 5d8c3ce8..7590e2c4 100644 --- a/tests/test_measure.py +++ b/tests/test_measure.py @@ -11,6 +11,7 @@ from imops._configs import measure_configs from imops.backend import Backend, Scipy from imops.measure import center_of_mass, label +from imops.utils import make_immutable np.random.seed(1337) @@ -195,6 +196,7 @@ def test_stress(connectivity, ndim): if ndim == 4 or np.random.binomial(1, 0.2) else np.random.randint(0, 5, size=np.random.randint(32, 64, size=ndim)) ) + sk_labeled, sk_num_components = sk_label(inp, connectivity=connectivity, return_num=True) labeled, num_components = label(inp, connectivity=connectivity, return_num=True) @@ -273,13 +275,15 @@ def test_center_of_mass(num_threads, backend, dtype): if dtype == 'bool' or 'u' in dtype else (32 * np.random.randn(*shape)).astype(dtype) + 4 ) + if np.random.binomial(1, 0.5): + make_immutable(inp) out = center_of_mass(inp, num_threads=num_threads, backend=backend) desired_out = scipy_center_of_mass(inp) assert isinstance(out, tuple) assert isinstance(desired_out, tuple) - allclose(out, desired_out, err_msg=(inp, inp.shape), rtol=1e-4) + allclose(out, desired_out, err_msg=(inp, inp.shape), rtol=1e-3) def test_scipy_warning(num_threads, backend, dtype): @@ -337,6 +341,10 @@ def test_labeled_center_of_mass(backend, dtype, label_dtype): if label_dtype != 'bool' else np.random.choice(np.array([[False], [True], [False, True], [True, False]], dtype=object)) ) + if np.random.binomial(1, 0.5): + make_immutable(inp) + if np.random.binomial(1, 0.5): + make_immutable(labels) out = center_of_mass(inp, labels, index, num_threads=1, backend=backend) desired_out = scipy_center_of_mass(inp, labels, index) @@ -345,4 +353,4 @@ def test_labeled_center_of_mass(backend, dtype, label_dtype): assert isinstance(x, tuple) assert isinstance(y, tuple) - allclose(out, desired_out, err_msg=(inp, inp.shape), rtol=1e-5) + allclose(out, desired_out, err_msg=(inp, inp.shape), rtol=1e-3) diff --git a/tests/test_morphology.py b/tests/test_morphology.py index bd87cd8a..071f9a8b 100644 --- a/tests/test_morphology.py +++ b/tests/test_morphology.py @@ -14,6 +14,7 @@ from imops.backend import Backend, Scipy from imops.morphology import binary_closing, binary_dilation, binary_erosion, binary_opening, distance_transform_edt from imops.pad import restore_crop +from imops.utils import make_immutable np.random.seed(1337) @@ -163,9 +164,15 @@ def take_by_coords(array, coords): else: inp = np.random.binomial(1, 0.5, shape).astype(bool) + if np.random.binomial(1, 0.5): + make_immutable(inp) + footprint_shape = footprint_shape_modifier(np.random.randint(1, 4, size=inp.ndim)) footprint = np.random.binomial(1, 0.5, footprint_shape) if np.random.binomial(1, 0.5) else None + if footprint is not None and np.random.binomial(1, 0.5): + make_immutable(footprint) + if backend == Scipy() and boxed: with pytest.raises(ValueError): imops_op(inp, footprint, backend=backend, boxed=boxed) @@ -236,6 +243,9 @@ def test_stress_edt(backend, num_threads, return_distances, return_indices, samp shape = np.random.randint(64, 128, size=size) x = np.random.randint(5, size=shape) + if np.random.binomial(1, 0.5): + make_immutable(x) + out = distance_transform_edt( x, sampling=sampling, diff --git a/tests/test_numeric.py b/tests/test_numeric.py index d7970b12..34b77a7e 100644 --- a/tests/test_numeric.py +++ b/tests/test_numeric.py @@ -8,6 +8,7 @@ from imops._configs import numeric_configs from imops.backend import Backend from imops.numeric import _STR_TYPES, copy, fill_, full, pointwise_add +from imops.utils import make_immutable np.random.seed(1337) @@ -221,9 +222,10 @@ def test_stress_pointwise_add_inplace(backend, num_threads, dtype): def test_stress_fill_(backend, num_threads, dtype): def sample_value(dtype): - x = dtype.type(32 * np.random.randn(1)) + x = 32 * np.random.randn(1) if isinstance(x, np.ndarray): x = x[0] + x = dtype.type(x) dice = np.random.randint(3) @@ -254,9 +256,10 @@ def sample_value(dtype): def test_stress_full(backend, num_threads, dtype): def sample_value(dtype): - x = dtype.type(32 * np.random.randn(1)) + x = 32 * np.random.randn(1) if isinstance(x, np.ndarray): x = x[0] + x = dtype.type(x) dice = np.random.randint(3) @@ -288,6 +291,9 @@ def test_stress_copy(backend, num_threads, dtype): shape = np.random.randint(32, 64, size=np.random.randint(1, 5)) nums = (32 * np.random.randn(*shape)).astype(dtype) + if np.random.binomial(1, 0.5): + make_immutable(nums) + old_nums = np.copy(nums) copy_nums = copy( nums, diff --git a/tests/test_radon.py b/tests/test_radon.py index b7504a89..9cb6d790 100644 --- a/tests/test_radon.py +++ b/tests/test_radon.py @@ -8,6 +8,7 @@ from imops.backend import Backend from imops.radon import inverse_radon, radon from imops.testing import fill_outside, sample_ct, sk_iradon, sk_radon +from imops.utils import make_immutable np.random.seed(1337) @@ -50,14 +51,21 @@ def test_inverse_radon(backend): image = sample_ct(slices, size) sinogram = sk_radon(image) + ref = sk_iradon(sinogram) + ref_0 = sk_iradon(sinogram[[0]]) + + if np.random.binomial(1, 0.5): + make_immutable(sinogram) + inv = inverse_radon(sinogram, axes=(1, 2), backend=backend) - almost_eq(sk_iradon(sinogram), inv, 3) + + almost_eq(ref, inv, 3) almost_eq(inv[:2], inverse_radon(sinogram[:2], axes=(1, 2), backend=backend)) almost_eq(inv[[0]], inverse_radon(sinogram[[0]], axes=(1, 2), backend=backend)) almost_eq(inv[0], inverse_radon(sinogram[0])) - almost_eq(sk_iradon(sinogram[[0]]), inverse_radon(sinogram[[0]], axes=(1, 2), backend=backend), 3) + almost_eq(ref_0, inverse_radon(sinogram[[0]], axes=(1, 2), backend=backend), 3) almost_eq(inv, np.stack(list(map(partial(inverse_radon, backend=backend), sinogram))), 3) @@ -66,13 +74,19 @@ def test_radon(backend): for size in [40, 47, 64]: image = sample_ct(slices, size) + ref = sk_radon(image) + ref_0 = sk_radon(image[[0]]) + + if np.random.binomial(1, 0.5): + make_immutable(image) + sinogram = radon(image, axes=(1, 2), backend=backend) - almost_eq(sk_radon(image), sinogram) + almost_eq(ref, sinogram) almost_eq(sinogram, radon(fill_outside(image, -1000), axes=(1, 2), backend=backend)) almost_eq(sinogram[:2], radon(image[:2], axes=(1, 2), backend=backend)) almost_eq(sinogram[[0]], radon(image[[0]], axes=(1, 2), backend=backend)) almost_eq(sinogram[0], radon(image[0], backend=backend)) - almost_eq(sk_radon(image[[0]]), radon(image[[0]], axes=(1, 2), backend=backend), 3) + almost_eq(ref_0, radon(image[[0]], axes=(1, 2), backend=backend), 3) almost_eq(sinogram, np.stack(list(map(partial(radon, backend=backend), image))), 3) diff --git a/tests/test_utils.py b/tests/test_utils.py index e31f7b57..45adfe6d 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -11,6 +11,9 @@ build_slices, check_len, imops_num_threads, + isin, + make_immutable, + make_mutable, normalize_num_threads, set_num_threads, ) @@ -21,6 +24,21 @@ MANY_THREADS = 42069 +@pytest.fixture(params=[np.int16, np.int32, np.int64]) +def dtype(request): + return request.param + + +@pytest.fixture(params=[bool, np.float32, np.uint8]) +def bad_dtype(request): + return request.param + + +@pytest.fixture(params=[1, 2, 3]) +def num_threads(request): + return request.param + + def test_check_len(): assert check_len([]) is None assert check_len([0]) is None @@ -110,3 +128,40 @@ def test_broadcast_axis(): with pytest.raises(ValueError): broadcast_axis([0, 1], 2, *arrays) + + +def test_isin(num_threads, dtype): + for _ in range(16): + shape = np.random.randint(32, 64, size=np.random.randint(1, 5)) + elements = np.random.randint(0, 5, size=shape, dtype=dtype) + test_elements = np.random.randint(0, 5, size=np.random.randint(1, 10)) + + if np.random.binomial(1, 0.5): + make_immutable(elements) + if np.random.binomial(1, 0.5): + make_immutable(test_elements) + + assert (np.isin(elements, test_elements) == isin(elements, test_elements, num_threads)).all() + + +def test_bad_dtype(num_threads, bad_dtype): + shape = np.random.randint(32, 64, size=np.random.randint(1, 5)) + elements = np.random.randint(0, 5, size=shape).astype(bad_dtype) + test_elements = np.random.randint(0, 5, size=np.random.randint(1, 10)) + + with pytest.raises(ValueError): + isin(elements, test_elements, num_threads) + + +def test_make_immutable(): + x = np.ones(3) + make_immutable(x) + with pytest.raises(ValueError): + x[0] = 1337 + + +def test_make_mutable(): + x = np.ones(3) + make_immutable(x) + make_mutable(x) + x[0] = 1337 diff --git a/tests/test_zoom.py b/tests/test_zoom.py index 353bcd40..83cb155e 100644 --- a/tests/test_zoom.py +++ b/tests/test_zoom.py @@ -9,7 +9,7 @@ from imops._configs import zoom_configs from imops.backend import Backend -from imops.utils import ZOOM_SRC_DIM, get_c_contiguous_permutaion, inverse_permutation +from imops.utils import ZOOM_SRC_DIM, get_c_contiguous_permutaion, inverse_permutation, make_immutable from imops.zoom import _zoom, zoom, zoom_to_shape @@ -250,6 +250,8 @@ def test_stress(backend, order, dtype): scale = np.random.uniform(0.5, 2, size=inp.ndim if np.random.binomial(1, 0.5) else 1) if len(scale) == 1: scale = scale[0] + if np.random.binomial(1, 0.5): + make_immutable(inp) without_borders = np.index_exp[:-1, :-1, :-1, :-1, :-1][: inp.ndim] desired_out = scipy_zoom(inp, scale, order=order)[without_borders]