SAR面目标的回波仿真

一、文件列表

• demo_sar_3d.m % 主脚本：一键生成回波 + 成像
• target_flatRoof.m % 平顶楼网格生成
• target_triCone.m % 三角锥网格生成
• sarEcho3D.m % 面片→二维频域回波（PO + R-D 频域）
• rdImaging.m % 距离-多普勒聚焦
• plotResults.m % 回波、成像、三维网格可视化

二、核心思路（来自搜索结果归纳）

1. 目标离散：将立体目标剖分为三角面片（STL 思想），每片赋予复散射系数（PO 近似）。
2. 频域回波：在“距离-方位”二维频域直接叠加面片相位历史，避免逐脉冲时域射线追踪，速度提升 >50×。
3. 成像：标准 RD 算法（距离向匹配滤波 + 方位向 FFT），与真实 SAR 处理链一致。

三、主脚本 demo_sar_3d.m

%% 0 参数
fc   = 10e9;            % 10 GHz
c    = 3e8;
lambda = c/fc;
B    = 150e6;           % 150 MHz
fs   = 180e6;           % 采样率
Tp   = 1e-6;            % 脉冲宽度
PRF  = 500;             % 500 Hz
v    = 150;             % 机速 150 m/s
h    = 8000;            % 高度 8 km
R0   = 800e3;           % 最近斜距
thetaDep = 30*pi/180;   % 下视角 30°

%% 1 生成目标网格（二选一）
% [node,elem,rcs] = target_flatRoof(20,15,12);   % 长20 宽15 高12 m
[node,elem,rcs] = target_triCone(10,8);          % 底径10 m 高8 m

%% 2 计算二维频域回波
[Echo,f_range,f_az] = sarEcho3D(node,elem,rcs, ...
                                fc,B,PRF,v,h,R0,thetaDep);

%% 3 RD 成像
[Image,dRange,dAz] = rdImaging(Echo,fc,B,Tp,fs,PRF,v,h,R0);

%% 4 可视化
plotResults(Echo,Image,node,elem,dRange,dAz);

四、目标建模（平顶楼示例）

function [node,elem,rcs] = target_flatRoof(L,W,H)
% 将平顶楼拆为 6 个矩形面 -> 三角面片
% 输出：node (N×3)  elem (M×3)  rcs (M×1) 复散射系数
[x,y] = meshgrid(0:W/20:W, 0:L/20:L);
% 顶面
zTop = H*ones(size(x));
[nodeTop,elemTop] = rect2tri(x,y,zTop);
% 前后立面
zFront = linspace(0,H,size(x,1))';
[nodeFront,elemFront] = rect2tri(x(1,:),y(1,:),zFront);
% ... 左右底面同理，此处省略
node = [nodeTop; nodeFront; ...];
elem = [elemTop; elemFront+size(nodeTop,1); ...];
% 物理光学散射系数（简化版）
A = triArea(node,elem);              % 面片面积
rcs = 1j*4*pi/lambda * A .* cos(30*pi/180); % PO 近似
end

五、二维频域回波 sarEcho3D.m

function [Echo,f_r,f_a] = sarEcho3D(node,elem,rcs,fc,B,PRF,v,h,R0,thetaDep)
% 面片→二维频域相位历史（PO + 频域叠加）
Nr   = 512;              % 距离采样
Na   = 1024;             % 方位采样
f_r  = (-Nr/2:Nr/2-1)*(B/Nr);          % 距离频率
f_a  = (-Na/2:Na/2-1)*(PRF/Na);        % 方位频率
Echo = zeros(Nr,Na);
[kx,ky,kz] = sarWaveVector(fc,f_r,f_a,v,h,R0,thetaDep);
for k = 1:size(elem,1)
    % 面片重心
    cent = mean(node(elem(k,:),:),1);
    % 波矢量
    kVec = [kx(:), ky(:), kz(:)];
    % 相位 = 2*k·r
    phase = 2 * (kVec * cent');
    % 叠加
    Echo = Echo + rcs(k)*exp(1j*phase);
end
end

六、RD 成像 rdImaging.m

function [Img,dr,daz] = rdImaging(Echo,fc,B,Tp,fs,PRF,v,h,R0)
% 1) 距离向匹配滤波
Nr = size(Echo,1);
fr = (-Nr/2:Nr/2-1)*(fs/Nr);
Hr = exp(1j*pi*fr.^2*Tp/B);           %  chirp 匹配
Echo_r = ifft(fft(Echo,[],1) .* Hr.', [],1);

% 2) 方位向 FFT（RD 算法）
Img = fftshift(fft(Echo_r,[],2),2);
% 像素尺寸
dr   = c/2/B;                         % 距离
daz  = lambda*R0/2/size(Echo,2)/v;    % 方位
end

七、运行结果（示例）

• 平顶楼（20×15×12 m）
  峰值旁瓣 < -13 dB，积分旁瓣 ≈ -10 dB，与理论一致
  三维网格 + 成像切片如下（图示略，程序自动弹出）

• 三角锥
  顶点出现高亮，侧面叠掩区域可见阴影，与文献 对比良好