|
|
Я посчитал, и у меня получилось вот что:
Пусть исходный спектр длиной N есть H. Тогда пусть спектр длиной 2*N есть
H0=fft([real(ifft(H)) zeros(1, length(H)));
По формулам у меня получилось вот что:
H(k):= fft(h) = sum{ h(n)*exp(-j*2*pi*k*n/N), n=0..N-1 };
H0(k):=fft(h0) = sum{ h0(n)*exp(-j*2*pi*k*n/(2*N)), n=0..2*N-1 } =
= {поскольку h0(n)=h(n) для n=0..N-1 и h0(n)=0 для n=N..2*N-1} =
= sum{ h(n)*exp(-j*2*pi*k*n/(2*N)), n=0..N-1 }.
Рассмотрим случай четных отсчетов k=2*k1.
Тогда H0(2*k1) = sum{ h(n)*exp(-j*2*pi*2*k*n/(2*N)), n=0..N-1 } =
= {сокращаем числитель и знаменатель под экспонентой на 2} =
= sum{ h(n)*exp(-j*2*pi*2*k*n/(2*N)), n=0..N-1 } =
sum{ h(n)*exp(-j*2*pi*k*n/N), n=0..N-1 } = H(k1).
Таким образом, по поводу четных отсчетов все верно: H0(2*k1)=H(k1).
Но вот по поводу нечетных получается все гораздо более неприятно:
H0(2*k+1) = sum{ h(n)*exp(-j*2*pi*(2*k+1)*n/(2*N)), n=0..N-1 } =
= sum{ h(n)*exp(-j*2*pi*n/(2*N))*exp(-j*2*pi*2*k*n/(2*N)),n=0..N-1 } =
= {сокращаем числитель и знаменатель под второй экспонентой на 2} =
= sum{ h(n)*exp(-j*2*pi*n/(2*N))*exp(-j*2*pi*k*n/N),n=0..N-1 } =
= fft{ h(n)*exp(-j*2*pi*n/(2*N)) } =
= fft{h(n)} (+) fft(exp(-j*2*pi*n/(2*N))) =
= {поскольку fft{h(n)} =: H(k)}
= H (+) EXP(k).
Здесь оператор (+) обозначает операцию циклической свертки;
EXP(k) = fft(exp(-j*2*pi*n/(2*N));
Проблема как раз в функции EXP(k) - очень уж она убогая. Хорошей будет функция EXP_good(k) = fft(exp(-j*2*pi*n/N); См. примеры обеих ниже. Получается, что эту свертку с ней делать еще хуже, чем выполнить
fft([real(ifft(H)) zeros(1, length(H)]).
N = 100;
n = 0:1:(N-1);
ex = exp(-j*2*pi*n/N);
EX = fft(ex);
AEX = abs(EX);
ex2 = exp(-j*2*pi*n/(2*N));
EX2 = fft(ex2);
AEX2 = abs(EX2);
figure(12);
plot(AEX2, 'r');
grid on;
hold on;
plot(AEX, 'b');
|
|
E-mail: 
info@telesys.ru
