匹配滤波器的实现-时域卷积法和频域法

clear all;
close all;
clc;
%%******************* 线性调频与脉冲压缩 ********************%%
set(0,'defaultfigurecolor','w'); %设置MATLAB的图形窗口背景为白色
%% 线性调频信号参数设置
T = 1e-6;      %时宽
B = 200e6;     %带宽
Fs = 4*B;      %采样率

K = B/T;                          %调谐频率
N = round(T/(1/Fs) );             %采样点数，四舍五入
%% 产生线性调频信号
t = linspace( -T/2 ,T/2 ,N);
st = ( abs(t) < T/2 ) .* exp( 1j * pi * K * t.^2 ); 
phase = pi * K * t.^2;          %信号相位
f= K * t;                       %信号频率
%% 频谱
freq = linspace(-Fs/2,Fs/2,N);  %频域采样
Sf = fftshift( fft(st) );
%% 时域匹配滤波
% conj()函数用于计算复数的共轭值；
% fliplr()函数将数组从左到右翻转.
ht = conj( fliplr(st) );    %时域匹配滤波为发射信号时间反褶再取共轭：h(t)=s*(-t)
s1 = conv(st,ht);           %线性调频信号经过匹配滤波器后的输出(时域卷积)
N1 = N+N-1 ;                %线性卷积后信号长度变为 N1+N2-1
t1 = linspace( -T/2 , T/2 , N1);
%% 频域匹配滤波1 (复制发射脉冲进行时间反褶并取共轭，计算补零DFT)
N2 = 2*N;                    %循环卷积长度 （N2应当≥N+N-1，其中弃置区位于长度大于N+N-1的部分）
t2 = linspace( -T/2 , T/2 , N2);
Hf1 = fftshift(fft(ht,N2));  %频域匹配滤波器
Sf2 = fftshift(fft(st,N2));  %频域信号
S2 = Sf2 .* Hf1;             %频域乘积
s2 = ifft(S2);
%% 频域匹配滤波2（复制发射脉冲补零后计算FFT，再取共轭）
Hf2 = fftshift(conj(fft(st,N2)));
S3  = Sf2.* Hf2;
s3  = fftshift(ifft(S3));
%% 频域匹配滤波3（直接生成频域滤波器）
f_h = (-N2/2:N2/2-1)*Fs/N2;
Hf3 = exp(1i*pi*f_h.^2/K);       
S4 = Sf2.*Hf3;
s4 = ifft(S4);
%% 绘图
% 时域
figure,plot( t*1e6, real(st) ),xlabel('t /us'),ylabel('幅度'),title('Chirp信号 实部');
figure,plot( t*1e6, imag(st) ),xlabel('t /us'),ylabel('幅度'),title('Chirp信号 虚部');
figure,plot( t*1e6, f*1e-6 ),xlabel('t /us'),ylabel('频率 /MHz'),title('Chirp信号 频率');
figure,plot( t*1e6, phase ),xlabel('t /us'),ylabel('相位 /rad'),title('Chirp信号 相位');
% 频域
figure,plot( freq*1e-6,abs(Sf) ),xlabel('f /MHz'),ylabel('幅度谱'),title('Chirp信号 幅度谱');
figure,plot( freq*1e-6,-pi*freq.^2/K+pi/4 ),xlabel('f /MHz'),ylabel('相位谱'),title('Chirp信号 相位谱');
% 时域匹配滤波
figure,plot( t1*1e6 , abs(s1) ),xlabel('t /us'),ylabel('幅度谱'),title('时域匹配滤波');
% 频域匹配滤波1
figure,plot( t2*1e6 , abs(s2) ),xlabel('t /us'),ylabel('幅度谱'),title('频域匹配滤波1');
% 频域匹配滤波2
figure,plot( t2*1e6 , abs(s3) ),xlabel('t /us'),ylabel('幅度谱'),title('频域匹配滤波2');
% 频域匹配滤波3
figure,plot( t2*1e6 , abs(s4) ),xlabel('t /us'),ylabel('幅度谱'),title('频域匹配滤波3');