1. 使用Cython扩展numpy

下面的代码将图像与滤镜进行二维离散卷积.它是有效的Python和有效的Cython代码.我将其称为Python版本的convolve_py.py和Cython版本的convolve1.pyx

%%writefile convolve_py.py import numpy as np def naive_convolve(f, g):     # f is an image and is indexed by (v, w)     # g is a filter kernel and is indexed by (s, t),     #   it needs odd dimensions     # h is the output image and is indexed by (x, y),     #   it is not cropped     if g.shape[0] % 2 != 1 or g.shape[1] % 2 != 1:         raise ValueError("Only odd dimensions on filter supported")     # smid and tmid are number of pixels between the center pixel     # and the edge, ie for a 5x5 filter they will be 2.     #     # The output size is calculated by adding smid, tmid to each     # side of the dimensions of the input image.     vmax = f.shape[0]     wmax = f.shape[1]     smax = g.shape[0]     tmax = g.shape[1]     smid = smax // 2     tmid = tmax // 2     xmax = vmax + 2*smid     ymax = wmax + 2*tmid     # Allocate result image.     h = np.zeros([xmax, ymax], dtype=f.dtype)     # Do convolution     for x in range(xmax):         for y in range(ymax):             # Calculate pixel value for h at (x,y). Sum one component             # for each pixel (s, t) of the filter g.             s_from = max(smid - x, -smid)             s_to = min((xmax - x) - smid, smid + 1)             t_from = max(tmid - y, -tmid)             t_to = min((ymax - y) - tmid, tmid + 1)             value = 0             for s in range(s_from, s_to):                 for t in range(t_from, t_to):                     v = x - smid + s                     w = y - tmid + t                     value += g[smid - s, tmid - t] * f[v, w]             h[x, y] = value     return h 

Overwriting convolve_py.py 

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 

from convolve_py import naive_convolve 

naive_convolve(f,g) 

array([[      0,       0,       1, ...,    2056,    1477,     792],        [      0,     109,     329, ...,    8858,    6227,    3275],        [    900,    2127,    3684, ...,   23106,   16050,    8349],        ...,        [1850400, 3730389, 5639970, ..., 6230334, 4183464, 2106687],        [1329300, 2678435, 4047407, ..., 4445402, 2983649, 1501849],        [ 712800, 1435572, 2168317, ..., 2369524, 1589761,  799920]]) 

%timeit -n10 naive_convolve(f,g) 

436 ms ± 12.4 ms per loop (mean ± std. dev. of 7 runs, 10 loops each) 

使用cython编译带numpy的代码需要在setup.py中指定include_dirs=[numpy.get_include()]

1.1. 第一版迭代--使用cython编译以提高性能

代码不用改,直接编译以提高性能

%%writefile convolve1.pyx  import numpy as np def naive_convolve_1(f, g):     # f is an image and is indexed by (v, w)     # g is a filter kernel and is indexed by (s, t),     #   it needs odd dimensions     # h is the output image and is indexed by (x, y),     #   it is not cropped     if g.shape[0] % 2 != 1 or g.shape[1] % 2 != 1:         raise ValueError("Only odd dimensions on filter supported")     # smid and tmid are number of pixels between the center pixel     # and the edge, ie for a 5x5 filter they will be 2.     #     # The output size is calculated by adding smid, tmid to each     # side of the dimensions of the input image.     vmax = f.shape[0]     wmax = f.shape[1]     smax = g.shape[0]     tmax = g.shape[1]     smid = smax // 2     tmid = tmax // 2     xmax = vmax + 2*smid     ymax = wmax + 2*tmid     # Allocate result image.     h = np.zeros([xmax, ymax], dtype=f.dtype)     # Do convolution     for x in range(xmax):         for y in range(ymax):             # Calculate pixel value for h at (x,y). Sum one component             # for each pixel (s, t) of the filter g.             s_from = max(smid - x, -smid)             s_to = min((xmax - x) - smid, smid + 1)             t_from = max(tmid - y, -tmid)             t_to = min((ymax - y) - tmid, tmid + 1)             value = 0             for s in range(s_from, s_to):                 for t in range(t_from, t_to):                     v = x - smid + s                     w = y - tmid + t                     value += g[smid - s, tmid - t] * f[v, w]             h[x, y] = value     return h 

Overwriting convolve1.pyx 

Overwriting convolve_py.py 0

Overwriting convolve_py.py 1

Overwriting convolve_py.py 2

Overwriting convolve_py.py 3

Overwriting convolve_py.py 4

Overwriting convolve_py.py 5

array([[      0,       0,       1, ...,    2056,    1477,     792],        [      0,     109,     329, ...,    8858,    6227,    3275],        [    900,    2127,    3684, ...,   23106,   16050,    8349],        ...,        [1850400, 3730389, 5639970, ..., 6230334, 4183464, 2106687],        [1329300, 2678435, 4047407, ..., 4445402, 2983649, 1501849],        [ 712800, 1435572, 2168317, ..., 2369524, 1589761,  799920]]) 

Overwriting convolve_py.py 7

Overwriting convolve_py.py 8

第一版什么也不改就可以提高1/4的性能

1.2. 第二版迭代--静态化参数

将函数的参数以及中间变量都申明为静态类型以提高运行效率

Overwriting convolve_py.py 9

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 0

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 1

Overwriting convolve_py.py 1

Overwriting convolve_py.py 2

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 4

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 5

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 6

array([[      0,       0,       1, ...,    2056,    1477,     792],        [      0,     109,     329, ...,    8858,    6227,    3275],        [    900,    2127,    3684, ...,   23106,   16050,    8349],        ...,        [1850400, 3730389, 5639970, ..., 6230334, 4183464, 2106687],        [1329300, 2678435, 4047407, ..., 4445402, 2983649, 1501849],        [ 712800, 1435572, 2168317, ..., 2369524, 1589761,  799920]]) 

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 8

import numpy as np N = 100 f = np.arange(N*N, dtype=np.int).reshape((N,N)) g = np.arange(81, dtype=np.int).reshape((9, 9)) 9

1.3. 第三版迭代--“缓冲”语法

提高np数组的效率,我们用一个特殊的“缓冲”语法来做到这一点，它必须告诉数据类型（第一个参数）和维数（“ndim”仅关键字参数，如果不提供，则假定一维

from convolve_py import naive_convolve 0

from convolve_py import naive_convolve 1

from convolve_py import naive_convolve 2

Overwriting convolve_py.py 1

Overwriting convolve_py.py 2

from convolve_py import naive_convolve 5

from convolve_py import naive_convolve 6

from convolve_py import naive_convolve 7

array([[      0,       0,       1, ...,    2056,    1477,     792],        [      0,     109,     329, ...,    8858,    6227,    3275],        [    900,    2127,    3684, ...,   23106,   16050,    8349],        ...,        [1850400, 3730389, 5639970, ..., 6230334, 4183464, 2106687],        [1329300, 2678435, 4047407, ..., 4445402, 2983649, 1501849],        [ 712800, 1435572, 2168317, ..., 2369524, 1589761,  799920]]) 

from convolve_py import naive_convolve 9

naive_convolve(f,g) 0

提高了150倍的性能