FFT: FFTW Matlab FFT2 mystery

Question

Я унаследовал старый код fortran с подпрограммой fft, и я не могу отследить источник этой программы. Единственное, что я знаю, это вызов ff2prp() и вызов fft2() для выполнения 2D прямого и обратного DFT. Чтобы узнать, что делает код, я беру DFT 2D-матрицы (матрицы) 4x4, и результаты сильно отличаются от результатов Matlab и FFTW.

Вопрос: Может ли кто-нибудь сказать, что делает код, глядя на выход. Вход и выход оба реальные массивы

Входной массив

0.20000  0.30000  1.00000  1.20000 
0.00000 12.00000  5.00000  1.30000 
0.30000  0.30000  1.00000  1.40000 
0.00000  0.00000  0.00000  1.50000 

После прямого преобразования Фурье с fft2() Fortran подпрограмма

0.16875 -0.01875 -0.05000  0.05625 
0.00000 12.00000  5.00000  1.30000 
0.30000  0.30000  1.00000  1.40000 
0.00000  0.00000  0.00000  1.50000 

Matlab выход выполнения ДКП: dct2 (вход)

6.3750 -0.8429 -3.4250 -2.4922 
    2.4620 0.6181 -2.6356 -0.9887 
    -4.2750 -0.9798 4.2250 2.2730 
    -4.8352 -1.2387 5.0695 3.4819 

Вывод кода C++ с использованием библиотеки FFTW. DCT от FFTW

(6.3750, 0.00) (-0.8429, 0.00) (-3.4250, 0.00) (-2.4922, 0.00) 
(2.4620, 0.00) (0.6181, 0.00) (-2.6356, 0.00) (-0.9887, 0.00) 
(-4.2750, 0.00) (-0.9798, 0.00) (4.2250, 0.00) (2.2730, 0.00) 
(-4.8352, 0.00) (-1.2387, 0.00) (5.0695, 0.00) (3.4819, 0.00) 

вперед FFT с Matlab - fft2 (вход)

25.5000 + 0.0000i -6.5000 - 7.2000i -10.5000 + 0.0000i -6.5000 + 7.2000i 
    -0.3000 -16.8000i -12.3000 + 4.8000i 0.1000 + 6.8000i 12.1000 + 5.2000i 
-14.1000 + 0.0000i 3.5000 +11.2000i 9.1000 + 0.0000i 3.5000 -11.2000i 
    -0.3000 +16.8000i 12.1000 - 5.2000i 0.1000 - 6.8000i -12.3000 - 4.8000i 

вперед FFT с FFTW

(25.50, 0.00) (-6.50, -7.20) (-10.50, 0.00) (-6.50, 7.20) 
(-0.30, -16.80) (-12.30, 4.80) (0.10, 6.80) (12.10, 5.20) 
(-14.10, 0.00) (3.50, 11.20) (9.10, 0.00) (3.50, -11.20) 
(-0.30, 16.80) (12.10, -5.20) (0.10, -6.80) (-12.30, -4.80) 

Как вы тростниковый видите выход Matlab и FFTW согласуется с друг друга, но не с выходом кода fortran. Я хотел бы использовать FFTW, но результаты отличаются от FFT. Я не могу понять, что делает FFT программа fortran. Может кто-нибудь сказать, глядя на выход.

Answer 1

Насколько я могу судить, fft2 кажется, рассчитали 1D FFT первой строки (оставляя все 3 другие неизменными), с результатом масштабируются 1/16 и упаковано в формате r0, r2, r1, i1.

Другими словами, выход может быть построен в Matlab с помощью:

input = [0.2 0.3 1 1.2;0 12 5 1.3;0.3 0.3 1 1.4;0 0 0 1.5]; 
N = size(input,2); 
A = fft(input(1,:))/16; 
B = reshape([real(A);imag(A)],1,2*N); 
B(2) = B(N+1); 
output = [B(1:N);A(2:size(input,1),:)]; 

Если у вас есть основания полагать, что fft2 должны быть вычисления 2D БПФ, то, возможно, какая-то проблема в том, как вы передать ваши данные к этой рутине, которая вызывает неправильные результаты. Кроме того, дополнительные тестовые примеры (или то, как вы вызываете ff2prp) для входов разных размеров, могут дать более полное представление о выборе коэффициента масштабирования (например, это 1/N^2 или 1/4N или что-то еще).