考虑不确定性的低空无人机运行风险评估方法*

王文涛，甘旭升*，吴亚荣，任谨慎，王翠香

（空军工程大学，a. 空管领航学院；b. 基础部，陕西 西安 710051）

0 引言

随着航空技术与信息技术日趋成熟，无人机系统技术发展迅猛，同时，凭借其低成本、高机动、操作便捷等特点，应用日益多样化［1-2］。而低空空域（特别是400～500 m 高度以下）作为“低、慢、小”无人机运行场所，能够满足无人机多样化应用（侦察监视、物流快递、救援救灾、农业耕种等）需求［3-6］。但低空无人机数目多、体积小、执行任务种类多样，现有的空中交通管理技术与手段难以适用于低空目标。因此，如何保证在低空复杂态势下，实现无人机安全和高效的运行，已成为全球工业界和学术界的共同挑战［7］。

目前，为解决低空无人机安全高效运行问题，各国都普遍在无人机运行服务系统（UAV service system，USS）系统和飞行情报管理系统（flight information management system，FIMS）支撑下，构建无人机空中交通管理（unmanned aircraft system traffic management，UTM）体系，并将 UTM 作为管控无人机的主要手段。美国联邦航空管理局（federal aviation administration，FAA）通过 FIMS，在 USS、ATM 和国家空域系统（national airspace system，NAS）之间的协调机制，并通过UTM 维持空中态势感知，维持空域和交通运行的监管与运营，同时给空域用户提供使用权限［8］。

风险评估主要是指利用科学的计算方法和模型，对潜在的或实际事故造成的影响后果进行定量估计，这是决定无人飞行系统能否安全融入国家空域系统的关键环节。为提升UTM 概念中的实时安全评估能力，推进无人机系统逐步融入NAS，近年来已经有一些学者对无人机低空运行风险进行评估。Giulio Avanzini 等［9］提出了一种适用于在人口稠密地区的无人航空系统风险评估程序。该方法根据无人机系统可靠性、人口密度以及撞击点的位置来评估发生故障后对地撞击风险。Stefano Primatesta 等［10］利用风险图对无人机运行风险进行描述，该方法可以对无人机运行环境的风险程度进行量化。其中风险图通过概率方法生成，包括人口密度、遮蔽因子、禁飞区和障碍物等层级，无人机每个运行位置都有其对应的量化风险值。P. Gonçalves等［11］提出一种基于Petri 网无人机安全评估模型，通过Petri 网显示已识别的无人机频率，以及无人机在出现故障后的反应能力，对无人机运行中操作授权过程进行优化。Achim Washington 等［12］基于贝叶斯方法，充分考虑系统安全评估过程中固有的不确定性，对风险评估与决策框架进行定义，并以无人机安全评估系统为研究案例，验证了所提方法的合理性。Xuejun Zhang 等［13］考虑并分析了无人机地面撞击操作和空中碰撞的2 个关键危害，通过计算满足不同无人机类别的安全目标水平所需的系统可靠性，来构建多因素影响下的风险模型，从而实现地面安全风险评估。REN X H 等［14］为评估低空运行的无人机的空域风险，定义了无人机在城市运行中的第三方风险——无人系统对地面无关人员造成的风险。并且明确了风险的来源和对象。考虑到无人机坠毁、噪声、机载摄像头和地面环境因素，构建了无人机城市物流风险指数评价模型。

从当前研究现状来看，低空无人机运行风险评估研究已有一些进展，但往往忽略复杂低空环境参数具有不确定性特征，以确定的数值进行仿真计算，往往与实际运行情况相差过大。基于以上考虑，本文提出一种考虑无人机坠落不确定性的风险评估方法。根据无人机下降过程中的受力模型，推导出无人机的坠落区域，并考虑无人机运行过程中的不确定性，计算出无人机坠落时的事故危害程度。

1 事故发生危害性建模

与有人机相比，无人机在性能参数、任务场景等方面都有显著差异，因此在评估无人机坠落地面撞击事件造成的危害程度时，需要结合UAV 自身特点，对事故发生危害性进行建模评估。如图1 所示，从UAV 系统来看，与地面撞击事件危害程度相关的因素包括无人机自身属性、坠落影响面积等。从任务场景角度看，无人机运行区域内地面人口分布、属性与遮蔽参数是地面撞击事故需要重点考虑的因素。

图1 UAV 地面撞击事故安全风险要素Fig.1 Safety risk elements of UAV ground impact accident

当前，主要用UAV 地面撞击事故死亡人数P 作为衡量事故等级的标准［15］，表示为

式中：P( f ) 为被撞击人员死亡率；Nexp为暴露在UAV 地面撞击事故中的地面人数，可由UAV 任务场景下的人口分布密度ρ 与撞击事故影响区域面积Aexp进行计算：

无人机系统由于运行故障，坠落形式通常为抛体运动，在下落过程中受多因素影响，其运动方程［16］为

式中：m，G 分别为无人机质量与重力；x，z 分别为下落过程中的水平与垂直方向的距离；t 为下落时间；Dx，Dz分别为无人机下落过程中所受的空气阻力。

式中：CD为空气阻力系数；ρa为空气密度；g 为重力加速度；Ax，Az分别为无人机侧向与垂向迎风面积。

式中：vx0，vz0分别为无人机水平和垂直方向上的初始速度；x0，z0分别为无人机水平和垂直方向上的初始位置。

无人机坠落时间为

式中：h0为无人机坠落时的运行高度。

将式（6）代入式（5）中，可得无人机在水平初速度v0和运行高度h0坠落条件下的坠落撞击点位置。实际无人机运行过程中的高度精度和水平精度分别记为 δh0，δxz。

无人机撞击致死率表示为

式中：ps为遮蔽参数，表示无人机运行环境下人群的暴露程度，取值范围为（0，1］，参数取值越大，同等撞击条件下致死率越低，各类型区域遮蔽参数设置如表1 所示。α 定义为ps= 0.5 时，致死率达到50%需要的撞击能量；β 定义为ps= 0 时，发生死亡所需的撞击能量临界值。Eimp指撞击事故中无人机动中无人机水平方向与垂直方向的碰撞速度vx，vy可由式（5），（6）确定，表示为

表1 遮蔽参数Table 1 Masking parameters

2 构建三维风险图

由于无人机低空运行区域内人口分布、地面遮蔽效果、飞行高度等参数均不相同，每个运行位置点存在的风险值也不一致，根据式（1）～（5），可对三维空间内运行风险进行量化，每一个地理坐标上的位置都对应一个风险值［17］。构建 n × m × l 的三维风险矩阵R，其中，第k(k = 1，2，…，l)层高度所代表的二维风险矩阵为

式 中 ：r(pi，j，k)(i = 1，2，…，n；j = 1，2，…，m；k =1，2，…，l)表示三维低空空域内位置处于 pi，j，k时的风险值。根据上文分析，可用UAV 地面撞击事故死亡人数P 来代替。

3 事故发生风险评估

无人机在低空运行空域所受多类型不确定因素影响，例如，各类导航、制导与控制（guidance，navigation，and control，GNS）误差，地面撞击事故死亡人数P 同样也难以用确定的数值进行描述，特别是无人机坠落实际撞击点将随机分布在标准撞击点周围，其概率密度函数（probability density function，PDF）依赖于撞击误差的统计分布。因此地面撞击区域的统计影响区域（statistical impact footprint，SIF）可通过 Monte-Carlo 方法来确定，图 2 展示了无人机以水平初速度v0和高度h0飞行时，抛体坠落运动地面撞击点二维示意图。

图2 无人机抛体坠落运动地面撞击点二维示意图Fig.2 Two-dimensional schematic diagram of the ground impact point of a UAV projectile falling motion

三维空域内无人机坠落撞击点分布也有类似的规律，通过Monte-Carlo 模拟可以确定双变量概率密度函数：

式中：Δx = x - μx，Δy = y - μy，其中 μx，μy分别为轨迹方向与垂直轨迹方向的位置平均值；ρxy为x 与y之间的相关系数；σx＞ 0，σy＞ 0 分别为轨迹方向与垂直轨迹方向的位置方差。

无人机坠落撞击点Monte-Carlo 模拟分布如图3所示，坠地撞击点三维分布图如图4 所示。其中红色曲线表示坠落撞击点分别在x 轴和y 轴的分布，是立体三维分布在二维面上的投影。

图3 无人机坠落撞击点Monte-Carlo 模拟分布Fig.3 Monte-Carlo simulation distribution of the impact point of a UAV fall

图4 坠地撞击点三维分布图Fig.4 Three-dimensional distribution map of the impact points of a fall

不同概率密度函数下的撞击影响点分布与影响区域面积如图5 所示，按照撞击区域3σ 准则，为充分考虑无人机对地面不同区域的撞击概率不同，准确衡量无人机在空域内某个运行位置的风险值，将撞击点分布图区分为 Sσ、S2σ与 S3σ3 种类别，若无人机坠落影响区域经过栅格化后，分为n 个基本计算单位，则风险值表示为

图5 不同撞击概率下影响面积示意图Fig.5 Schematic diagram of affected area under different impact probabilities

式中：a，b，c 分别为 Sσ、S2σ- Sσ与 S3σ- S2σ3 种区域面积内覆盖的栅格个数，满足a + b + c = n。

4 案例分析

低空空域环境内地面障碍物类型复杂多样、人口密度高，且处于动态变化之中。选取108°54′44″N～108°55'13″N，34°14'28″E～34°14'41″E 作为无人机运行区域，如图6 所示。该地区包含建筑、空地、绿化等多类型地形，同时人口密集，是典型的人口聚集区，适合用来衡量无人机运行风险。考虑到实际运行评估过程中的可操作性，确保无人机在高于地面障碍物的高度上运行，其不干扰民航飞机的运行，选择该区域上空30 m 空域作为无人机运行高度进行风险评估。

图6 无人机运行区域Fig.6 UAV operating area

选取大疆经纬M210-RTK 型无人机作为风险评估对象，该型号无人机的参数如表2 所示。

表2 大疆经纬M210-RTK 型无人机参数Table 2 Parameters of DJI Jingwei M210-RTK UAV

人口密度是衡量地面人员分布的重要指标，是无人机运行风险评估中的一个关键变量。由于区域内每日人口空间分布具有较强的规律性与动态性，因此本研究选取 2020-04-20T12：00：00 固定时刻的城市热力图，作为衡量人口分布的依据，结果如图7 所示，将运行区域栅格化，形成80 × 41 个网格，网格中不同颜色代表不同的人口分布密度。

图7 人口密度分布图Fig.7 Population density distribution map

遮蔽参数是衡量地面不同设施对人员的保护效果的重要指标，如图8 所示，运行区域栅格化，形成80 × 41 个网格，根据表1 中定义，网格中不同颜色代表不同的遮蔽参数数值。

图8 遮蔽参数分布图Fig.8 Distribution of shadowing parameters

如图9 所示，A 点为无人机任务起始点，经过A-B-C-D 路线最终返回A 点，运行高度设置在所有地面设施高度之上，为简化评估流程，假设无人机在任务全程保持运行速度为v0= 20 m/s，运行高度h0= 30 m 不变。风险评估所需的参数如表3 所示。

图9 无人机任务场景Fig.9 UAV mission scenario

表3 所需参数取值Table 3 Values of required parameters

为衡量无人机运行过程中的危害性，需要对离散化的运行区域进行逐一分析，对于所选取的区域，规模为：740 m×401 m，则每个栅格网络的规模大约为：10 m×10 m。现在对初始点状态的无人机进行风险评估，其余路径点类似。

初始点无人机速度为v0= 20 m/s，高度h0= 30 m，以坠地撞击点纵向分布为例，根据统计信息，其均值服从均值为42.51 m，标准差为3.17 m 的正态分布，坠地撞击地纵向分布特征如图10 所示。

图10 坠地撞击点纵向分布特征Fig.10 Characteristics of longitudinal distribution of impact points on the ground

按照事故发生危害性评估流程，可计算出初始点的撞击点分布图如图11 所示，计算得出：Sσ= 19.91 m2，S2σ= 79.63 m2，S3σ= 179.17 m2，根据式（11）定义，可计算出风险值为14.81。同理，整个无人机执行任务路径点均可按照此方式实时计算风险值。事故发生危害性等级标准如表4 所示。对照标准，无人机在初始点的危害等级为轻微危害。

表4 事故发生危害性等级标准Table 4 Accident hazardous level standard

图11 初始点撞击点分布图Fig.11 Distribution of initial point impact points

无人机在运行过程中，通过计算事故危害性数值，结合事故发生危害性等级标准，实时判断无人机所处的安全水平，整个运行的实时评估结果如图12 所示。为更精细地显示无人机整个运行过程中的风险等级，绘制相应的实时风险等级安全谱，如图13 所示。从安全谱中可以直观地看出，无人机运行风险在中度危害及以下，只有极个别时刻处于中度危害等级，无严重危害等级，整体处于比较安全的水平。

图12 无人机实时风险评估结果Fig.12 UAV real-time risk assessment results

图13 无人机实时运行安全谱Fig.13 UAV real-time operation safety spectrum

由此可见，基于提出的风险评估框架能够实时评估无人机所处的安全水平，为无人机在复杂低空运行提供理论以及与标准。

为验证本文所提风险评估的优越性，将近年来应用在无人机风险评估领域的安全性分析方法的优缺点进行对比，具体信息如表5 所示。本文所提的风险评估方法能够无人机坠落的不确定性因素，符合无人机实际飞行情况，更加准确地衡量无人机风险水平。

表5 风险评估方法对比Table 5 Comparison of risk assessment methods

5 结论

作为UTM 中安全运行能力的体现，风险评估是无人飞行系统能否安全融入国家空域系统的关键环节。本文提出考虑无人机坠落过程不确定性的风险评估方法，可将该功能嵌入到UTM 服务中，能够实时评估无人机的安全水平，用于风险缓解以及突发事件管理。

（1）提出基于无人机坠落动力学模型与Monte-Carlo 方法的事故危害性建模，能够充分考虑坠落过程中的不确定性，结合事故危害性等级标准，确定无人机运行过程中的风险值。

（2）以108°54'44″N～108°55'13″N，34°14′28″E～34°14'41″E 区域为例，选取大疆经纬 M210-RTK 型无人机作为评估对象，验证了风险评估方法的可行性。