纯方位无源定位：无人机圆形与锥形编队建模与调整

简介：2022年全国大学生数学建模竞赛本科组B题河北省一等奖论文，聚焦无人机纯方位无源定位与编队位置调整，面向准备数模竞赛、学习定位建模及MATLAB实现的高校学生。内容围绕圆形编队三问展开：由方位信息与余弦定理反解实际位置，借“三圆共点”讨论增加发射信号无人机后的可定位性，并以0.01米步长遍历完成偏差位置调整，同时延伸到锥形编队的定位与调整方案。资源含1个PDF文件，压缩包约2.63MB，完整呈现摘要、问题重述、模型假设、符号说明、模型建立求解及MATLAB编程思路，可作为赛题复盘、论文写作和算法复现的参考。已有8041人学习，适合需要理解纯方位无源定位建模流程、枚举搜索与数值计算细节的读者查阅。

1. 从十架无人机的圆阵说起：纯方位无源定位要解决什么

十架无人机在夜空里排成一个圆，一架钉在圆心，九架沿半径一百米的圆周均匀铺开，相邻夹角四十度。地面站只报方位角，不报距离——这就是纯方位无源定位要处理的局面，也是全国大学生数学建模竞赛本科组 B 题整套方案的起点。难点不在飞控，在几何：被动接收信号的那架无人机，怎么靠三个夹角把自己算回本来该在的位置；发射信号的机器又要开几台才够用。这套方案面向两类人：要交数学建模课设、毕设的学生，和需要在编队控制里落一套可复现定位流程的从业者。下面按问题拆解、模型建立、代码落地、避坑、验收的顺序往下推。

2. 赛题拆解与建模约定：圆形与锥形两套编队各自要交什么

2.1 三个小问的交付边界

问题一的三小问不是三个独立题目，是一条递进的链条，任何一环没交代清楚，后面就接不上。第一小问要求建立被动接收信号无人机的定位模型，这里有个容易被含糊过去的点：题目说的“定位”其实是两个方向，一是根据接收到的方位信息反解出自己当前实际在哪儿，二是算出自己应该到哪儿。很多方案只做前者，答辩一问“你的目标位置怎么来的”就卡住了。所以我把定位模型拆成“实际位置反解”和“理想位置求解”两个子模型，分别对应两组不同的已知条件。

第二小问是一个数量下界问题：0 号和 1 号已经发射信号，还需要额外增加几架发射信号无人机，才能让其余无人机的定位唯一。结论是一架就够。这个结论的证明方式不是纯几何推导，而是“枚举—匹配”：把所有可能的发射机组合和接收机组合铺开算一遍，看哪一组的角度信息与实测偏差最小。第三小问才是完整的调整方案，给定所有无人机带偏差的初始位置，用 0 号加圆周上最多三架发射信号，让其余无人机一轮轮调回理想位置。

问题二把队形从圆形换成锥形，直线方向上相邻两架间距相等，相邻三架构成等边三角形。定位原理不变，变的是拓扑：圆形编队里任意两架发射机加圆心就能覆盖全周，锥形编队里会出现一个盲区——某架无人机恰好落在三架已定位无人机构成的三角形边的延长线上时，只能收到一个有效方位信息，位置仍然解不出来。这个盲区的处理方式决定了第二轮、第三轮怎么排，是问题二真正的采分点。

提示：把“定位实际位置”和“定位目标位置”分开写进模型假设，后面所有公式的符号体系才不会打架。

2.2 坐标系、符号与五条假设怎么定

坐标系选择看着是小事，实际决定了后面公式的复杂度。我一般把 0 号无人机（圆心那架）放在原点，把 1 号固定在 x 轴正方向上，圆周九架按逆时针每四十度排开。这样任意一架圆周无人机的理想坐标就是 $(r\cos 40(k-1)^\circ,\ r\sin 40(k-1)^\circ)$，不用做任何旋转平移，余弦定理里出现的边长也能直接用编号差算出来。

符号体系必须在动笔推导之前定死，中途改符号是返工的头号来源。下面这套是我实际用的版本：

模型假设只留五条，多了是给自己挖坑：无人机飞行不受外界干扰；接收到的方向信息准确无误；每次位置调整耗时不计，可以理解为“知道该往哪动”的瞬间就已经到位；第一小问第三问中 0 号与 1 号在飞行过程中保持初始相对位置关系不变；问题二假定初始时刻有两架无人机位置无偏差。第三条尤其关键，它把动态过程压成了静态迭代，否则整道题要变成最优控制问题，本科组的篇幅根本装不下。

编号差算边长的公式是整套代码里复用率最高的一行：$a = \sqrt{2r^2\left(1-\cos(40^\circ\cdot|i-j|)\right)}$。它的几何含义是等腰三角形顶角为 $40^\circ|i-j|$、腰长为 $r$ 时的底边。写代码时别直接抄 $\sqrt{2r^2-2r^2\cos\theta}$ 这种半展开形式，括号一多就容易漏因子，血泪经验。

2.3 MATLAB 在方案里承担的三类分工

这套方案用 MATLAB 落地，但四个程序不是同一类东西，混着写会很难维护。我把它们分成三层：第一层是正演，由编号和半径算出理想坐标与理想方位角，纯代数，没有迭代；第二层是反演，由方位角反解坐标，用solve做符号求解再double转数值，或者直接用fsolve做数值求解；第三层是搜索，包括枚举匹配和步长遍历，用双重for循环加最小值追踪。

符号求解那层有个细节值得单说。三个余弦定理方程联立，解出来通常是多组根，其中一部分是复数或者负距离，属于数学上存在但物理上无意义的分支。处理方式是取实部并筛掉负值，同时用“接收机应当落在圆形编队内部或附近”这个先验做一次合理性判断。不做这一步，程序在某些编号组合下会吐出明显离谱的坐标，你还以为模型错了。

搜索层则是另一回事。枚举 56 种情况、步长遍历 360 个方向，这些循环量都不大，MATLAB 里裸写for完全够用，不需要向量化。真正要留意的是循环里别反复调用input——把输入提到循环外，或者改写成函数参数，否则调试一次要敲几十遍回车，非常劝退。

注意：三层程序共用同一套编号规则和半径定义，把这两个量写成脚本顶部的常量，别在每个文件里各写一遍。

3. 被动定位模型：余弦定理反解与三圆共点的两条求解路径

3.1 实际位置反解：三条边三个未知量的封闭方程

接收机当前实际位置的求解思路，本质上是“三边未知、三角已知”的三角形拼接。设两架发射机为 $i$、$j$，接收机为 $k$。发射机之间距离 $a_1$ 已知，三个夹角 $\alpha_1$（接收机对 $i$、$j$ 的张角）、$\alpha_2$、$\alpha_3$ 已知，未知的是 $k_0$（接收机到圆心距离）、$a_2$、$a_3$ 三条边。三个三角形各写一次余弦定理，正好三个方程三个未知数，方程组的封闭性是有保证的。

写成代码就是三个余弦定理表达式，用符号求解器一次解出，再挑物理上合理的那组根：

% 由两架发射机编号和半径，先算出它们之间的距离 a1 r = 100; % 圆形编队半径，单位 m FY0i = 0; FY0j = 1; % 发射信号无人机编号 a1 = sqrt(2*r^2*(1 - cosd(abs(FY0i - FY0j)*40))); % 顶角 40*|i-j| 的等腰三角形底边 % 三个夹角由接收机测得，作为模型输入 alpha1 = 30; alpha2 = 35; alpha3 = 40; % 单位：度 syms k0 a2 a3 % 三个三角形各写一次余弦定理，构成封闭方程组 eq1 = (k0^2 + a3^2 - r^2) / (2*k0*a3) - cosd(alpha1); eq2 = (k0^2 + a2^2 - r^2) / (2*k0*a2) - cosd(alpha2); eq3 = (a2^2 + a3^2 - a1^2) / (2*a2*a3) - cosd(alpha3); sol = solve([eq1, eq2, eq3], [k0, a2, a3]); K0 = double(sol.k0); % 可能有多组根，需筛掉复数与负值 K0 = K0(imag(K0) == 0 & real(K0) > 0); k0 = real(K0(1)); % 取物理上合理的一组

三个方程对应三个三角形：圆心—发射机 $i$—接收机、圆心—发射机 $j$—接收机、发射机 $i$—发射机 $j$—接收机。cosd用角度制，别和cos混用，这是最常见的沉默错误——算出来不报错，但结果全错。求完 $k_0$ 之后，极角由 $k_0$ 与已知边的夹角关系推出，再转成直角坐标 $x=k_0\cos\theta$、$y=k_0\sin\theta$。这一步得到的坐标就是“实际位置”，它和理想坐标的差值，正是后面调参要消掉的东西。

3.2 理想位置：定角对定弦的“眼睛”与“葫芦”轨迹

理想位置的求解走的是另一条路，也是整道题最漂亮的一段几何。一个动点到两个定点的夹角是定值，这个动点的轨迹是上下对称的两段圆弧合成的封闭曲线。夹角为钝角时，拼出来的形状像一只眼睛；夹角为锐角时，像一只葫芦。这两种形态不是修辞，它们直接决定了轨迹方程里半径和圆心的符号，写错一个，整个交点就偏到编队外面去。

三架位置无偏差的无人机提供三条这样的轨迹，三条轨迹的公共交点就是接收机的理想位置，也就是“三圆共点”。求解时不必真的去画圆求交，继续用余弦定理把三个方位角用已知边和已知角表示即可，公式链条是：先由 $a_1$ 与 $\alpha$ 求接收机到发射机的距离，再求各边对应的角，最后合成极角。这条链条里的每一步都有明确的几何对应，出问题时能定位到具体是哪一步错了。

钝角情形的轨迹方程推导可以给出一个可验证的形式。设两定点距离为 $2c$，动点对它们的张角为钝角 $\varphi$，则轨迹满足一个配平方后的圆方程，圆心落在两定点连线的中垂线上，半径由 $\sin\varphi$ 决定。锐角情形同理，只是圆心跑到连线另一侧，两段圆弧的拼接方式变化，视觉上就成了葫芦。这两种情形都要在代码里显式分支处理，不能靠一个通式蒙混过去。

实际写方案时，我会把“轨迹形状”当成一个自检手段：程序算出来的交点如果落在两发射机连线的外侧、且离编队圆心很远，那多半是锐角/钝角分支选反了，回头查 $\cos$ 那一项的符号，比逐行打印变量快得多。

3.3 优角情形与代码落地

有一类位置关系必须单独说明：当接收机落在某个特殊方位，使得其中一个夹角变成优角时，这个优角与它和另一个角相加构成周角。因为 $\cos(360^\circ-\varphi)=\cos\varphi$，两个角的余弦值完全一致，所以无论是优角还是它的补角，代入余弦定理的结果都不变。这意味着定位模型的正确性不受无人机相对位置关系影响，不需要为优角单独写一套公式。

这个结论看着是省事，实际是救命。因为如果不做这层说明，评审很容易追问“你的模型在夹角大于 180 度时还成立吗”，而现场推一遍余弦的周期性既费时间又容易紧张出错。把它写进模型建立的正文里，一句话带过，反而显得考虑周全。

正演部分——由编号和半径算出理想坐标与理想方位角——代码量不大，但它是所有校验的基准，值得单独成一个函数：

function [xy, ang] = idealCircle(r, num) % 输入：r 编队半径，num 无人机编号（0 为圆心机） % 输出：xy 理想直角坐标，ang 该机在极坐标下的极角（度） if num == 0 xy = [0, 0]; ang = 0; return; end theta = 40 * (num - 1); % 圆周九架，相邻夹角 40 度 xy = [r*cosd(theta), r*sind(theta)]; ang = theta; end

num取值 1 到 9 时按逆时针生成九个等分点，与表 1 的理想坐标一一对应；num=0单独返回原点。这个函数在后续三个程序里都会被调用，写成独立文件而不是复制粘贴，是为了保证“理想位置”这个基准永远只有一个版本。理想坐标表如下，$r$ 取 100 时直接代入即可：

系数是 $40^\circ$ 整数倍的余弦与正弦，$0.77$ 对应 $\cos 40^\circ$，$0.64$ 对应 $\sin 40^\circ$，$0.17$ 对应 $\cos 80^\circ$，$0.98$ 对应 $\sin 80^\circ$。把这张表打印出来贴在推导稿旁边，任何时候算出的坐标与它差得离谱，都能立刻发现是公式抄错了还是编号搞反了。

4. 发射机数量与位置调整：枚举匹配加 0.01 米步长遍历

4.1 为什么一架就够：从单方位信息到三圆共点

第二小问的答案是一架，但“一架”这个结论必须给出下界证明，否则只是猜。论证路径是：不增加任何发射机时，除 0 号和 1 号外其余任何一架接收机只能收到一个方位信息，它可能的落点是整条“眼睛”或“葫芦”轨迹上的任意一点，位置不唯一，定位失败。增加一架之后，接收机同时拿到三个方位信息，对应三条轨迹，三圆共点给出唯一交点，定位成功。既然一架就能唯一确定，增加两架、三架只会让信息冗余，不会有额外收益。

这里的逻辑严密性来自“至少三个方位信息”这个必要条件。两个方位信息对应两条轨迹，两条封闭曲线一般有两个交点，仍然不唯一；只有第三条轨迹进来做筛选，才能收缩到单点。把这个论证写清楚，比直接写“枚举发现一架就够”要有说服力得多，答辩时也更抗追问。

4.2 56 种情况枚举与最小平方和匹配

枚举的规模是这样的：0 号和 1 号固定发射，从剩余八架里选一架作为额外发射机，再把八架里除自己以外的一架当作接收机，组合数是 $8\times7=56$。每一种情况都能算出一组理想的三个方位角，全部列成表备用。实际运行时，接收机给出实测的三个角，方案遍历这张表，取使三个角差值的平方和最小的那一行，对应的额外发射机编号就是答案：

% 预置：alpha1/alpha2/alpha3 为 8x8 的理想方位角表 % 输入：beta1/beta2/beta3 为接收机实测的三个方位角 FY0k = 8; % 接收信号无人机编号 best = inf; bestIdx = 0; for i = 2:9 % 遍历可作为额外发射机的八架 if i == FY0k continue; % 自己不能给自己发信号 end % 三个角差值的平方和，作为匹配代价 cost = (beta1 - alpha1(i, FY0k))^2 ... + (beta2 - alpha2(i, FY0k))^2 ... + (beta3 - alpha3(i, FY0k))^2; if cost < best best = cost; bestIdx = i; % 记录当前最优的额外发射机编号 end end

alpha1(i, FY0k)这类索引的含义是“以 $i$ 为额外发射机、以FY0k为接收机时的理想方位角”，行是发射机、列是接收机，和结果表的方向一致，别转置。cost用平方和而不是绝对值和，是因为三个角的量纲相同、权重相当，平方和对大偏差更敏感，匹配更稳。

这套匹配有个隐含前提：实测角度必须足够接近某一组理想角度。如果输入的角度本身来自一个明显偏离编队的位置，最小代价可能也很大，此时应当先判断best是否超过阈值，超了就不要硬给答案，而是提示“当前偏差过大，建议先粗调”。这条判断在原方案里没写，但实际调试时几乎必然遇到，补上会让程序健壮很多。

举两个可以手工验算的例子：半径取 10，接收机为 8 号、实测角为某组特定值时，程序输出额外发射机为 9 号；接收机为 3 号、另一组实测角时，输出为 5 号。手工验算的价值在于验证索引方向，一旦行列弄反，这两个例子会立刻给出错误的编号，比盲调快得多。

4.3 圆形编队的南北分向步进调整

第三小问的调整方案，本质是在极坐标下做一次贪心搜索。接收机并不知道自己该往哪个方向动，于是设定步长 0.01 米，围绕当前位置向四周全方位遍历，每走一步就用定位模型重新算一次自己的坐标，与理想坐标比较，误差进入 0.01 米以内就停下。之所以按南北分向描述，是因为极坐标下极径的变化对应沿半径方向的移动，极角的微调对应切向移动，两个方向组合起来就能覆盖整个邻域。

单步移动后的极径与极角更新，仍然靠余弦定理完成：当前位置 $B$ 已知、移动一个微小步长后的位置 $C$ 已知、$B$ 处的方位角已知，那么 $BC$ 长度与 $C$ 点的新极角都能解出来。写成代码就是一层while套一层for，外层控制误差是否收敛，内层枚举 360 个方向：

njl = 0.01; % 收敛阈值，单位 m r = 100; % 编队半径 theta = [0 40 80 120 160 200 240 280 320]; % 九架理想极角 rowNow = 98.5; % 当前极径（带偏差的初始值） thetaNow = 3.0; % 当前极角（带偏差的初始值） jl = inf; while jl > njl bestStep = inf; bestR = rowNow; bestT = thetaNow; for d = 0:359 % 全方位枚举移动方向 % 沿方向 d 走一个步长后的候选位置 rTry = rowNow + njl*cosd(d); tTry = thetaNow + njl*sind(d); % 与 2 号机理想位置 (r, 40 度) 的偏差 dx = rTry*cosd(tTry) - r*cosd(40); dy = rTry*sind(tTry) - r*sind(40); err = sqrt(dx^2 + dy^2); if err < bestStep bestStep = err; bestR = rTry; bestT = tTry; end end rowNow = bestR; thetaNow = bestT; % 采纳本轮最优方向 jl = bestStep; % 更新收敛判据 end

外层while的终止条件是本轮最优偏差小于 0.01 米，内层for遍历 360 个方向并记录最优点。rowNow/thetaNow是极坐标表示，转直角坐标再做差，是为了避免极角在 0 度和 360 度附近回绕时算出巨大的虚假误差——这个坑不踩一次很难想到。

按原方案，0 号和 1 号在飞行过程中保持初始相对位置不变，因此始终是位置无偏差的基准；2 号先被定位并调整到位，随后 0、1、2 三架都进入“无偏差元组”，可以用它们去定位 3 号，如此逐架推进。由于调整时间不计，整个过程可以理解为一次串行扫描，每架机只在自己被定位的那一轮动一次。

4.4 锥形编队的轮次推进与延长线盲区

锥形编队沿用同一套定位原理，但要先解决“基准从哪来”的问题。如果只有零架或一架位置无偏差，根本无法启动定位链条，所以必须假设初始时刻有两架是准的，方案里取 1 号和 5 号，并假定它们之间的间距保持不变。有了这两架，2 号就能通过两次“单角测量”被锁定：以 1、2 为发射机、5 号为接收机得一个角，以 2、5 为发射机、1 号为接收机得另一个角，三角形里已知两角，第三角随之确定，2 号的位置就解出来了，随后按 4.3 的步长遍历把它推到理想位置。

第一轮之后，1、2、5 三架构成定位基准，可以继续定位其他无人机。但这里出现一个几何盲区：如果某架无人机恰好落在 1、2、5 构成的三角形某条边的延长线上，它只能收到一个有效方位信息，位置仍然不唯一，无法定位。按原方案的推导，第二轮里 4、7、9、11、13、14 号就卡在这个盲区上，只能先定位 3、6、8、10、12、15 号并调整到位。

第三轮的处理是这套方案里最见巧思的一步：重新选取 2、6、8 作为发射信号无人机，此时剩余未定位的无人机都不在这三架构成的三角形边的延长线上，全部可获得唯一解。这一步的价值不只是“能定位”，而是把原本可能需要四五轮的推进压缩到三轮，效率提升明显。写方案时我会把每一轮的基准组合、可定位编号、盲区编号列成一张小表，评审一眼就能看出你是在按图论式的覆盖思路推进，而不是碰运气换组合。

锥形编队的理想坐标与圆形不同，按等边三角形铺排：原点一架，纵轴方向每隔一个间距布一架，两侧按 $50\sqrt{3}$、$100\sqrt{3}$ 这类倍数展开。坐标表在方案里要给出，但不必逐行解释每个数字的来源，说清“相邻三架构成等边三角形、直线方向等距”这两条生成规则就够了，剩下的交给代码算。

5. 避坑与排查：这套模型最容易翻车的五个地方

5.1 几何与方程层面的三个坑

第一个坑是角度制与弧度制混用。现象是程序不报错，但解出来的坐标离编队十万八千里，或者solve直接返回空解。原因在于 MATLAB 里cos吃弧度、cosd吃角度，题目的方位角全是度，一旦某个式子里混进一个cos，整条余弦定理链条就废了。解决办法是全程统一用cosd/sind/acosd，并在脚本开头写一行注释声明“本文件所有角度均为度”，改代码时先看这行。我自己的习惯是搜一遍全文的cos(，出现一个改一个。

第二个坑是锐角/钝角分支选反。现象是轨迹交点落在两发射机连线的外侧远处，或者三圆共点解出的位置明显偏离编队圆周。原因是“眼睛”和“葫芦”两种轨迹的圆心在连线两侧，公式里 $\sin\varphi$ 的符号一翻，圆就跑到另一边去了。解决办法是在推导阶段就把两种情形的轨迹方程分别写出来并画图确认，代码里用if显式分支，不要试图写一个通式通吃。判断依据可以取“三个方位角中钝角的个数”，大于等于两个时走“眼睛”分支，否则走“葫芦”分支。

第三个坑是三圆共点退化成两圆相交。现象是解出来的位置有两个候选，程序随便取了一个，某些编号下结果对、换一组就错。原因是有两条轨迹几乎重合，交点收缩不干净。解决办法是检查三条轨迹对应的发射机是否共线或近似共线，如果是，换一组发射机重新组合；同时在代码里加一条判据：两个候选点到圆心的距离差超过阈值时，取更接近编队半径的那个。

5.2 数值与流程层面的两个坑

第四个坑是符号求解返回多组根，取根策略粗暴。现象是同一组输入，改一下初值或换个编号，坐标就跳到完全不同的象限。原因是三个余弦方程联立后往往有两到四组解，其中包含复数根和负距离根。解决办法是先按虚部为零筛一遍，再按距离为正筛一遍，最后用“接收机应位于圆周附近”这个先验做第三次筛选。三次筛选之后通常只剩一组，如果还剩多组，说明输入的三个角本身不自洽，该回头查测量数据。

第五个坑是步长遍历陷入局部循环。现象是while跑了几万次还不收敛，误差在 0.02 米附近来回震荡。原因是 0.01 米的步长与 0.01 米的阈值太接近，走到最优点附近时每一步都会略微越过，然后又走回来。解决办法是把收敛阈值放宽到步长的 1.5 到 2 倍，或者改成“连续两轮最优偏差不再下降就退出”。这类数值细节在论文里可以不写，但代码里必须有，否则演示时当场卡住，非常尴尬。

提示：所有涉及while的循环都加一个最大迭代次数保护，跑满就退出并打印当前误差，比无限等待强。

6. 验收与答辩：用三张表把这套方案讲圆

方案做完不等于能交，验收环节要能自证。我一般准备三张表，按顺序摆在答辩稿里。第一张是理想位置对照表，列出十架无人机的理论坐标，用来证明坐标系和编号规则没写错；第二张是枚举匹配表，56 种组合的理想方位角全列出来，随便挑一行让程序复算，结果必须一致；第三张是收敛记录表，记下每架无人机从初始偏差调到 0.01 米以内所用的轮次和最终误差。三张表覆盖了“模型对不对、匹配准不准、迭代收不收”三个层面，评委问哪个都有据可查。

验收标准我按下面这套卡：

闭环测试是这套验收里最有用的一条：拿理想坐标正演出三个方位角，再喂给反解程序，看能不能还原出原来的坐标。这一步能同时抓出角度制错误、索引转置、取根策略错误三类问题，比单独测每个函数高效得多。锥形编队的验收则重点看盲区编号对不对——如果程序报告第二轮就能定位 4 号，那基本可以断定延长线的判据写漏了。

答辩有个具体技巧：不要从头讲公式，先讲那张盲区表。锥形编队第二轮为什么定不了 4、7、9、11、13、14 号，第三轮换成 2、6、8 为什么就全通了，这一问一答里既有几何直觉又有覆盖策略，讲清楚了两分钟，评委的印象分基本就拿到了。公式细节留到追问环节再展开，主动讲反而容易在符号上被绕住。

从那以后我每次做这类几何建模题，都强制走一遍“正演—反演闭环 + 盲区枚举”这两步再动笔写论文，哪怕时间紧也先跑通最小闭环，因为模型再漂亮，跑不通就是零。希望帮到你。

本文还有配套的精品资源，点击获取