【三】FMCW雷达距离-多普勒图的计算

原理

为什么不能直接通过脉冲（慢时间）之间物体距离的变化直接得到速度呢？这是因为目标在一个脉冲时间内移动的距离通常小于雷达的距离分辨率。也就是说在 Range-FFT 的结果中，目标的峰值会一直停留在同一个距离索引中，肉眼看不出它在移动。

虽然距离峰值在坐标轴上没有跳变，但目标微小的位移会引起反射波相位的剧烈变化。既然相位在旋转，旋转的快慢就代表了速度的大小。根据物理定义，频率是相位随时间的变化率 $f = \frac{1}{2\pi} \cdot \frac{d\phi}{dt}$ 。而目标运动产生的相位旋转频率就是多普勒频率。

如果场景中有多个目标（比如两个人在不同速度行走，或者一个人在呼吸的同时还在挥手），在同一个距离索引里会叠加多个不同旋转频率的相位信号。慢时间 FFT 能够将这些重叠在一起的旋转相位分离出来。每一个多普勒频率峰值就对应一个独立的速度。

一般流程是：先沿着快时间轴做FFT运算，区分不同距离的目标；再沿着慢时间轴做FFT运算，区分相同距离不同速度的目标。

移除杂波和静止背景

data = mat - np.average(mat)

i_ant = 0 # 第一个天线
num_ant = 1 # 天线数量
mti_alpha = 0.8
mti_history = np.zeros((config.num_chirps, config.chirp.num_samples, num_ant))

data_mti = data - mti_history[:, :, i_ant]
mti_history[:, :, i_ant] = data * mti_alpha + mti_history[:, :, i_ant] * (1 - mti_alpha)

杂波就是移除接收信号的平均值。

移除静止背景，就需要采集历史数据（长期不动的背景），然后将当前数据减去历史数据。

历史数据的更新是将新收集的数据一定程度上与历史数据进行混合，得到更新的历史数据。

mti_alpha是MTI滤波器的平滑系数，取值范围为0到1，较大的值会更快地适应变化，但可能会引入更多的噪声，较小的值会更慢地适应变化，但可以更好地抑制噪声。

2D FFT得到距离-多普勒图（RD图）

构造2D-FFT函数。

def fft_2d_spectrum(mat, range_window, doppler_window):
    [num_chirps, num_samples] = np.shape(mat)

    avgs = np.average(mat, 1).reshape(num_chirps, 1)
    mat = mat - avgs

    # prepare for range FFT
    mat = np.multiply(mat, range_window)
    zp1 = np.pad(mat, ((0, 0), (0, num_samples)), 'constant')
    range_fft = np.fft.fft(zp1, axis=1) / num_samples
    range_fft = 2 * range_fft[:, range(int(num_samples))]


    # prepare for doppler FFT
    range_fft_windowed = np.multiply(range_fft, doppler_window)
    zp2 = np.pad(range_fft_windowed, ((0, num_chirps), (0, 0)), 'constant')
    doppler_fft = np.fft.fft(zp2, axis=0) / num_chirps
    doppler_fft = np.fft.fftshift(doppler_fft, (0,)) # 将FFT结果移到中心点
    return doppler_fft

步骤：

再取一次平均，抑制平均值。
np.multiply(mat, range_window)是在做FFT之前进行加窗。
首先对快时间域做FFT变换（axis是指”沿着哪个轴做变换“，axis=1表示沿着列做FFT变换，即对每一行做FFT变换），与之前距离估计类似，也是补零平滑FFT图、取前半段内容，补偿幅度。
其次对慢时间域做FFT变换，同样需要补零平滑FFT图。但需要注意，经过距离变换后range_fft不再是实数，对复数序列的 FFT，频谱不再满足共轭镜像对称，正负多普勒分别代表接近/远离的速度方向，所以必须保留全谱并用fftshift把零速度移到中间。

坐标转换

根据距离计算那篇文章，我们知道

\Delta d\text{（单格分辨率）}= \frac{c}{4B_s}\text{(单位：m)}

而对于速度轴的转换，公式如下：

\Delta v=\frac{c}{2\text{PRT}\cdot f_c\cdot N_\text{FFT}}=\frac{\lambda}{2 \text{PRF}\cdot N_\text{FFT}}=\frac{\lambda\cdot \text{PRF}}{2N_\text{FFT}}​

其中 $N_\text{FFT}$ 表示沿着行做FFT的点数，也就是行数。因此完成坐标转换

bandwidth_hz = abs(config.chirp.end_frequency_Hz - config.chirp.start_frequency_Hz)
range_bin_m = 3e8 / (4 * bandwidth_hz)
range_bins_m = np.arange(config.chirp.num_samples) * range_bin_m

fc_hz = 0.5 * (config.chirp.start_frequency_Hz + config.chirp.end_frequency_Hz)
wavelength_m = 3e8 / fc_hz
prf_hz = 1.0 / config.chirp_repetition_time_s
doppler_bins = doppler_fft.shape[0]
doppler_bin_hz = prf_hz / doppler_bins
velocity_bin_mps = 0.5 * wavelength_m * doppler_bin_hz
velocity_bins_mps = (np.arange(doppler_bins) - doppler_bins / 2) * velocity_bin_mps

doppler_mag = 20 * np.log10(np.abs(doppler_fft) + 1e-6) # 将单位转换成dB

最后画图，即可得到下图效果：

rdmap

注：要得到更干净的图，考虑调整set_clim中的vmin和vmax。

参考文献

基于英飞凌官方例程中DopplerAlgo.py的代码。

← Back to blog

Table of contents

【三】FMCW雷达距离-多普勒图的计算

原理

移除杂波和静止背景

2D FFT得到距离-多普勒图（RD图）

坐标转换

参考文献