9

これは「プログラミング」の問題ではありません。しかし、それはこのコミュニティで広く知られ、理解されていることだと確信しています。

画像 x と、はるかに小さい画像 y があり、それらの FFT を乗算して 2 つを畳み込む必要があります。しかし、それらは同じサイズではないため、周波数領域の乗算を行う方法がわかりません。

x (次元 4096 x 4096 の整数行列) の (2 次元) FFT を取得すると、周波数領域表現 X (複素数の行列であり、次元は 2048 x 2048 だと思います) が得られます。 ）。

同様に、(次元 64 x 64 の整数行列である y の 2 次元 FFT を取得します。これにより、周波数領域表現 Y (複素数の行列でもあり、次元は 32 であると思います) が得られます。 ×32）。

Numerical Recipes で fourn 関数を使用しているため、入力行列 x と y を 1 次元配列に折りたたむ必要があり、これが離散フーリエ変換 X と Y に置き換えられます。画像の二次元問題、私は一次元配列で作業しています。

まったく同じ次元の x と y の 2 つの画像を畳み込みたいとします。それはすべて非常に簡単です：

 X = FFT(x)

 Y = FFT(y)

 Z = X * Y (term by term multiplication)

 Convolution of x and y = IFFT(Z)

しかし、X と Y の長さが異なる場合、どのように乗算を行うのでしょうか?

1 つの可能性は、x と同じ次元になるように y をパディングすることです。しかし、これは恐ろしく非効率的です。別の可能性は、Y をパディングして X と同じ次元にすることです。しかし、これが周波数空間で何を意味するのかわかりません。

この質問の別の言い方は次のとおりです。FFT を使用して非常に異なる次元の 2 つの画像を畳み込みたい場合、それらのスペクトル (周波数領域表現) の乗算を行うことができます。その乗算を行うにはどうすればよいですか?

ありがとう、

〜マイケル。

4

2 に答える 2

8

入力画像サイズ（行列x）に一致するように、小さい方の配列（畳み込みカーネル、この場合はy）をゼロでパディングするのが標準的なアプローチです。空間領域で畳み込みを行う場合、これは非常に非効率的ですが、FFT を乗算する場合は必要であり、パディングされた配列の FFT を計算するコストはそれほど悪くありません。

于 2009-09-03T22:38:42.757 に答える
4

2 つの周波数間隔が同じである必要があると考えるのは正しいことです。1D の例を見てみましょう (私は Matlab 構文を使用しています):

N = 4096;
M = 64;

x = randn(N, 1);
y = hann(M, 'symmetric');

zLinear = conv(x,y);
zCircular = ifft( fft(x) .* fft(y,N) );

disp(max(abs(zLinear(65:4096) - zCircular(65:4096))));

2 つの手法の差は ~2e-14 なので、丸め誤差です。線形畳み込みと循環畳み込みの違いにより、最初の 64 サンプルをスキップする必要があることに注意してください。

zCircular の計算では、fft(y,N) に注意してください。これは、fft を計算する前に最大 N までゼロで y 信号をパディングするための Matlab 構文です。これはメモリ使用効率が悪いと考えられるかもしれませんが、速度を比較してください。

線形畳み込み: 各 64 の 4096 の乗算/加算 = 262144 の乗算/加算

巡回畳み込み: 4096 の 2 つの FFT + 2 * 4096 要素の 1 つの複素乗算 + 1 つの逆 FFT
= 3 * 4096 * log2(4096) + 4096 * 6 = 172032 (複素数の乗算に 6 つの操作を想定)

基本的に、FFT の NlogN 速度は、3 つ必要な場合でも、M が非常に短い場合を除き、N * M 畳み込み操作よりも優れています。

EDIT 2Dケースの速度推定を追加

2D データの場合、速度の利点が拡大されることを付け加えておく価値があります。2D FFT は N * N * log2(N * N) 演算を必要とするため、N = 4096 の 3 つの FFT + 複素数 N^2 配列の乗算は 1.3e10 演算です。しかし、直接畳み込みは N^2 * M^2 = 6.9e10 演算であり、約 50 倍遅くなります。

于 2009-09-03T22:39:32.753 に答える