python
提供了多种卷积方案,相比之下,定义在ndimage
中的卷积函数,在功能上比numpy
和signal
中的卷积要稍微复杂一些,这点仅从输入参数的多少就可略窥一二
numpy.convolve(a, v, mode='full') scipy.ndimage.convolve1d(input, weights, axis=-1, output=None, mode='reflect', cval=0.0, origin=0) scipy.signal.convolve(in1, in2, mode='full', method='auto') scipy.ndimage.convolve(input, weights, output=None, mode='reflect', cval=0.0, origin=0)
前两者为1维卷积函数,且ndimage可对多维数组沿着单个坐标轴进行卷积操作,后两者为多维卷积。
numpy和signal中的卷积函数,其mode都有三种,用以调节卷积后的边缘特性,如果输入的两个卷积对象的维度分别是N NN和M MM,则这三种模式的输出结果为
full
: 输出维度N + M − 1 N+M-1N+M−1,其最后一点的信号完全不交叠,故而边缘效应明显。same
:输出维度max ( M , N ) \max(M,N)max(M,N),边缘效应仍然可见valid
:输出维度∣ M − N ∣ |M-N|∣M−N∣,只返回完全交叠的区域,相当于把存在边缘效应的点都率除掉了
ndimage
中的convolve
针对边缘效应,对图像进行扩展,而其mode
决定的就是扩展之后的填充格式,设待滤波数组为a b c d
,则在不同的模式下,对边缘进行如下填充
reflect
d c b aa b c dd c b aconstant
k k k ka b c dk k k knearest
a a a aa b c dd d d dmirror
d c ba b c dc b awrap
a b c da b c da b c d其中,k
通过参数cval
设定。
这五种修改边界的方法,在scipy.ndimage
的函数中十分普遍,尤其是涉及到卷积的滤波函数,堪称标配。
接下来针对这些不同的卷积函数,做一下性能测试,用5 × 5的卷积模板,对1000 × 1000的矩阵进行卷积计算,来看一下不同实现方案的卷积,其速度如何
import numpy as np import scipy.signal as ss import scipy.ndimage as sn from timeit import timeit A = np.random.rand(1000,1000) B = np.random.rand(5,5) timeit(lambda : ss.convolve(A, B), number=10) # 0.418 timeit(lambda : sn.convolve(A, B), number=10) # 0.126
相比之下,ndimage
中的卷积显然是更高效的。
接下来测试一下一维卷积的表现
A = np.random.rand(10000) B = np.random.rand(15) timeit(lambda : np.convolve(A, B), number=1000) # 0.15256029999727616 timeit(lambda : ss.convolve(A, B), number=1000) # 0.1231262000001152 timeit(lambda : sn.convolve(A, B), number=1000) # 0.09218210000108229 timeit(lambda : sn.convolve1d(A, B), number=1000) # 0.03915820000111125
相比之下,convolve1d
不愧是写明了1d
的卷积函数,速度最快,而numpy
中提供的函数速度最慢。
卷积操作经常被作用在图像滤波以及边缘提取上,例如,通过类似下面这样的矩阵,可以将图像的纵向的边缘提取出来。
下面做一个简单的测试
from scipy.misc import ascent import matplotlib.pyplot as plt img = ascent() temp = np.zeros([3,3]) temp[:,0] = -1 temp[:,2] = 1 edge = sn.convolve(img, temp) fig = plt.figure() ax = fig.add_subplot(121) ax.imshow(img) ax = fig.add_subplot(122) ax.imshow(edge) plt.show()
【本文由:香港云服务器 http://www.558idc.com/ne.html 复制请保留原URL】