手把手教你用MATLAB实现麦克风阵列的Delay-and-Sum波束形成（附完整代码）

麦克风阵列波束形成实战：从Delay-and-Sum原理到MATLAB工程实现全解析

在智能语音交互设备、会议系统乃至专业录音设备中，我们常常惊叹于它们能在嘈杂环境中“听清”特定方向的声音。这背后的一项核心技术，便是麦克风阵列波束形成。对于许多初次接触阵列信号处理的工程师和开发者而言，理论公式与算法论文往往令人望而生畏，而一个清晰、可运行、能直接上手的工程实现方案，才是打通从理论到应用的关键桥梁。本文将彻底抛开复杂的数学推导，聚焦于最经典、最基础的Delay-and-Sum (DAS) 波束形成算法，手把手带你用MATLAB从零构建一个完整的波束形成器。我们将深入探讨时域与频域两种延迟补偿方案的实现细节、分数延迟处理的工程权衡，并提供可直接复用的代码模块。无论你是希望快速验证算法效果的学生，还是需要在产品中集成基础声源增强功能的开发者，这篇文章都将为你提供一条清晰的实践路径。

1. 波束形成核心思想与Delay-and-Sum原理

想象一下，在一个喧闹的会议室里，你希望听清坐在长桌另一端同事的发言。如果你只有一只耳朵，这可能会很困难。但如果你有两只耳朵，并且可以微微侧头，将“听觉焦点”对准那位同事，他的声音就会变得相对清晰。麦克风阵列波束形成，本质上就是为机器赋予这种“侧耳倾听”甚至“聚焦收音”的能力。

一个由多个麦克风按特定几何形状（如线性、圆形）排列而成的系统，构成了一个麦克风阵列。声波从声源传播到阵列中不同位置的麦克风，会产生微小的时间差（Time Difference of Arrival, TDOA）。波束形成算法通过巧妙地处理这些时间差，可以增强来自特定方向的声音，同时抑制来自其他方向的干扰噪声。

Delay-and-Sum (DAS) 是众多波束形成算法中最直观、最基础的一种。它的名字直白地揭示了其两个核心步骤：

1. Delay (延迟)：根据目标声源的方向，计算每个麦克风通道信号相对于某个参考点（通常是阵列中心或第一个麦克风）需要补偿的时间延迟。
2. Sum (求和)：将所有通道经过延迟补偿后的信号进行加权平均（求和）。

其根本目的是让来自目标方向的信号在叠加时同相相加，实现相长干涉，从而幅度增强；而来自其他方向的信号则可能异相抵消，实现相消干涉，从而被抑制。

提示：DAS算法属于固定波束形成范畴，其权重（此处的延迟量）是预先根据目标方向计算好的，不随输入信号变化。这使得它计算量小、实现简单、鲁棒性强，但对抗复杂干扰（如移动声源、强混响）的能力有限。

在工程实现中，一个无法回避的挑战是分数延迟问题。我们采集到的音频信号是离散的数字信号，采样间隔固定（例如，采样率Fs=16kHz时，间隔为1/16000秒）。计算出的理想延迟时间，很可能不是采样间隔的整数倍。例如，需要延迟0.375个采样点。我们无法直接“移动”信号0.375个点。这就引出了实现DAS的两种主流思路：

接下来的章节，我们将分别深入这两种方法的MATLAB实现细节。

2. 工程准备：阵列几何与延迟计算

在编写任何一行信号处理代码之前，我们必须先明确系统的“物理配置”。这包括麦克风阵列的几何布局、目标声源的方向假设，以及如何将这些物理量转换为信号处理所需的采样点延迟。

2.1 定义阵列与声源

我们以一个简单的均匀线性阵列为例。假设有5个麦克风，沿x轴等间距排列，间距为4厘米。声速在空气中通常取343米/秒。

%% 1. 系统参数设置
clear; close all; clc;

% 基本参数
c = 343;               % 声速 (m/s)
Fs = 16000;            % 采样率 (Hz)
num_mics = 5;          % 麦克风数量
d = 0.04;              % 麦克风间距 (米)

% 定义均匀线性阵列位置 (单位：米)
% 假设阵列沿x轴排列，阵列中心在原点
mic_positions = zeros(3, num_mics);
for i = 1:num_mics
    mic_positions(1, i) = (i - (num_mics+1)/2) * d; % x坐标
    % y和z坐标均为0，表示在同一水平线上
end
disp('麦克风阵列位置 (x, y, z 坐标，单位：米):');
disp(mic_positions');

2.2 从方向到时间延迟

目标方向通常用方位角和俯仰角表示。我们假设目标声源位于阵列的远场，即声波到达阵列时可视为平面波。这意味着对于所有麦克风，声波的方向是平行的，延迟差仅由麦克风在波前法线方向上的投影决定。

计算延迟的核心是导向矢量。对于远场平面波，从方向单位向量到每个麦克风的时间延迟差可以通过简单的点积运算得到。

%% 2. 计算目标方向的延迟 (以采样点为单位)
% 假设目标声源方向：方位角 azimuth = 30度，俯仰角 elevation = 0度（水平方向）
azimuth_deg = 30;   % 方位角，从正x轴逆时针旋转
elevation_deg = 0;  % 俯仰角，从xy平面向上为正

% 将角度转换为弧度
azimuth = deg2rad(azimuth_deg);
elevation = deg2rad(elevation_deg);

% 计算目标方向的单位向量 (笛卡尔坐标系)
% 注意：这里的方向向量指向声源，而延迟计算需要的是波传播的方向（相反方向）。
% 通常约定：信号从(az, el)方向传来，其传播方向向量为：
% u = -[cos(el)*cos(az); cos(el)*sin(az); sin(el)]
u = -[cos(elevation)*cos(azimuth);
      cos(elevation)*sin(azimuth);
      sin(elevation)];

% 计算每个麦克风相对于原点的时延（秒）
% 时延 = (麦克风位置向量 点乘 波传播方向单位向量) / 声速
tau_seconds = (mic_positions' * u) / c; % 列向量，每个元素对应一个麦克风

% 将时延转换为采样点延迟（可能是小数）
delay_samples = tau_seconds * Fs;

disp('每个麦克风的理论延迟（秒）:');
disp(tau_seconds');
disp('每个麦克风的理论延迟（采样点数，含小数）:');
disp(delay_samples');

运行这段代码，你会得到类似下面的输出。注意，延迟值有正有负，负延迟意味着该麦克风接收到的信号比参考点（原点）更早，在补偿时我们需要对其进行正向延迟（即让它的信号“等一等”），所以后续处理中通常使用 -delay_samples 作为补偿量。

麦克风阵列位置 (x, y, z 坐标，单位：米):
   -0.0800         0         0
   -0.0400         0         0
         0