张宏鸣 谭紫薇 韩文霆 朱珊娜 张姝茵 葛晨宇
(1.西北农林科技大学信息工程学院, 陕西杨凌 712100; 2.西北农林科技大学机械与电子工程学院, 陕西杨凌 712100)
农情遥感监测是以遥感技术为主对农业生产进行动态监测的过程。在整个植被时期监测作物是精准农业的先决条件,在精准农业中占有重要地位[1-2]。随着遥感技术不断发展,遥感监测作物株高成为可能。玉米在解决粮食安全、饲料保障、发展国民经济以及缓解能源危机等方面具有重要作用[3]。玉米植株高度可以间接地反映生物量积累,从而估算玉米产量,是进行生产调控的重要参考因素之一。因此在生长周期内,估算玉米植株高度具有重要意义。
传统的玉米株高测量方法采用刻度尺人工测量,速度慢、工作量大且准确率较低,已不能满足农业生产需要。无人机遥感系统具有运载便利,灵活性高,作业周期短,影像数据分辨率高等一系列优点[4],为大范围玉米株高信息的快速、准确、动态监测提供重要的技术手段,有效弥补了地面调查的部分缺陷[5]。
文献[6-7]利用作物的3D模型对作物的高度进行计算。文献[8-9]利用高光谱数据建立植株高度估算模型。文献[10-12]通过将航空数字影像转换为作物表面模型(Digital surface model,DSM)提取作物高度。文献[13]利用数字图像处理技术处理作物图像,从图像中取得株高等长势信息。文献[14-15]利用地面激光雷达实现作物冠层高度测量。文献[16]运用机器视觉和图像处理技术间接计算出作物株高,但是这种方法仅依靠二维图像特征对植物生长信息进行测量,植株生长中的弯曲变形对测量误差影响较大,降低误差较困难。
图像分割及识别技术在农业工程领域中得到广泛运用。文献[17-18]利用颜色特征对绿色植物图像进行分割。文献[19-20]通过图像处理和支持向量机(Support vector machine,SVM)进行植物的分割及识别。文献[21]采用直方图阈值化和形态学滤波方法实现绿色植物与土壤背景的分割。文献[22]提出一种基于植物整体外观特征提取的植物自动识别方案,但是当图像存在复杂背景时,可能会将背景分割到植物对象区域中。
目前关于无人机数字正射影像(Digital orthophoto map,DOM)与DSM结合进行玉米株高提取的研究较少。本文通过无人机获取大面积、高精度的夏玉米DOM和DSM,采用K-means算法、遗传神经网络算法和骨架算法提取玉米株高,与实地测量株高进行对比和精度评价。
实验区位于内蒙古自治区鄂尔多斯市达拉特旗昭君镇,其地势南高北低,南部为丘陵山区,北部为肥沃的黄河冲积平原。实验区域的地理位置为东经109°36′24″~109°36′28″,北纬40°25′58″~40°26′2″,常年干旱少雨,气候为温带大陆性气候。实验区为半径60 m的圆形区域,总面积约为9 498.5 m2。实验区玉米行距58 cm,株距25 cm。将其划分为5个不同水分实验区域,如图1a所示,各时期的实际灌溉量和降雨量分别通过安装在喷灌机上的流量计(MIK-2000H型)及和实验地相邻的标准气象站采集。在每个区域内选择3个地势高低有差异的6 m×6 m的实验规划区域,每个区域的对角线上平均选取3株玉米作为标记植株,在每个标记植株附近随机采集4株玉米的高度,每个样本区取5株玉米高度求平均值,以该平均值作为地面采集玉米的高度值,从而获取不同的玉米植株高度。为精确统计株高,采用2 m×2 m作为实验样本区,共计45个样本区。研究区域位置及样本区分布如图1a所示。标记植株的选取方法及实验样本区如图1b所示。部分实验样本区如图2所示。
图1 研究区域Fig.1 Study area
图2 实验样本区Fig.2 Experimental sample regions
采用大疆精灵4Pro型无人机进行遥感影像的采集。大疆精灵4Pro型无人机具有飞行稳定、续航时间长、防止快速移动过程中产生拖影、障碍感知等优点,可以稳定获取无畸变失真、可拼接的遥感影像。该无人机系统主要由飞行器、稳定云台、影像传感器等组成。其中飞行器是影像传感器及其稳定云台的搭载平台,是获取农业遥感数据的基础;稳定云台使得影像传感器在飞行过程中保持相对地面稳定的状态,从而避免了遥感影像的几何畸变,同时也保证了影像采集过程中成像角度的相对稳定。无人机及传感器如图3所示,其主要技术参数如表1所示。
图3 无人机和影像传感器Fig.3 UAV and image sensor
参数数值/类型轴距/mm350起飞质量/g1388续航时间/min30最大遥控距离/mFCC:7000影像传感器1英寸2000万像素CMOS
在无人机遥感影像获取时,为了避免云朵遮挡,选取太阳光辐射强度稳定的正午12:00,天空晴朗无云的天气情况进行采集。采用Pix4Dmapper软件快速生成专业的、精确的DOM及DSM数据。Pix4Dmapper处理数据大致流程如下:①导入影像(格式为TIF)和位置与姿态系统(Position and orientation system,POS)数据。②导入地面控制点(Ground control point,GCP)文件。③根据不同要求设置参数。④“一键式”全自动处理进行点云提取和立体模型建立,获得DOM和DSM[23]。技术路线如图4所示。
图4 技术路线图Fig.4 Technical roadmap
由Pix4Dmapper拼合得到的DSM和DOM数据的精度与无人机飞行高度、影像数量、影像重叠度、影像比例等因素有关。飞行高度越小,模型表面越光滑,具有更多的细节内容;影像数量多则产生的模型表面会更加光滑清晰,细节更加突出;重叠度越高,相邻影像间的差异越小,同名点的匹配也越容易,相对定向精度越高,建成的模型质量越好[24],而影像重叠度不够则会缺失模型信息;影像比例越大生成点越多,模型表现细节越丰富。
参照有关无人机遥感研究的设置参数[12,25-29],综合分析实验区域情况,最终选用飞行高度50 m拍摄的854幅影像作为数据源,地面分辨率为1.25 cm。基于精度要求,样区范围共布设5个地面控制点,使用实时动态定位(Real-time kinematic,RTK)进行测量,可用于空三运算和精度检测。同时用这些点来检测影像集合定位精度,保证校正影像在作物株高提取中的基本应用需求。GCP坐标系选择CGCS2000/3-degree Gauss- Kruger CM111E (egm96),拼接时Pix4Dmapper自动获取相机型号FC63108.85472x3648(84ddd3d74c736564626ec d8e10c57f19)(RGB)。其他参数设置如表2所示。
表2 Pix4Dmapper拼接主要参数Tab.2 Main parameters of splicing in Pix4Dmapper
无人机影像实验数据于2018年6月12日至2018年7月8日在内蒙古自治区鄂尔多斯市达拉特旗昭君镇进行采集,共5期数据,分别记作D0、D1、D2、D3和D4,对应时期为t0~t4。在t0时获取的数据记作D0,此时实验田为播种后至出苗前的裸土,获取其无人机可见光影像,结合GCP生成实验田的DSM,记作Ddsm0,可以得到实验田高精度的高低起伏变化情况,作为后期H数据提取的地表基准面(因为此时地面近乎为裸土,无植株高度变化,不进行株高提取实验,仅作为后期辅助计算的地表基准面)。在ti(i=1,2,3,4)时期使用与t0时相同的GCP,生成后期各生长阶段的Ddsmi(i=1,2,3,4),与Ddsm0作差可以得到对应ti时期玉米的植株高度为
Hi=Ddsmi-Ddsm0(i=1,2,3,4)
(1)
研究区域的玉米DSM如图5a所示,基于DSM的株高提取原理如图5b所示。在玉米生长过程中,实地测量值以顶端完全展开的叶子为基准测量的高度作为玉米的高度。
图5 DSM及株高提取原理Fig.5 DSM and plant height extraction principles
从图5a可看出高度为-4.581~4.775 m,这是因为两期数据的喷灌机位置不同。t0时,喷灌机位于图5a标记位置1处,此时的Ddsm0为喷灌机高度,后期Ddsm1较喷灌机Ddsm0低,则Ddsm1与Ddsm0作差产生负值;t=1时,喷灌机位于图5a标记位置2处,此时的Ddsm1为喷灌机高度,后期Ddsm1较Ddsm0高,则Ddsm1与Ddsm0作差产生正值。根据统计,99.8%的高度数据在正常范围内,可以进行研究。
采用无人机DOM与DSM结合的方法进行玉米株高提取。主要步骤包括:①获取玉米高清可见光航拍影像。②完成航拍影像的拼接以及预处理等,生成所需的DSM及DOM,选择播种前的DSM作为地表基准面,用之后测量的各生长时期的DSM与其相减得到不同时期的玉米DSM。③对DOM进行图像增强等多种预处理,通过K-means算法[30]、遗传神经网络算法[31]和骨架算法[32]得到玉米植株区域。④提取的玉米区域进行几何配准后经过遥感影像处理软件生成掩膜。⑤运用得到的掩膜与DSM套和得到影像上的玉米高度。⑥影像高度与实际株高进行对比及精度评价,得出精确度较好的株高模型。方法流程如图6所示。
图6 方法流程图Fig.6 Flow chart of method
2.1.1遥感影像预处理
无人机拍摄的原始影像仅仅记录了实验区部分区域,为了便于数据分析,需要对原始影像进行拼接,得到实验区的整体影像。本文采用Pix4Dmapper软件进行影像拼接。导入无人机可见光航拍影像、POS数据,结合5个控制点对影像进行几何校正,通过全自动空三加密,生成DOM及DSM。
2.1.2图像增强
由于光线原因会造成图像局部过亮或过暗,对图像进行拉伸,使之覆盖较大的取值区间,提高图像的对比度,便于后期提取玉米区域。图像增强公式为
J=iimadjust(I,[llow.in;hhigh.in],[llow.out;hhigh.out])
(2)
将图像I中的亮度值映射到图像J中的新值,即将llow.in至hhigh.in之间的值映射到llow.out至hhigh.out之间的值。llow.in以下与hhigh.in以上的值则被剪切掉,它们都可以使用空矩阵,默认值是[0 1]。图7a原图增强后效果如图7b所示。
2.1.3色彩空间转换
通过对RGB、HSV和YCbCr色彩空间模型[29]进行色彩空间转换并对比,如图7c、7d所示。
2.1.4OTSU阈值分割
最大类间方差法简称OTSU,可以根据图像的灰度特性,将图分为前景和背景两部分。
对于图像I(x,y),前景(即目标)和背景的分割阈值记作T,属于前景的像素点数占整幅图像的比例记为ω0,平均灰度为u0;背景像素点数占整幅图像的比例为ω1,平均灰度为u1;整幅图像的平均灰度记为u,类间方差记为g。则有
u=ω0u0+ω1u1
(3)
g=ω0(u0-u)2+ω1(u1-u)2
(4)
联立式(3)、(4)可得
g=ω0ω1(u0-u1)2
(5)
当方差g最大时,可以认为此时前景和背景差异最大,此时的灰度T为最佳阈值。
通过对比图7e、7f可发现,原始图像进行分割后会产生大量噪声,增强后能够有效解决此类问题。对比图7f、7g、7h可以发现,HSV色彩空间进行图像分割会丢失大量玉米区域信息,增强后的图像进行图像分割则存在大量非玉米植株区域,而通过YCbCr色彩空间进行图像分割减少了大量非玉米植株区域,同时保证了较完整的玉米区域信息。因此,采用YCbCr色彩空间进行玉米株高的骨架提取。
图7 图像预处理Fig.7 Image preprocessing
2.2.1骨架提取算法
骨架包含了图像特征的最有效数字化信息,能够对图像进行有效描述。它是由一些细(或者比较细)的弧线和曲线集合构成的一个能够保持原始形状相连性的表示[33]。文献[34]提出利用火烧稻草的模型进行骨架提取,即“草火法”。假设一个模型的内部被稻草填满,火从边缘的每一点以相同的速度燃烧,则中轴点即为稻草边界上的火源同时向内燃烧的相遇点。因为火烧的速度相同,所以在同一个时间相对的两个边缘点的前进距离是相同的,于是他们停止的位置同样是这两个边缘点的对称中心,也就是骨架点,如图8所示。
图8 草火法示意图Fig.8 Diagram of grass fire
一般说来,骨架必须保持3个特性:①连续性,即连通结构必须细化成连通线结构。②中心性,即骨架与图像具有结构同一性。③最小宽度为l。
骨架提取是指根据不同的定义和算法提取原始形状骨架的方法。本文主要应用Zhang-Suen骨架算法[35]结合形态学处理进行玉米骨架的提取。
假设已知目标标记为1,背景点标记为0。定义边界点为1,其8连通邻域至少有1个点标记为0。根据算法,需要对边界点进行以下处理:
(1)在以8连通邻域为中心的边界点上,中心点为P1,顺时针方向的相邻点记为P2,P3,…,P9。其中P2位于中心点P1之上(图9)。首先,选择满足要求的点
(6)
N(P1)为非零相邻点的个数,S(P1)为P2~P9~P2序列从0到1的变化个数。经过对所有边界点的检查,所有标记点都被删除。
图9 8邻域图Fig.9 Eight neighborhood diagram
(2)满足
(7)
经检查,标识点已删除。
循环上述两步骤,直到都没有像素被标记为删除为止,输出的结果即为二值图像细化后的骨架。
2.2.2玉米骨架掩膜与高度提取
2.2.2.1提取骨架制作掩膜
由图10c、10f可以看出,单像素宽骨架制作掩膜不连贯,会丢失大量玉米骨架区域信息,通过合适的结构元素对其进行形态学处理得到2像素宽骨架,可以得到较完整的玉米骨架区域。因为掩膜不能与DSM叠加显示,所以为了方便观察,使用制作掩膜的面替代掩膜进行显示。
图10 制作的掩膜Fig.10 Making mask
2.2.2.2玉米高度提取
通过K-means算法(图11a)、遗传神经网络算法(图11c)和骨架算法(图11e)分别对无人机DOM中的玉米区域进行提取,生成掩膜,与DSM套和(图11b与K-means算法套和,图11d与遗传神经网络算法套和,图11f与骨架算法套和),获取玉米高度信息。3种方法提取效果如图11所示,图中从左到右依次为D1、D2、D3和D4时期,白色部分表示玉米植株区域。
图11 玉米区域提取与DSM结合Fig.11 Maize area extraction combining with DSM
2.3.1方法对比
通过K-means算法、遗传神经网络算法及骨架算法提取玉米区域与DOM叠加(图12),可以看到本文的骨架算法在提取完整玉米植株区域的同时避免了较低叶片和大范围阴影区域对株高提取的干扰。而K-means会将大量较低叶片和阴影区域保留下来,遗传神经网络方法则会丢失部分玉米植株区域,对株高提取干扰严重。
2.3.2精度验证
(1)提取算法精度验证
为验证本文方法的准确性,用精确度S对3种提取算法效果进行评价,越大表示提取算法精确度越高。评估方法在文献[36]提出的方法基础上改进得到。计算公式为
(8)
式中S——提取算法的精确度
T——预测为玉米的玉米样本
F——预测为玉米的非玉米样本
(2)提取高度精度验证
为验证本文方法的准确性,用平均绝对误差(Mean absolute error,MAE)和均方根误差(Root mean squared error,RMSE)两个指标对模型精度进行验证,值越小,表示提取株高与实测株高越接近;用决定系数R2(Coefficient of determination)对模型拟合优度进行评价,越大表示提取株高与实测株高拟合效果越好。
图12 掩膜+DOM提取效果对比(部分)Fig.12 Comparison of mask and DOM extraction effect(part)
从图12中红框区域可以看出,K-means算法及遗传神经网络算法会将大量较低叶片和阴影区域保留下来,经过掩膜和DSM结合提取的高度会受到一定影响;从图12中黄框区域可以看出,遗传神经网络算法则会丢失部分玉米植株区域,对于玉米株高提取也会产生一定的误差。而骨架算法提取的是玉米植株的中心骨架,在提取完整玉米植株区域的同时避免了较低叶片和大范围阴影区域对株高提取的干扰,同时骨架算法在提取玉米骨架的过程中同样存在误提取的情况,如图12中蓝框区域,在苗期出现误把杂草或者石子当作玉米区域而保留下的情况,在后期生长过程中这种情况均有所减少。
对4期共180幅DOM的灰度直方图进行分析,选择合适的阈值区分玉米和其他区域(包括阴影和较低叶片),通过提取算法的精确度S对3种方法进行对比,结果如表3所示。通过K-means算法提取玉米区域的精确度最低,通过遗传神经网络算法提取玉米区域的精确度较高,通过骨架算法提取玉米区域的精确度最高。3种方法的S平均值分别为0.57、0.68、0.83。结合图12可知,D1时期阴影和叶片较小,遗传神经网络算法和骨架算法提取效果差距不大,K-means算法对土壤区域识别能力较差,对提取效果影响较大;在D2~D4期间,叶片逐渐长大、展开,产生大量阴影区域,K-means算法和遗传神经网络算法误提取大量阴影区域和较低叶片,提取效果较差;相比而言,骨架算法在D1~D4期间的提取效果较稳定,误提取的土壤区域更少,设定合适的阈值区分玉米和其他区域后,通过骨架算法可以有效减少阴影和较低叶片等区域对提取玉米区域的影响,精确度更高。
表3 玉米区域提取精度Tab.3 Maize area extraction accuracy evaluation
通过决定系数(R2)、均方根误差(RMSE)和平均绝对误差(MAE)对3种方法提取的株高进行对比,如表4所示。K-means算法提取株高的R2为0.853,RMSE为15.886 cm,MAE为13.743 cm;遗传神经网络算法提取株高的R2为0.877,RMSE为14.519 cm,MAE为11.884 cm;骨架算法提取株高的R2为0.923,RMSE为11.493 cm,MAE为8.927 cm。由此可知,通过骨架算法提取株高,效果更好。
表4 精度评价Tab.4 Evaluation of precision
(1)利用K-means算法、遗传神经网络算法和骨架算法分别对玉米区域提取,生成掩膜,并结合DSM提取株高。3种方法的R2分别为0.853、0.877、0.923,RMSE分别为15.886、14.519、11.493 cm,MAE分别为13.743、11.884、8.927 cm。表明基于玉米生长阶段的无人机DOM提取玉米区域,效果较好,结合DSM能为玉米株高提取研究提供借鉴,且骨架算法提取精度较高,可为大面积的田间株高测量提供新的技术手段。
(2)利用无人机影像进行玉米区域提取的精度受天气影响较为明显。在无人机获取数据期间,如有风则会导致实验材料的高度改变;如阴天则会影响无人机遥感影像的拼接,大量阴影则不利于玉米植株区域的提取。这些都会明显影响数据的质量,因此获取数据时尽量选择晴朗无风的天气,且在正午时分、无外界因素影响的情况下进行拍摄。
(3)土地、石子、阴影、较低叶片等区域对提取玉米区域影响较大,骨架算法可以解决多数问题,且提取效果较稳定。但对于部分阴影仍需手工进行一定干预才能达到较好的效果,可以通过结合纹理特征的方法进行进一步的研究。