でFFTベースのサブピクセルシフト(変換)アルゴリズムを実装しようとしていPython
ます。フーリエシフト定理により、アレイをサブピクセル量で次のように変換できます。1.順方向FFTアレイ2.フーリエ空間での線形位相ランプによるアレイの乗算3.逆方向FFTアレイ
このアルゴリズムは、numpy / scipyを使用してPythonで簡単に実装できますが、256 ** 2アレイの場合、シフトごとに非常に低速(〜10msec)になります。scipy.weave.inlineを使用してPythonから直接cコードを呼び出すことで、これを高速化しようとしています。
しかし、複雑な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_complex配列「i」にコピーする必要があることがわかります。そして最後に、結果「i」を出力のnumpy配列「b」にコピーする必要があります。fftwでnumpy配列「a」と「b」を直接使用する方がはるかに効率的ですが、これを機能させることはできません。
fftwが複雑なnumpy配列を直接使用する方法について誰かが考えていますscipy.weave.inline
か?
ありがとう