Cプログラムによる画像の高速フーリエ変換FFT

Cプログラムにる画像処理について教えてください。逆高速フーリエ変換（IFFT)によって周波数領域から実空間領域に変換する際、周波数帯域を（低周波数領域に）制限してIFFTを実行できるという文献があるのですが可能でしょうか？通常の逆離散フーリエ変換（IDFT)では可能でした。
以下のコードの修正で可能でしょうか？
例）　５１２＊５１２の空間周波数領域データを212*212（仮に）の低周波数領域に帯域制限してIFFTで画像に戻したいのです。
以下の以下の一次元FFTコードを用いてX方向とY方向にIFFTしたいと考えています。

void FFT(int ir, int nx, float *xr, float *xi, float *si, float *co, unsigned short *brv)
// 1次元フーリエ変換
// int ir; 順変換(1)と逆変換(-1)
// int nx; 1次元FFTのデータ数
// float *xr; 実部のデータ xr[nx]
// float *xi; 虚部のデータ xi[nx]
// float *si; FFT用のサインデータ si[nx/2]
// float *co; FFT用のコサインデータ co[nx/2]
// unsigned short *brv; FFT用の入れ替えデータ brv[nx]
{
int i, j, n1, n2=nx, j3, j4, k, l, ll, d=1, g;
float a, b, c, s;

for(l = 1; l <= nx/2; l *= 2, d += d) {
g = 0;
ll = n2;
n2 /= 2;
for(k = 1; k <= n2; k++) {
n1 = k-ll;
c = co[g];
s = -ir*si[g];
g += d;
for(j = ll; j <= nx; j += ll) {
j3 = j+n1-1;
j4 = j3+n2;
a = xr[j3]-xr[j4];
b = xi[j3]-xi[j4];
xr[j3] += xr[j4]; xi[j3] += xi[j4];
xr[j4] = c*a+s*b; xi[j4] = c*b-s*a;
}
}
}

通常のDFTでの帯域制限はこのように行い可能でした。
void fourier1d(int ir, float *fr, float *fi, int nx)
{
int i, j, n = 1;
float *gr, *gi;
double u, x;

gr = (float *)malloc((unsigned long)nx*sizeof(float));
gi = (float *)malloc((unsigned long)nx*sizeof(float));
for(i = 0 ; i < nx ; i++) {
u = i-nx/2;
gr[i] = gi[i] = 0;
for(j = 150 ; j < 362 ; j++) {
x = j-nx/2;
gr[i] += (float)( fr[j]*cos(2*PI*u*x/nx)+ir*fi[j]*sin(2*PI*u*x/nx));
gi[i] += (float)(-ir*fr[j]*sin(2*PI*u*x/nx)+fi[j]*cos(2*PI*u*x/nx));
}
}
if(ir == -1) n = 212; // 逆変換はデータ数で割る
for(i = 0 ; i < nx ; i++) {
fr[i] = gr[i]/n;
fi[i] = gi[i]/n;
}
free(gr);
free(gi);
}

No.2 回答者： Tacosan 回答日時： 2011/05/10 11:30

文章は逆フーリエ変換だけどコードは (逆じゃない) フーリエ変換だ. どっちが正しいんだろう.

いずれにしても, 「FFT ないし IFFT の途中で帯域制限」はできないことはないだろうけどする必要があるかどうかが疑問. #1 でも言われてるけどおそらくかなり面倒なので, よほど狭い帯域に制限するのでなければ普通にやった方がお手軽かつ高速な可能性もあります.

# でもって, そこまで狭帯域に制限するなら DFT の方が速かったりする.

No.1 回答者： ICE_FALCON 回答日時： 2011/05/09 23:47

あんまりソース見ていないですけど、

(1)帯域制限をしてから(2)IFFTするということですか？
それなら簡単にできます。
FFTとDFTは演算量が違うだけで結果は同じですから。
圧縮では普通に使われる技術ですね（DCTだったりもしますが）。

それとも、
IFFTを実行しつつ帯域制限もして（制限した帯域の）無駄な演算も行わないということでしょうか？
やったこと無いですが、めんどうそうですね、
制限する帯域によっては効率よくできそうですが、普通のIFFT演算ではなくなりそうです。

この回答へのお礼

IFFTを実行しつつ帯域制限もして（制限した帯域の）無駄な演算も行わないということです。帯域制限位相限定相関法（BLPOC)というものに使用します。IFFTを用いてやりたいのですがコード作成が難しくてできません。

お礼日時：2011/05/10 01:31