Matlab, fft2, баг?

[fillsv 0]

Столкнулся с необходимостью делать двумерное фурье набора картинок. Это можно делать двумя путями:

1) через цикл прогнать все картинки по одной.

2) загнать из в трехмерную матрицу и обработать одной командой.

Забавность заключается в том, что при выполнении двумерного фурье двуми, казалось бы, идентичными способами, получаются слегка разные результаты. Ниже пример этого в коде:

Создаем исходный массив:

A = rand(100,100,10);

Обработка через цикл:

for ii = 1:10
    f1(:,:,ii) = fft2(A(:,:, ii));
end

Обработка единой матрицы:

f2 = fft2(A);

Сравнение результатов:

std(f1(:)-f2(:))

ans =

   1.6847e-14

С практической точки зрения, отличие двух результатов, конечно же, ничтожное. Однако таки интересно, почему оно есть? В хелпе по этому поводу написано следующее:

If the dimensionality of X is greater than 2, the fft2 function returns the 2-D DFT for each higher dimensional slice of X. For example, if size(X) = [100 100 3], then fft2 computes the DFT of X(:,:,1), X(:,:,2) and X(:,:,3).

Таки это баг или все в пределах нормы?

[peselnik 0]

Это нормально. Если посмотреть код функции fft2, то она вызывает либо fftn, либо fft.

В первом случае вызывается fftn, потому что на входе матрица. Во втором случае вызывается fft, потому что на входе не матрица.

 

Обе функции fftn и fft - встроенные в MATLAB и обращаются к библиотеке FFTW.

 

О вычислительной разнице между одномерным и многомерным преобразованием написано в документации FFTW:

 

http://www.fftw.org/fftw3_doc/Multi_002ddi...Transforms.html

 

We should remind the user that the separable product of 1d transforms along each dimension, as computed by FFTW, is not always the same thing as the usual multi-dimensional transform. A multi-dimensional R2HC (or HC2R) transform is not identical to the multi-dimensional DFT, requiring some post-processing to combine the requisite real and imaginary parts, as was described in The Halfcomplex-format DFT.

[fillsv 0]

Спасибо! Это можно на досуге покопать поподробнее)