是否可以对访问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(以及其他)提供过滤图像的工具:

In [51]: import scipy.ndimage

从1维数组创建样本图像.重塑创建视图而不是复制:

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().前者创建原始二维数组的一维视图,后者创建一个展平的副本:

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])

Answer 2

假设arr代表输入数组，并且thresh是与每个输入元素进行比较的阈值。现在，我们可以根据给定阈值对输入数组进行阈值设置，从而得到一个掩码/布尔数组。然后，我们可以执行 2D 卷积并1s从中减去True阈值数组的值。

因此，实现看起来像这样 -

from scipy.signal import convolve2d

# Get thresholded mask as int array & set first, last cols and rows as 0s
mask = (arr > thresh).astype(int)
mask[[0,-1]] = 0
mask[:,[0,-1]] = 0

# Perform 2D convolution and subtract 1s corresponding to True elems in mask
out = convolve2d(mask,np.ones((3,3),dtype=int),'same') - mask