MATLAB实战：从零推导合成孔径雷达(SAR)后向投影(BP)算法核心公式与成像全流程

1. 从“听回声”到“画地图”：理解SAR与BP算法的本质

大家好，我是老张，在雷达信号处理这个行当里摸爬滚打了十几年，从最初的懵懂到后来带着团队做项目，踩过的坑比走过的路还多。今天，我想和你聊聊合成孔径雷达（SAR）成像，特别是那个听起来很玄乎、但原理却异常直观的后向投影算法。咱们不整那些虚头巴脑的理论堆砌，就从一个你绝对能听懂的场景开始：听回声，画地图。

想象一下，你站在一个漆黑的大山洞里，手里只有一个手电筒和一个能精确记录回声时间的设备。你想知道山洞的墙壁长什么样。你会怎么做？很简单，你朝一个方向打一束光（相当于雷达发射电磁波），然后仔细听回声。回声回来的时间越晚，说明那个反射点离你越远。如果你只在一个位置听，你只能知道前方一条线上物体的远近，这就像最原始的雷达，只能测距，得不到“图像”。

那怎么得到一幅二维的“地图”呢？你开始沿着山洞的一侧慢慢走动，每走一小步，就朝同一个方向打一束光并记录回声。你走的这条路径，就形成了一个“合成孔径”。当你把所有位置听到的回声，按照它们各自的空间位置（你当时的位置、回声的时间）巧妙地“投影”回你心中设想的那张地图网格上时，一幅清晰的山洞墙壁图像就逐渐浮现出来了。这个“把不同位置听到的回声，反向投影回去构建图像”的过程，就是后向投影算法的核心思想。

SAR干的就是这个事，只不过它是在飞机或卫星上，对地面进行照射。它发射的不是简单的光脉冲，而是线性调频信号，这种信号像鸟叫一样，频率由低到高变化，拥有巨大的带宽，从而能获得极高的距离分辨率。然后，它接收地面散射回来的、经过延迟和调制的回波。BP算法的任务，就是从这一大堆看似混乱的回波数据中，反推出地面上每一个小格子（像素）的反射强度，最终形成一幅高分辨率雷达图像。

所以，如果你手头有一组SAR原始回波数据，感觉像是一团乱麻，别慌。接下来，我就带你用MATLAB，从最基础的公式推导开始，一步步把这团“乱麻”织成一幅清晰的“锦绣”。我们会紧扣“公式-代码”的对应关系，你不仅会知道每一步要做什么，更会明白为什么这么做，以及如何在MATLAB里把它实现出来。我们的目标很明确：推导出核心公式，写出每一行都能理解的代码，最终对一个点目标实现完美成像。

2. 基石：深入理解LFM信号与回波模型

万事开头难，而理解我们发射的信号，就是这最关键的第一步。SAR之所以能看得清，一大半功劳要归于它发射的线性调频信号。

2.1 LFM信号：雷达的“鸟叫声”

为什么用LFM信号？你可以把它想象成一把声音洪亮、音调变化丰富的“哨子”。普通的单频脉冲就像一声短暂的“嘀”，能量分散，距离分辨能力差。而LFM信号在持续时间内，频率线性地从f_start扫到f_end，就像一声悠长的“嘀~~~~~~~~~~”，它把能量分散在了更宽的频带上。接收时，我们通过一种叫“脉冲压缩”的技术，又能把这声长“哨”压缩成一个极其尖锐的脉冲，从而在距离维度上区分开两个非常近的目标。

它的数学表达式，是咱们所有推导的起点。在MATLAB里，我们首先得把它造出来。公式看起来有点唬人：

s_t(t) = rect(t / T) * exp(1j * 2 * pi * (f_c * t + 0.5 * K_r * t.^2))

别怕，我们拆开看：

• rect(t / T)：这是一个矩形窗，表示信号只在时间长度T（脉冲宽度）内存在。在代码里，这就是一个从0到T的时间向量。
• exp(1j * ...)：这是复数表达形式，1j在MATLAB里就是虚数单位i。它描述了信号的相位变化。
• f_c * t：这是载频带来的相位变化，很平稳。
• 0.5 * K_r * t.^2：这是关键！ K_r是调频率，单位是Hz/s。正是这个t的平方项，让信号的频率随着时间t线性增加（因为瞬时频率是相位对时间的导数）。

我们来写代码生成它。这里要注意，MATLAB是离散采样的，所以我们的t是一个时间向量。

% 参数设置
T = 10e-6;         % 脉冲宽度 10微秒
B = 50e6;          % 信号带宽 50MHz
Fs = 100e6;        % 采样率 100MHz
f_c = 5.3e9;       % 载频 5.3GHz (C波段)

% 计算调频率K_r: K_r = B / T
K_r = B / T;

% 生成快时间轴（一个脉冲内的采样点）
t_fast = (0:1/Fs:T-1/Fs).'; % 列向量，方便后续处理
% 注意：这里T可能不是Fs的整数倍，实际工程中会调整。我们简化处理。

% 生成LFM发射信号
s_tx = exp(1j * 2 * pi * (f_c * t_fast + 0.5 * K_r * t_fast.^2));
% 添加矩形窗（实际上指数部分已经隐含在时间内，这里显式加一下）
s_tx = s_tx .* (abs(t_fast - T/2) <= T/2); % 矩形窗

% 画图看看时域和频域
figure;
subplot(2,1,1);
plot(t_fast*1e6, real(s_tx));
xlabel('时间 (\mus)'); ylabel('幅度');
title('LFM信号实部（时域）');
grid on;

subplot(2,1,2);
freq_axis = linspace(-Fs/2, Fs/2, length(s_tx));
S_tx = fftshift(fft(s_tx));
plot(freq_axis/1e6, abs(S_tx));
xlabel('频率 (MHz)'); ylabel('幅度');
title('LFM信号频谱');
grid on;

运行这段代码，你会看到一个时域上振幅恒定但相位变化越来越“密”的信号，以及一个频域上近似矩形的频谱。这就对了，它就是我们后续所有处理的“源头活水”。

2.2 回波信号：带了“距离烙印”的LFM