1、技术|Technology68 风能 Wind Energy 面向无人机自主巡检的风电机组姿态计算和航迹规划方法*文|魏昂昂，姚中原，张宇，王恩民，李志轩，司泽随着无人机和图像处理技术的发展,基于旋翼无人机的风电机组无人化自主巡检技术获得了较大的发展。在进行风电机组的无人机自主巡检过程中,需要保证无人机在按照规划的航迹飞行时不与风电机组发生碰撞。因此,风电机组姿态的准确获取,以及航迹规划的有效优化,是保证风电机组无人自主巡检高效、准确、安全的关键因素。目前,无人机巡检较普遍地应用于电力行业中,巡检范围包括风力/光伏发电领域、变电站、输电线路。在当前的实际应用中,无人机巡检仍主要依赖于人工控制。
2、为实现巡检的自动化和无人化水平,越来越多的学者开始介入无人机自动航迹规划技术研究。例如,徐晨1基于激光雷达建立电力线路的巡检航迹。徐维2基于机载激光雷达对电力走廊进行了场景重构与三维测量。赵猛等3利用三维激光雷达进行了无人机场景的建模。孙英君等4利用机载激光雷达对分布式光伏板进行了建模分析。此外,还有研究人员通过双目摄像头或深度相机与激光雷达的配合使用进行了无人机巡检中的航迹规划。然而,由于机载激光雷达实现商业化与工业化的时间并不长,以此为基础的各项计算机技术尚未完全成熟,且存在成本高、种类杂、精度不一等问题。此外,广泛使用的多旋翼无人机电池续航能力一直较弱,在必备的相机载荷外再增加雷达载荷,
3、会影响无人机的巡检时长和巡检效率。针对目前技术方案存在的载荷较大、航时不足、成本较高等问题,本文采用单个相机载荷采集风电机组及周边环境的图像信息,并生成风电机组点云数据,通过无人机上的机载相机快速获取物体表面大量采样点三维空间坐标,重构对象的三维模型；再通过模型分割、点云投影、二值化等算法计算获取风电机组的偏航角、叶轮转角、位置等参数；最后结合风电机组原有经纬度坐标信息等,确定风电机组的巡检航行轨迹。面向无人机自主巡检的风电机组姿态计算和航迹规划本文所述的面向无人机自主巡检的风电机组姿态计算和航迹规划方法主要包含以下步骤（图 1）：首先通过无人机上的机载相机快速获取物体表面大量采样点三维空间坐
4、标,重构对象的三维模型；再通过模型分割等算法确定风电机组的偏航角及转角等姿态信息；最后结合风电机组原有经纬度坐标信息等确定风电机组的巡检航行轨迹。一、三维重建三 维 重 建 的 过 程 主 要 可 以 分 为 两 个 步 骤：SFM（Structure From Motion）、MVS（Multi View Stereo）。其中,SFM 是使用一系列图片重建出拍照场景三维结构的*中国华能集团科技项目（HNKJ20-H44）1:徐展.一种激光雷达导航的全自主智能无人机巡线系统 J.浙江电力，2017，36（6）：4447.2:徐维.基于机载激光雷达的电力线信息重建技术研究 J.激光杂志，2020
5、，41（12）：114117.3:赵猛，姚凯，杨国栋，等.基于三维激光雷达的无人机场景建模 J.现代电子技术，2020，43（14）：4247.4:孙英君，崔信国，亓伟，等.激光雷达数据的建筑物光伏潜力评估 J.测绘科学，2019，44（7）：147151.图1 无人机巡检的风电机组姿态计算和航迹规划算法三维重建点云模型分割偏航计算算法转角计算巡检路径规划Technology|技术2023年第01期 69 （a）三维重建图片拍摄路径图3 以风电机组为中心的部分点云数据图2 基于点云数据的风电机组三维建模（b）基于SFM生成的拍摄场景稠密点云方法,该系列图片通过无人机的机载相机拍摄获取。图片摄取
6、方式如图 2（a）所示,在风电机组正上方按照红色轨迹绕行一周,该圆形轨迹直径和高度参考塔筒高度和叶片长度确定。其输入是相同物体不同视角的有序图片,相邻图片之间有一定的重叠区域,输出是物体的三维结构和各个拍照点相机的内外参数。典型的 SFM 系统工作流程可以分为三个阶段：特征检测和提取、特征匹配和几何验证、物体结构和相机运动重建。通过这一流程可以得到各图像中互相匹配的特征点在三维空间中的位置,生成稀疏点云。MVS 是在 SFM 输出结果的基础上计算每张照片中每个像素的深度和法线信息,融合深度图和法线图生成拍摄场景的稠密点云,如图 2（b）所示。二、点云模型分割点云分割的目的是得到风电机组各部分模
7、型,以便后续进一步进行处理。对于该任务,为了增加计算效率,减少无用计算,需要对点云进行关键区域提取,保留以风电机组为中心的部分点云数据。根据施工信息,可以获得风电机组基座经纬度坐标和海拔高度,结合风电机组塔筒高度可以粗略得到风电机组塔筒顶端中心的空间位置（cx,cy,cz）,结合叶片长度 l 即可确定风轮区域。由于塔筒顶端中心位置不是精确的实际位置,需要预留一定的范围 delta,则风轮区域是以（cx,cy,cz）为中心,l+delta 为边长的矩形区域,提取重建点云中对应区域即得到风电机组区域的点云,如图 3 所示。在实现点云模型分割前,需要获得塔筒中心精确的实际位置。本文对塔筒中心位置的计
8、算采用常用于立体模型虚拟中心位置计算的投影法,具体计算方式如下：首先,对风电机组点云自下而上间隔 1m 进行切片,将切片区域点投影到平面上,并与之前的结果进行累加计算,所有段计算完成后的截面累加效果如图 4（a）所示。由图可知,由于塔筒区域的垂直投影重叠度较高,累加后对应的塔筒位置为一个明显的圆环。接着,对累加图进行阈值化处理,阈值 th 为累加图所有点中最大累积数量的 0.2倍,去除累加图中所有小于 th 的点,得到如图 4（b）所示的二值图 G。通过施工信息获取塔筒最大直径 max_pillar 和最小直径 min_pillar,生成如图 5 所示的环形卷积核 F。用该卷积核对二值图进行卷
9、积运算,公式如下：R(i,j)=mnF(m,n)G(i m,j n)（1）技术|Technology70 风能 Wind Energy 计算结果 R 如图 6 所示。将卷积结果 R 中响应最强的点转换为通用坐标,即可得到风电机组塔筒的精确中心点位置（x,y）。根据塔筒精确中心点位置、塔筒直径,即可分离出塔筒,如图 7（a）所示。根据塔筒最高处坐标和机舱长度,在剩余点云中去除以塔筒最高点坐标为中心、机舱长度为半径的球形区域,得到如图 7（b）所示的风电机组叶片点云。三、偏航计算算法风电机组偏航角由叶片所在平面法向量确定,叶片所在平面由空间点拟合得到。建模点云存在较多噪声,且密度、分布位置不均匀。
10、为保证叶片平面的精度,叶片平面拟合分为两步：粗平面拟合和平面精调,二者均以最小二乘法平面拟合为基础,使用不同区域的点对平面进行拟合。空间平面拟合：记待拟合空间点集合为(x1,y1,z1),(x2,y2,z2),(xi,yi,zi),(xn,yn,zn),记矩阵 A1、B 如下：（2）（3）（a）累加图图5 环形卷积核F图6 二值图G的卷积计算结果图4 风电机组点云间隔1m水平切片的垂直投影（b）二值图GTechnology|技术2023年第01期 71 （a）塔筒点云图8 叶片平面拟合图7 风电机组点云分离（b）叶片点云通过以下公式进行空间直线拟合得到计算结果 R：（4）平面参数为：（5）空间
11、平面方程为：（6）平面法向量为：（7）首先,使用全部叶片点进行上述平面拟合,得到叶片点云沿平面法向量投影到平面的三视图（图 8）,该平面即为粗拟合平面 P1。由于所有叶片点并非都在同一平面上,且点沿平面法向量方向分布较散,最小二乘法拟合平面时会尽可能兼顾所有空间点,因此,拟合平面并不能很好地技术|Technology72 风能 Wind Energy 贴合风电机组叶片的朝向,需要进一步调整以获得更加精确的拟合平面。为了得到更精确、更贴合实际叶片位置的平面,使用叶片迎风面最外侧点进行二次拟合。为了确定这些点的位置,将叶片点云投影至上述粗拟合平面,在投影的平面图形中使用形态学方法寻找叶片投影的骨架
12、点,由骨架点反向投影至点云即可得到迎风面最外侧点。空间点到平面投影过程如下所示：记待投影空间点集合为(x1,y1,z1),(x2,y2,z2),(xi,yi,zi),(xn,yn,zn)。投影平面法向量为 V=(A,B,C),本例中 C=-1,平面上一点为(x0,y0,z0),平面参数 D 为：D=(A x0+B y0+C z0)（8）空间点(xi,yi,zi)在该平面上的投影点(Xi,Yi,Zi)：t=(A xi+B yi+C zi+D)/(A2+B2+C2)（9）Xi=xi AtYi=yi BtZi=ziC t （10）（a）叶片点云在粗拟合平面P1上的投影（a）拟合平面投影点图9 叶片平
13、面确定方法中的平面精调过程图10 叶片点云的粗影投影点和二值化投影点（b）使用二值图骨架提取方法提取的叶片骨架（b）二值化投影点Technology|技术2023年第01期 73 （a）提取前图12 叶片骨架提取图11 精调平面P2（b）提取后平面精调过程：图 9（a）为叶片点云在粗拟合平面P1 上的投影,图 9（b）是使用二值图骨架提取方法提取的叶片骨架。将叶片投影骨架点从粗拟合平面 P1 反投影,找到投影骨架点对应叶片点,如图 10 所示。对这些点进行平面拟合得到精调平面 P2,其三视图如图 11 所示。相比于 P1 平面,P2 平面与叶片拟合程度更好,记 P2法向量为 V2,则此时的风电
14、机组偏航角如下：yaw=arccos(V21V202+V212+V222)（11）俯仰角如下：pitch=arccos(V202+V212V202+V212+V222)（12）四、转角计算将叶片点沿 P2 法向量投影至 P2 平面,并用二值图骨架提取方法提取骨架,结果如图 12 所示。对二值骨架使用霍夫变换进行直线检测,检测到的直线和其端点如图 13 所示。确定转角的前提是找到轮毂中心点,即旋转中心。当前条件下,检测到的多条直线不一定会相交于同一个点,待求解的是一个超定方程,因此,使用最小二乘法来求解中心点,使得解出的中心点到各条直线的距离尽可能近,技术|Technology74 风能 Win
15、d Energy 从而保证中心点的准确性。求解方式如下：设有 n 条二维直线,其中,第 i 条直线的起点为(xi,yi),方向为(ni,vi),记矩阵 G1为一个 2n 行 2 列的矩阵,G2为一个 2n 行 n 列的矩阵,其中,G1为 n 个 2 行 2 列的单位阵叠加,G2矩阵的第 k 列的 2k 到 2(k+1)行元素分别为-ni、-vi,G1、G2横向叠加得到矩阵 G,则中心点 m为一个 2 行 1 列矩阵,矩阵两元素分别为中心点横纵坐标,中心点 m 计算方式为：（13）中心点为(m00,m10)。确定中心点后,保留每条检测线（线段）距离中心点最远的端点,将图 13 中靠近轮毂中心区域
16、的检测线的端点都去除,仅保留检测线最外侧的端点,连接中心点和该最外侧端点用于计算转角,如图 14 所示。由于风电机组的 3 支叶片是相同的,以每支叶片为参考可以计算出 3 组转角,目前以其中相对最小的转角作为计算结果。在上述条件下即计算每组由中心点-直线远端端点组成的线段与竖直方向的夹角,这些直线远端端点均是图 14 中保留的距离中心点最远的端点,其绝对值最小的角度为当前风电机组的转角。五、巡检路径规划通过无人机自带的自动测量单元将无人机的 GPS 信息与相机的图像信息进行融合,实现了对风电机组塔筒、风电机组塔高、飞机叶片姿态的准确估计。在完成计算后,根据这些数据建立风电机组初始模型。如图 15 所示,规划出 9 条自动化飞行的路径（每支叶片 3 条飞行路径）,构成一个以叶片贯穿线为中心的圆柱体包络线。巡检航线需要覆盖轮毂和风电机组的 3 支叶片,轮毂点是整个航线的中心点。风电机组有 3 支叶片,任意 2 支叶片之间呈 120,每支叶片的材料、外观形状、动力特性都完全相同,因此,可以针对抽象出的单支叶片进行航线设计,再根据风电机组角度生成所有 3 支叶片的飞行航线。单支叶片参考坐标系如
