一种基于梯度复数场的多旋翼无人直升机航线规划方法

1.本发明属于无人直升机航线规划领域，涉及到使用无人直升机性能和动、静态威胁等信息规划出最优飞行航线方法，具体是一种基于梯度复数场的多旋翼无人直升机航线规划方法。


背景技术：

2.到目前为止，国内外在综合复杂动态和静态环境下航线规划技术方面的研究明显增多，但距飞行航线实施规划上存在较大的差距。复杂多变环境的建模和实时无人直升机航线优化目前还没有得到很好的解决，如何能应用目标特征来规划安全性能最大的飞行航线，规划出满足任务要求、安全性能等约束的较优飞行航线是亟待解决的问题，此问题的解决对军用无人机、民用无人机都具有极其重要的现实意义。


技术实现要素：

3.(1)构建综合威胁模型t
t
，数学描述如下：
4.首先利用公式(1)得到地形威胁模型v1：
[0005][0006]
上述公式中，v1(x，y)表示横坐标x及纵坐标y处的地形威胁模型值，单位为公里；的是一种(-x
2-(y+1)2)次方指数函数；
[0007]
其次利用公式(2)～(7)得到x轴坐标和y轴坐标处人群威胁模型v2：
[0008][0009]
mu＝round(rand(3，2)-0.5)*19+30
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
[0010]
sigma＝round(rand(3，2)*40)/10+30
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)
[0011]
gm＝gmdistribution.fit(x，3，
′
options
′
，statset(
′
display
′
，
′
off
′
))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(5)
[0012]
ii＝ezsurf(@(x，y)pdf(gm，[x，y])，[-15 15]，[-15 10])
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0013]
v2＝ii.cdata
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)
[0014]
公式(2)～(7)中，[]为矩阵符号，x为人群分布三维坐标矩阵，mvnrnd、mu、sigma分别为产生人群高斯分布数据的函数及其对应的均值、协方差；gmdistribution.fit、gm、statset分别为三维混合高斯分布参数拟合函数及其创建的结构体对象和创建统计选项结构函数，options、display、off分别表示选项、显示和关闭，ezsurf、ii、ii.cdata分别表示创建一个的表面图及其句柄和三维表面数据矩阵，pdf为求概率密度的函数，@为构造函数表达式的符号；
[0015]
最后利用公式(8)得到融合静态地形威胁和动态人群威胁的综合威胁模型v：
[0016]
v＝v1+v2
ꢀꢀꢀꢀ
(8)
[0017]
公式(8)中，v为综合威胁模型；
[0018]
(2)计算威胁模型梯度值及其对应的复数场pf：
[0019]
首先利用公式(9)～(10)得到x轴坐标梯度分量px、y轴坐标梯度分量py和复数场pg：
[0020]
[px，py]＝gradient(v，0.2，0.2)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)
[0021]
pg＝px+py*i
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(10)
[0022]
公式(9)～(10)中，gradient为求梯度的函数，px、py、pg分别为梯度横、纵坐标分量及其对应的梯度复数场，i为序数表示符；
[0023]
其次利用公式(11)～(12)得到x轴坐标梯度分量px、y轴坐标梯度分量py和梯度矢量复数场pf：
[0024]
iq＝quiver(1：60，1：60，v)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)
[0025]
pf＝iq.udata+iq.vdata*j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)
[0026]
公式(9)～(10)中，quiver是matlab中绘制二维矢量场的函数，iq.udata、iq.vdata、pf分别为梯度矢量场横、纵坐标分量及其对应的梯度矢量复数场，h为序数表示符；
[0027]
(3)遴选适飞区域，构建飞行航迹点集fpc：
[0028]
利用公式(13)～(16)遴选适飞区域x、y轴坐标集iqc及梯度矢量场横、纵坐标分量比值集iqvm，并随机选出飞行航迹点集fpc：
[0029][0030][0031]
id＝floor(rand(1，90)*length(iqc)+1)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(15)
[0032]
fpc＝iqc(id)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(16)
[0033]
公式(13)～(16)中，iqc、iqvm分别表示遴选适飞区域平面坐标集及其对应位置的梯度矢量场横、纵坐标分量比值集iqvm，inf表示无穷大，abs、floor、length、rand分别为求绝对值、向下取整数、求最大数组维度的长度和产生数值范围在0～1之间随机矩阵的函数，id、fpc分别为随机选择的飞行航迹点在iqc中的序号和飞行航迹点集；
[0034]
(4)执行航迹点优选函数，生成初始最优飞行航线：
[0035]
首先利用公式(17)～(20)建立目标区域引导模型ze及其梯度x、y轴坐标轴分量dx、dy，构造飞行航线航向引导复数场pfz：
[0036]
[mx，my]＝meshgrid(0∶60)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(17)
[0037]
ze＝-((mx-60)2+(my-60)2)/2
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(18)
[0038]
[dx，dy]＝gradient(ze，1，1)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(19)
[0039]
pfz＝dx+dy*i
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(20)
[0040]
公式(17)～(20)中，mx、my分别表示x、y轴坐标离散矩阵，meshgrid为生成二维网格坐标函数，dx、dy、pfz分别表示目标区域引导模型梯度x、y轴坐标轴分量和飞行航线航向引导复数场；
[0041]
利用公式(21)～(26)定义目标区域引导函数fun及其对x、y轴坐标得导数funx、
funy，并初始化航迹点初值坐标x0、y0：
[0042]
fum＝inline(
′‑
ze
′
，’x
′
，
′y′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(21)
[0043]
dfunx＝inline(
′
x
′
，
′
x
′
，
′y′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(22)
[0044]
dfuny＝inline(
′y′
，
′
x
′
，
′y′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(23)
[0045]
indx＝floor(rand(1)*length(iqc))+1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(24)
[0046]
x0＝iqc(indx，1)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(25)
[0047]
，y0＝iqc(indx，2)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(26)
[0048]
公式(21)～(26)中，fun、funx、funy分别表示目标区域引导函数和其对x、y的导数，inline为定义函数函数，indx、x0、y0分别随机生成的初始化航迹点在iqc中的序号及其横、纵坐标轴；
[0049]
利用公式(27)～(28)利用改进的梯度下降算法tidu1对适飞区域内的飞行航迹点集进行选取，组成优化的航迹点集合fp：
[0050]
[xtd，ytd，n，point]＝tidu1(fun，dfunx，dfuny，x0，y0，epsilon1，lambda1，iqc)
ꢀꢀꢀꢀꢀ
(27)
[0051]
fp＝[[xb yb zb]；[point v(point)]；[xe ye ze]]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(28)
[0052]
公式(27)～(28)中，xtd、ytd、n、point分别表示选出的最优航迹点的横、纵坐标及航迹点个数、航迹点平面坐标集合，tidu1为改进的根据目标引导场从适飞区域选择航迹点的梯度下降算法，fp为选出航线点组成的飞行航线，epsilon1、lambda1分别为精度和步长/更新率，xb、yb、zb分别为航线起始点的x坐标、y坐标、z坐标，xe、ye、ze分别为航线终止点的x坐标、y坐标、z坐标；
[0053]
(5)采用分维b-spline方法改进最优飞行航线，实现航线可飞行性优化：
[0054]
利用公式(29)～(32)分别应用b-spline方法实现航线点集的y坐标和z坐标分维优化，生成优化后的航迹点集pb：
[0055]
xxgh＝0：.25：60
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(29)
[0056]
yygh＝spline(xgh，ygh，xxgh)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(30)
[0057]
zzgh＝spline(xgh，zgh，xxgh)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(31)
[0058]
pb＝[xxgh，yygh，zzgh]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(32)
[0059]
公式(26)～(27)中，xxgh、yygh、zzgh分别为航线点优化后的x坐标、y坐标、z坐标spline样条曲线优化函数，pbx为航线优化后航迹点序列。
[0060]
本发明与现有技术相比，其优点在于：
[0061]
a.实现复杂动态和静态环境下全局最优航线的快速规划。采用遴选适飞区域和改进的最速下降方法相结合，为避免无人直升机坠落的问题提供了一种有效的方法。
[0062]
b.利用分维航线实现二维样条曲线的三维航线优化。利用离散的x坐标通过二维spline函数实现y、z坐标光滑优化处理以达到较快速且具有可操作性的特点，在动态和静态环境下三维航线快速生成，实现动、静态威胁快速回避方面具有重要的作用。
具体实施方式
[0063]
下面结合实施例对本发明技术方案做进一步说明。实施例中有1架无人直升机飞行，规划空域左下角坐标为[0km，0km，0km]，右上角坐标为[60km，60km，60km]，在(11km，
16km)、(9km，9km)、(15km，1km)、(11km，5km)、(7km，12km)、(12km，12km)处共有6个地形威胁，在(17km，32km)、(30km，19km)、(40km，14km)处共有3个人群威胁；航程的长度大于50km，续航时间大于45minh，巡航速度大于60km/h。出发点和目标点的坐标分别为(1km，1km，5m)、(60km，60km，5m)；
[0064]
(1)构建综合威胁模型t
t
，数学描述如下：
[0065]
首先利用公式(1)得到地形威胁模型v1：
[0066][0067]
上述公式中，v1(x，y)表示横坐标x及纵坐标y处的地形威胁模型值，单位为公里；的是一种(-x
2-(y+1)2)次方指数函数；
[0068]
其次利用公式(2)～(7)得到x轴坐标和y轴坐标处人群威胁模型v2：
[0069][0070]
mu＝round(rand(3，2)-0.5)*19+30
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
[0071]
sigma＝round(rand(3，2)*40)/10+30
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)
[0072]
gm＝gmdistribution.fit(x，3，
′
options
′
，statset(
′
display
′
，
′
off
′
))
ꢀꢀꢀꢀꢀꢀꢀꢀ
(5)
[0073]
ii＝ezsurf(@(x，y)pdf((gm，[x，y])，[-15 15]，[-15 10])
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0074]
v2＝ii.cdata
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)
[0075]
公式(2)～(7)中，[]为矩阵符号，x为人群分布三维坐标矩阵，mvnrnd、mu、sigma分别为产生人群高斯分布数据的函数及其对应的均值、协方差；gmdistribution.fit、gm、statset分别为三维混合高斯分布参数拟合函数及其创建的结构体对象和创建统计选项结构函数，options、display、off分别表示选项、显示和关闭，ezsurf、h、ii.cdata分别表示创建一个的表面图及其句柄和三维表面数据矩阵，pdf为求概率密度的函数，@为构造函数表达式的符号；
[0076]
最后利用公式(8)得到融合静态地形威胁和动态人群威胁的综合威胁模型v：
[0077]
v＝v1+v2
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(8)
[0078]
公式(8)中，v为综合威胁模型；
[0079]
(2)计算威胁模型梯度值及其对应的复数场pf：
[0080]
首先利用公式(9)～(10)得到x轴坐标梯度分量px、y轴坐标梯度分量py和复数场pg：
[0081]
[px，py]＝gradient(v，0.2，0.2)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)
[0082]
pg＝px+py*i
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(10)
[0083]
公式(9)～(10)中，gradient为求梯度的函数，px、py、pg分别为梯度横、纵坐标分量及其对应的梯度复数场，i为序数表示符；
[0084]
其次利用公式(11)～(12)得到x轴坐标梯度分量px、y轴坐标梯度分量py和梯度矢量复数场pf：
[0085]
iq＝quiver(1：60，1：60，v)
ꢀꢀꢀꢀ
(11)
[0086]
pf＝iq.udata+iq.vdata*j
ꢀꢀꢀꢀ
(12)
[0087]
公式(9)～(10)中，quiver是matlab中绘制二维矢量场的函数，iq.udata、iq.vdata、pf分别为梯度矢量场横、纵坐标分量及其对应的梯度矢量复数场，j为序数表示符；
[0088]
(3)遴选适飞区域，构建飞行航迹点集fpc：
[0089]
利用公式(13)～(16)遴选适飞区域x、y轴坐标集iqc及梯度矢量场横、纵坐标分量比值集iqvm，并随机选出飞行航迹点集fpc：
[0090][0091][0092]
id＝floor(rand(1，90)*length(iqc)+1)(15)
[0093]
fpc＝iqc(id)
ꢀꢀꢀꢀ
(16)
[0094]
公式(13)～(16)中，iqc、iqvm分别表示遴选适飞区域平面坐标集及其对应位置的梯度矢量场横、纵坐标分量比值集iqvm，inf表示无穷大，abs、floor、length、rand分别为求绝对值、向下取整数、求最大数组维度的长度和产生数值范围在0～1之间随机矩阵的函数，id、fpc分别为随机选择的飞行航迹点在iqc中的序号和飞行航迹点集；
[0095]
(4)执行航迹点优选函数，生成初始最优飞行航线：
[0096]
首先利用公式(17)～(20)建立目标区域引导模型ze及其梯度x、y轴坐标轴分量dx、dy，构造飞行航线航向引导复数场pfz：
[0097]
[mx，my]＝meshgrid(0：60)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(17)
[0098]
ze＝-((mx-60)2+(my-60)2)/2
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(18)
[0099]
[dx，dy]＝gradient(ze，1，1)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(19)
[0100]
pfz＝dx+dx*i
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(20)
[0101]
公式(17)～(20)中，mx、my分别表示x、y轴坐标离散矩阵，meshgrid为生成二维网格坐标函数，dx、dy、pfz分别表示目标区域引导模型梯度x、y轴坐标轴分量和飞行航线航向引导复数场；
[0102]
利用公式(21)～(26)定义目标区域引导函数fun及其对x、y轴坐标得导数funx、funy，并初始化航迹点初值坐标x0、y0：
[0103]
fun＝inline(
′‑‑
ze
′′
x
′
，
′y′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(21)
[0104]
dfunx＝inline(
′
x
′
，
′
x
′
，
′y′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(22)
[0105]
dfuny＝inline(
′y′
，
′
x
′
，
′y′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(23)
[0106]
indx＝floor(rand(1)* length(iqc))+1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(24)
[0107]
x0＝iqc(indx，1)
ꢀꢀꢀꢀ
(25)
[0108]
，y0＝iqc(indx，2)
ꢀꢀꢀꢀ
(26)
[0109]
公式(21)～(26)中，fun、funx、funy分别表示目标区域引导函数和其对x、y的导数，inline为定义函数函数，indx、x0、y0分别随机生成的初始化航迹点在iqc中的序号及其横、纵坐标轴；
[0110]
利用公式(27)～(28)利用改进的梯度下降算法tidu1对适飞区域内的飞行航迹点集进行选取，组成优化的航迹点集合fp：
[0111]
[xtd，ytd，n，point]＝tidu1(fun，dfunx，dfuny，x0，y0，50，5，iqc)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(27)
[0112]
fp＝[[xb yb zb]；[point v(point)]；[xeye ze]]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(28)
[0113]
公式(27)～(28)中，xtd、ytd、n、point分别表示选出的最优航迹点的横、纵坐标及航迹点个数、航迹点平面坐标集合，tidu1为改进的根据目标引导场从适飞区域选择航迹点的梯度下降算法，fp为选出航线点组成的飞行航线，epsilon1、lambda1分别为精度和步长/更新率，xb、yb、zb分别为航线起始点的x坐标、y坐标、z坐标，xe、ye、ze分别为航线终止点的x坐标、y坐标、z坐标；
[0114]
(5)采用分维b-spline方法改进最优飞行航线，实现航线可飞行性优化：
[0115]
利用公式(29)～(32)分别应用b-spline方法实现航线点集的y坐标和z坐标分维优化，生成优化后的航迹点集pb：
[0116]
xxgh＝0：.25：60
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(29)
[0117]
yygh＝spline(xgh，ygh，xxgh)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(30)
[0118]
zzgh＝spline(xgh，zgh，xxgh)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(31)
[0119]
pb＝[xxgh，yygh，zzgh]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(32)
[0120]
公式(26)～(27)中，xxgh、yygh、zzgh分别为航线点优化后的x坐标、y坐标、z坐标spline为样条曲线优化函数，pbx为航线优化后航迹点序列。
[0121]
实验结果为：规划的航程为83.9852公里，满足航程大于50km的约束，规划时间为0.999秒，采用本方法规划的航线，能够满足在线规划的需求(5秒)，航线规划时间83.99min满足续航时间大于45min的需求，同时航线选择了安全的地势较低的区域进行飞行，能够实现在线飞行航线规划，为紧急事故救援提供了安全保障。

技术特征：
1.一种基于梯度复数场的多旋翼无人直升机航线规划方法，其特征采用以下实现步骤：(1)构建综合威胁模型t
t
，数学描述如下：首先利用公式(1)得到地形威胁模型v1：上述公式中，v1(x，y)表示横坐标x及纵坐标y处的地形威胁模型值，单位为公里；的是一种(-x
2-(y+1)2)次方指数函数；其次利用公式(2)～(7)得到x轴坐标和y轴坐标处人群威胁模型v2：mu＝round(rand(3，2)-0.5)*19+30
ꢀꢀꢀꢀꢀꢀꢀꢀ
(3)sigma＝round(rand(3，2)*40)/10+30
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)gm＝gmdistribution.fit(x，3，
′
options
′
，statset(
′
display
′
，
′
off
′
))
ꢀꢀꢀꢀꢀꢀ
(5)ii＝ezsurf(@(x，y)pdf(gm，[x，y])，[-15 15]，[-15 10])
ꢀꢀꢀꢀꢀꢀꢀ
(6)v2＝ii.cdata
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)公式(2)～(7)中，[]为矩阵符号，x为人群分布三维坐标矩阵，mvnrnd、mu、sigma分别为产生人群高斯分布数据的函数及其对应的均值、协方差；gmdistribution.fit、gm、statset分别为三维混合高斯分布参数拟合函数及其创建的结构体对象和创建统计选项结构函数，options、display、off分别表示选项、显示和关闭，ezsurf、ii、ii.cdata分别表示创建一个的表面图及其句柄和三维表面数据矩阵，pdf为求概率密度的函数，@为构造函数表达式的符号；最后利用公式(8)得到融合静态地形威胁和动态人群威胁的综合威胁模型v：v＝v1+v2
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(8)公式(8)中，v为综合威胁模型；(2)计算威胁模型梯度值及其对应的复数场pf：首先利用公式(9)～(10)得到x轴坐标梯度分量px、y轴坐标梯度分量py和复数场pg：[px，py]＝gradient(v，0.2，0.2)
ꢀꢀꢀꢀꢀꢀꢀ
(9)pg＝px+py*i
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(10)公式(9)～(10)中，gradient为求梯度的函数，px、py、pg分别为梯度横、纵坐标分量及其对应的梯度复数场，i为序数表示符；其次利用公式(11)～(12)得到x轴坐标梯度分量px、y轴坐标梯度分量py和梯度矢量复数场pf：iq＝quiver(1：60，1：60，v)
ꢀꢀꢀꢀꢀꢀꢀ
(11)pf＝iq.udata+iq.vdata*j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)公式(9)～(10)中，quiver是matlab中绘制二维矢量场的函数，iq.udata、iq.vdata、p分别为梯度矢量场横、纵坐标分量及其对应的梯度矢量复数场，j为序数表示符；(3)遴选适飞区域，构建飞行航迹点集fpc：
利用公式(13)～(16)遴选适飞区域x、y轴坐标集iqc及梯度矢量场横、纵坐标分量比值集iqvm，并随机选出飞行航迹点集fpc：集iqvm，并随机选出飞行航迹点集fpc：id＝floor(rand(1，90)*length(iqc)+1)
ꢀꢀꢀꢀꢀꢀꢀꢀ
(15)fpc＝iqc(id)公式(13)～(16)中，iqc、iqvm分别表示遴选适飞区域平面坐标集及其对应位置的梯度矢量场横、纵坐标分量比值集iqvm，inf表示无穷大，abs、floor、length、rand分别为求绝对值、向下取整数、求最大数组维度的长度和产生数值范围在0～1之间随机矩阵的函数，id、fpc分别为随机选择的飞行航迹点在iqc中的序号和飞行航迹点集；(4)执行航迹点优选函数，生成初始最优飞行航线：首先利用公式(17)～(20)建立目标区域引导模型ze及其梯度x、y轴坐标轴分量dx、dy，构造飞行航线航向引导复数场pfz：[mx，my]＝meshgrid(0：60)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(17)ze＝-((mx-60)2+(my-60)2)/2
ꢀꢀꢀꢀꢀꢀꢀꢀ
(18)[dx，dy]＝gradient(ze，1，1)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(19)pfz＝dx+dy*i
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(20)公式(17)～(20)中，mx、my分别表示x、y轴坐标离散矩阵，meshgrid为生成二维网格坐标函数，dx、dy、pfz分别表示目标区域引导模型梯度x、y轴坐标轴分量和飞行航线航向引导复数场；利用公式(21)～(26)定义目标区域引导函数fun及其对x、y轴坐标得导数funx、funy，并初始化航迹点初值坐标x0、y0：fun＝inline(
′‑
ze
′
，
′
x
′
，
′
y
′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(21)dfunx＝inline(
′
x
′
，
′
x
′
，
′
y
′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(22)dfuny＝inline(
′
y
′
，
′
x
′
，
′
y
′
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(23)indx＝floor(rand(1)*length(iqc))+1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(24)x0＝iqc(indx，1)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(25)y0＝iqc(indx，2)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(26)公式(21)～(26)中，fun、funx、funy分别表示目标区域引导函数和其对x、y的导数，inline为定义函数函数，indx、x0、y0分别随机生成的初始化航迹点在iqc中的序号及其横、纵坐标轴；利用公式(27)～(28)利用改进的梯度下降算法tidu1对适飞区域内的飞行航迹点集进行选取，组成优化的航迹点集合fp：[xtd，ytd，n，point]＝tidu1(fun，dfunx，dfuny，x0，y0，epsilon1，lambda1，iqc)
ꢀꢀꢀꢀꢀꢀ
(27)fp＝[xb yb zb]；[point v(point)]；[xe ye ze]]
ꢀꢀꢀꢀꢀꢀꢀ
(28)
公式(27)～(28)中，xtd、ytd、n、point分别表示选出的最优航迹点的横、纵坐标及航迹点个数、航迹点平面坐标集合，tidu1为改进的根据目标引导场从适飞区域选择航迹点的梯度下降算法，fp为选出航线点组成的飞行航线，epsilon1、lambda1分别为精度和步长/更新率，xb、yb、zb分别为航线起始点的x坐标、y坐标、z坐标，xe、ye、ze分别为航线终止点的x坐标、y坐标、z坐标；(5)采用分维b-spline方法改进最优飞行航线，实现航线可飞行性优化：利用公式(29)～(32)分别应用b-spline方法实现航线点集的y坐标和z坐标分维优化，生成优化后的航迹点集pb：xxgh＝0：.25：60
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(29)yygh＝spline(xgh，ygh，xxgh)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(30)zzgh＝spling(xgh，zgh，xxgh)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(31)pb＝[xxgh，yygh，zzgh]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(32)公式(26)～(27)中，xxgh、yygh、zzgh分别为航线点优化后的x坐标、y坐标、z坐标spline为样条曲线优化函数，pbx为航线优化后航迹点序列。

技术总结
本发明涉及一种基于梯度复数场的多旋翼无人直升机航线规划方法，其步骤为：构建综合威胁模型，计算威胁模型梯度值及其对应的复数场，遴选适飞区域，构建飞行航迹点集，执行航迹点优选函数，生成初始最优飞行航线，采用分维B-Spline方法改进最优飞行航线，实现航线可飞行性优化。本发明是考虑到综合形威胁空间布局特征、飞行安全性和飞机性能等的航线选优问题。其算例结果表明，利用复数场优选算法在复杂威胁环境下能实现最优航线的规划，该方法能够给无人直升机提供有效飞行服务同时提高了飞行的安全性和可靠性。飞行的安全性和可靠性。


技术研发人员：刘丽峰 牛宗伟 张强
受保护的技术使用者：山东理工大学
技术研发日：2023.06.20
技术公布日：2023/8/6