调频连续波(FMCW)雷达目标轨迹提取

大白菜菜籽

已于 2022-10-16 11:25:33 修改

阅读量5.8k

点赞数 22

分类专栏：雷达 FMCW 文章标签： matlab 其他

于 2022-04-27 17:32:12 首次发布

本文链接：https://blog.csdn.net/weixin_42686221/article/details/124451239

版权

雷达同时被 2 个专栏收录

10 篇文章 59 订阅

订阅专栏

FMCW

2 篇文章 3 订阅

订阅专栏

调频连续波雷达目标轨迹提取

本篇文章内容主要是基于调频连续波雷达采集的人体运动回波，提取回波中的距离、角度信息，并对目标运动轨迹进行绘制。
主要涉及到的信号处理的相关知识有：脉冲压缩、加窗、MTI/MTD、多天线测角（DOA估计）等

调频连续波雷达简介

调频连续波雷达是一种发射的信号为调频连续波的雷达。调频连续波信号就是频率随时间变化，而且占空比为1的信号。
但是一直发射信号也要考虑一些其他的因素，所以有些调频连续波雷达还是可以设置成脉冲的形式，比如TI的XWR系列。
我们只需要关心雷达采集到的三维回波矩阵（快时间、慢时间、通道）就可以了，所以对于这个问题就不再这里深究了。
调频连续波的两种主要形式：三角波、锯齿波。

下面给出三角波示意图（图片来自百度）

图片来自百度
关于三角波与锯齿波在信号处理中的差距，对本文的目的影响不大，此处便不再赘述，有感兴趣者可自查。
该链接可参考：https://www.zhihu.com/question/401476711

发射信号参数配置

%% 雷达为24GHz的三角波调频连续波雷达，天线为一发两收。
C          = 3e8;                 %光速              m/s
Fc         = 24.125e9;            %中心频率          Hz
B          = 250e6;               %信号带宽          Hz
Tr         = 2e-3;                %重复周期          s
PRF        = 1/Tr;                %重复频率          Hz
Fs         = 500e3;               %快时间维采样频率  Hz
Ts         = 1/Fs;                %快时间维采样周期  s 
D          = 0.0145;              %天线间隔

回波信号形式

回波矩阵

本文的雷达只有两个接收天线，所以通道维的个数=2，分为channelA 和 channelB。
快时间维是每个信号每个周期的采样点数，采样点数=采样时间 X 采样频率=1000。
慢时间就是发射的脉冲个数（可以自己设定，根据每帧多少个脉冲，发射多少帧这样，不过不重要），本文脉冲数为	13838。

回波处理

 在得到回波矩阵之后就可以，对数据进行处理了。首先，为了方便处理，首先将三角波矩阵转变成了锯齿波。
（很明显只需要取快时间前500就可以了，但是这个时候，相当于采样率没变，采样点数少了一半，在计算距离分辨率的时候，我们再进行分析）

距离向脉冲压缩

为了得到目标距离信息，需要进行脉冲压缩。脉冲压缩主要有两种方法：匹配滤波 和 Dechirp+FFT。
前者针对LFM信号，后者针对混频后的中频信号。本文调频连续波雷达采用混频得到中频信号属于第二种形式。
所以只需要对回波快时间维进行FFT即可得到回波的一维距离像。
考虑到有环境背景噪声的影响，本文采取了加窗操作。

下面展示加窗前后距离像对比。为了提高频率分辨率，快时间维进行了2048点的FFT，并对回波位置进行了截取。
加窗前后对比

下面给出距离向处理代码

echoes=load('~');
channelA=echoes.echoes.channelA;
channelB=echoes.echoes.channelB;
channelB=imag(channelB)+1j*real(channelB);   %%由于雷达数据格式问题，B通道实部虚部需要互调
[Nr,Na]=size(channelB);
%三角波转锯齿波
channelAA=channelA(1:500,:);
channelBB=channelB(1:500,:);
clear channelB channelA

%%距离向FFT
w=hamming(500)*ones(1,13838);  %% 窗函数。
figure
channelAA_fw=fftshift(fft(channelAA.*w,2048,1),1);
dataA_FFT=channelAA_fw(minbin_range:maxbin_range,:);
channelBB_fw=fftshift(fft(channelBB.*w,2048,1),1);
dataB_FFT=channelBB_fw(minbin_range:maxbin_range,:);
imagesc(abs(dataA_FFT))
title("通道A加窗后")

得到距离信息之后需要计算实际距离，其实只需要知道距离分辨率就可以了。
在采样率为 ${F_s}$ 的情况下进行 ${N}$ 点的采样，可得频率分辨率为
$Δf=\frac{F_s}{N}$ 根据调频率 ${K=B/\frac{T_r}{2}}$ ，就可以计算时间分辨率 ${Δt=\frac{Δf}{K}}$ ，进而可以计算距离分辨率 ${ΔR=\frac{CΔt}{2}=\frac{C\frac{\frac{F_s}{N}}{B/\frac{T_r}{2}}}{2}=\frac{CT_rF_s}{4BN}}$ 需要注意：由于是三角波，所以调频率的计算只用半个周期。

MTI/MTD

MTI本质就是个频域滤波器：可以抑制固定目标和慢速杂波。进行双脉冲对消，直接慢时间维差分就可以了。
MTD是一个带通滤波器组，一般使用慢时间FFT来实现，也就意味着在多普勒域对目标进行锁定，提高检测性能

下面给出示例代码

%%MTI/MTD 脉冲对消
%两脉冲对消
figure
dataA_FFT_mti=diff(dataA_FFT,1,2);
dataB_FFT_mti=diff(dataB_FFT,1,2);
%%MTD处理：方位向FFT，然后在距离-多普勒图中搜索最大值
figure
%方位向分帧数，帧数决定方位FFT个数，影响速度分辨率，这里Nfft_A取了64
framenum=floor((Na-1)/Nfft_A);
dataA_MTD=zeros(size(dataA_FFT_mti,1),framenum);
dataB_MTD=zeros(size(dataA_FFT_mti,1),framenum);
for k=1:framenum
    dataA_mtd=fftshift(fft(dataA_FFT_mti(:,(k-1)*Nfft_A+1:k*Nfft_A),Nfft_A,2),2);
    dataB_mtd=fftshift(fft(dataB_FFT_mti(:,(k-1)*Nfft_A+1:k*Nfft_A),Nfft_A,2),2);
    dataA_RV=dataA_mtd(:,minbin_v:maxbin_v);
    dataB_RV=dataB_mtd(:,minbin_v:maxbin_v);
    [maxrow,maxcol]=find(abs(dataA_RV)==max(max(abs(dataA_RV))));
    dataA_MTD(:,k)=dataA_RV(:,maxcol);
    dataB_MTD(:,k)=dataB_RV(:,maxcol);
end
imagesc(abs(dataA_MTD))

下面给出MTI/MTD处理过程图
在这里插入图片描述
结果显示MTI+MTD后得到的距离像已经较为理想，下一步可以进行角度提取。

角度提取

对角度信息的提取主要有天线维比相：长短基线法、天线维FFT，超分辨DOA估计：MUSIC、ESPRIT，等具体理论知识自查哈...
由于本文只针对单目标，且只有两个平面天线，所以可以简单的进行双天线比相，得到平面角。

双通道比相具体步骤：
获得channelA距离向最大值位置，将该位置上对应的channelA和channelB的数据进行共轭相乘后取相角：
$\begin{array}{l} \Delta \varphi {\rm{ = }}angle(channel{A_{\max }}.*channelB_{\max }^*)\end{array}$ 得到的相角即为两个通道的相位角度差，根据几何关系可以得到来波方向： $\begin{array}{l}\theta = \arcsin (\frac{{\lambda \Delta \varphi }}{{2\pi d}}) \end{array}$ 下面给出代码示例

%%比相测角
target_theta=[];
target_range=[];
dataA_MTD_dB=20*log10(abs(dataA_MTD)/max(max(abs(dataA_MTD))));
for k=1:framenum
    [maxrow,maxValue]=find(abs(dataA_MTD(:,k))==max(abs(dataA_MTD(:,k))));
    if dataA_MTD_dB(maxrow,k)>-25
        target_range=[target_range rangRes*maxrow];
        target_theta=[target_theta angleEstimate(dataA_MTD(maxrow,k),dataB_MTD(maxrow,k))];
    end
end
%%%%   angleEstimate函数   %%%%%
%%   phi=angle(a.*conj(b));
%%   theta=asin(lamba* phi/2/pi/D)/pi*180;
%平滑滤波
target_theta=smooth(target_theta);
target_range=smooth(target_range);

轨迹绘制

得到目标距离、角度信息之后就可以绘制目标轨迹了，下面分别给出绘制动态极坐标以及直角坐标图，并保存为了gif。

%图像显示
filename='track.gif';
h=figure;
for i=1:length(target_range)
    polarplot(target_theta(i)*pi/180,target_range(i),'bo')
    thetalim([-90,90]);
    hold on
%     pause(0.001)
    drawnow 
      % Capture the plot as an image 
      frame = getframe(h); 
      im = frame2im(frame); 
      [imind,cm] = rgb2ind(im,256); 
     if i == 1 
          imwrite(imind,cm,filename,'gif', 'Loopcount',inf,'DelayTime',0.05); 
      else 
          imwrite(imind,cm,filename,'gif','WriteMode','append','DelayTime',0.05); 
     end
end
%直角坐标系
X=target_range.*sind(target_theta);
Y=target_range.*cosd(target_theta);
b=figure;
filename='track_xy.gif';
for i=1:length(target_range)
    plot(X(i),Y(i),'bo')
    xlabel('水平距离（m）')
    ylabel('垂直距离（m）')
    axis([-4 4 0 8])
    hold on
%     pause(0.001)
        drawnow 
      % Capture the plot as an image 
      frame = getframe(b); 
      im = frame2im(frame); 
      [imind,cm] = rgb2ind(im,256); 
     if i == 1 
          imwrite(imind,cm,filename,'gif', 'Loopcount',inf,'DelayTime',0.05); 
      else 
          imwrite(imind,cm,filename,'gif','WriteMode','append','DelayTime',0.05); 
     end
end

结果展示

本文到这里就结束了一个雷达目标估计提取的小实验处理的全部流程代码都已经给出可以供大家学习参考！😄

大白菜菜籽

关注

22
点赞
踩
116

收藏

觉得还不错? 一键收藏
15
评论
调频连续波(FMCW)雷达目标轨迹提取

调频连续波雷达目标轨迹提取本篇文章内容主要是基于调频连续波雷达采集的人体运动回波，提取回波中的距离、角度信息，并对目标运动轨迹进行绘制。主要涉及到的信号处理的相关知识有：脉冲压缩、加窗、MTI/MTD、多天线测角（DOA估计）等调频连续波雷达简介调频连续波雷达是一种发射的信号为调频连续波的雷达。调频连续波信号就是频率随时间变化，而且占空比为1的信号。但是一直发射信号也要考虑一些其他的因素，所以有些调频连续波雷达还是可以设置成脉冲的形式，比如TI的XWR系列。我们只需要关心雷达采集到的三维回波矩
复制链接

扫一扫