Доброго времени суток. Разбираюсь с FFT и IFFT. Уже какой день парюсь над обратным преобразованием Фурье. Добрые люди, может объясните в чем может быть проблема?
Звуковой сигнал 2 кГц, который имитирует программка справа, поступает на микрофон. Данные заносятся в массив и рисуется входной сигнал - верхняя диаграмма.
Этот массив с данными раскладываю в FFT, получаю массив fft и рисую диаграмму - второй график сверху. Проверял на разных частотах, все здесь правильно.
И наконец, я делаю обратное FFT, взяв массив fft. В итоге, 3-й график, вроде и частота сигнала та же, что и на входе, но почему-то блин амплитуда сигнала падает. И видно,что на стыке двух спектров (основного и зеркального) амплитуда минимальна.
Делал эксперимент, если в спектре четко только одна гармоника, то при IFFT все нормально, амплитуда сигнала равномерна по всей длине, как на входном графике.
Как только в спектре есть две гармоники и более, то сразу же график как на картинке. Как будто перемножение что ли где-то, или что может быть?
Буду благодарен знатокам!
Код:
Код:
procedure TFFTBase.InitMem(WindowSize: Integer);
begin
fWindowSize:=WindowSize;
tTlbSize:=(2*fWindowSize+1)*SizeOf(Single);
fSinTbl:=AllocMem(tTlbSize);
fCosTbl:=AllocMem(tTlbSize);
end;
procedure TFFTBase.InitSinCosTbl;
var i :integer;
begin
for i :=1 to 2*fWindowSize do begin
fSinTbl[i] := (-1)*Sin(PI/i);
fCosTbl[i] := Cos(PI/i);
end;
end;
NV2 :=N shr 1;
NM1 :=N-1;
J :=1;
if Inverse then
for i :=0 to N-1 do
TComplexArray(F^)[i].Im :=-TComplexArray(F^)[i].Im;
for I :=1 to NM1 do begin
if I<J then begin
T :=TComplexArray(F^)[J-1];
TComplexArray(F^)[J-1] :=TComplexArray(F^)[I-1];
TComplexArray(F^)[I-1] :=T;
end;
K :=NV2;
while K < J do begin
J :=J - K;
K :=K shr 1;
end;
J :=J + K;
end;
for L :=1 to M do begin
LE :=2 shl (L-1);
LE1 :=LE shr 1;
U.Re :=1.0; U.Im :=0.0;
W.Re :=fCosTbl[LE1];
W.Im :=fSinTbl[LE1];
for J :=1 to LE1 do begin
I :=J;
while I <= N do begin
IP :=I + LE1;
T.Re :=TComplexArray(F^)[IP-1].Re * U.Re - TComplexArray(F^)[IP-1].Im * U.Im;
T.Im :=TComplexArray(F^)[IP-1].Re * U.Im + TComplexArray(F^)[IP-1].Im * U.Re;
TComplexArray(F^)[IP-1].Re :=TComplexArray(F^)[I-1].Re - T.Re;
TComplexArray(F^)[IP-1].Im :=TComplexArray(F^)[I-1].Im - T.Im;
TComplexArray(F^)[I-1].Re:=TComplexArray(F^)[I-1].Re+T.Re;
TComplexArray(F^)[I-1].Im:=TComplexArray(F^)[I-1].Im+T.Im;
Inc(I,LE);
end;
Uo :=U;
U.Re :=(Uo.Re * W.Re) - (Uo.Im * W.Im);
U.Im :=(Uo.Re * W.Im) + (Uo.Im * W.Re);
end;
end;
ImDummy :=1/Sqrt(N);
if Inverse then ImDummy :=-ImDummy;
ReDummy :=Abs(ImDummy);
for I :=1 to N do
begin
TComplexArray(F^)[i-1].Re :=TComplexArray(F^)[i-1].Re*ReDummy;
TComplexArray(F^)[i-1].Im :=TComplexArray(F^)[i-1].Im*ImDummy;
end;
DelMem;
end;