有效
一种非合作无人机入侵合作无人机的冲突风险评估方法
羊钊、王艳、李娜、朱仁伟、杨晓方、胡锦标、张颖、王兵
南京航空航天大学
摘要
本发明公开了一种非合作无人机入侵合作无人机的冲突风险评估方法,包括:采集多种类型的无人机真实轨迹数据,生成轨迹时间序列切片;根据非合作无人机与合作无人机在未来某时刻的预测信息,以及设定的保护区范围进行冲突识别;对有冲突的点对组合进行模拟仿真,生成蒙特卡洛仿真样本;评估非合作无人机入侵合作无人机的危险行为发生的可能性和严重程度,确定风险评估指标并进行评级;对非合作无人机入侵合作无人机的场景进行冲突风险评估及等级计算。本发明方法可对非合作无人机与单架无人机的冲突场景进行风险评估并量化冲突风险等级,从而提供更精准的报警信息。
1.一种非合作无人机入侵合作无人机的冲突风险评估方法,其特征在于,步骤如下:1)采集多种类型的无人机真实轨迹数据,并进行轨迹数据处理与分割,生成轨迹时间序列切片;2)根据非合作无人机与合作无人机在未来某时刻的预测信息,以及设定的保护区范围进行冲突识别;3)基于步骤2)中的冲突识别结果,对有冲突的点对组合进行模拟仿真,生成蒙特卡洛仿真样本;4)基于步骤3)中生成的仿真样本,评估非合作无人机入侵合作无人机的危险行为发生的可能性和严重程度,确定风险评估指标并进行评级;5)基于步骤4)中确定的风险评估指标,对非合作无人机入侵合作无人机的场景进行冲突风险评估及等级计算;所述步骤2)的具体过程如下:21)判断未来某时刻非合作无人机的预测位置与合作无人机的冲突保护区的相对位置关系;211)定义合作无人机的冲突保护区和碰撞保护区;以未来某时刻合作无人机的预测位置作为保护区水平及垂直方向的中点,考虑冲突及碰撞两个过程,设置冲突保护区和碰撞保护区形状为圆柱体,碰撞保护区外扩形成冲突保护区;碰撞保护区的尺寸计算公式如下:式中,h为碰撞保护区圆柱体的高度,r为碰撞保护区圆柱体的半径,h u 为无人机外形尺寸的高度,r u 为无人机水平尺寸最大值的二分之一,G y 与G h 分别为无人机在GPS模式下垂直和水平方向的悬停精度;冲突保护区的圆柱体高度D ver 和半径D hor 在碰撞保护区的圆柱体高度和半径的基础上按一定比例扩大;212)判断无人机的冲突间隔;判断作为质点的无人机B在第t时刻的三维位置坐标预测值在水平、垂直方向上是否均位于合作无人机A的冲突保护区内,若是,则下面公式(3)、(4)同时成立,无人机A与无人机B发生冲突,具体判断公式如下:|P(z) i-coo -P(z) i-non |<D ver (4)式中,P(x) i-coo 、P(y) i-coo 、P(z) i-coo 分别为合作无人机在水平方向的x坐标、水平方向的y坐标和垂直方向的z坐标,P(x) i-non 、P(y) i-non 、P(z) i-non 分别为非合作无人机在水平方向上的x坐标、水平方向的y坐标和垂直方向上的z坐标;若公式(3)、(4)均不成立,则无人机A与无人机B不发生冲突;若公式(3)、(4)中仅一个成立,则进入步骤22);22)基于步骤21)的判断结果,若公式(3)、(4)非同时满足,则进一步结合非合作无人机在未来时刻的预测飞行意图进行冲突识别;其中,将未来某时刻至下一时刻时间段内的飞行意图表征为三个方面,包括水平飞行行为、纵向飞行行为和加速度情况;水平飞行行为具体为两相邻时刻之间的航向角区间,以90°为间隔划分为四个区间;纵向飞行行为具体为两相邻时刻之间的俯仰角区间,包括下降,爬升及平飞三个区间;加速度情况由两相邻时刻之间加速、匀速和减速三种情况组成;若第t时刻下非合作无人机B的三维位置坐标预测值仅在水平或垂直方向位于合作无人机A的冲突保护区内,且另一个方向两者间的距离小于无人机单个步长下以最大速度运行的飞行距离,则判断非合作无人机B在对应方向的飞行行为区间与合作无人机A位于相同区间的概率是否超过预定阈值,若超过,则表示第t时刻非合作无人机B与合作无人机A预测轨迹具有冲突风险;若未超过,则合作无人机A与非合作无人机B不发生冲突;23)基于步骤21)、步骤22)中的结果得到第t时刻非合作无人机与合作无人机冲突行为的预测标签;同时,非合作无人机在第t时刻的真实位置与合作无人机真实轨迹坐标的水平及垂直方向的相对距离小于/等于设置的冲突保护区距离时,得到冲突事件的真实标签为是,否则为否;24)进行冲突识别结果有效性评价;241)冲突识别精度计算;冲突识别精度为预测正确的样本数量占总样本数的比例Acc,具体计算公式如下:式中,y l 、y l '分别为样本l的真实标签和预测标签,n为预测样本数量;242)冲突识别查准率和查全率计算;在非合作无人机入侵场景的冲突识别中,查准率和查全率具体计算公式如下:式中,R为识别查全率,P为识别查准率;243)冲突识别综合评价指标计算;通过平衡因子α表征识别任务对查准率P和查全率R的侧重程度,α>1时,侧重于查全率R;α<1时,侧重于查准率P;将P和R综合为一个评价指标F α ,具体计算公式如下:设置综合评价指标的阈值为F α * ,当F α >F α * 时,进入步骤3);否则,返回步骤1)。
2.根据权利要求1所述的非合作无人机入侵合作无人机的冲突风险评估方法,其特征在于,所述步骤1)的具体过程如下:采集多种类型的无人机真实轨迹数据,轨迹数据中的各轨迹点包括时间戳、身份ID、纬度、经度、飞行相对高度、速度及姿态角信息,将WGS-84坐标系下的经度、纬度、高度转化为笛卡尔坐标系中的(x,y,z),对无人机真实轨迹数据进行等间隔采样,采用滑动时间窗口方法对无人机的轨迹序列数据进行切片划分,并进行归一化处理,生成各类型无人机的轨迹时间序列切片。
3.根据权利要求1所述的非合作无人机入侵合作无人机的冲突风险评估方法,其特征在于,所述步骤3)的具体过程如下:31)建立无人机运动学方程;建立无人机运动学方程,以保证每个相邻时刻内生成的随机运动轨迹的物理可行性与真实性,具体公式如下:式中,a为加速度矢量,V为速度矢量,T为推力矢量,K为经验阻力系数,m为无人机的质量,g为加速度;其中,K与无人机姿态的变化有关,在不同的倾斜角度 下,K与轴向阻力系数K a 、径向阻力系数K s 的关系为:设定无人机阻力项轴对称,结合无人机在普通档和运动挡下的相关性能参数进行求解,以确定初始推力和阻力,具体公式如下:式中, 和 分别为无人机在运动档和普通档下的最大倾斜角度,V 1 和V 2 分别为无人机在运动档和普通档下的最大水平速度;32)生成指定速度矢量时所需的推力矢量和姿态倾角;设定在每个时间步长初始时发生姿态和速度变化,无人机的推力矢量与姿态矢量在同一个方向上,且阻力矢量的方向与速度矢量的方向相反,单位姿态矢量 与单位速度矢量 具体计算公式如下:式中,ψ与θ分别为航向角与俯仰角, 与 分别为绕x轴和绕y轴方向的旋转角度;随机生成的无人机速度和姿态满足运动轨迹的物理可行性需求,具体公式如下:式中,D为阻力矢量,依据限制条件 在加速度矢量a=0的情况下,结合式(9)和式(14),通过迭代改变 与 直至推力矢量T与单位姿态矢量 垂直时,生成给定速度矢量时所需的推力矢量和姿态倾角,使随机生成的速度和姿态满足运动轨迹的物理可行性需求;33)生成蒙特卡洛仿真样本;每个蒙特卡洛仿真样本的飞行轨迹在每个时间步长随机生成速度、航向角和仰角,速度、航向角和仰角的初始生成遵循正态分布,分别生成不同时刻下无人机的蒙特卡洛仿真样本,具体生成方式如下:式中,|V|表示速度大小;|V target |表示平均速度的大小;|V max |表示最大速度的大小;ψ与θ分别为航向角与俯仰角;V target 、ψ target 、θ target 分别为速度、航向角和俯仰角的均值,σ V 、σ ψ 、σ θ 分别为速度、航向角和俯仰角的方差。
4.根据权利要求1所述的非合作无人机入侵合作无人机的冲突风险评估方法,其特征在于,所述步骤4)的具体过程如下:41)评估非合作无人机入侵合作无人机危险行为发生可能性;411)计算无人机预计入侵概率;预计入侵概率表示预测时刻非合作无人机入侵合作无人机冲突保护区和碰撞保护区的行为发生概率,通过计算位于合作无人机冲突保护区或碰撞保护区内蒙特卡洛仿真样本的百分比得到,具体计算公式如下:式中,μ 1 为预计入侵概率,f为蒙特卡洛仿真样本个数,s k ∈I表示所生成的第k个蒙特卡洛仿真样本位置位于合作无人机保护区I c 内,其中,c=1时为冲突保护区,c=2时为碰撞保护区;412)计算无人机预计冲突时间;设定预测步长为i时,合作无人机的预测位置为(x i ,y i ,z i ),非合作无人机的蒙特卡洛仿真样本分布位置为 其中, j为生成的仿真样本个数;通过计算合作无人机至非合作无人机各样本点的距离,得到样本距离集合 为预测步长i时合作无人机预测点与非合作无人机第p个仿真样本之间的距离,具体计算公式如下:合作无人机与非合作无人机在预测步长i时的最小距离 具体计算公式如下:设定在合作无人机相对静止,不执行任何机动的情况下,以非合作无人机按最大飞行速度径直朝向合作无人机方向的最坏意图假设进行预计冲突时间的计算,具体计算公式如下:式中,μ 2 为预计冲突时间, 为合作无人机与非合作无人机在预测步长i时的最小距离, 为非合作无人机的最大飞行速度;413)计算无人机预计冲突意图;设定非合作无人机在预测步长i→i+1时段内的水平意图区间为I h ,垂直方向意图区间I v ,速度变化情况I a 在加速时I a =1,匀速时I a =0,减速时I a =-1;合作无人机在步长为i时的预测位置为(x i ,y i ,z i );通过结合意图区间的相对位置及速度变化情况判断冲突意图等级;当(x i ,y i )∈I h ,z i ∈I v 且I a =1,非合作无人机在下一时刻将进入合作无人机在水平及垂直方向的位置区间,且加速运行时,表示非合作无人机具有入侵意图,发生事故的风险等级最高;当 且 I a =-1,非合作无人机在下一时刻水平及垂直方向的位置区间均与合作无人机隔离,且减速运行时,表示非合作无人机进行无人机间的安全间隔保持,发生事故的风险等级最低;414)计算无人机距离规避指数;利用冲突主体之间的位置和运动方向来进行计算无人机距离规避指数,具体计算公式如下:式中,μ 4 为距离规避指数,X为预测时刻合作无人机与非合作无人机之间的距离,D mar 为合作无人机安全执行规避机动所需的最小距离,表示最小避让范围;415)计算无人机方位规避指数;利用速度障碍法的几何定义判断合作无人机和非合作无人机是否处于潜在的碰撞路线上,以预测时刻非合作无人机的位置为圆心,冲突保护区水平安全间隔为半径作圆;以预测时刻合作无人机的位置作圆的两条切线段,形成的圆锥区域为碰撞锥;若合作无人机的速度V h 与非合作无人机的速度V r 的矢量和V s 位于碰撞锥内,则存在冲突的可能,计算方位规避指数μ 5 ,具体如下:式中,ω为合作无人机探测角度的二分之一,V s 表示合作无人机与非合作无人机的速度矢量和,α 1 、α 2 分别为V s 与主机探测角度的两个边界所形成的夹角;416)确定非合作无人机入侵冲突区域概率、非合作无人机入侵碰撞区域概率、非合作无人机与合作无人机预计冲突时间、非合作无人机预计冲突意图、合作无人机距离规避指数、合作无人机水平方位规避指数、合作无人机垂直方位规避指数作为危险行为发生的可能性的评估指标;42)评估非合作无人机入侵合作无人机危险行为发生严重程度;421)获取无人机坠毁区域;无人机坠毁区域分为垂直碰撞区域和水平碰撞区域,垂直碰撞区域具体计算公式如下:Re v =π(r p +r uav ) 2 (22)式中,Re v 为垂直碰撞区域,r p 为人体平均半径,r uav 为无人机最大尺寸半径;水平碰撞区域的计算公式如下:Re c =2(r p +r uav )·d+π(r p +r uav ) 2 (23)式中,Re c 为水平碰撞区域,d表示坠落无人机在下降到人体的高度时行进的水平距离,d=H p /tanγ;H p 为人体平均高度,γ为速度矢量与水平面或人体表面撞击形成的角度;422)估计无人机地面碰撞动能;计算无人机碰撞点速度V imp ,通过最大飞行速度V x 和相对于最大飞行高度的自由落体速度V y 得到,具体如下:通过V x 与V y ,求解无人机碰撞角度γ,具体如下:根据无人机在撞击时刻的速度大小V imp 和最大起飞质量M,可求解坠毁无人机在撞击点的动能E c ,具体如下:423)计算无人机地面碰撞死亡概率;考虑遮蔽因子下的死亡概率P f 具体计算公式如下:式中,P s 表示遮蔽因子,E h 表示在遮蔽因子为6的情况下,死亡概率为50%时对应的冲击能量,E d 表示导致死亡的最低冲击能量阈值,E c 为坠毁无人机在撞击点的动能,μ表示校正因子,用于对小于E d 阈值或位于其附近的低动能值进行改进估计,具体如下:424)计算每飞行小时地面伤亡人数;基于公式(22)-公式(28),得到无人机每飞行小时地面伤亡人数计算公式如下:N p =P b ·Re c ·D p ·P f (29)式中,N P 为每飞行小时地面伤亡人数,P b 为无人机故障概率,由于考虑入侵场景下的冲突或碰撞事件可能导致无人机发生事故及故障的行为,因此,设定P b =μ 1 ,Re c 为水平碰撞区域面积,D p 为无人机运行区域的人口密度,P f 为死亡概率;425)确定考虑遮蔽因子下的死亡概率和每飞行小时地面伤亡人数作为危险行为发生的严重程度的评估指标。
5.根据权利要求1所述的非合作无人机入侵合作无人机的冲突风险评估方法,其特征在于,所述步骤5)的具体过程如下:51)设置无人机性能参数和环境参数,包括无人机最大起飞质量、最大运行速度、最大尺寸半径、探测距离、探测角度、商业区人口密度、商业区遮蔽系数、站立人体平均高度、人体平均半径、重力加速度;52)分别对冲突及碰撞预测仿真示例情况进行预测步长5s、10s、15s、20s下的风险计算及预警等级确定,包括入侵冲突及碰撞保护区域的概率、预计冲突时间、预计冲突意图、距离规避指数、水平及垂直方向的方位规避指数、死亡概率、每飞行小时地面伤亡人数。






