需求的提出：
在利用numpy进行图像频域滤波，其中一步需要获取fshift矩阵中心点到矩阵其余各点的直线距离。书本上基本用双层for循环实现，这样对于一个解释性语言来说简直就是噩梦！这让我非常不爽！！！
寻求解决：
我上网搜寻方案，浏览器翻了5页也没找到，只能放弃寻找，自己想办法。
分析问题：
问题：
求一个二维矩阵中某一点到其余各点距离的同形矩阵。
分析：
在二维矩阵中每个点都可以用（x ，y）表示，一共两个数，能不能分开计算呢？好像可以！
二维空间的欧式距离的计算公式就是把x和y分开计算的。我们可以再把矩阵的 行 和 列 拆成 两个 一维向量序列 把中心点拆成 x y，利 用numpy的广播机制，实现快速计算，如下演示：
(x[:] - x0)^2 (y[:] - y0)^2
再利用numpy 的维度扩增函数，对两条向量进行维度扩充，再利用numpy 广播加法 就能得到一个完整的二维矩阵，进而对矩阵进行开方，得到某一点对其余点的完整距离矩阵。这样就能大幅提高计算速度，避免使用双层for循环。
验证方法：
根据上述思路，初步实现设想代码，如下：
import numpy as np
h,w = 6, 6 # src mat shape
tmp = [np.arange(h), np.arange(w)]
center = -0.5 + np.array([h, w])/2
dis_mat = (np.expand_dims((tmp[0] - center[0])**2, axis=1)
+ np.expand_dims((tmp[1] - center[1])**2,axis=0))**0.5
print("tmp:\n",tmp, end='\n\n')
print("center:\n", center, end='\n\n')
print("dis_mat:\n", dis_mat)
实验结果：
Example.png
由上图可以看出，我们得到了以原矩阵正中心到其他点距离的同形矩阵。
解决问题：
例：
数字图像处理的频域滤波中常用的ButterWorth 高频滤波实验。
环境：
|AMD-R7-4800u | Win10-21H2 | WSL2 | Ubuntu-20.04LTS | 中使用Vscode for jupyternotebook 运行代码。必要库：numpy、skimage、pyplot、time 。
验证方式：
通过在函数执行部分嵌入记录代码，记录操作完成需要的时间并打印。
原书ButterworthHFilter函数代码：
def butterworthPassFilter(img, d, n):
import numpy as np
from math import sqrt
start = time.time()
f = np.fft.fft2(img)
fshift = np.fft.fftshift(f)
d_mat = np.zeros(img.shape)
center = tuple(map(lambda x:(x-1)/2, fshift.shape))
for i in range(d_mat.shape[0]):
for j in range(d_mat.shape[1]):
dis = sqrt((center[0]-i)**2 + (center[1]-j)**2)
d_mat[i, j] = 1/(1+(d/(dis+1e-6))**(2*n))
new_img = np.abs(np.fft.ifft2(np.fft.ifftshift(fshift*d_mat)))
end = time.time()
print(end - start)
return new_img
本人改良ButterworthHFilter函数代码：
def butterworthPassFilter(img, d, n):
import numpy as np
start = time.time()
f = np.fft.fft2(img)
fshift = np.fft.fftshift(f)
tmp = [np.arange(fshift.shape[0]),np.arange(fshift.shape[1])]
center = -0.5 + np.array([fshift.shape[0], fshift.shape[1]])/2.0
dis = (np.expand_dims((tmp[0] - center[0])**2, axis=1)
+ np.expand_dims((tmp[1] - center[1])**2, axis=0))**0.5
d_mat = 1/(1+(d/(dis+1e-6))**(2*n)) # Butterworth 高通滤波
new_img = np.abs(np.fft.ifft2(np.fft.ifftshift(fshift*d_mat)))
end = time.time()
print(end - start)
return new_img
调用代码：
from matplotlib import pyplot as plt
from skimage import io, color
import time
img = io.imread('../medium/images/9.png')
img = color.rgb2gray(img)
plt.subplot(221)
plt.axis('off')
plt.title('Src')
plt.imshow(img, cmap='gray')
plt.subplot(222)
plt.axis('off')
plt.title('Butter 10 1')
butter = butterworthPassFilter(img, 10, 1)
plt.imshow(butter, cmap='gray')
plt.subplot(223)
plt.axis('off')
plt.title('Butter 30 1')
butter = butterworthPassFilter(img, 30, 1)
plt.imshow(butter, cmap='gray')
plt.subplot(224)
plt.axis('off')
plt.title('Butter 10 5')
butter = butterworthPassFilter(img,10, 5)
plt.imshow(butter, cmap='gray')
plt.show()
实验结果
原书双层for循环：
losers.png
本人改良版：
superfast.png
结果分析：
在与原书代码运行时间相比，运行速度平均提升可超5倍多！！！
由此可见利用numpy进行矩阵运算优化的效果非常好！
注意：
文章为本人在简书平台原创，盗用必究！
