LFM信号加噪、时频分析、滤波

前言

线性调频（Linear Frequency Modulation,LFM）信号具有很大的时宽带宽积，可获得很大的脉冲压缩比，是雷达系统和声呐系统广泛采用的一种信号形式。
LFM信号的数学表达式为：

其中，A为信号幅度，t为时间，T为脉冲持续时间（脉冲宽度），fc为载波频率，K为信号的线性调频率，rect()为矩形窗函数，数学表达式如下：

假设信号带宽B=10MHz，脉冲宽度T=100us，载波频率fc=100MHz，采样频率为fs=2B=20MHz。线性调频率K=B/T.

设计要求如下：
1.用Matlab对该参数的LFM信号进行时频分析，得出结论；
2.对信号噪声干扰（单一频率噪声或高斯白噪声），对加噪后的信号时频分析；
3.设计合适的滤波器，对含噪的LFM信号进行滤波，比较滤波前后的效果。

1.LFM信号的时频分析

线性调频（LFM）信号是指瞬时频率随时间成线性变化的信号，也称为Chirp信号。表达式中han’you含有时间t的项对t求导得到K*t，即瞬时频率。在Matlab中，生成了2000点长的LFM信号，绘制实部、虚部、相位、瞬时频率如下图：

LFM信号的时域频域波形如下：

2.加入噪声干扰

为了比较不同种类的噪声对信号的影响，分别加入频率为6MHz的噪声和高斯白噪声（加噪后SNR=10dB）。
加入单一频率噪声的时域图频谱图如下：

加入高斯白噪声的时域图、频谱图如下：

3.含噪LFM信号滤波

3.1 滤除单一频率噪声

利用Matlab里的Filter Design工具箱，设计一个IIRLPF，其低通频率为5500Hz，截止频率为5700Hz，阻带衰减为40dB。

将含有单一频率的信号通过该滤波器，观察输出的时域图和频域图。6000Hz频率噪声明显滤除了，但是信号存在部分失真。

3.2 滤除信号中的高斯白噪声

由于高斯白噪声的频率在LFM信号的任何频率都覆盖，所以不能用滤波器的方法滤波。比较常用的方法是自适应信号滤波，例如LMS滤波、RLS滤波等。原理在此不过多赘述。本文采用RLS滤波方法。
在Matlab编写的RLS滤波程序如下：

function [e,w]=my_rls(lambda,M,u,d,delta)
% recursive least squares,rls.
% Call:
% [e,w]=rls(lambda,M,u,d,delta)
%
% Input arguments:
% lambda = constant, (0,1]
% M = filter length, dim 1x1
% u = input signal, dim Nx1
% d = desired signal, dim Nx1
% delta = constant for initializaton, suggest 1e-7.
%
% Output arguments:
% e = estimation error, dim Nx1
% w = final filter coefficients, dim Mx1
% Step1:initialize
w=zeros(M,1);
P=eye(M)/delta;
u=u(:);
d=d(:);
% input signal length
N=length(u);
% error vector
e=d.';
% Step2: Loop, RLS
for n=M:N
    uvec=u(n:-1:n-M+1);
    e(n)=d(n)-w'*uvec;
    k=lambda^(-1)*P*uvec/(1+lambda^(-1)*uvec'*P*uvec);
    P=lambda^(-1)*P-lambda^(-1)*k*uvec'*P;
    w=w+k*conj(e(n));
end

滤波前、后的LFM信号时域图如下：

滤波后的LFM信号频谱如下：

（仿真分析留给各位大佬来完成！）

代码

clear;
clc;
close all;
%% 参数设置
B=10e6;    %信号带宽10MHz
tao=10e-5;   %脉宽100us
fs=2*B;    %采样频率
Ts=1/fs;   %采样周期
K=B/tao;     %线性调频率
fc =1e8;   %载波频率
%% 产生LFM信号
t=-tao/2:1/fs:tao/2-1/fs;
LFM=exp(j*2*pi*(fc*t+0.5*K*t.^2));  %发射信号
theta =  pi*K*t.^2;      %信号弧度
ft = K*t;                 %信号频率

figure();
subplot(2,2,1);
plot(real(LFM));
xlabel('信号点数n');ylabel('幅度');
title('LFM信号实部');
subplot(2,2,2);
plot(imag(LFM));
xlabel('信号点数n');ylabel('幅度');
title('LFM信号虚部');
subplot(2,2,3);
plot(theta);
xlabel('信号点数n');
title('信号相位');
subplot(2,2,4);
plot(ft);
xlabel('信号点数n');
title('信号频率 Hz');

X=fftshift(fft(LFM));
f=linspace(0,fs,length(t))-fs/2;        %设置频率变量
figure();
subplot(211);
plot(t,real(LFM));
xlabel('时间/t');ylabel('幅度');
title('LFM信号时域');
subplot(212);
plot(f,abs(X));
xlabel('频率/f');ylabel('幅度');
title('LFM信号频谱');

%% 加入单一频率噪声
fn1=6e6;        %噪声频率6MHz
n1=0.1*cos(2*pi*fn1*t);
S1=LFM+n1;
figure();
subplot(211);
plot(t,S1);
xlabel('时间/t');ylabel('幅度');
title('加噪后信号时域');
subplot(212);
plot(f,abs(fftshift(fft(S1))));
xlabel('频率/f');ylabel('幅度');
title('加噪后信号频谱');

%% 设计带通滤波器
y=filter(IIRBPF2,S1);       %用Matlab自带的filter design工具箱，设计低通滤波器,代码见IIRBPF2.m
[b,a]=tf(IIRBPF2);        %把低通滤波器IIRBPF2转化成传输函数,系数为b,a
figure();
freqz(b,a);     %作IIRBPF2 幅频特性曲线图

figure();
subplot(211);
plot(t,y);
xlabel('时间/t');ylabel('幅度');
title('去噪后信号时域图');
subplot(212)
plot(f,abs(fftshift(fft(y))));  %去噪后的频谱
xlabel('频率/f');ylabel('幅度');
title('去噪后信号频谱');
%% 加高斯白噪声
S2=awgn(LFM,10,'measured');    %加入高斯白噪声，信噪比为10dB
sn=S2-LFM;                     %噪声

figure();
subplot(211);
plot(t,S2);
xlabel('时间/t');ylabel('幅度');
title('LFM+高斯白噪声的时域波形');
subplot(212);
plot(f,abs(fftshift(fft(S2))));
xlabel('时间/t');ylabel('幅度');
title('LFM+高斯白噪声频谱');

%% RLS算法滤波
mu2=0.002;
M=2;
espon=1e-5;
delta=1e-7;
lambda=0.99;
[en,w2]=my_rls(lambda,M,sn,S2,delta);%RLS算法子程序
er=en-LFM;            %er为误差信号，滤波器输出-输入
figure();
subplot(311);
plot(t,S2);
xlabel('时间/t');ylabel('幅度');
title('LFM+高斯白噪声的时域波形');
subplot(312);
plot(t,en);
xlabel('时间/t');ylabel('幅度');
title('滤波后信号时域波形');
subplot(313);
plot(t,er);
xlabel('时间/t');ylabel('幅度');
title('误差信号');

figure;
plot(f,abs(fftshift(fft(en))));
xlabel('频率/f');ylabel('幅度');
title('滤除高斯白噪声后信号频谱');