直接在scipy.weave.inlin

Question

我试图在Python中实现一个基于FFT的亚像素移位（翻译）算法。Fourier移位定理允许通过以下方式将数组转换为亚像素量： 1前向FFT阵列 2傅里叶空间中线性相位斜坡倍增阵 三。逆FFT阵列

这个算法很容易在python中使用numpy/scipy来实现，但是对于256**2数组来说，它的每移位速度非常慢（大约10毫秒）。我正试图通过使用scipy.weave.inline. 在

但是，在将复杂的numpy数组传递给FFTW时遇到了麻烦。c代码看起来像：

    #include <fftw3.h>
    #include <stdlib.h>

    #define INVERSE +1
    #define FORWARD -1


    fftw_complex *i, *o;
    int n, m;
    fftw_plan pf, pi;
    #line 22 "test_scipy_weave.py" 

    i = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * xdim*ydim);
    o = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * xdim*ydim);

    pf = fftw_plan_dft_2d(xdim, ydim, i, o, -1, FFTW_PATIENT);
    pi = fftw_plan_dft_2d(xdim, ydim, o, i,  1, FFTW_PATIENT);

    # Copy data to fftw_complex array. How to use python arrays directly
    for (n=0; n<xdim;n++){
        for (m=0; m<ydim; m++){
            i[n*xdim+m][0]=a[n*xdim+m].real();
            i[n*xdim+m][1]=a[n*xdim+m].imag();
        }
    }

    fftw_execute(pf);

    /* Mult by linear phase ramp here */

    fftw_execute(pi);

    for (n=0; n<xdim;n++){
        for (m=0; m<ydim; m++){
            b[n*xdim+m] = std::complex<double>([in*xdim+m][0], i[n*xdim+m][1]);
        }
    }

    fftw_destroy_plan(p);

所以您可以看到我必须将存储在numpy数组“a”中的数据复制到fftw峎复杂数组“I”中。最后，我必须将结果“I”复制到输出numpy数组“b”中。直接在fftw中使用numpy数组“a”和“b”会更有效，但我无法使其工作。在

有人知道如何让fftw直接在scipy.weave.inline中使用复杂的numpy数组吗？在

谢谢

Answer 1

根据fftw manual，您可以在fftw.h之前导入{}，这将保证{}与本机C数据类型相对应。我很确定numpy数据类型也保证（或者在实践中很可能）与本机C数据类型兼容。在

在这种情况下，您可以访问指向数组数据的指针a.data_as(ctypes.c_void_p)。不幸的是，ctypes不能识别复杂的类型，但是希望转换成一个空指针就可以了。在

执行此操作时，必须注意数组a以C-连续的方式存储，在创建数组时由参数order='C'指定。在