Цитата(Alexeis @ 11.11.2008, 21:07  ) | | Какого еще зеркального эффекта? |
С верху как должно быть а снизу получается в итоге:
Применяемый мною код БПФ:
| Код |
unit fft;
interface
uses Math, Ap, Sysutils;
procedure FastFourierTransform(var a : TReal1DArray; nn : Integer;
InverseFFT : Boolean);
implementation
uses Unit1;
(*************************************************************************
Быстрое преобразование Фурье
Алгоритм проводит быстрое преобразование Фурье комплексной
функции, заданной nn отсчетами на действительной оси.
В зависимости от переданных параметров, может выполняться
как прямое, так и обратное преобразование.
Входные параметры:
nn - Число значений функции. Должно быть степенью
двойки. Алгоритм не проверяет правильность
переданного значения.
a - array [0 .. 2*nn-1] of Real
Значения функции. I-ому значению соответствуют
элементы a[2*I] (вещественная часть)
и a[2*I+1] (мнимая часть).
InverseFFT
- направление преобразования.
True, если обратное, False, если прямое.
Выходные параметры:
a - результат преобразования. Подробнее см.
описание на сайте.
*************************************************************************)
procedure FastFourierTransform(var a : TReal1DArray; nn : Integer;
InverseFFT : Boolean);
var
Jj, n ,Mmax ,m , j, istep, i, isign: Integer;
wtemp, wr, wpr, wpi, wi, theta, tempr, tempi: Extended;
begin
if InverseFFT then isign := -1
else isign := 1;
n := 2 * nn;
j := 1;
I := 1;
mmax := 2;
while I <= n do // Реверс
begin
if j > i then
begin
tempr := a[j - 1];
tempi := a[j];
a[j - 1] := a[i - 1];
a[j] := a[i];
a[i - 1] := tempr;
a[i] := tempi;
end;
m := nn;
while (m >= 2) and (j > m) do
begin
j := j - m;
m := m div 2;
end;
j := j + m;
I := I + 2;
end; // Реверс
while n > mmax do // FFT
begin
istep := 2 * mmax; // 4, 8, 16..2 * NN
theta := 2 * Pi * isign / mmax; // Прямое 2 * Pi * (-1) / (2, 4, 8, 16..2 * NN) Обратное 2 * Pi * 1 / (2, 4, 8, 16..2 * NN)
wpr := -2 * sqr(sin(0.5 * theta)); // Cos(X) - 1
wpi := sin(theta); // Sin(X)
wr := 1;
wi := 0;
M := 1;
while M <= mmax do
begin
for Jj := 0 to (n - m) div istep do
begin
i := m + Jj * istep;
j := i + mmax;
tempr := wr * a[j - 1] - wi * a[j];
tempi := wr * a[j] + wi * a[j - 1];
a[j - 1] := a[i - 1] - tempr;
a[j] := a[i] - tempi;
a[i - 1] := a[i - 1] + tempr;
a[i] := a[i] + tempi;
end;
wtemp := wr;
wr := wr * wpr - wi * wpi + wr;
wi := wi * wpr + wtemp * wpi + wi;
M := M + 2;
end;
mmax := istep;
end; //FFT
if InverseFFT then // Обратное FFT
for I := 0 to N - 1 do
a[I] := a[I] / nn;
end;
end.
|
|
