python - 是否可以对访问 numpy 数组中不同元素的函数进行矢量化？

Question

假设我希望在 python 中实现以下代码

此函数将图像作为一维数组，并迭代数组中的各个元素(输入图像中的像素)，这会影响输出数组，该输出数组也是表示为一维数组的图像

例子:输入图像(红色)中的单个像素会影响(橙色)中的 8 个周围像素

C 中的基本实现是

/* C version 
 * Given an input image create an 
 * output image that is shaped by individual pixels
 * of the input image
 */

int image[width * height]; //image retrieved elsewhere
int output [width * height]; //output image
int y = 0, x = 0;
for( y = 1; y < height-1 ; ++ y) {
    for(x = 1; x < width-1; ++ x) {
        if (image[y * width + x] > condition) {
            /* pixel affects the surrounding 8 pixels in the output image */

            output[(y-1) * width + x - 1]++; /* upper left  */
            output[(y-1) * width + x    ]++; /* above       */
            output[(y-1) * width + x + 1]++; /* upper right */
            output[y * width + x + 1    ]++; /* right       */
            output[y * width + x - 1    ]++; /* left        */
            output[(y+1) * width + x - 1]++; /* lower left  */
            output[(y+1) * width + x    ]++; /* below       */
            output[(y+1) * width + x + 1]++; /* lower right */


        }
    }
}

python 中的天真方法是使用完全相同的元素明智访问，如下所示

#Python version
input  = blah # formed elsewhere  
output = np.zeros(width * height)
for y in xrange(1, height-1):
    for x in xrange(1, width-1):
        if input[y * width + x] > condition:
            output[(y-1) * width + x - 1]+= 1; # upper left  
            output[(y-1) * width + x    ]+= 1; # above       
            output[(y-1) * width + x + 1]+= 1; # upper right 
            output[y * width + x + 1    ]+= 1; # right       
            output[y * width + x - 1    ]+= 1; # left        
            output[(y+1) * width + x - 1]+= 1; # lower left  
            output[(y+1) * width + x    ]+= 1; # below       
            output[(y+1) * width + x + 1]+= 1; # lower right 

有没有更好的方法来实现这个？是否可以矢量化此函数？

Answer 1

如果我正确理解了这个问题，那么这个方法可以颠倒过来:如果一个像素在其邻域中有与条件匹配的像素，则每次匹配时将其递增 1。对所有像素执行此操作。 Scipy(以及其他)为 filtering images 提供工具:

In [51]: import scipy.ndimage

从一维数组创建示例图像。 Reshape 创建 View 而不是复制:

In [62]: I1d
Out[62]: 
array([  0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,
         0,   0,   0,   0,   0, 129,   0, 129, 129,   0,   0,   0,   0,
         0,   0,   0,   0,   0,   0,   0,   0,   0,   0, 129,   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, 129])

In [63]: height
Out[63]: 8

In [64]: width
Out[64]: 8

In [65]: I = I1d.reshape((height, width))

In [66]: I
Out[66]: 
array([[  0,   0,   0,   0,   0,   0,   0,   0],
       [  0,   0,   0,   0,   0,   0,   0,   0],
       [  0,   0, 129,   0, 129, 129,   0,   0],
       [  0,   0,   0,   0,   0,   0,   0,   0],
       [  0,   0,   0,   0, 129,   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, 129]])

使用卷积创建一个图像，该图像从超过条件(此处为 128)的像素的二进制掩码中保存原始像素中每个像素的增量:

In [67]: scipy.ndimage.filters.convolve(
    (I > 128).astype(np.int),  # conditioned binary image
    weights=np.array([[1, 1, 1],  # each match weighted as 1
                      [1, 0, 1],
                      [1, 1, 1]]),
    mode='constant', cval=0)  # Use zeros as constant fill values on edges
Out[67]: 
array([[0, 0, 0, 0, 0, 0, 0, 0],
       [0, 1, 1, 2, 2, 2, 1, 0],
       [0, 1, 0, 2, 1, 1, 1, 0],
       [0, 1, 1, 3, 3, 3, 1, 0],
       [0, 0, 0, 1, 0, 1, 0, 0],
       [0, 0, 0, 1, 1, 1, 0, 0],
       [0, 0, 0, 0, 0, 0, 1, 1],
       [0, 0, 0, 0, 0, 0, 1, 0]])

In [68]: conv = _

如果最终目标是添加原始和增量:

In [69]: I + conv
Out[69]: 
array([[  0,   0,   0,   0,   0,   0,   0,   0],
       [  0,   1,   1,   2,   2,   2,   1,   0],
       [  0,   1, 129,   2, 130, 130,   1,   0],
       [  0,   1,   1,   3,   3,   3,   1,   0],
       [  0,   0,   0,   1, 129,   1,   0,   0],
       [  0,   0,   0,   1,   1,   1,   0,   0],
       [  0,   0,   0,   0,   0,   0,   1,   1],
       [  0,   0,   0,   0,   0,   0,   1, 129]])

要同时输出一维数组，请使用 ravel()或 flatten() .前者创建一个原始二维数组的一维 View ，后者创建一个扁平化的副本:

In [70]: conv.ravel()
Out[70]: 
array([0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 2, 2, 2, 1, 0, 0, 1, 0, 2, 1, 1, 1,
       0, 0, 1, 1, 3, 3, 3, 1, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 1,
       0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 1, 0])