scipy / scipy

SciPy library main repository
https://scipy.org
BSD 3-Clause "New" or "Revised" License
12.71k stars 5.08k forks source link

BUG: scipy.ndimage.zoom produce wrong output with zoom factor of 1 #20999

Open arsalanfiroozi opened 3 weeks ago

arsalanfiroozi commented 3 weeks ago

Describe your issue.

Hey, I am using scipy.ndimage.zoom for 3d images. Based on the input parameters, I might need to zoom or not. So, I expect the function to act as an identity function if the zoom factor is one. You can see the code and the output. mask

Reproducing Code Example

import numpy as np
from scipy.ndimage import zoom
import matplotlib.pyplot as plt

size = 100
mask = np.zeros((100, 100, 100))
center = np.multiply(size,(0.5,0.5,0.5))
print(center)
radius = size/8
for x in range(size):
    for y in range(size):
        for z in range(size):
            if (x - center[0]) ** 2 + (y - center[1]) ** 2 + (z - center[2]) ** 2 < radius ** 2:
                mask[x, y, z] = 1

zoomed_mask = zoom(mask, (1,1,1))

# Plot the original and zoomed masks in 3d
fig = plt.figure()
ax = fig.add_subplot(121, projection='3d')
ax.voxels(mask, edgecolor='k')
ax.set_title('Original Mask')
ax = fig.add_subplot(122, projection='3d')
ax.voxels(zoomed_mask, edgecolor='k')
ax.set_title('Zoomed Mask')
plt.savefig('mask.png')

Error message

There is no error message. See the output image.

SciPy/NumPy/Python version and system information

Build Dependencies:
  blas:
    detection method: pkgconfig
    found: true
    include directory: /usr/local/include
    lib directory: /usr/local/lib
    name: openblas
    openblas configuration: USE_64BITINT=0 DYNAMIC_ARCH=1 DYNAMIC_OLDER= NO_CBLAS=
      NO_LAPACK= NO_LAPACKE= NO_AFFINITY=1 USE_OPENMP= ZEN MAX_THREADS=64
    pc file directory: /usr/local/lib/pkgconfig
    version: 0.3.27
  lapack:
    detection method: pkgconfig
    found: true
    include directory: /usr/local/include
    lib directory: /usr/local/lib
    name: openblas
    openblas configuration: USE_64BITINT=0 DYNAMIC_ARCH=1 DYNAMIC_OLDER= NO_CBLAS=
      NO_LAPACK= NO_LAPACKE= NO_AFFINITY=1 USE_OPENMP= ZEN MAX_THREADS=64
    pc file directory: /usr/local/lib/pkgconfig
    version: 0.3.27
  pybind11:
    detection method: config-tool
    include directory: unknown
    name: pybind11
    version: 2.12.0
Compilers:
  c:
    commands: cc
    linker: ld.bfd
    name: gcc
    version: 10.2.1
  c++:
    commands: c++
    linker: ld.bfd
    name: gcc
    version: 10.2.1
  cython:
    commands: cython
    linker: cython
    name: cython
    version: 3.0.10
  fortran:
    commands: gfortran
    linker: ld.bfd
    name: gcc
    version: 10.2.1
  pythran:
    include directory: ../../tmp/pip-build-env-mnl4e8vy/overlay/lib/python3.10/site-packages/pythran
    version: 0.15.0
Machine Information:
  build:
    cpu: x86_64
    endian: little
    family: x86_64
    system: linux
  cross-compiled: false
  host:
    cpu: x86_64
    endian: little
    family: x86_64
    system: linux
Python Information:
  path: /opt/python/cp310-cp310/bin/python
  version: '3.10'
nickodell commented 2 weeks ago

More information

Hello, could you include the version info for SciPy/NumPy/Python?

You can get them with this:

import sys, scipy, numpy; print(scipy.__version__, numpy.__version__, sys.version_info);

Workaround

I looked into why this plot looks like this. When I plot a single 2D layer of this 3D array, I get the following result for mask:

plt.imshow(mask[50])
plt.colorbar()
plt.show()

Untitled

And the following result for zoomed mask.

plt.imshow(zoomed_mask[50])
plt.colorbar()
plt.show()

Untitled

Which look visually identical.

Looking at the arrays, mask has zeros where the sphere is.

>>> print(mask[50])
[[0. 0. 0. ... 0. 0. 0.]
 [0. 0. 0. ... 0. 0. 0.]
 [0. 0. 0. ... 0. 0. 0.]
 ...
 [0. 0. 0. ... 0. 0. 0.]
 [0. 0. 0. ... 0. 0. 0.]
 [0. 0. 0. ... 0. 0. 0.]]

But zoomed_mask has values very close to zero.

>>> print(zoomed_mask[50])
[[-2.73008327e-64  4.48669571e-63  2.98858092e-63 ...  7.35583558e-62
   4.81237224e-63  5.64220157e-63]
 [ 3.47812810e-64  9.82065113e-63  1.14177549e-62 ... -5.28495731e-61
  -5.63448647e-62  1.48921256e-62]
 [-9.62627161e-63  1.21406027e-62  6.31602054e-62 ... -1.36521778e-60
  -2.06049781e-61  1.28043091e-61]
 ...
 [-6.97146241e-62 -1.37985552e-61  5.45084572e-61 ... -5.55510792e-60
  -5.19043953e-61 -6.48359415e-62]
 [-2.95606196e-62 -7.19974141e-62 -2.14201518e-61 ... -4.31188291e-61
   3.37309667e-61  2.80366236e-61]
 [-4.08125980e-63 -1.64422845e-62 -1.88761873e-62 ...  2.78968899e-61
   6.02296151e-63 -1.62186346e-62]]

However, according to the documentation for matplotlib, voxels interprets non-zero values as filled values:

A 3D array of values, with truthy values indicating which voxels to fill

So this gives us a work-around, which is to avoid filling in near-zero values.

# Plot the original and zoomed masks in 3d
fig = plt.figure()
ax = fig.add_subplot(121, projection='3d')
ax.voxels(mask, edgecolor='k')
ax.set_title('Original Mask')
ax = fig.add_subplot(122, projection='3d')
ax.voxels(~np.isclose(zoomed_mask, 0), edgecolor='k')
ax.set_title('Zoomed Mask')
plt.show()

Untitled