（转）基于毫米波雷达的舱内乘员检测原理

https://zhuanlan.zhihu.com/p/272797271?utm_source=wechat_session

1.简介

随着家用车的普及，汽车走进越来越多的家庭，其中儿童被困车内的事故时有发生，甚至造成死亡，下图是由于儿童被困车内致死事故的统计。

儿童被困车内致死人数

为此，ENCAP Roadmap 2025提出，将儿童舱内检测技术作为汽车的安全性评定的加分项，后续将有法规逐步实施。

ENCAP Roadmap 2025

随着安全意识的提高，不少车型开始搭载舱内乘员监测系统。通常采用毫米波雷达的方案，将其布置与B柱上部或者顶棚，检测舱内的生命体。

毫米波雷达方案优势对比

本人从事相关领域，以TI的Demo 为例，将其原理及实现进行简单的介绍。

2.舱内检测原理

舱内检测的原理主要包括以下步骤：

前端（Front end）：依据天线和波形等配置，实现毫米波雷达FMCW的信号发送和接收；
模数转换（ADC）：对中频信号进行采样，并将 ADC输出样本缓冲在ADC输出缓冲器中，以供DSP部分访问；
DSP :舱内检测的算法主要在DSP中运行，该算法主要包括Rang FFT、Clutter Removal、Heatmap Generation、Feature Extraction、Zone Detect 等步骤，下文中将详细介绍
ARM ：该ARM MCU可执行应用程序代码，包括进一步的信号处理操作和其他更高级别的功能。

舱内乘员检测的基本原理图

2.1 距离FFT Range Processing

该部分主要将ADC采样后的中频信号数据进行1维加窗和1D FFT,获得目标的距离信息。同时，该过程在有效线性调频时间和帧结束之间的空闲时间（“帧间”周期）内进行。

2.2 去除杂波 Clutter Removal

在杂乱的环境中，尤其是在室内，检测具有较小RCS（雷达横截面）的物体（例如行人和生命形式的物体）具有较大的难度。一些具有较大RCS的杂波（例如建筑结构和家具）可能会主导接收到的信号。因此，具有小RCS的物体可能在杂波的强烈干扰下被掩盖，并且变得难以检测。当雷达传感器处于静止状态时，杂波和所有静止物体的零多普勒频率将被用于改善其余雷达信号处理的性能。

2.3 生成热力图 Range-Azimuth Heatmap Generation

每帧使用去除杂波的1D FFT输出数据，结合基于协方差的高分辨率测角方法DOA，生成“距离-方位角”的热力图。为避免数值不稳定性问题，沿协方差估计矩阵的对角线添加一个较小的修正量。“热力图”实际上是一个速度图，其幅值越大，表示非静态对象的运动越多。

DoA角度估计

$X_{n,p,k}$ :表示经过Range Processing和Clutter Removal两个步骤后的输出， $n,p,k$ 分别表示range bin index， chirp index， virtual receive antenna index。

对于第 $n$ 个range bin，空间协方差的估计值为：

$hat{R_{n}}=frac{1}{N_c}sum_{k=1}^{N_c}X_{n,k}X_{n,k}^{H}$

其中： $X_{n,k}$ 为第n个range bin和第k个chirp的空间向量，它是由数据立方体通过在虚拟天线上堆叠样本而形成的，N_c表示天线个数。

例如：

假设： $X_{n,k}=[X_{n,k,1},X_{n,k,2},...,X_{n,k,N_{a}}]^T$

为了减少计算量，可以使用小于 $N_{c}$ 的值 $N_{a}$ 来估计空间协方差。
当入射雷达射频（RF）信号来自方位角 $heta$ 时，转向矢量 $alpha( heta)$ 表示每个虚拟接收天线上的相位差。

对于线性一维虚拟接收机的天线阵列，定义方位角 $heta$ 的转向矢量 $alpha( heta)$ 为：

$alpha( heta)=[1,exp(j2pi dsin( heta)),...,exp(j2pi (N_a-1)dsin( heta))]$

其中，d是归一化后的天线间距。

第n个rangebin的角度谱可由下式计算：

$P_{n}( heta)=frac{1}{alpha( heta)^Hhat{R_{n}}^{-1}alpha( heta)}$

相应的波束成形向量为：

$W_{n}( heta)=frac{hat{R_{n}}^{-1}alpha( heta)}{alpha( heta)^Hhat{R_{n}}^{-1}alpha( heta)}$

将每个rangebin的DoA频谱表示为矩阵形式，其中第n行表示为第n个rangebin的DoA频谱，称为“距离-方位热力图”

——————————————————————————————————————

Diagonally-Loading Method

在某些情况下，协方差估计是奇异的或接近奇异的，这可能会在进行矩阵求逆时引起数值稳定性问题。

通常沿协方差估计矩阵的对角线添加一个修正量来避免数值稳定性问题

$P_{n}( heta)=frac{1}{alpha( heta)^H[hat{R_{n}}+alpha_nI]^{-1}alpha( heta)}$

可以使用如下的噪声方差估计来确定对角因子 $alpha_n$ :

$alpha_n=eta overline P_{n}$

对于常数 $eta$ ，其中 $P_{n}=frac{1}{N_alpha}trace{R_n}$ 为协方差估计矩阵对角项的平均值。

2.4 热力图调整 Range-Azimuth heatmap conditioning

为去除杂波和DoA问题，需要对热图进行进一步处理，以抑制ARC，并平滑因杂波去除中帧变化而导致的掉帧。

Arc Removal

for each range row:
	if (range_row_min_val > Nr)
		row_range = range_row_max_val - range_row_min_val
		scale_factor = 1.0 + ((range_row_min_val - Nr) / row_range
		
	for each cell in the range row:
	//Translate the row to zero and rescale
	rescaled_val = (cell - range_row_min_val) * scale_factor
	//Translate the scaled row to the noise floor
	cell = rescaled_val + Nr

Heatmap Smoothing

for each cell in the heatmap:
	cell_thresh = previous_cell * threshold
	
	if (current_cell < cell_thresh)
	current_cell = cell_thresh

对热力图进行去底噪、平滑之后，即可进行下一步的操作。

2.5 特征提取 Feature Extraction

Notations

- $N_z$ : 位置的个数；

- $r_n$ : 第 $n$ 个range bin所对应的距离值, $n=1,2,..,N_r$ ；

- $heta_m$ :第 $m$ 个方位角索引对应的角度值， $m=1,2,..,N_{az}$ ；

- $S[t]$ ：第 $t$ 帧处 $N_r×N_{az}$ 的热力图矩阵；

——————————————————————————————————————

定义区域 Definition of Zones

在范围和方位角域中定义区域，其取值与雷达传感器的安装位置和姿态有关。可将每个乘客座椅的位置定义为目标区域。其定义方式如下：

$Z_i={(n,m):r_l^j≤r_n<r_u^j, heta_l^i≤ heta_m< heta_u^i}$

- $Z_i$ 是第 $i$ 个区域（Range-Azimuth Heatmap中的矩形边界）内所距离-方位角网格的集合，其由四个参数 $r_l^j、r_u^j、 heta_l^i、 heta_u^i$ 定义。

——————————————————————————————————————

特征 Features

为区分不同目标区域是否有乘员，需要利用不同的区域特征进行判断。通常目标区域是有乘员对于该区域的功率能量影响较大。

平均功率能量 the average zone-power ，表征第 $t$ 帧的某个区域能量，对于第 $i$ 个区域，第 $t$ 帧时，该区域的平均能量为：
$Q_i[t]=frac{1}{|Z_i|}sum S_{n,m}[t]$

通过平均功率能量这一指标可推导出三个具有代表意义的特征：移动平均功率能量、平均功率比、区域功率相关系数。

移动平均功率能量Moving-averaged zone power ，表征前 $L$ 帧的平均功率，对于第 $i$ 个区域，第 $t$ 帧时，其移动平均功率能量为：

$overline Q_i[t]=frac{1}{L}sum Q_{i}[l]$

其中， $l=t-L+1$ ， $L$ 是移动平均值窗的长度。

2.平均功率比 Moving-averaged power ratio ,表征某个区域能量占总能量的比值。对于第 $i$ 个区域，第 $t$ 帧时，其平均功率比为：

$q_i[t]=frac{ overline Q_i[t]}{sum_{i=1}^{N_z} overline Q_i[t]}$

3.区域功率相关系数 Correlation coefficient of zone power,表示不同目标区域间的相关性,对于区域 $i$ 和区域 $j$ ,其相关系数为：

$ho_{i,j}[t]=frac{1}{L-1} sum_{l=t-L+1}^{t}frac{Q_i[l]-overline Q_i[l]}{sigma_i[t]}$

对于三个目标位置的检测场景，即 $N_z=3$ 时，其特征向量可表示为：

$X[t]=[overline Q_1[t],overline Q_2[t],overline Q_3[t],q_1[t],q_2[t],q_3[t], ho_{1,2}[t], ho_{2,3}[t], ho_{3,1}[t]]$

2.6 区域检测算法 Zone Detection Algorithms

对于 $N_z$ 个目标区域，其可能的乘员分布有 $2^{N_z}$ 种，将特征向量和离线生成的训练系数相乘，将得到一系列的决策结果，其中“1”代表有人，“0”表示空。

对于 $N_z=3$ ，乘员分布的所有状态如下：

$(O_1[t],O_2[t],O_3[t])in{(0,0,0),(1,0,0),(0,1,0),(0,0,1),(1,1,0),(1,0,1),(0,1,1),(1,1,1)}$

对于每种可能的乘员分布状态，定义 $N+1$ 维的决策权值向量：

$W_{(o_1,o_2,...,o_{N_z})}=[W_{(o_1,o_2,...,o_{N_z})}^0,W_{(o_1,o_2,...,o_{N_z})}^1,...W_{(o_1,o_2,...,o_{N_z})}^N]^T$

其中， $N$ 是特征向量 $X[t]$ 的元素个数。

***

目标空间上有乘员的概率定义为：

$hat p_{(o_1,o_2,...,o_{N_z})}[t]=g(W_{(o_1,o_2,...,o_{N_z}）}^{ ilde x[t]})=gsum_{i=0}^{N}W_{(o_1,o_2,...,o_{N_z}）}i^{x[t]}$

其中： $ilde X[t]=[1;X[t]]=[X_0[t],X_1[t],...,X_N[t]]^T$ ，

$g(cdot)$ 为非线性函数，逻辑回归函数, $g(z)=frac{1}{1+e_{-z}}$ 。

***

按照下述估计目标空间上是有有乘员：

$(hat O_1[t],hat O_2[t],...,hat O_{N_z}[t])=arg_{（O_1,O_2,..,0_{N_Z}）}in{{0,1}^{N_z}}maxhat p_{(o_1,o_2,...,o_{N_z})}[t]$

因为 $g(cdot)$ 为增函数，上式可改写为：

$(hat O_1[t],hat O_2[t],...,hat O_{N_z}[t])=arg_{（O_1,O_2,..,0_{N_Z}）in{{0,1}^{N_z}}}maxW_{(o_1,o_2,...,o_{N_z}）}^{T} ilde X[t]$

***

进一步，区域检测可通过矩阵矢量乘法实现：

将矩阵 $W_{(o_1,o_2,...,o_{N_z})}$ 写为 $2^{N_z}×（1+N）$ 的矩阵 $W$ :

$W=[egin{matrix}W_{0,0,...,0}^{T}\W_{1,0,...,0}^{T}\...\W_{1,1,...,1}^{T}end{matrix}]$

因此，可通过下式得到检测结果：

$(hat O_1[t],hat O_2[t],...,hat O_{N_z}[t])=de2bi[argmaxg(W ilde x[t])]$

其中： $g(cdot)$ 为做分类的非线性函数， $arg max(cdot)$ 是返回最大值的索引， $de2bi[cdot]$ 是将十进制转化为二进制。

***

当 $N_z=3$ 时，其索引与二进制关系如下：

index
0	（0,0,0）
1	（1,0,0）
2	（0,1,0）
3	（1,1,0）
4	（0,0,1）
5	（1,0,1）
6	（0,1,1）
7	（1,1,1）

3 总结

通过上述操作即可实现舱内乘员的检测，首先通过1DFFT及DOA生成“距离-方位”的热力图，在得到热力图后，提取其特征，并采集数据进行离线的训练优化，得到权重矩阵，依据权重矩阵及区域检测算法，得到乘员分布情况。

上述参考TI的官方指南，如有需要可私信
参考视频：车内乘员检测TI:https://edu.21ic.com/video/3813