简介MATLAB仿真技术在半波振子天线子阵方向图研究中具有重要应用。此份zip资源包面向天线设计与电磁仿真初学者、通信工程专业学生及相关工程人员提供完整的半波振子子阵方向图仿真实现帮助理解振子天线辐射特性和子阵方向图生成原理。压缩包共2个文件包含1个MATLAB脚本文件和1个仿真结果图片整体大小仅45KB脚本涵盖天线参数设置、子阵构建、仿真计算与方向图绘制等关键环节图片直观展示了半波振子子阵的辐射强度分布。目前已有613人学习下载。通过运行该脚本读者可快速获得可视化方向图并能通过调整振子长度、馈电点位置、子阵排列等参数深入探究天线性能为实际阵列设计与教学科研提供可复用的仿真参考也可为赋形波束、宽角扫描等复杂辐射模式研究提供基础。1. 半波振子子阵方向图仿真单元方向图和阵因子到底谁说了算平时画半波振子的方向图读到的几乎都是那个“8”字形截面一旦把它拼成子阵整个方向图就会多出若干瓣和零陷主瓣也被压窄。真正决定主瓣方向和零点位置的往往不是单元方向图而是子阵阵因子单根振子的单元方向图只相当于给阵因子加了一个“罩子”。banbozhenzi.zip 里的 banbozhenzi.m做的就是把这套方向图乘积定理落成 MATLAB 代码最后生成类似 3.2figure.jpg 的极坐标子阵方向图。它适合三类人正在看阵列天线教材的学生、要做波束扫描的射频工程师、以及想快速对照全波仿真结果的系统设计人员。这个压缩包里的程序不长核心逻辑就是“单元方向图 × 阵因子”但每一行参数设置都直接反映天线工程习惯。下面从单元建模开始一步步拆开这套 matlab 仿真。2. 半波振子的单元方向图建模从解析式到 MATLAB 代码2.1 为什么半波振子单元方向图可以直接公式化半波振子是最典型的线天线振子长度等于半个工作波长电流沿振子呈近似正弦分布、末端电流趋于零。它的远场方向性有解析表达式不需要像微带贴片那样非得靠全波求解器才能得到稳定结果。若振子沿 z 轴放置球坐标系下的归一化电场方向图写作[ E_{unit}(\theta) \frac{\cos\left(\frac{\pi}{2}\cos\theta\right)}{\sin\theta} ]这个式子的特点是水平面θ 90°方向值最大归一化后为 1沿振子轴线方向θ 0° 和 180°辐射为零所以方向图是典型的“8”字形。做子阵仿真时单元方向图通常提前算好并归一化后续直接和阵因子相乘。工程里常有人直接把这段公式抄进 MATLAB结果画出来的方向图在轴向出现一个很大的跳变原因是分母 sinθ 在 θ 接近 0 或 π 时趋近于零分子也趋近于零代码里直接相除会得到 NaN。这不是物理问题而是数值处理问题需要在计算网格上特殊处理。2.2 参数初始化把频率、波长、阵元数写清楚打开 banbozhenzi.m你会发现第一段几乎都在做参数声明。这类脚本我习惯按下面的顺序组织c 2.99792458e8; % 光速单位 m/s freq 2.4e9; % 工作频率单位 Hz这里以 2.4 GHz 为例 lambda c / freq; % 波长约 0.125 m N 8; % 子阵单元个数 d lambda / 2; % 单元间距取半波长 theta linspace(0, pi, 721); % 极角网格从 0 到 pi分辨率 0.25°参数说明N 是参与合成的振子数量它直接决定阵因子波束宽度d 是相邻振子相位中心的距离通常取半波长取太大会在可见区出现栅瓣这个后面第五章专门讲。theta 网格建议取 721 或 3601 点点数太少会导致 -3 dB 波束宽度读数粗糙点数太多则没必要因为解析公式算起来很快瓶颈只在绘图。2.3 单元方向图计算的数值处理单元方向图的计算代码可以这样写theta_safe theta; % 将轴向sin theta 接近0置为 NaN先避开除零 theta_safe(abs(sin(theta)) 1e-6) NaN; E_unit cos(pi/2 * cos(theta_safe)) ./ sin(theta_safe); % 轴向辐射物理上为 0这里把之前置 NaN 的位置补回 0 E_unit(abs(sin(theta)) 1e-6) 0; % 归一化让最大值等于 1 E_unit E_unit / max(abs(E_unit));逻辑说明第一处替换是为了让 MATLAB 在计算 NaN 时不会报警告第二处替换是利用半波振子轴方向电场为零的物理结论把轴向值强制设为 0。这里的阈值 1e-6 不是固定标准如果 theta 网格点足够密取 1e-8 也一样关键是保证 sinθ 不为零。有人在网上抄的代码喜欢写E_unit(isnan(E_unit)) 1对半波振子来说这是错的。轴向零点是半波振子的固有特性不是数值异常改成 1 会把后续子阵方向图的零点全部抬高副瓣电平直接失真。3. 方向图乘积定理子阵方向图的合成逻辑3.1 阵因子和单元方向图之间的边界在哪里子阵方向图不是把每根振子的方向图简单相加而是先算阵因子再和单元方向图相乘。方向图乘积定理的适用条件是单元之间互耦可忽略、阵列处于远场观察区域、各单元方向图相同。它的表达式很简洁[ F_{total}(\theta) F_{unit}(\theta) \times AF(\theta) ]AF 就是阵因子它只跟阵列几何布局、单元间距、激励幅度和相位有关。对于沿 z 轴等间距排布的 N 元线阵阵因子可写成[ AF(\theta) \sum_{n0}^{N-1} w_n e^{j n k d \cos\theta} ]其中 k 2π/λ 是波数w_n 是第 n 个单元的复激励。这个公式的物理含义是每个单元到远场观察点的路径差不同产生的相位差叠加后形成干涉图样。均匀激励时 w_n 1阵因子有解析闭合形式但用循环求和的方式写更直观也方便以后加幅度锥削。3.2 均匀线阵的法向阵因子实现先看不扫描、无幅度加权的情况MATLAB 代码可以这样写psi 2 * pi * d / lambda * cos(theta); % 相邻单元相对观察方向的空间相位差 AF zeros(size(theta)); for n 0:N-1 AF AF exp(1j * n * psi); end AF abs(AF) / N; % 除以 N使最大值为 1代码逻辑循环变量 n 从 0 开始对应第 0 个阵元取参考相位psi 是一个随 theta 变化的一维数组exp(1jnpsi) 表示第 n 个单元相对参考单元的相位偏移。最后取模并除以 N 完成归一化。如果想写成更精简的闭合形式可以用abs(sin(N*psi/2) ./ (N*sin(psi/2)))但注意 psi 接近 0 和 2π 的整数倍时会出现 0/0需要额外做数值保护。工程上我更推荐循环写法因为后续加相位扫描、加幅度锥削时改动更小。均匀激励的阵列副瓣电平理论上固定在 -13.26 dB和单元数无关但主瓣宽度会随 N 增加而变窄。这是判断代码是否写对的第一条标准。3.3 合成子阵方向图并转换成分贝值把单元方向图和阵因子乘起来就得到完整的子阵方向图E_total E_unit .* AF; E_total_db 20 * log10(abs(E_total) eps);后面加上 eps 是为了避免 0 dB 以下的极小值在取对数时出现 -Inf影响后续绘图。要注意的是阵列做波束扫描以后单元方向图 E_unit 基本不变变的只有 AF所以如果看到扫描后方向图整体被压下去那不是代码错了而是单元方向图的“罩子”作用。不同单元数下的均匀线阵方向图特征如下单元数 N间距 d第一副瓣电平3 dB 波束宽度约4λ/2-13.26 dB25.4°8λ/2-13.26 dB12.7°16λ/2-13.26 dB6.4°这个表在初步验证 matlab 仿真结果时非常有用。看到 -13.26 dB 附近的副瓣说明阵因子代码基本正确如果副瓣变成了 -6 dB 左右大概率是归一化时少除了 N。4. 复现 banbozhenzi.m 的实际流程从工作区到 3.2figure.jpg4.1 脚本结构按四段式组织这一类半波振子仿真脚本结构上通常分成四段第一段清空工作区并声明参数第二段计算单元方向图第三段生成阵因子第四段合成并绘图。banbozhenzi.m 虽然文件名简短但跑通后输出的 3.2figure.jpg 正是第四段绘图的产物。我重写这类脚本时会先把骨架立出来clear; close all; clc; % 参数声明 freq 2.4e9; lambda 3e8 / freq; N 8; d lambda / 2; theta linspace(0, pi, 721); % 单元方向图 theta_safe theta; theta_safe(abs(sin(theta)) 1e-6) NaN; E_unit cos(pi/2 * cos(theta_safe)) ./ sin(theta_safe); E_unit(abs(sin(theta)) 1e-6) 0; E_unit E_unit / max(abs(E_unit)); % 阵因子 AF zeros(size(theta)); for n 0:N-1 AF AF exp(1j * n * 2*pi*d/lambda * cos(theta)); end AF abs(AF) / N; % 合成方向图 E_total E_unit .* AF; E_total_db 20 * log10(abs(E_total) eps);这段代码没有任何工具箱依赖只要有 MATLAB 基础环境就能运行。有些教材会把 theta 范围写成 -90° 到 90°那是把坐标原点放在阵面法向这里用 0 到 π是因为单元方向图的解析式天然以振子轴线为基准混用两套坐标容易把方向图画反。4.2 极坐标绘图与坐标轴参数生成 3.2figure.jpg 这一步推荐用 polarplot 而不是 plot因为天线方向图在极坐标下能直接读出主瓣指向和零陷夹角figure(Name, 半波振子子阵方向图); polarplot(theta, E_total_db, LineWidth, 1.5); rlim([-40 0]); % 动态范围 40 dB低于 -40 dB 的旁瓣压到圆外 thetalim([0 180]); % 只显示 0 到 180 度对应半空间 grid on;rlim 参数是这里最关键的一步。方向图全动态范围可能超过 60 dB直接绘图会让主瓣占满整个极坐标半径副瓣细节全部被压扁。工程上一般取 30 至 40 dB 动态范围既能看清主瓣和第一副瓣又能避免噪声底抬起来干扰判断。thetalim 是否要做限制取决于你要看全空间还是只看前半空间。半波振子在轴线方向本来就存在零点所以 0° 附近自然会凹下去不用刻意放大。4.3 从方向图数据里读取副瓣和波束宽度运行脚本后除了看图还应该直接从数组里算指标。手动把鼠标放到图上读数不同版本的 MATLAB 图形窗口读数习惯不同效率太低。可以这样自动提取E_total_db E_total_db - max(E_total_db); % 重新归一化到 0 dB % 主瓣位置 [~, idx_max] max(E_total_db); theta_deg rad2deg(theta); theta_main theta_deg(idx_max); % 3 dB 波束宽度 idx_hp find(E_total_db -3); HPBW theta_deg(max(idx_hp)) - theta_deg(min(idx_hp)); % 第一副瓣电平把主瓣周围 ±5 度剔除后取最大值 mainlobe_mask abs(theta_deg - theta_main) 5; pattern_for_sll E_total_db; pattern_for_sll(mainlobe_mask) -inf; SLL max(pattern_for_sll);说明HPBW 计算使用的是角度网格步长步长越大结果越粗糙如果 theta 用 721 点步长 0.25°波束宽度读数精度已经足够。SLL 的剔除宽度取 ±5° 是经验值只适用于窄波束情况波束扫描到 60° 以后主瓣宽度明显加大需要按主瓣宽度的一半动态设置剔除区间。5. 参数如何影响扫描、栅瓣和副瓣调参实战5.1 阵元间距 d 与栅瓣边界把 d 从 λ/2 加大到 λ甚至 1.5λ可见区会出现栅瓣也就是第二个幅度与主瓣相同的主瓣。栅瓣出现的条件是[ \sin\theta_{gl} \frac{m\lambda}{d}, \quad m \pm 1, \pm 2, \dots ]当这个方程在 [-1, 1] 范围内有解时栅瓣就进入方向图。下表总结了不同间距下的表现阵元间距 d栅瓣情况可见区内主瓣数适用场景λ/4无栅瓣1单元数多、尺寸敏感λ/2无栅瓣1常规窄波束设计λ栅瓣在端射方向2不适合常规设计1.5λ栅瓣在 ±41.8°3仅特殊分集场景仿真时可以直接改 d 参数重跑脚本看到第二主瓣出现时不用怀疑代码那是间距过大导致的物理现象。需要注意的是均匀线阵在 d λ 时栅瓣恰好出现在 θ 90° 和 θ 90° 对称位置d 继续加大栅瓣会进入前半空间并与主瓣争夺能量。5.2 用相位梯度实现波束扫描实际子阵设计里经常要求主瓣指向某个特定方向 θ_scan这时需要对每个单元施加相位补偿theta_scan deg2rad(30); % 期望主瓣指向 30° % 计算相邻单元之间的扫描相位差 phase_step 2 * pi * d / lambda * cos(theta_scan); w exp(-1j * phase_step * (0:N-1)); % 共轭相位使主瓣向指定方向移动 AF zeros(size(theta)); for n 1:N AF AF w(n) * exp(1j * 2*pi*d/lambda * cos(theta) * (n-1)); end AF abs(AF) / max(abs(AF));相位为什么要取共轭因为阵因子求和里第 n 个单元在 θ_scan 方向本来就自带一个正相位exp(j n k d cosθ_scan)想让它们在该方向同相叠加就给每个单元乘一个大小相等、符号相反的相位。这样在 θ_scan 方向所有单元贡献同相相加而其他方向相位对消。扫描以后有两个边界条件要同时检查一是主瓣扫描到 60° 以上时单元方向图增益明显下降这是阵元的“罩子效应”不是阵列故障二是 d 必须满足[ \frac{d}{\lambda} \le \frac{1}{1 |\cos\theta_{scan}|} ]例如扫描到 60°cosθ 0.5此时 d/λ 必须小于 0.667超过这个值就可能看到栅瓣进入可见区。这是调 d 参数时必须对照的判据。5.3 幅度锥削对副瓣的抑制均匀分布的副瓣电平固定 -13.26 dB想进一步压副瓣就得做幅度锥削。最常用的是道尔夫-切比雪夫加权它能在给定副瓣电平下让主瓣宽度最窄N 8; SLL_desired -30; % 目标副瓣电平 w_amp chebwin(N, -SLL_desired); % 输入正值例如 30 % 把幅度加权乘进阵因子 AF zeros(size(theta)); for n 1:N AF AF w_amp(n) * exp(1j * 2*pi*d/lambda * cos(theta) * (n-1)); end AF abs(AF) / max(abs(AF));chebwin 需要 Signal Processing Toolbox如果没有该工具箱可以用汉明窗做近似代价是副瓣电平会均匀下降到约 -40 dB但主瓣会展宽 1.3 倍左右。做 matlab 仿真时副瓣和波束宽度之间的取舍关系比具体窗函数更重要代码写对的前提下副瓣压得越低主瓣必然越宽这不是数组问题。6. 子阵方向图的快速校验与全波对拍方法6.1 用 3 dB 波束宽度反向验证方向性仿真结果对不对不能只看主瓣是不是在 0°。均匀线阵的法向波束宽度有理论近似公式[ HPBW \approx 0.886 \frac{\lambda}{N d} ]以 N 8、d λ/2 为例HPBW ≈ 0.886 / 4 rad ≈ 12.7°。代码里算出的 HPBW 如果明显偏离这个值优先检查 theta 网格是否包含了完整的可见区其次检查阵因子归一化是否正确。扫描到 θ_scan 后波束宽度按 1/cosθ_scan 展宽如果扫描到 60° 还不变宽说明相位加权代码写错了。6.2 与 HFSS、CST 全波结果对拍方向图乘积定理忽略了单元互耦而全波仿真会包含互耦、边缘截断效应和实际馈电结构的影响。工程上最常见的做法是先跑 HFSS 得到一个单元的辐射方向图和 S 参数再用 MATLAB 做阵因子合成。如果两者在副瓣区域差异超过 2 dB通常是互耦导致的单元方向图畸变在起作用。修正办法不是把全波结果直接当作单元方向图而是提取阵列中每个单元的有源方向图再重新合成这是阵列天线仿真中更贴近真实的一步。6.3 把脚本结果导出成对比图多次调参后最好把不同扫描角度的方向图叠在一张图里保存。用 R2020a 之后的版本可以这样导出figure(Name, 多波束对比); hold on; for scan_deg [0 15 30 45] % 此处重新计算 AF 并合成 E_total_db polarplot(theta, E_total_db, DisplayName, [num2str(scan_deg) °]); end rlim([-40 0]); legend(show, Location, southoutside); exportgraphics(gcf, subarray_scan_compare.png, Resolution, 300);exportgraphics 比 saveas 更稳定不会出现白色边距裁不掉的问题。这样保存出来的图可以直接贴进仿真报告审阅人一眼就能看出波束扫描过程中主瓣宽度和副瓣电平的变化趋势。本文还有配套的精品资源点击获取 半波振子子阵方向图MATLAB仿真：方向图乘积定理与参数解析
简介MATLAB仿真技术在半波振子天线子阵方向图研究中具有重要应用。此份zip资源包面向天线设计与电磁仿真初学者、通信工程专业学生及相关工程人员提供完整的半波振子子阵方向图仿真实现帮助理解振子天线辐射特性和子阵方向图生成原理。压缩包共2个文件包含1个MATLAB脚本文件和1个仿真结果图片整体大小仅45KB脚本涵盖天线参数设置、子阵构建、仿真计算与方向图绘制等关键环节图片直观展示了半波振子子阵的辐射强度分布。目前已有613人学习下载。通过运行该脚本读者可快速获得可视化方向图并能通过调整振子长度、馈电点位置、子阵排列等参数深入探究天线性能为实际阵列设计与教学科研提供可复用的仿真参考也可为赋形波束、宽角扫描等复杂辐射模式研究提供基础。1. 半波振子子阵方向图仿真单元方向图和阵因子到底谁说了算平时画半波振子的方向图读到的几乎都是那个“8”字形截面一旦把它拼成子阵整个方向图就会多出若干瓣和零陷主瓣也被压窄。真正决定主瓣方向和零点位置的往往不是单元方向图而是子阵阵因子单根振子的单元方向图只相当于给阵因子加了一个“罩子”。banbozhenzi.zip 里的 banbozhenzi.m做的就是把这套方向图乘积定理落成 MATLAB 代码最后生成类似 3.2figure.jpg 的极坐标子阵方向图。它适合三类人正在看阵列天线教材的学生、要做波束扫描的射频工程师、以及想快速对照全波仿真结果的系统设计人员。这个压缩包里的程序不长核心逻辑就是“单元方向图 × 阵因子”但每一行参数设置都直接反映天线工程习惯。下面从单元建模开始一步步拆开这套 matlab 仿真。2. 半波振子的单元方向图建模从解析式到 MATLAB 代码2.1 为什么半波振子单元方向图可以直接公式化半波振子是最典型的线天线振子长度等于半个工作波长电流沿振子呈近似正弦分布、末端电流趋于零。它的远场方向性有解析表达式不需要像微带贴片那样非得靠全波求解器才能得到稳定结果。若振子沿 z 轴放置球坐标系下的归一化电场方向图写作[ E_{unit}(\theta) \frac{\cos\left(\frac{\pi}{2}\cos\theta\right)}{\sin\theta} ]这个式子的特点是水平面θ 90°方向值最大归一化后为 1沿振子轴线方向θ 0° 和 180°辐射为零所以方向图是典型的“8”字形。做子阵仿真时单元方向图通常提前算好并归一化后续直接和阵因子相乘。工程里常有人直接把这段公式抄进 MATLAB结果画出来的方向图在轴向出现一个很大的跳变原因是分母 sinθ 在 θ 接近 0 或 π 时趋近于零分子也趋近于零代码里直接相除会得到 NaN。这不是物理问题而是数值处理问题需要在计算网格上特殊处理。2.2 参数初始化把频率、波长、阵元数写清楚打开 banbozhenzi.m你会发现第一段几乎都在做参数声明。这类脚本我习惯按下面的顺序组织c 2.99792458e8; % 光速单位 m/s freq 2.4e9; % 工作频率单位 Hz这里以 2.4 GHz 为例 lambda c / freq; % 波长约 0.125 m N 8; % 子阵单元个数 d lambda / 2; % 单元间距取半波长 theta linspace(0, pi, 721); % 极角网格从 0 到 pi分辨率 0.25°参数说明N 是参与合成的振子数量它直接决定阵因子波束宽度d 是相邻振子相位中心的距离通常取半波长取太大会在可见区出现栅瓣这个后面第五章专门讲。theta 网格建议取 721 或 3601 点点数太少会导致 -3 dB 波束宽度读数粗糙点数太多则没必要因为解析公式算起来很快瓶颈只在绘图。2.3 单元方向图计算的数值处理单元方向图的计算代码可以这样写theta_safe theta; % 将轴向sin theta 接近0置为 NaN先避开除零 theta_safe(abs(sin(theta)) 1e-6) NaN; E_unit cos(pi/2 * cos(theta_safe)) ./ sin(theta_safe); % 轴向辐射物理上为 0这里把之前置 NaN 的位置补回 0 E_unit(abs(sin(theta)) 1e-6) 0; % 归一化让最大值等于 1 E_unit E_unit / max(abs(E_unit));逻辑说明第一处替换是为了让 MATLAB 在计算 NaN 时不会报警告第二处替换是利用半波振子轴方向电场为零的物理结论把轴向值强制设为 0。这里的阈值 1e-6 不是固定标准如果 theta 网格点足够密取 1e-8 也一样关键是保证 sinθ 不为零。有人在网上抄的代码喜欢写E_unit(isnan(E_unit)) 1对半波振子来说这是错的。轴向零点是半波振子的固有特性不是数值异常改成 1 会把后续子阵方向图的零点全部抬高副瓣电平直接失真。3. 方向图乘积定理子阵方向图的合成逻辑3.1 阵因子和单元方向图之间的边界在哪里子阵方向图不是把每根振子的方向图简单相加而是先算阵因子再和单元方向图相乘。方向图乘积定理的适用条件是单元之间互耦可忽略、阵列处于远场观察区域、各单元方向图相同。它的表达式很简洁[ F_{total}(\theta) F_{unit}(\theta) \times AF(\theta) ]AF 就是阵因子它只跟阵列几何布局、单元间距、激励幅度和相位有关。对于沿 z 轴等间距排布的 N 元线阵阵因子可写成[ AF(\theta) \sum_{n0}^{N-1} w_n e^{j n k d \cos\theta} ]其中 k 2π/λ 是波数w_n 是第 n 个单元的复激励。这个公式的物理含义是每个单元到远场观察点的路径差不同产生的相位差叠加后形成干涉图样。均匀激励时 w_n 1阵因子有解析闭合形式但用循环求和的方式写更直观也方便以后加幅度锥削。3.2 均匀线阵的法向阵因子实现先看不扫描、无幅度加权的情况MATLAB 代码可以这样写psi 2 * pi * d / lambda * cos(theta); % 相邻单元相对观察方向的空间相位差 AF zeros(size(theta)); for n 0:N-1 AF AF exp(1j * n * psi); end AF abs(AF) / N; % 除以 N使最大值为 1代码逻辑循环变量 n 从 0 开始对应第 0 个阵元取参考相位psi 是一个随 theta 变化的一维数组exp(1jnpsi) 表示第 n 个单元相对参考单元的相位偏移。最后取模并除以 N 完成归一化。如果想写成更精简的闭合形式可以用abs(sin(N*psi/2) ./ (N*sin(psi/2)))但注意 psi 接近 0 和 2π 的整数倍时会出现 0/0需要额外做数值保护。工程上我更推荐循环写法因为后续加相位扫描、加幅度锥削时改动更小。均匀激励的阵列副瓣电平理论上固定在 -13.26 dB和单元数无关但主瓣宽度会随 N 增加而变窄。这是判断代码是否写对的第一条标准。3.3 合成子阵方向图并转换成分贝值把单元方向图和阵因子乘起来就得到完整的子阵方向图E_total E_unit .* AF; E_total_db 20 * log10(abs(E_total) eps);后面加上 eps 是为了避免 0 dB 以下的极小值在取对数时出现 -Inf影响后续绘图。要注意的是阵列做波束扫描以后单元方向图 E_unit 基本不变变的只有 AF所以如果看到扫描后方向图整体被压下去那不是代码错了而是单元方向图的“罩子”作用。不同单元数下的均匀线阵方向图特征如下单元数 N间距 d第一副瓣电平3 dB 波束宽度约4λ/2-13.26 dB25.4°8λ/2-13.26 dB12.7°16λ/2-13.26 dB6.4°这个表在初步验证 matlab 仿真结果时非常有用。看到 -13.26 dB 附近的副瓣说明阵因子代码基本正确如果副瓣变成了 -6 dB 左右大概率是归一化时少除了 N。4. 复现 banbozhenzi.m 的实际流程从工作区到 3.2figure.jpg4.1 脚本结构按四段式组织这一类半波振子仿真脚本结构上通常分成四段第一段清空工作区并声明参数第二段计算单元方向图第三段生成阵因子第四段合成并绘图。banbozhenzi.m 虽然文件名简短但跑通后输出的 3.2figure.jpg 正是第四段绘图的产物。我重写这类脚本时会先把骨架立出来clear; close all; clc; % 参数声明 freq 2.4e9; lambda 3e8 / freq; N 8; d lambda / 2; theta linspace(0, pi, 721); % 单元方向图 theta_safe theta; theta_safe(abs(sin(theta)) 1e-6) NaN; E_unit cos(pi/2 * cos(theta_safe)) ./ sin(theta_safe); E_unit(abs(sin(theta)) 1e-6) 0; E_unit E_unit / max(abs(E_unit)); % 阵因子 AF zeros(size(theta)); for n 0:N-1 AF AF exp(1j * n * 2*pi*d/lambda * cos(theta)); end AF abs(AF) / N; % 合成方向图 E_total E_unit .* AF; E_total_db 20 * log10(abs(E_total) eps);这段代码没有任何工具箱依赖只要有 MATLAB 基础环境就能运行。有些教材会把 theta 范围写成 -90° 到 90°那是把坐标原点放在阵面法向这里用 0 到 π是因为单元方向图的解析式天然以振子轴线为基准混用两套坐标容易把方向图画反。4.2 极坐标绘图与坐标轴参数生成 3.2figure.jpg 这一步推荐用 polarplot 而不是 plot因为天线方向图在极坐标下能直接读出主瓣指向和零陷夹角figure(Name, 半波振子子阵方向图); polarplot(theta, E_total_db, LineWidth, 1.5); rlim([-40 0]); % 动态范围 40 dB低于 -40 dB 的旁瓣压到圆外 thetalim([0 180]); % 只显示 0 到 180 度对应半空间 grid on;rlim 参数是这里最关键的一步。方向图全动态范围可能超过 60 dB直接绘图会让主瓣占满整个极坐标半径副瓣细节全部被压扁。工程上一般取 30 至 40 dB 动态范围既能看清主瓣和第一副瓣又能避免噪声底抬起来干扰判断。thetalim 是否要做限制取决于你要看全空间还是只看前半空间。半波振子在轴线方向本来就存在零点所以 0° 附近自然会凹下去不用刻意放大。4.3 从方向图数据里读取副瓣和波束宽度运行脚本后除了看图还应该直接从数组里算指标。手动把鼠标放到图上读数不同版本的 MATLAB 图形窗口读数习惯不同效率太低。可以这样自动提取E_total_db E_total_db - max(E_total_db); % 重新归一化到 0 dB % 主瓣位置 [~, idx_max] max(E_total_db); theta_deg rad2deg(theta); theta_main theta_deg(idx_max); % 3 dB 波束宽度 idx_hp find(E_total_db -3); HPBW theta_deg(max(idx_hp)) - theta_deg(min(idx_hp)); % 第一副瓣电平把主瓣周围 ±5 度剔除后取最大值 mainlobe_mask abs(theta_deg - theta_main) 5; pattern_for_sll E_total_db; pattern_for_sll(mainlobe_mask) -inf; SLL max(pattern_for_sll);说明HPBW 计算使用的是角度网格步长步长越大结果越粗糙如果 theta 用 721 点步长 0.25°波束宽度读数精度已经足够。SLL 的剔除宽度取 ±5° 是经验值只适用于窄波束情况波束扫描到 60° 以后主瓣宽度明显加大需要按主瓣宽度的一半动态设置剔除区间。5. 参数如何影响扫描、栅瓣和副瓣调参实战5.1 阵元间距 d 与栅瓣边界把 d 从 λ/2 加大到 λ甚至 1.5λ可见区会出现栅瓣也就是第二个幅度与主瓣相同的主瓣。栅瓣出现的条件是[ \sin\theta_{gl} \frac{m\lambda}{d}, \quad m \pm 1, \pm 2, \dots ]当这个方程在 [-1, 1] 范围内有解时栅瓣就进入方向图。下表总结了不同间距下的表现阵元间距 d栅瓣情况可见区内主瓣数适用场景λ/4无栅瓣1单元数多、尺寸敏感λ/2无栅瓣1常规窄波束设计λ栅瓣在端射方向2不适合常规设计1.5λ栅瓣在 ±41.8°3仅特殊分集场景仿真时可以直接改 d 参数重跑脚本看到第二主瓣出现时不用怀疑代码那是间距过大导致的物理现象。需要注意的是均匀线阵在 d λ 时栅瓣恰好出现在 θ 90° 和 θ 90° 对称位置d 继续加大栅瓣会进入前半空间并与主瓣争夺能量。5.2 用相位梯度实现波束扫描实际子阵设计里经常要求主瓣指向某个特定方向 θ_scan这时需要对每个单元施加相位补偿theta_scan deg2rad(30); % 期望主瓣指向 30° % 计算相邻单元之间的扫描相位差 phase_step 2 * pi * d / lambda * cos(theta_scan); w exp(-1j * phase_step * (0:N-1)); % 共轭相位使主瓣向指定方向移动 AF zeros(size(theta)); for n 1:N AF AF w(n) * exp(1j * 2*pi*d/lambda * cos(theta) * (n-1)); end AF abs(AF) / max(abs(AF));相位为什么要取共轭因为阵因子求和里第 n 个单元在 θ_scan 方向本来就自带一个正相位exp(j n k d cosθ_scan)想让它们在该方向同相叠加就给每个单元乘一个大小相等、符号相反的相位。这样在 θ_scan 方向所有单元贡献同相相加而其他方向相位对消。扫描以后有两个边界条件要同时检查一是主瓣扫描到 60° 以上时单元方向图增益明显下降这是阵元的“罩子效应”不是阵列故障二是 d 必须满足[ \frac{d}{\lambda} \le \frac{1}{1 |\cos\theta_{scan}|} ]例如扫描到 60°cosθ 0.5此时 d/λ 必须小于 0.667超过这个值就可能看到栅瓣进入可见区。这是调 d 参数时必须对照的判据。5.3 幅度锥削对副瓣的抑制均匀分布的副瓣电平固定 -13.26 dB想进一步压副瓣就得做幅度锥削。最常用的是道尔夫-切比雪夫加权它能在给定副瓣电平下让主瓣宽度最窄N 8; SLL_desired -30; % 目标副瓣电平 w_amp chebwin(N, -SLL_desired); % 输入正值例如 30 % 把幅度加权乘进阵因子 AF zeros(size(theta)); for n 1:N AF AF w_amp(n) * exp(1j * 2*pi*d/lambda * cos(theta) * (n-1)); end AF abs(AF) / max(abs(AF));chebwin 需要 Signal Processing Toolbox如果没有该工具箱可以用汉明窗做近似代价是副瓣电平会均匀下降到约 -40 dB但主瓣会展宽 1.3 倍左右。做 matlab 仿真时副瓣和波束宽度之间的取舍关系比具体窗函数更重要代码写对的前提下副瓣压得越低主瓣必然越宽这不是数组问题。6. 子阵方向图的快速校验与全波对拍方法6.1 用 3 dB 波束宽度反向验证方向性仿真结果对不对不能只看主瓣是不是在 0°。均匀线阵的法向波束宽度有理论近似公式[ HPBW \approx 0.886 \frac{\lambda}{N d} ]以 N 8、d λ/2 为例HPBW ≈ 0.886 / 4 rad ≈ 12.7°。代码里算出的 HPBW 如果明显偏离这个值优先检查 theta 网格是否包含了完整的可见区其次检查阵因子归一化是否正确。扫描到 θ_scan 后波束宽度按 1/cosθ_scan 展宽如果扫描到 60° 还不变宽说明相位加权代码写错了。6.2 与 HFSS、CST 全波结果对拍方向图乘积定理忽略了单元互耦而全波仿真会包含互耦、边缘截断效应和实际馈电结构的影响。工程上最常见的做法是先跑 HFSS 得到一个单元的辐射方向图和 S 参数再用 MATLAB 做阵因子合成。如果两者在副瓣区域差异超过 2 dB通常是互耦导致的单元方向图畸变在起作用。修正办法不是把全波结果直接当作单元方向图而是提取阵列中每个单元的有源方向图再重新合成这是阵列天线仿真中更贴近真实的一步。6.3 把脚本结果导出成对比图多次调参后最好把不同扫描角度的方向图叠在一张图里保存。用 R2020a 之后的版本可以这样导出figure(Name, 多波束对比); hold on; for scan_deg [0 15 30 45] % 此处重新计算 AF 并合成 E_total_db polarplot(theta, E_total_db, DisplayName, [num2str(scan_deg) °]); end rlim([-40 0]); legend(show, Location, southoutside); exportgraphics(gcf, subarray_scan_compare.png, Resolution, 300);exportgraphics 比 saveas 更稳定不会出现白色边距裁不掉的问题。这样保存出来的图可以直接贴进仿真报告审阅人一眼就能看出波束扫描过程中主瓣宽度和副瓣电平的变化趋势。本文还有配套的精品资源点击获取