最小化 sinc 谱瓣

信息处理 fft 离散信号 卷积 窗函数
2022-01-30 05:34:39

由于某些原因,我想最小化矩形信号的 FFT 变换的次级波瓣(显示在下图的左上角)

在此处输入图像描述 所以我使用了窗口函数来做到这一点(左下角)。但是,如果我想直接在信号的复变换上进行(右下),则次级瓣的衰减值与时间情况不同。

有下面的测试代码,我想知道它有什么问题:

传递给编辑

     Fs = 1000; L = 2000; t = (0:L-1)*1/Fs;  
% the signal
y = zeros(1,L);  y(500:1500) = 1;     % rectangle signal

figure; 
subplot(2,2,1);
plot(Fs*t,y); title('Temporal signal'); xlabel('time (milliseconds)')

NFFT = 6*2^nextpow2(L);
Y1 = 2*abs(fft(y,NFFT)); Y1 = Y1/max(Y1);
f = 0:Fs/NFFT:Fs/2;

% the FFT of the signal
subplot(2,2,2);
semilogy(f, Y1(1:NFFT/2+1), '.-'); axis tight; grid on;
title('Signal''s spectrum'); xlabel('Frequency (Hz)'); ylabel('|Y(f)|');

w = window(@blackman, L);
y1 = y.*w';
Y2 = 2*abs(fft(y1,NFFT)); Y2 = Y2/max(Y2);

% the FFT of the windowed signal, in time
subplot(2,2,3);
semilogy(f, Y2(1:NFFT/2+1), '.-'); axis tight; grid on;
title('FFT windowing - temporal product'); xlabel('Frequency (Hz)'); ylabel('|Y(f)|')

% the FFT of the windowed signal, in complex
Y3 = 2*abs(fftshift(conv(Y1, fft(w,NFFT), 'same'))); Y3 = Y3/max(Y3);
subplot(2,2,4);
semilogy(f, Y3(1:NFFT/2+1), '.-');  axis tight; grid on;
title('FFT windowing - spectral convolution'); xlabel('Frequency (Hz)'); ylabel('|Y(f)|')

编辑
我已经根据@Matt 的评论更改了代码,但最后一张图仍然存在问题(FFT 卷积上的噪声)

NFFT = 6*2^nextpow2(L);
Y1 = fft(y,NFFT);
YY1 = 2*abs(Y1); YY1 = YY1/max(YY1);
f = 0:Fs/NFFT:Fs/2;

% the FFT of the signal
subplot(2,2,2);
semilogy(f, YY1(1:NFFT/2+1), '.-'); axis tight; grid on;
title('Signal''s spectrum'); xlabel('Frequency (Hz)'); ylabel('|Y(f)|');

w = window(@blackman, L);
y1 = y.*w';
YY2 = 2*abs(fft(y1,NFFT)); YY2 = YY2/max(YY2);

% the FFT of the windowed signal, in time
subplot(2,2,3);
semilogy(f, YY2(1:NFFT/2+1), '.-'); axis tight; grid on;
title('FFT windowing - temporal product'); xlabel('Frequency (Hz)'); ylabel('|Y(f)|')

% the FFT of the windowed signal, in complex
YY3 = 2*abs(cconv(Y1, fft(w,NFFT)', NFFT)); YY3 = YY3/max(YY3);
subplot(2,2,4);
semilogy(f, YY3(1:NFFT/2+1), '.-');  axis tight; grid on;
title('FFT windowing - spectral convolution'); xlabel('Frequency (Hz)'); ylabel('|Y(f)|')'

在此处输入图像描述

1个回答

两个时域序列的乘法对应于两个序列的 DFT 的循环(循环)卷积,而不是代码中实现的线性卷积。因此,如果x[n]y[n]是两个长度N时域序列,和X[k]Y[k]是它们各自的 DFT,以下对应关系成立:

DFTN{x[n]y[n]}[k]=1Nm=0N1X[m]Y[(km)N],0kN1

在哪里

(km)N=km+lN

与整数l选择使得结果在范围内[0,N1].

另请查看相关问题的答案。

编辑:您也不能Y1在卷积中使用,因为它不是 DFT,而是 DFT 的大小fft(y)您应该与循环卷积fft(w)