Python 调用 C/C++实现卷积
scipy与numpy中很多函数用c实现。
本篇以最简单的convolve为例做一个相应的动态库。
1.[转载] Python 调用 C/C++（基础篇）
下面转载自Jerry Jho在知乎“如何实现 C/C++ 与 Python 的通信？”问题下的回复。
示例简单，说明清晰，在我电脑上测试无误。
这种做法称为Python扩展。比如说，我们有一个功能强大的C函数：期望在Python里这样使用：int great_function(int a) { return a + 1; }考虑最简单的情况。我们把功能强大的函数放入C文件 great_module.c 中。>>> from great_module import great_function >>> great_function(2) 3除了功能强大的函数great_function外，这个文件中还有以下部分：#include <Python.h> int great_function(int a) { return a + 1; } static PyObject * _great_function(PyObject *self, PyObject *args) { int _a; int res; if (!PyArg_ParseTuple(args, "i", &_a)) return NULL; res = great_function(_a); return PyLong_FromLong(res); } static PyMethodDef GreateModuleMethods[] = { { "great_function", _great_function, METH_VARARGS, "" }, {NULL, NULL, 0, NULL} }; PyMODINIT_FUNC initgreat_module(void) { (void) Py_InitModule("great_module", GreateModuleMethods); }
- 包裹函数_great_function。它负责将Python的参数转化为C的参数（PyArg_ParseTuple），调用实际的great_function，并处理great_function的返回值，最终返回给Python环境。
- 导出表GreateModuleMethods。它负责告诉Python这个模块里有哪些函数可以被Python调用。导出表的名字可以随便起，每一项有4个参数：第一个参数是提供给Python环境的函数名称，第二个参数是_great_function，即包裹函数。第三个参数的含义是参数变长，第四个参数是一个说明性的字符串。导出表总是以{NULL, NULL, 0, NULL}结束。
- 导出函数initgreat_module。这个的名字不是任取的，是你的module名称添加前缀init。导出函数中将模块名称与导出表进行连接。
在Windows下面，在Visual Studio命令提示符下编译这个文件的命令是
cl /LD great_module.c /o great_module.pyd -IC:\Python27\include C:\Python27\libs\python27.lib/LD 即生成动态链接库。编译成功后在当前目录可以得到 great_module.pyd（实际上是dll）。这个pyd可以在Python环境下直接当作module使用。
在Linux下面，则用gcc编译：
gcc -fPIC -shared great_module.c -o great_module.so -I/usr/include/python2.7/ -lpython2.7在当前目录下得到great_module.so，同理可以在Python中直接使用。
测试结果如下
2.[转载] Python调用C/C++（使用SWIG）
下面转载自Jerry Jho在知乎“如何实现 C/C++ 与 Python 的通信？”问题下的回复。
示例简单，说明清晰，在我电脑上测试无误。
用C/C++对脚本语言的功能扩展是非常常见的事情，Python也不例外。除了SWIG，市面上还有若干用于Python扩展的工具包，比较知名的还有Boost.Python、SIP等，此外，Cython由于可以直接集成C/C++代码，并方便的生成Python模块，故也可以完成扩展Python的任务。答主在这里选用SWIG的一个重要原因是，它不仅可以用于Python，也可以用于其他语言。如今SWIG已经支持C/C++的好基友Java，主流脚本语言Python、Perl、Ruby、PHP、JavaScript、tcl、Lua，还有Go、C#，以及R。SWIG是基于配置的，也就是说，原则上一套配置改变不同的编译方法就能适用各种语言（当然，这是理想情况了……）SWIG的安装方便，有Windows的预编译包，解压即用，绿色健康。主流Linux通常集成swig的包，也可以下载源代码自己编译，SWIG非常小巧，通常安装不会出什么问题。还记得本文第2节的那个great_function吗？有了SWIG，事情就会变得如此简单：
/* great_module.i */ %module great_module %{ int great_function(int a) { return a + 1; } %} int great_function(int a);换句话说，SWIG自动完成了诸如Python类型转换、module初始化、导出代码表生成的诸多工作。
对于C++，SWIG也可以应对。例如以下代码有C++类的定义：
//great_class.h #ifndef GREAT_CLASS #define GREAT_CLASS class Great { private: int s; public: void setWall (int _s) {s = _s;}; int getWall () {return s;}; }; #endif // GREAT_CLASS
对应的SWIG配置文件
/* great_class.i */ %module great_class %{ #include "great_class.h" %} %include "great_class.h"
这里不再重新敲一遍class的定义了，直接使用SWIG的%include指令
SWIG编译时要加-c++这个选项，生成的扩展名为cxx
Windows下编译：swig -c++ -python great_class.i
cl /LD great_class_wrap.cxx /o _great_class.pyd -IC:\Python27\include C:\Python27\libs\python27.lib
Linux，使用C++的编译器
在Python交互模式下测试：g++ -fPIC -shared great_class_wrap.cxx -o _great_class.so -I/usr/include/python2.7/ -lpython2.7
也就是说C++的class会直接映射到Python classSWIG非常强大，对于Python接口而言，简单类型，甚至指针，都无需人工干涉即可自动转换，而复杂类型，尤其是自定义类型，SWIG提供了typemap供转换。而一旦使用了typemap，配置文件将不再在各个语言当中通用。>>> import great_class >>> c = great_class.Great() >>> c.setWall(5) >>> c.getWall() 5
3.Python调用C/C++（convolve函数，使用SWIG）
由上可见，使用swig相比第一种方法要简单很多。
需要实现的函数在Python中的代码是
def conv1d(a,b): m=len(a) n=len(b) t = m + n - 1 r = [0] * t for i in range(t): for j in range(min(m,i+1)): if (i - j < n): r[i] += a[j] * b[i-j] return r
使用swig时遇到几个问题：
1. 输入的Python数组无法被c/c++识别
这里采用
%typemap(in) int[ANY](int temp[$1_dim0])
将参数中的PySequence转化为c中的数组。
$1_dim0指的是数组传入大小，如下例中为20，PyObject_Length($input)为输入数据中PyList的长度。
通过一个循环，将每个数据用PySequence_GetItem($input,i)获取并传递给 $1作为转化后的参数。
2. 输出的c数组指针无法被Python识别
这里采用
%typemap(out) int* conv1d
将返回值c数组转化为PyList。
通过一个循环，将每个数据用PyList_SetItem($input,i)赋值并传递给 $result作为转化后的返回值。
3.将c数组转化为PyList时，不知道长度。
这里采用
%typemap(in,numinputs=0,noblock=1) size_t *len
将长度传递给上一步中的typemap。
numinputs=0 指不需要数据输入。
Typemaps支持一个名为noblock的特殊属性，其中可以使用{...}分隔符，但是分隔符实际上不会生成到代码中。 效果类似于使用“”或％{％}分隔符，但代码通过预处理器运行。
%module t %include <stdint.i> %typemap(in,numinputs=0,noblock=1) size_t *len { size_t templen; $1 = &templen; } %typemap(in) int[ANY](int temp[$1_dim0]) { int i; if (!PySequence_Check($input)) { PyErr_SetString(PyExc_TypeError,"Expecting a sequence"); return NULL; } if (PyObject_Length($input) > $1_dim0) { PyErr_SetString(PyExc_ValueError,"Expecting a sequence with less elements"); return NULL; } int j=PyObject_Length($input); for (i =0; i < j; i++) { PyObject *o = PySequence_GetItem($input,i); if (!PyInt_Check(o)) { Py_XDECREF(o); PyErr_SetString(PyExc_ValueError,"Expecting a sequence of ints"); return NULL; } temp[i] = PyInt_AsLong(o); Py_DECREF(o); } $1 = &temp[0]; } %typemap(out) int* conv1d { int i; $result = PyList_New(templen); for (i = 0; i < templen; i++) { PyObject *o = PyInt_FromLong((long)$1[i]); PyList_SetItem($result,i,o); } } %inline %{ int *conv1d(int a[20],int b[20],int n, int m, size_t *len) { int j,i; const int t = n+m-1; int *r = (int *) malloc((t)*sizeof(int)); for (i = 0; i < t; i++) for (r[i]=0,j = 0; j < m && j <= i; j++) if (i - j < n) r[i] += a[j] * b[i - j]; *len = t; return r; } %}
4.Python调用C/C++（convolve函数，使用Cython）
以下介绍引用自用Cython来提高Python代码速度 [一]
什么是Cython？
Cython并不是C，也不是Python，他是一个将介于C和Python之间的语言结构编译成C，然后在Python中运行的接口。它保持了大多数Python的语法结构，支持所有Python的命令，并通过一定方法加入C的结构，从而绕过Python的重重转译和类型检查，来实现快速的程序。
大多数Python整合包都会包括Cython。你也可以通过pip来安装：
pip install Cython
与Python不一样，你并不能直接使用Cython的源代码。需要通过编译来完成从Python调用C程序的目的。此外，没有错，如果你有用C直接编写的程序，可以用Cython来调用纯C程序！
Cython同样支持C++
Cython自带一些C和C++的标准库
ipython支持直接在命令行编译和使用cython程序！然而我的实践发现可能出现奇怪bug。
Cython的参考书非常少，最常用的一本是Cython A GUIDE FOR PYTHON PROGRAMMERS。然而……此书非常不详尽，绝对不是noob friendly…… 建议有什么问题还是stackoverflow求助大牛。
什么时候用Cython?
切记：Cython并不是所有python速度问题的良药。它牺牲了flexibility来换取速度，你必须要考虑所有细节问题，包括很多python新手甚至老手忽略了的、Python帮你暗地里处理好的问题。
-
如果你需要大量使用python的类，请慎重考虑是否使用。
-
如果并没有大量循环，或者最内层的循环有一些复杂python的类，请慎重使用。
-
如果不需要重复使用脚本，请慎重使用。
-
如果提升只有不到20%，请慎重使用。
-
如果这个函数或者类会被很多人debug和改进，请慎重使用。除非很多人熟悉Cython，不然maintain程序会非常麻烦。
-
如果不能确定性能瓶颈是你要改写为cython的函数，请慎重使用。
-
如果你不会从200%的速度提升中获得快感，请慎重使用。
根据上面的需求，我们制作了如下pyx文件。
import numpy as np cimport cython cimport numpy as np @cython.boundscheck(False) @cython.wraparound(False) def _cython_convolve(np.ndarray[np.int_t, ndim=1] a, np.ndarray[np.int_t, ndim=1] b): cdef int m = a.shape[0] cdef int n = b.shape[0] cdef int i,j,t = m + n - 1 cdef np.ndarray[np.int_t, ndim=1] r = np.zeros(t, dtype=np.int) for i in range(t): for j in range(min(m,i+1)): if (i - j < n): r[i] += a[j] * b[i-j] return r def cython_convolve(a, b): return _cython_convolve(np.array(a), np.array(b))
以及如下的setup.py文件
# setup.py from distutils.core import setup, Extension from Cython.Build import cythonize import numpy setup(ext_modules = cythonize(Extension( 'cython_convolve', sources=['cython_convolve.pyx'], language='c', include_dirs=[numpy.get_include()], library_dirs=[], libraries=[], extra_compile_args=[], extra_link_args=[] )))
在终端运行
python setup.py build_ext --inplace
生成的cython_convolve.so即可用于调用。
5.在Python中的使用效果
自定义的conv1d最慢，使用swig实现的最快。使用cython实现的比numpy自带的略慢。
import scipy.signal.signaltools as sg import numpy as np import timeit import t import cython_convolve a = [[1,2,3,4,1,2,3,4,1,2,3,4,1,2,3,4]] b = [[1,2,3,4,1,2,3,4,1,2,3,4,1,2,3,4]] def conv1d(a,b): m=len(a) n=len(b) t = m + n - 1 r = [0] * t for i in range(t): for j in range(min(m,i+1)): if (i - j < n): r[i] += a[j] * b[i-j] return r def npconv(): np.convolve(a[0],b[0]) def sgconv(): sg.convolve2d(a, b) def newconv(): conv1d(a[0], b[0]) def tconv(): t.conv1d(a[0],b[0],16,16) def cconv(): cython_convolve.cython_convolve(a[0], b[0]) if __name__=='__main__': print sg.convolve2d(a, b) print np.convolve(a[0], b[0]) print conv1d(a[0], b[0]) print t.conv1d(a[0],b[0],16,16) print cython_convolve.cython_convolve(a[0], b[0]) print timeit.timeit('sgconv()',setup='from __main__ import sgconv',number=100000) print timeit.timeit('npconv()',setup='from __main__ import npconv',number=100000) print timeit.timeit('newconv()',setup='from __main__ import newconv',number=100000) print timeit.timeit('tconv()', setup='from __main__ import tconv', number=100000) print timeit.timeit('cconv()', setup='from __main__ import cconv', number=100000) """ [[ 1 4 10 20 27 32 36 40 53 60 62 60 79 88 88 80 103 108 94 60 77 80 68 40 51 52 42 20 25 24 16]] [ 1 4 10 20 27 32 36 40 53 60 62 60 79 88 88 80 103 108 94 60 77 80 68 40 51 52 42 20 25 24 16] [1, 4, 10, 20, 27, 32, 36, 40, 53, 60, 62, 60, 79, 88, 88, 80, 103, 108, 94, 60, 77, 80, 68, 40, 51, 52, 42, 20, 25, 24, 16] [1, 4, 10, 20, 27, 32, 36, 40, 53, 60, 62, 60, 79, 88, 88, 80, 103, 108, 94, 60, 77, 80, 68, 40, 51, 52, 42, 20, 25, 24, 16] [ 1 4 10 20 27 32 36 40 53 60 62 60 79 88 88 80 103 108 94 60 77 80 68 40 51 52 42 20 25 24 16] 1.34557986259 0.578319072723 5.11206579208 0.260015964508 0.591851968765 """
