我需要一种快速的方法来遍历Python中图像/堆栈的像素

Question

我创建了一个3D中值滤镜，它可以正常工作，并且如下所示：

def Median_Filter_3D(image,kernel):

window = np.zeros(shape=(kernel,kernel,kernel), dtype = np.uint8)
n = (kernel-1)/2    #Deals with Image border
imgout = np.empty_like(image)
w,h,l = image.shape()

每个像素的%%启动循环

for y in np.arange(0,(w-n*2),1):
    for x in np.arange(0,(h-n*2),1):
        for z in np.arange(0,(l-n*2),1):
            window[:,:,:] = image[x:x+kernel,y:y+kernel,z:z+kernel]
            med = np.median(window)
            imgout[x+n,y+n,z+n] = med 
return(imgout)

因此，在每个像素处，它都会创建一个大小为kernelxkernelxkernel的窗口，找到该窗口中像素的中值，然后用新的中值替换该像素的值。

我的问题是，它非常缓慢，我需要处理数千个大图像。必须有一种更快的方法来遍历所有这些像素，并且仍然能够获得相同的结果。

提前致谢！！

Answer 1

首先，在python中循环3D矩阵是一个非常非常非常糟糕的主意。为了循环播放大型3D矩阵，最好使用Cython或C / C ++ / Fortran并创建python扩展。但是，对于这种特殊情况，scipy已经包含n维数组的中值滤波器的实现：

>>> from scipy.ndimage import median_filter
>>> median_filter(my_large_3d_array, radious)

简而言之，没有更快的方法来遍历python中的体素（也许numpy迭代器会有所帮助，但不会显着提高性能）。如果您需要在python中执行更复杂的3D内容，则应考虑在Cython的循环接口中进行编程，或者使用诸如Dask之类的分块库进行编程，该库对数组的块实现并行操作。

Python的问题是，如果for循环极慢，特别是如果它们是嵌套的并且有大数组的话。因此，没有用于在数组上获得有效迭代的标准pythonic方法。通常，提高速度的方法是通过矢量化操作和numpy-ticks，但是这些都是非常特定于问题的，没有通用技巧，您将在SO中学到很多numpy技巧。

作为一种通用方法，如果确实需要遍历数组，则可以使用Cython编写代码。Cython是Python的类似C的扩展。您可以使用Python语法编写代码，但是使用int或来指定变量类型（例如在C中），float然后该代码会自动编译为C，并且可以从python调用。以下是一个简单的示例：

示例Python循环函数：

import numpy as np

def iter_A(A):
    B = np.empty(A.shape, dtype=np.float64)

    for i in range(A.shape[0]):
        for j in range(A.shape[1]):
            B[i, j] = A[i, j] * 2
    return B

我知道上面的代码有点多余，可以写成B = A * 2，但是它的目的只是为了说明python循环非常慢。

Cython版本的功能：

import numpy as np
cimport numpy as np

def iter_A_cy(double[:, ::1] A):
    cdef Py_ssize_t H = A.shape[0], W = A.shape[1]
    cdef double[:, ::1] B = np.empty((H, W), dtype=np.float64)
    cdef Py_ssize_t i, j

    for i in range(H):
        for j in range(W):
            B[i, j] = A[i, j] * 2

    return np.asarray(B)

测试两种实现的速度：

>>> import numpy as np
>>> A = np.random.randn(1000, 1000)
>>> %timeit iter_A(A)
1 loop, best of 3: 399 ms per loop
>>> %timeit iter_A_cy(A)
100 loops, best of 3: 2.11 ms per loop

注意：您不能照原样运行Cython函数。您需要将其放在单独的文件中，然后首先进行编译（或%%cython在IPython Notebook中使用magic）。

它显示原始python版本需要400ms迭代整个数组，而仅 2ms适用于Cython版本（x200speedup）。