简介面向雷达数据处理与信号分析场景的MATLAB工具包聚焦DCA1000采集数据的读取与雷达信号峰值检测两大任务。readDCA1000.m负责打开二进制文件、解码IQ采样并分离I/Q分量将原始数据转换为可处理的时间序列findAllPeak.m则在预处理基础上执行寻峰算法通过幅度阈值或导数特征识别目标峰位输出峰值位置与幅度。两个脚本协同可覆盖从底层数据解析到目标特征提取的完整链路适合雷达工程初学者快速上手实测数据处理也为研究人员验证算法提供可直接修改的参考代码。压缩包仅含2个m文件、约2KB体量轻但功能聚焦目前已有363人学习下载适合正在学习雷达信号处理、需要DCA1000数据读取或峰值检测现成实现的读者。1. findAllPeak 到底解决雷达目标检测里的哪个痛点做完距离维 FFT 和速度维 FFT 之后你手里是一张 RD 图目标在这个二维谱图上表现为一个个凸起的峰。很多刚接触雷达数据处理的工程师会随手写max()拿到全图最大值结果发现只输出一个目标可雷达前方明明可能站着两个人、停着两辆车多目标场景下全被这一个全局最大点吃掉了。findAllPeak 这类峰值检测算法解决的就是这个痛点它不是让你找“最大”而是把满足幅值、间距、局部凸起条件的所有峰值全部提取出来作为后续点迹聚合、目标跟踪的输入。这篇文章面向从“跑通 FFT”走向“真正做雷达数据处理”的工程师也适合已经用 CFAR 拉出候选点但被重复峰困扰的老手。读完你可以自己实现一个可用的 findAllPeak并且知道参数怎么调、坑在哪。2. 雷达数据处理链条里findAllPeak 是点迹提取的最后一道闸2.1 从 ADC 到点迹findAllPeak 在信号处理流程中的位置常见毫米波雷达数据处理的链路是这样的ADC 采样得到一串复数 IQ 数据按 chirp 排成矩阵先沿着快时间维度做距离维 FFT把回波时延转成距离再沿着慢时间维度做多普勒维 FFT把相位变化转成速度得到的就是 Range-Doppler 图也就是 RD 图。如果是 4D 毫米波雷达会在后续做角度维 FFT把多个接收天线的相位差转成方位角和俯仰角。之后一般做恒虚警检测CFAR把超过自适应门限的单元标成候选目标最后才是峰值检测与点迹输出。findAllPeak 在整条链里的位置非常靠后它工作在 CFAR 之后或者在没有任何检测算法时直接作用在谱图上。CFAR 解决的问题是“哪里可能超过噪声门限”findAllPeak 解决的问题是“这个超门限连通域里到底哪个点代表目标”。原因是雷达的点扩展函数有主瓣也有旁瓣一个目标往往会点亮相邻好几个距离-多普勒单元CFAR 会把这些单元全部判成过门限如果不做峰值提取一个目标会被输出成十几个“点迹”后端的关联滤波器会直接疯掉。所以峰值检测是把检测结果从“区域”压缩成“点”的最后一道闸。2.2 峰值检测的数学定义局部极大值判定从一个一维频谱讲起。给定一段幅度序列x[i]所谓峰值第一层含义是局部极大这个点的幅度同时大于它的左邻点和右邻点也就是满足x[i] x[i-1]且x[i] x[i1]。如果序列里存在平坦段就要改成x[i] x[i-1]且x[i] x[i1]或选用带方向的比较策略避免在平顶上输出一堆重复点。这是 findAllPeak 最基本的筛选条件。只靠局部极大显然不够因为噪声也会产生大量毛刺每个毛刺都是局部极大。所以第二层含义是幅度门限峰的幅度必须大于某个绝对高度height。第三层是凸起门限峰与周围最低点的差值要超过threshold用来剔除落到噪声基底上的小起伏。第四层是相邻峰最小间距distance两个峰之间至少隔多少个采样点不满足的只保留幅度更高的那个。把这几层组合起来才是一个真正可用的 findAllPeak。严格说它不像 CFAR 那样是统计检测器而是一个确定性提取器把已经确认存在的目标从谱峰邻域里“挑”出一个代表点。2.3 为什么 max 不够用多目标与漏检有人会问既然目标在谱图上是一个峰那取最大值不就行了吗问题是雷达场景天然是多目标。同一距离门上可能有两个径向速度不同的目标距离维 FFT 后它们重叠在一起但多普勒维会分开全图取 max 只能留下速度更高的那一个另一个直接丢了。同一速度上可能有两个不同距离的目标距离维 FFT 分开后max 又只能留距离最近那个。在安检、周界、矿山等场景里这种漏检是致命的。还有一个更隐蔽的问题目标的真实回波不一定是全局最大。近处有一个强目标远处有一个弱目标弱目标幅度低但依然超过检测门限max()只能看到强目标而 findAllPeak 配合 CFAR 门限后可以把强目标主瓣之外的所有显著峰都列出来。雷达距离方程告诉我们回波功率与距离的四次方成反比强弱目标动态范围可以相差几十分贝所以“全局最大”在雷达数据处理里几乎没有意义。真正有意义的只有一个问题哪些峰是可信的目标而不是旁边峰拖出来的尾巴。findAllPeak 的存在就是把这些尾巴切掉留下干净的目标位置。3. 用 Python 写一个 findAllPeak核心实现与参数调优3.1 一维峰值检测最小实现基于 NumPy常见做法是先把信号转成浮点数组然后一次性用 NumPy 的向量化比较找出所有局部极大再逐个应用门限和去重。下面是一个可以直接拿去跑的版本import numpy as np def find_all_peak(x, heightNone, distance1, thresholdNone): 在一维频谱 x 中找出所有局部极大值点。 参数 ---- x : np.ndarray 输入频谱幅度可为距离维 FFT 输出或 RD 图某一切片 height : float 峰的最低绝对幅度低于该值直接丢弃 distance : int 两个峰之间的最小样本数小于该距离的峰只保留幅值高的 threshold : float 峰相对两侧基线的凸起高度用于抑制噪声毛刺 返回 ---- peaks : np.ndarray 每个峰在 x 中的索引 peak_heights : np.ndarray 每个峰对应的幅度值 x np.asarray(x, dtypenp.float64) if len(x) 3: return np.array([], dtypeint), np.array([], dtypefloat) # 候选峰严格大于左右相邻点边界点不参与 left np.concatenate(([np.inf], x[:-1])) right np.concatenate((x[1:], [np.inf])) is_local_max (x left) (x right) peak_idx np.nonzero(is_local_max)[0] # 绝对幅度门限 if height is not None: peak_idx peak_idx[x[peak_idx] height] # 凸起高度门限 if threshold is not None: keep [] for i in peak_idx: # 取左右相邻点的较大值作为基线 base max(x[i - 1], x[i 1]) if i not in (0, len(x) - 1) else x[i] if x[i] - base threshold: keep.append(i) peak_idx np.array(keep, dtypeint) # 按幅度降序排列再按照 distance 做非极大值抑制 if len(peak_idx) 0: return np.array([], dtypeint), np.array([], dtypefloat) order np.argsort(x[peak_idx])[::-1] sorted_idx peak_idx[order] sorted_val x[sorted_idx] selected [] for idx, val in zip(sorted_idx, sorted_val): if all(abs(idx - s) distance for s in selected): selected.append(idx) selected np.array(sorted(selected), dtypeint) return selected, x[selected]这段代码的逻辑很直白先用向量化比较筛出所有严格局部极大再用height去掉低幅度杂波用threshold去掉相对周围没有明显凸起的平缓隆起最后把候选峰按幅度从高到低排序逐个检查是否与已选峰距离过近近则丢弃。最后一步非常重要因为窗函数旁瓣和频谱泄漏会让同一个目标在相隔一两个采样点的地方出现两个假峰按幅度排序后保留最高者就能把目标分裂压回去。distance的单位是采样点不是米或者米/秒。如果距离维 FFT 的 bin 间隔是 0.15 米distance3就表示两个峰至少要相隔 0.45 米才被当成不同目标。这个参数需要根据距离分辨率、过采样率来定不能拍脑袋。3.2 参数设计height、distance、threshold 的取值依据参数含义雷达场景推荐值设太大设太小height峰绝对幅度下限噪声平均幅度的 46 倍或 CFAR 门限值弱目标漏检大量噪声峰混入distance峰间最小样本间隔23 个 FFT bin相邻真实目标被合并目标分裂成多峰threshold相对两侧基线的凸起量36 dB 对应的线性倍数低信噪比目标被丢掉噪声毛刺被当成峰关于height很多从实验室走向现场的人会栽跟头直接用幅度绝对值换一个雷达芯片、换一版 ADC 增益设置参数全部失效。我一般建议输入 findAllPeak 之前先把谱图归一化比如除以整个热图的噪声平均功率或者让height接收“噪声倍数”而不是原始幅度。雷达距离方程决定回波动态范围很大近距离强目标可能比噪声高 40 dB远距离目标可能只高 6 dB固定绝对门限很难同时兼顾。所以更稳的做法是把height设成一个保守值比如噪声均值以上 6 dB宁多勿缺然后用 CFAR 或后处理去剔除假点。threshold是很多人容易忽略的参数。它的作用是判断一个点到底是不是真正的峰强目标主瓣很尖锐峰顶与两侧平台差很大噪声毛刺虽然也是局部极大但凸起量很小。用threshold可以滤掉大部分噪声峰而不用把height调得很高。我习惯把threshold设在 3 dB 左右也就是要求峰至少比左右相邻采样点高出 1.41 倍这个值在绝大多数调频连续波雷达里都能有效压制平台噪声。3.3 从一维扩展到二维 RD 图峰值检测不能直接分行找RD 图本质是一个二维矩阵目标在距离-多普勒平面是一个二维凸起不能简单对每一行调用一维 find_all_peak再对每一列调用一次。那样同一个目标会在行搜索和列搜索里各出现一次点迹直接翻倍。正确思路是使用二维局部最大值滤波器如果一个点的幅度在其周围n×n邻域内同时是最大值才把它当成二维峰。SciPy 提供了现成的maximum_filterfrom scipy.ndimage import maximum_filter def find_all_peak_2d(rd, heightNone, min_distance3, neighborhood5): rd : np.ndarray 距离-多普勒幅度图axis0 为多普勒axis1 为距离 # 邻域内最大值的掩码 local_max (rd maximum_filter(rd, sizeneighborhood, modenearest)) # 去掉边界导致的假极大值 pad neighborhood // 2 local_max[:pad, :] False local_max[-pad:, :] False local_max[:, :pad] False local_max[:, -pad:] False # 绝对幅度门限 if height is not None: local_max (rd height) ys, xs np.nonzero(local_max) vals rd[ys, xs] # 按照幅度排序再按照 min_distance 做非极大值抑制 order np.argsort(vals)[::-1] ys, xs, vals ys[order], xs[order], vals[order] selected_rows, selected_cols [], [] for y, x, v in zip(ys, xs, vals): too_close False for py, px in zip(selected_rows, selected_cols): if abs(y - py) min_distance and abs(x - px) min_distance: too_close True break if not too_close: selected_rows.append(y) selected_cols.append(x) return np.array(selected_rows), np.array(selected_cols), \ np.array([rd[y, x] for y, x in zip(selected_rows, selected_cols)])neighborhood的参数取值很关键。它表示一个目标在 RD 图上占据的“点扩展函数”尺寸。比如距离维用了汉宁窗主瓣会比矩形窗宽那么目标在距离方向会拖出两三个明显的旁瓣如果neighborhood太小旁瓣照样会超过邻域最大值产生假峰如果太大两个靠得很近的真实目标会被当成一个目标。常规做法是先画一张强目标的 RD 图切片数出主瓣上下左右大概几个格点再把neighborhood设为主瓣宽度的 1.5 倍左右。我一般先设为 5跑完看点迹数量如果发现目标分裂就调到 7 或 9这个参数本身不敏感但调差了效果很明显。4. 把 findAllPeak 接入真实雷达数据链以毫米波雷达 RD 图为例4.1 数据从哪里来AWR2243 读取与 RD 图生成以最常见的毫米波雷达数据处理为例很多板子和芯片比如 AWR2243通过以太网或 LVDS 输出原始 ADC 数据格式是多个发射天线、多个接收天线的复数矩阵。拿到数据后第一步是先把它按照 chirp 序号和采样点序号排成二维矩阵raw_data [n_chirp, n_sample]。这里是复数数据每一行代表一个 chirp 内的快时间采样同一行内的多个采样点按时间顺序排列。距离维 FFT 在快时间方向做速度维 FFT 在慢时间方向做。常见做法是先沿 axis1 做距离维 FFT再去掉距离门零频附近不需要的部分然后沿 axis0 做多普勒 FFT最后用fftshift把零多普勒移到中间import numpy as np # raw_data 形状: (n_chirp, n_sample)复数矩阵 range_fft np.fft.fft(raw_data, axis1)[:, :n_sample] # 距离维 # 可选对每根天线各自做完距离维后再做速度维 doppler_fft np.fft.fftshift( np.fft.fft(range_fft, axis0), axes0 ) rd_map np.abs(doppler_fft) # 幅度 RD 图这段代码里np.fft.fft(range_fft, axis0)是对每一列距离门做多普勒维 FFT。由于多普勒维有正负频率fftshift把负多普勒放到左边零多普勒放中间正多普勒放右边。后面做峰值检测时一个速度等于 0 的静止目标会出现在中间行附近这符合直觉。做完之后rd_map的形状是(n_chirp, n_sample)其中每一格代表一个距离-多普勒单元幅度代表该单元的回波强度。这个rd_map就是 findAllPeak 的输入。4.2 先 CFAR 后 findAllPeak 还是先峰后 CFAR这个问题的答案是必须先 CFAR 后 findAllPeak或者至少把 CFAR 门限作为height传入 find_all_peak_2d。原因是雷达背景噪声不平稳近距地面杂波强远距噪声低固定height在全图上效果很差。CFAR 在每个检测单元周围划定参考窗口当下背景噪声强就把门限抬高噪声弱就把门限降下来保证虚警率恒定。CFAR 的输出是一个与 RD 图同形的布尔矩阵每个True表示该单元超过局部门限。直接把rd_map中cfar_mask为 False 的位置清零再送进find_all_peak_2d能显著减少峰数量而且不会漏掉弱目标。示例def extract_peaks_with_cfar(rd_map, cfar_mask, min_distance3): # 只保留过门限位置的幅度其余置零 masked np.where(cfar_mask, rd_map, 0.0) # 二维局部最大避免一个目标输出多个峰 local_max (masked maximum_filter(masked, size5, modenearest)) ys, xs np.nonzero(local_max cfar_mask) # 按幅度排序 最小距离去重 vals masked[ys, xs] order np.argsort(vals)[::-1] ys, xs, vals ys[order], xs[order], vals[order] selected [] for y, x, v in zip(ys, xs, vals): if all(abs(y - py) min_distance or abs(x - px) min_distance for py, px in selected): selected.append((y, x)) return selected这里用了一个小技巧maximum_filter只会在局部极大位置保留原值其他位置会被邻域内更大的值替代然后和cfar_mask取交集这样即使 CFAR 标记出一个复杂连通域每个连通域也只输出一个峰。需要注意的是min_distance在二维里用“或”而不是“与”也就是多普勒方向或者距离方向只要有一个足够远就算不同目标如果同时很近才合并。你可以在实际项目里调整这个逻辑比如 4D 毫米波雷达还要把角度维一并纳入距离判断。4.3 从峰值坐标换算距离、速度、角度峰在 RD 图上的索引本身不算目标状态要换算成物理量。假设距离维 FFT 点数为N_range雷达带宽为B则距离分辨率是c / (2B)距离单元r_idx对应的距离是r_idx * c / (2B)。注意这里的r_idx是从 0 开始的而且距离维窗口可能因为滤波器延迟有一个固定偏置通常需要校准。速度维换算稍微复杂一点。调频连续波雷达一个 chirp 周期为T_chirp波长λ雷达最大不模糊速度是λ / (4 * T_chirp)多普勒维点数N_doppler对应速度分辨率λ / (2 * T_total)其中T_total是所有 chirp 的总时间。做完fftshift之后速度索引v_idx到速度的换算要先把索引减去N_doppler // 2再乘以速度分辨率。下面是一个换算示例量公式说明距离r r_idx * c / (2*B)距离分辨率不受 FFT 点数影响多普勒速度v (v_idx - N_doppler//2) * λ / (2*T_total)fftshift 后 v_idx 以 0 多普勒为中心角度theta arcsin(angle_fft_bin * λ / (2*pi*d))多个 RX 天线对距离-多普勒谱做角度 FFT峰值检测输出的坐标如果能提前换算成距离、速度再送进跟踪模块后端的处理会轻松很多。我习惯让 findAllPeak 的输出带上range、doppler、snr这三个字段而不是光给一个浮点索引。SNR 可以直接用峰的幅度减去同距离门噪声平均幅度单位 dB这个值后面做目标丢弃和质量评估很管用。5. findAllPeak 在雷达数据处理中的五个典型翻车现场与排查5.1 同一个目标被输出成两个相邻点迹现象RD 图上一辆车只有一片亮斑但点迹输出时同一个距离-多普勒位置出现两个相距 1 到 2 个 bin 的目标点迹。原因目标点扩展函数不是单点主瓣有一定宽度加上距离维 FFT 可能用了较少的点数频谱泄漏让主瓣顶部出现两个局部极大值CFAR 检测时这两个极大值都能过门限findAllPeak 若没有按距离去重就会当成两个峰输出。解决在 findAllPeak 的非极大值抑制里加大distance参数。常见做法是先看强度剖面图数出主瓣半宽在几个 bin再把distance设为主瓣宽度的 1.2 到 1.5 倍。条件允许的话最好在峰附近做质心聚合把主瓣内所有超过 3 dB 的格点合成一个峰避免误分裂。5.2 强目标旁瓣淹没了邻近弱目标现象一个强目标周围出现几个整齐的小峰幅度比真实目标低 10 dB 左右位置对称分布在强目标两侧且与强目标保持固定的间隔。原因这是窗函数旁瓣的直接体现。矩形窗第一旁瓣只有约 13 dB 抑制汉宁窗旁瓣抑制约 32 dB但旁瓣位置相对主瓣宽度是固定的强目标回波功率过大时旁瓣电平和弱目标真实回波相当findAllPeak 会把旁瓣峰值也选出来。解决优先在距离维和多普勒维加窗比如距离维加汉宁窗多普勒维再加汉宁窗把旁瓣压下去。然后调高threshold参数要求峰与两侧基线差值大于 6 dB这样旁瓣这种周围本来就有强主瓣背景的峰会被过滤掉。特别注意不要把height调太高那会连弱目标一起丢掉。5.3 静态杂波造成零多普勒附近的虚假峰现象RD 图中间行特别亮findAllPeak 在速度 0 附近输出几十个目标点迹真实运动目标反而淹没在里面。原因地面、天线直耦、静止建筑物反射的杂波集中在零多普勒附近幅度远高于运动目标回波。固定height和threshold拦不住这堵“直流墙”。解决在峰值检测之前先做一次多普勒高速滤波或 MTI 对消。最简单的两脉冲对消就是把多普勒维相邻两行相减doppler_fft[1:] doppler_fft[1:] - doppler_fft[:-1]这会直接削掉零频率附近的静止分量。做完之后再看 RD 图运动目标保持静止杂波被大幅压制这时 findAllPeak 的检测数量会恢复正常。如果对消后噪声抬升建议把height重新标定到对消后噪声均值上。5.4 二维搜索时行列分别找峰造成重复输出现象点迹数量突然翻倍每个目标有两个点迹一个距离坐标相同、多普勒坐标相差几个 bin或者反过来。原因这是最常见的 findAllPeak 使用错误——分别对每一行和每一列调用一维峰值搜索。行列交叉搜索会把同一个二维峰识别成两个一个来自行方向上的极大另一个来自列方向上的极大。解决改用二维局部最大值滤波器如上一章的find_all_peak_2d让一个峰只在其邻域内是最大值才被选中。如果工程上是 C 代码可以用 3×3 邻域遍历暴力实现对每个超过 CFAR 门限的格点检查它在上下左右和两个对角线方向是否全部大于邻居。这个 3×3 的最简版本其实已经能解决大多数重复点迹问题。5.5 阈值单位不一致导致参数不可移植现象同一套 findAllPeak 参数在 AB 雷达板上跑得好好的换了一台雷达 ADC 位数不同或模拟增益不同点迹数量暴涨或全部看不见。原因height使用了原始 ADC 幅度作为绝对值不同雷达的噪声电平、FFT 增益、窗函数增益都不一样绝对阈值没有物理意义。CFAR 门限是自适应的所以不受影响但 findAllPeak 里手写的绝对阈值就报废了。解决把 findAllPeak 的输入全部归一化到噪声功率。典型做法是先取 RD 图边缘区域无目标区域计算平均幅度noise_mean让rd_map rd_map / noise_mean再设置height为 613 dB 对应的倍数值。或者直接把 CFAR 输出的门限矩阵传给 findAllPeak用它替代height这样每个峰的检测标准都是“相对当下局部背景超多少 dB”换平台不需要改参数。6. 进阶用质心细化峰位输出点迹直连跟踪滤波6.1 质心细化突破频域栅格分辨率RD 图的 bin 间隔决定了距离/速度量化精度但目标真实参数并不刚好落在整数 bin 上。FFT 谱峰在离散采样后会呈现近似 sinc 形状峰顶附近的幅度加权中心靠近真实位置。可以在 findAllPeak 选出的峰周围取一个 3×3 邻域用幅度作为权值计算亚 bin 坐标def refine_peak(rd, y, x, half_win1): sub rd[y-half_win:yhalf_win1, x-half_win:xhalf_win1] w sub / np.sum(sub) ys np.arange(-half_win, half_win 1).reshape(-1, 1) xs np.arange(-half_win, half_win 1).reshape(1, -1) dy np.sum(ys * w) dx np.sum(xs * w) return y dy, x dx质心细化后距离坐标y dy不再是整数 bin。如果距离分辨率是 0.15 米质心偏移 0.3 bin 就对应 4.5 厘米的修正量在近距离测距应用里这个很重要。注意质心窗口不要开太大一般取峰左右各 1 个 bin 就够了窗口太大反而会把噪声旁瓣拉进来。6.2 输出结构化点迹直接接后端跟踪峰值检测的结果如果只输出一组(y, x)数组跟踪端还需要再去查原始数据这是设计缺陷。我习惯在峰值检测函数内部直接组装点迹结构体from dataclasses import dataclass dataclass class PointCloudTarget: range: float doppler: float azimuth: float snr: float peak_value: float每个目标携带从峰位换算出的物理距离、多普勒速度、角度、信噪比。这样做的好处是后端卡尔曼滤波不需要关心 FFT bin 和分辨率这些细节。4D 毫米波雷达只需要再多一步角度维 FFT然后把峰值对应的角度也填入azimuth字段整个方案从检测到跟踪就打通了。输出点迹时建议把 SNR 低于预置质量阈值的目标直接丢弃这比在峰值检测里硬调height更科学因为 SNR 在不同目标间已经归一化。6.3 性能预算峰值搜索的取舍最后说一个我在实际项目中养成的习惯永远不要让 findAllPeak 跑在全图每个点上。二维局部最大值滤波在超大 RD 图上代价不小尤其是 4D 毫米波雷达有多个接收天线数据量成倍增加。正确顺序是先跑一个轻量的 CFAR 或幅度门限筛出候选点再用小邻域做峰值判定这样计算量集中在少数候选区域。做嵌入式移植的时候还可以把峰值检测做成滑动窗比较器每来一行新数据就输出一次把延迟控制在单帧内。峰值搜索效果的好坏不取决于算法有多炫而在于参数是否能跟着场景走、输出的坐标是否能被后端直接用。这套思路帮我从“看波形毛刺”走到了“稳定输出空口点迹”希望帮到你。本文还有配套的精品资源点击获取