工具:matlabR2021b,vivado2018.3.
脉冲压缩的原理
脉冲压缩实际上就是对接收信号进行匹配滤波处理。根据发射的波形不同,脉冲压缩时选择不同的匹配滤波器系数。
数字脉冲压缩的实现方式有两种: 一是时域卷积法; 二是频域乘积法。依据傅里叶变换理论,时域卷积等效于频域乘法,因此两种方法实际上是等效的。时域处理方法直观、简单,但是运算量大,工程上一般采用频域法。
脉冲压缩的时域实现
时域脉冲压缩系统就是指在时域对信号进行卷积运算。
脉冲压缩的频域实现
脉冲压缩算法分别在频域上结论
频域脉冲压缩处理是指在频域实现时域线性卷积的过程,这个过程时域卷积和频域相乘时等效的。分别对回波信号和匹配滤波系数进行傅里叶变换后进行复乘运算。最后取傅里叶逆变换后得到结果。
脉冲压缩算法的matlab仿真
线性调频信号的产生
要获取较大的距离分辨率,可以增加发射信号的带宽。雷达的作用距离取决于信号的时宽,即要求信号要有大的时宽,而雷达的分辨精度取决于信号的带宽,即要求信号要有大的带宽。但是对于单载频脉冲信号,时宽和带宽的乘积近似等于1,所以同时得到大时宽和大带宽是矛盾的。
为了解决这一矛盾,必须采取同时具有大时宽和大带宽的复杂信号形式,最常用的就是线性调频信号(LFM),这种信号是在宽脉冲内附加载波线性调频,从而在大时宽的条件下扩展了带宽,通过脉冲压缩技术使宽脉冲变成窄脉冲,以获得高的距离分辨能力。
雷达的距离分辨率公式
其中c为光速,B为信号带宽。
关于线性调频信号
%%线性调频信号
T=10e-6; %脉冲持续时间10us
B=20e6; %线性调频信号的频带宽度20MHz
K=B/T; %调频斜率
Fs=5*B;Ts=1/Fs; %采样频率和采样间隔
N=512; %取点数,根据采样率取点数
t=linspace(-T/2,T/2,N); %持续时间范围-5us到5us
St=exp(j*pi*K*t.^2); %线性调频信号
时域和频域图谱
发射信号,创建了一个100us的时间长度,发射线性调频信号持续数据为10us。
% 生成 0 - 100us 的时间向量
t_extended = linspace(0, 100e-6, 4096 );
% 进行脉冲调制
pulsed_signal = zeros(1, length(t_extended));
pulsed_signal(1:length(St)) = St;
p_imag = imag(pulsed_signal);
接受回波信号,创建了一个100us的时间长度。分别在3KM,5KM和10KM模拟了目标,产生三个回波信号。并根据距离将回波的在100us回波的时间计算出来。
% 时间延迟,目标距离
R0 = 3e3; %目标距离
R1 = 5e3; %目标距离
R2 = 10e3; %目标距离
c = 3e8;
tr0 = 2*R0/c;
tr1 = 2*R1/c;
tr2 = 2*R2/c;
% 回波时间转化成起始点数
start_point0 = tr0/100e-6*4096;
start_point1 = tr1/100e-6*4096;
start_point2 = tr2/100e-6*4096;
start_point0 = ceil(start_point0);
start_point1 = ceil(start_point1);
start_point2 = ceil(start_point2);
%距离和时间的关系,光在这一来一回传输的时间
pulsed_signal2 = zeros(1, length(t_extended));
pulsed_signal2(start_point0+1:(start_point0+length(St))) = St;
pulsed_signal2(start_point1+1:(start_point1+length(St))) = St;
pulsed_signal2(start_point2+1:(start_point2+length(St))) = St;
进行频域脉冲压缩
对发射信号进行傅里叶变换后求共轭,得出匹配滤波器系数。求出回波信号的傅里叶变换后与匹配滤波器系数相乘。然后在进行逆傅里叶变换还原波形。在通过时间和距离的关系计算出目标距离。
%频域脉冲压缩算法
%求原始信号的傅里叶变换并算出匹配滤波器系数
Xs = fft(St,length(pulsed_signal)); % 原始信号的FFT
Xecho_0 = fft(pulsed_signal2,length(pulsed_signal2)); % 目标 1 回波信号的FFT
conj_s = conj(Xs);
Y_0 = conj(Xs).*Xecho_0; % 目标 1 的乘法器
% Y_0 = fftshift(Y_0);
y_0 = ifft(Y_0,length(pulsed_signal2)); % 目标 1 的 IFFT
r_0=t_extended*c/2; %回波时计算了来回的时间,这里将距离除以2获得准确的距离
得出脉冲压缩的结果为,可以看到我们设置了三个回波对应了三个目标,此时目标对应距离分别为3000m。5000m。10000m。
完整的matlab仿真代码
clc
clear
close all
%%线性调频信号
T=10e-6; %脉冲持续时间10us
B=20e6; %线性调频信号的频带宽度20MHz
K=B/T; %调频斜率
Fs=5*B;Ts=1/Fs; %采样频率和采样间隔
N=512; %取点数,根据采样率取点数
t=linspace(-T/2,T/2,N); %持续时间范围-5us到5us
St=exp(j*pi*K*t.^2); %线性调频信号
% 生成 0 - 100us 的时间向量
t_extended = linspace(0, 100e-6, 4096 );
% 进行脉冲调制
pulsed_signal = zeros(1, length(t_extended));
pulsed_signal(1:length(St)) = St;
p_imag = imag(pulsed_signal);
% 时间延迟,目标距离
R0 = 3e3; %目标距离
R1 = 5e3; %目标距离
R2 = 10e3; %目标距离
c = 3e8;
tr0 = 2*R0/c;
tr1 = 2*R1/c;
tr2 = 2*R2/c;
% 回波时间转化成起始点数
start_point0 = tr0/100e-6*4096;
start_point1 = tr1/100e-6*4096;
start_point2 = tr2/100e-6*4096;
start_point0 = ceil(start_point0);
start_point1 = ceil(start_point1);
start_point2 = ceil(start_point2);
%距离和时间的关系,光在这一来一回传输的时间
pulsed_signal2 = zeros(1, length(t_extended));
pulsed_signal2(start_point0+1:(start_point0+length(St))) = St;
pulsed_signal2(start_point1+1:(start_point1+length(St))) = St;
pulsed_signal2(start_point2+1:(start_point2+length(St))) = St;
% ========================================================
% 将回波信号存入coe
% 2. 提取实部和虚部
real_part = real(pulsed_signal2);
imag_part = imag(pulsed_signal2);
% 3. 将实部和虚部转换为 16 位整数
% 注意:需要将值缩放到合适的范围
real_int = int16(real_part * 2^14); % 进行缩放
imag_int = int16(imag_part * 2^14); % 进行缩放
% 4. 组合成 32 位整数
combined_int = bitshift(uint32(typecast(imag_int, 'uint16')), 16) + uint32(typecast(real_int, 'uint16'));
% 5. 转换为二进制字符串
hex_str = dec2hex(combined_int, 8);
num_points = length(hex_str);
% 保存到coe文件
fileID = fopen('pluse_echo.coe', 'w');
fprintf(fileID, 'memory_initialization_radix=16;\n');
fprintf(fileID, 'memory_initialization_vector=\n');
for i = 1:num_points
if i < num_points
fprintf(fileID, '%s,\n', hex_str(i,:));
else
fprintf(fileID, '%s;', hex_str(i,:));
end
end
fclose(fileID);
%频域脉冲压缩算法
%求原始信号的傅里叶变换并算出匹配滤波器系数
Xs = fft(St,length(pulsed_signal)); % 原始信号的FFT
Xecho_0 = fft(pulsed_signal2,length(pulsed_signal2)); % 目标 1 回波信号的FFT
conj_s = conj(Xs);
Y_0 = conj(Xs).*Xecho_0; % 目标 1 的乘法器
% Y_0 = fftshift(Y_0);
y_0 = ifft(Y_0,length(pulsed_signal2)); % 目标 1 的 IFFT
% ======================================================
% 画出原始信号傅里叶变换后的共轭存入coe
% 2. 提取实部和虚部
real_part = real(conj_s);
imag_part = imag(conj_s);
% 3. 将实部和虚部转换为 16 位整数
% 注意:需要将值缩放到合适的范围
real_int = int16(real_part * 100); % 进行缩放
imag_int = int16(imag_part * 100); % 进行缩放
% 4. 组合成 32 位整数
combined_int = bitshift(uint32(typecast(imag_int, 'uint16')), 16) + uint32(typecast(real_int, 'uint16'));
% 5. 转换为二进制字符串
hex_str = dec2hex(combined_int, 8);
num_points = length(hex_str);
% 保存到coe文件
fileID = fopen('conjxs.coe', 'w');
fprintf(fileID, 'memory_initialization_radix=16;\n');
fprintf(fileID, 'memory_initialization_vector=\n');
for i = 1:num_points
if i < num_points
fprintf(fileID, '%s,\n', hex_str(i,:));
else
fprintf(fileID, '%s;', hex_str(i,:));
end
end
fclose(fileID);
%%====================画图====================%%%
%线性调频信号
figure;
subplot(211)
plot(t*1e6,real(St));
xlabel('时间/us');
title('线性调频信号的实部');
grid on;axis tight;
subplot(212)
freq=linspace(-Fs/2,Fs/2,N);
plot(freq*1e-6,fftshift(abs(fft(St))));
xlabel('频率/MHz');
title('线性调频信号的幅频特性');
grid on;axis tight;
% 绘制时域图像
figure;
% subplot(211)
plot(t_extended*1e6, real(pulsed_signal));
xlabel('时间 (us)');
ylabel('幅度');
title('产生的线性调频信号时域图像');
% % 傅里叶变换求频谱
% spectrum = fftshift(fft(pulsed_signal));
% spectrum_abs = abs(spectrum);
% % 绘制频谱图
% subplot(212)
% %做FT变化,算平均功率谱,并画谱输出
% f = Fs * (0:length(spectrum) - 1) / length(spectrum);
% plot(f/1e6, abs(spectrum));
% xlabel('频率 (MHz)');
% ylabel('幅度');
% title('产生线性调频信号时域图像');
% 绘制时域图像
figure;
% subplot(211)
plot(t_extended*1e6, real(pulsed_signal2));
xlabel('时间 (us)');
ylabel('幅度');
title('回波线性调频信号时域图像');
% % 傅里叶变换求频谱
% spectrum_real = real(fft(pulsed_signal2));
% spectrum = fftshift(fft(pulsed_signal2));
% % 绘制频谱图
% subplot(212)
% f = Fs * (0:length(spectrum) - 1) / length(spectrum);
% plot(f/1e6, abs(spectrum));
% xlabel('频率 (MHz)');
% ylabel('幅度');
% title('回波线性调频信号的频谱图');
% ========= 脉冲压缩结果 ===========
figure;
subplot(211)
plot(t_extended*1e6, abs(y_0)/max(abs(y_0)));
xlabel('时间 (us)');
ylabel('幅度');
title('脉冲压缩后时域图像');
%求距离
subplot(212)
r_0=t_extended*c/2; %回波时计算了来回的时间,这里将距离除以2获得准确的距离
plot(r_0,(abs(y_0)/max(abs(y_0)))); % 目标 1 绘制 y 的幅度(以 dB 为单位)与距离 r 的关系
xlabel('距离/m');title('目标 1 脉冲压缩结果');
脉冲压缩算法的FPGA实现
为了在FPGA中仿真。将回波数据和匹配滤波器系数分别写成coe文件加载到FPGA的BRAM中。BRAM参数设置为。
逻辑实现的思路为。
回波数据从BRAM中读出。对回波数据进行FFT。完成后,从BRAM中读出匹配滤波器系数。通过回波的FFT和匹配滤波器的复乘完成后。对数据进行逆傅里叶变换,还原波形。
FFTIP设置为4096个点(结合matlab产生的波形数据)
回波数据进行FFT变换后,与匹配滤波器系数进行复乘运算。
可以看到上述的过程完成使用乘法器加法器即可完成。
FPGA测整体的仿真测试代码
`timescale 1ns / 1ps
//
// Company:
// Engineer:
//
// Create Date: 2024/08/15 15:45:55
// Design Name:
// Module Name: vtf_pluse
// Project Name:
// Target Devices:
// Tool Versions:
// Description:
//
// Dependencies:
//
// Revision:
// Revision 0.01 - File Created
// Additional Comments:
//
//
module vtf_pluse;
//bram echo
reg [11:0] addr_echo ;
wire [31:0] dout_echo ;
//system signal
reg clk ;
reg rst_n ;
reg data_gen ;
reg data_flog ;
reg douta_vld ;
reg [15:0] cnt ;
//fft
wire [31 : 0] s_axis_data_tdata ;
wire s_axis_data_tvalid ;
wire s_axis_data_tready ;
wire s_axis_data_tlast ;
wire [31 : 0] m_axis_data_tdata ;
wire [23 : 0] m_axis_data_tuser ;
wire m_axis_data_tvalid ;
wire m_axis_data_tready ;
wire m_axis_data_tlast ;
wire [15:0] echo_imag ;
wire [15:0] echo_real ;
wire [15:0] fft_imag ;
wire [15:0] fft_real ;
//输出脉冲压缩结果
wire [31:0] yk_data ; //
wire yk_valid ; //
wire [11:0] yk_cnt ; //
wire yk_ready ;
wire yk_tlast ;
wire [31 : 0] yk_data_tdata ;
wire [23 : 0] yk_data_tuser ;
wire yk_data_tvalid ;
wire yk_data_tready ;
wire yk_data_tlast ;
//=============================================
assign s_axis_data_tdata = dout_echo;
assign s_axis_data_tvalid = douta_vld;
assign s_axis_data_tlast = (cnt == 16'd4095);
assign m_axis_data_tready = 1'b1;
assign echo_imag = s_axis_data_tdata[31:16];
assign echo_real = s_axis_data_tdata[15:0];
assign fft_imag = m_axis_data_tdata[31:16];
assign fft_real = m_axis_data_tdata[15:0];
assign yk_tlast = {yk_cnt == 12'd4095};
assign yk_data_tready = 1'b1;
//======================================================
//回波读出
//======================================================
blk_mem_echo u_echo (
.clka (clk ), // input wire clka
.ena (1'b1 ), // input wire ena
.wea (1'b0 ), // input wire [0 : 0] wea
.addra (addr_echo ), // input wire [11 : 0] addra
.dina (32'd0 ), // input wire [31 : 0] dina
.douta (dout_echo )// output wire [31 : 0] douta
);
xfft_0 u_fft (
.aclk (clk ), // input wire aclk
.aresetn (rst_n ), // input wire aresetn
.s_axis_config_tdata (8'd1 ), // input wire [7 : 0] s_axis_config_tdata
.s_axis_config_tvalid (1'b1 ), // input wire s_axis_config_tvalid
.s_axis_config_tready ( ), // output wire s_axis_config_tready
.s_axis_data_tdata (s_axis_data_tdata ), // input wire [31 : 0] s_axis_data_tdata
.s_axis_data_tvalid (s_axis_data_tvalid ), // input wire s_axis_data_tvalid
.s_axis_data_tready (s_axis_data_tready ), // output wire s_axis_data_tready
.s_axis_data_tlast (s_axis_data_tlast ), // input wire s_axis_data_tlast
.m_axis_data_tdata (m_axis_data_tdata ), // output wire [31 : 0] m_axis_data_tdata
.m_axis_data_tuser (m_axis_data_tuser ), // output wire [23 : 0] m_axis_data_tuser
.m_axis_data_tvalid (m_axis_data_tvalid ), // output wire m_axis_data_tvalid
.m_axis_data_tready (m_axis_data_tready ), // input wire m_axis_data_tready
.m_axis_data_tlast (m_axis_data_tlast ), // output wire m_axis_data_tlast
.m_axis_status_tdata ( ), // output wire [7 : 0] m_axis_status_tdata
.m_axis_status_tvalid ( ), // output wire m_axis_status_tvalid
.m_axis_status_tready ( ), // input wire m_axis_status_tready
.event_frame_started ( ), // output wire event_frame_started
.event_tlast_unexpected ( ), // output wire event_tlast_unexpected
.event_tlast_missing ( ), // output wire event_tlast_missing
.event_status_channel_halt ( ), // output wire event_status_channel_halt
.event_data_in_channel_halt ( ), // output wire event_data_in_channel_halt
.event_data_out_channel_halt ( )// output wire event_data_out_channel_halt
);
//复乘
plural_mul u_plural_mul(
//输入sk回波FFT数据
.sk_data (m_axis_data_tdata[31:0] ), // (input ) (input )
.sk_cnt (m_axis_data_tuser[11:0] ), // (input ) (input )
.sk_valid (m_axis_data_tvalid ), // (input ) (input )
//输出脉冲压缩结果
.yk_data (yk_data[31:0] ), // (output) (output)
.yk_valid (yk_valid ), // (output) (output)
.yk_cnt (yk_cnt ), // (output) (output)
//
.clk (clk ), // (input ) (input )
.rst_n (rst_n ) // (input ) (input )
);
xfft_0 u_ifft (
.aclk (clk ), // input wire aclk
.aresetn (rst_n ), // input wire aresetn
.s_axis_config_tdata (8'd0 ), // input wire [7 : 0] s_axis_config_tdata
.s_axis_config_tvalid (1'b1 ), // input wire s_axis_config_tvalid
.s_axis_config_tready ( ), // output wire s_axis_config_tready
.s_axis_data_tdata (yk_data ), // input wire [31 : 0] s_axis_data_tdata
.s_axis_data_tvalid (yk_valid ), // input wire s_axis_data_tvalid
.s_axis_data_tready (yk_ready ), // output wire s_axis_data_tready
.s_axis_data_tlast (yk_tlast ), // input wire s_axis_data_tlast
.m_axis_data_tdata (yk_data_tdata ), // output wire [31 : 0] m_axis_data_tdata
.m_axis_data_tuser (yk_data_tuser ), // output wire [23 : 0] m_axis_data_tuser
.m_axis_data_tvalid (yk_data_tvalid ), // output wire m_axis_data_tvalid
.m_axis_data_tready (yk_data_tready ), // input wire m_axis_data_tready
.m_axis_data_tlast (yk_data_tlast ), // output wire m_axis_data_tlast
.m_axis_status_tdata ( ), // output wire [7 : 0] m_axis_status_tdata
.m_axis_status_tvalid ( ), // output wire m_axis_status_tvalid
.m_axis_status_tready ( ), // input wire m_axis_status_tready
.event_frame_started ( ), // output wire event_frame_started
.event_tlast_unexpected ( ), // output wire event_tlast_unexpected
.event_tlast_missing ( ), // output wire event_tlast_missing
.event_status_channel_halt ( ), // output wire event_status_channel_halt
.event_data_in_channel_halt ( ), // output wire event_data_in_channel_halt
.event_data_out_channel_halt ( )// output wire event_data_out_channel_halt
);
//================================================================
//================================================================
initial
begin
clk = 0;
rst_n=0;
data_gen =0;
#100;
rst_n =1;
#1000;
data_gen =1;
#100;
data_gen =0;
end
always @ (posedge clk or negedge rst_n)
begin
if(rst_n == 1'b0)begin
data_flog <= 1'b0;
end
else if(data_gen == 1'b1)begin
data_flog <= 1'b1;
end
else if(addr_echo >= 16'd4095)begin
data_flog <= 1'b0;
end
end
always @ (posedge clk or negedge rst_n)
begin
if(rst_n == 1'b0)begin
addr_echo <= 12'b0;
end
else if(data_flog == 1'b1)begin
addr_echo <= addr_echo+12'b1;
end
else begin
addr_echo <= 12'b0;
end
end
always @ (posedge clk or negedge rst_n)
begin
if(rst_n == 1'b0)begin
douta_vld <= 1'b0;
end
else begin
douta_vld <= data_flog;
end
end
always @ (posedge clk or negedge rst_n)
begin
if(rst_n == 1'b0)begin
cnt <= 16'd0;
end
else if(s_axis_data_tvalid == 1'b1 && s_axis_data_tready == 1'b1)begin
cnt <= cnt + 1'b1;
end
else begin
cnt <= cnt;
end
end
//================================================================
//================================================================
always #5 clk = ~clk;
endmodule
FPGA实现moldelsim仿真结果
可以看到脉冲压缩的结果高脉冲的点数在。
分别在820,1366,2731位置。计算结果和matlab结果一致。分别代表3KM,5KM以及10KM。