numpy数组中子矩阵的矢量化提取
Vectorized extraction of submatrices in numpy-array
我的目标是实现一个中值滤波器,这是一个函数,它用周围像素的中值替换(大部分)二维数组中的每个像素。可用于图像去噪。
我的实现从原始矩阵中提取子矩阵,其中包含像素本身及其相邻像素。这种提取目前是通过 for 循环完成的,正如您可以想象的那样,for 循环占用了大约 95% 的执行时间。
这是我当前的实现:
def median_blur(img, fsize=3):
img_intermediate = np.zeros((img.shape[0] + 2, img.shape[1] + 2), np.uint8) # First intermediate img, used for padding original image
img_intermediate[1:img.shape[0]+1, 1:img.shape[1]+1] = img
img_result = np.empty((*img.shape, fsize, fsize), np.uint8) # Will contain result, first receives kernel-submatrices
# Extract submatrices from image
for i in range(img.shape[0]):
for j in range(img.shape[1]):
img_result[i, j] = img_intermediate[i:i+fsize, j:j+fsize]
img_result = img_result.reshape(*img.shape, fsize**2) # Reshape Result-Matrix
img_result = np.median(img_result, axis=2) # Calculate median in one go
return img_result.astype('uint8')
如何使用矢量化运算提取这些子矩阵?
作为奖励,如果有计算机视觉经验的人正在阅读这篇文章:除了将中值滤波器应用于中间零填充矩阵之外,还有更好的方法来实现中值滤波器吗?
非常感谢。
这是一个向量化的解决方案。但是,您可以通过注意图像数组的内存顺序来提出更快的解决方案:
from numpy.lib.stride_tricks import as_strided
img_padded = np.pad(img, 1, mode='constant')
sub_shape = (fsize, fsize)
view_shape = tuple(np.subtract(img_padded.shape, sub_shape) + 1) + sub_shape
sub_img = as_strided(img_padded, view_shape, img_padded.strides * 2)
sub_img = sub_img.reshape(img.shape + (fsize**2,))
result = np.median(sub_img, axis=2).astype(int)
先使用 pad
函数进行填充,然后使用 as_strided
获得 sub_matrices,最后将 median
应用于您的步幅。
更新:使用@Divakar 在评论中建议的view_as_windows
:
from skimage.util.shape import view_as_windows
img_padded = np.pad(img, 1, mode='constant')
sub_shape = (fsize, fsize)
view_shape = view_as_windows(img_padded, sub_shape, 1)
sub_img = view_shape.reshape(img.shape + (fsize**2,))
result = np.median(sub_img, axis=2).astype(int)
view_as_windows
也提供了类似于strides的子矩阵。
示例图像和输出:
img:
[[ 1 2 3 4 5 6]
[ 7 8 9 10 11 12]
[13 14 15 16 17 18]
[19 20 21 22 23 24]]
median_filtered:
[[ 0 2 3 4 5 0]
[ 2 8 9 10 11 6]
[ 8 14 15 16 17 12]
[ 0 14 15 16 17 0]]
我的目标是实现一个中值滤波器,这是一个函数,它用周围像素的中值替换(大部分)二维数组中的每个像素。可用于图像去噪。
我的实现从原始矩阵中提取子矩阵,其中包含像素本身及其相邻像素。这种提取目前是通过 for 循环完成的,正如您可以想象的那样,for 循环占用了大约 95% 的执行时间。
这是我当前的实现:
def median_blur(img, fsize=3):
img_intermediate = np.zeros((img.shape[0] + 2, img.shape[1] + 2), np.uint8) # First intermediate img, used for padding original image
img_intermediate[1:img.shape[0]+1, 1:img.shape[1]+1] = img
img_result = np.empty((*img.shape, fsize, fsize), np.uint8) # Will contain result, first receives kernel-submatrices
# Extract submatrices from image
for i in range(img.shape[0]):
for j in range(img.shape[1]):
img_result[i, j] = img_intermediate[i:i+fsize, j:j+fsize]
img_result = img_result.reshape(*img.shape, fsize**2) # Reshape Result-Matrix
img_result = np.median(img_result, axis=2) # Calculate median in one go
return img_result.astype('uint8')
如何使用矢量化运算提取这些子矩阵?
作为奖励,如果有计算机视觉经验的人正在阅读这篇文章:除了将中值滤波器应用于中间零填充矩阵之外,还有更好的方法来实现中值滤波器吗?
非常感谢。
这是一个向量化的解决方案。但是,您可以通过注意图像数组的内存顺序来提出更快的解决方案:
from numpy.lib.stride_tricks import as_strided
img_padded = np.pad(img, 1, mode='constant')
sub_shape = (fsize, fsize)
view_shape = tuple(np.subtract(img_padded.shape, sub_shape) + 1) + sub_shape
sub_img = as_strided(img_padded, view_shape, img_padded.strides * 2)
sub_img = sub_img.reshape(img.shape + (fsize**2,))
result = np.median(sub_img, axis=2).astype(int)
先使用 pad
函数进行填充,然后使用 as_strided
获得 sub_matrices,最后将 median
应用于您的步幅。
更新:使用@Divakar 在评论中建议的view_as_windows
:
from skimage.util.shape import view_as_windows
img_padded = np.pad(img, 1, mode='constant')
sub_shape = (fsize, fsize)
view_shape = view_as_windows(img_padded, sub_shape, 1)
sub_img = view_shape.reshape(img.shape + (fsize**2,))
result = np.median(sub_img, axis=2).astype(int)
view_as_windows
也提供了类似于strides的子矩阵。
示例图像和输出:
img:
[[ 1 2 3 4 5 6]
[ 7 8 9 10 11 12]
[13 14 15 16 17 18]
[19 20 21 22 23 24]]
median_filtered:
[[ 0 2 3 4 5 0]
[ 2 8 9 10 11 6]
[ 8 14 15 16 17 12]
[ 0 14 15 16 17 0]]