龙 翔,赵庆展 ,王学文,马永建,江 萍
(1.石河子大学机械电气工程学院,新疆石河子 832000;2.兵团空间信息工程技术研究中心,新疆石河子 832000; 3.石河子大学信息科学与技术学院,新疆石河子 832000 )
【研究意义】新疆棉花产量占全国棉花总产量2/3以上[1]。对棉花的生长状况监测结果一直是水肥管理、产量预估的重要依据。近年来无人机载高光谱传感器大量应用于遥感应用领域,该方式可获得高时空分辨率的高光谱影像,提高遥感监测精度与效率[2,3]。与使用无人机搭载高光谱传感器获取冬小麦多期高光谱影像,并对冬小麦叶面积指数反演,基于比值植被指数的对数形式建立的反演模型,展现出较优的线性关系。高光谱影像中普遍存在混合像元,影响图像分类精度以及目标监测效果,制约定量高光谱遥感的发展与应用[4]。【前人研究进展】使用高光谱遥感技术在作物遥感监测应用领域开展了研究[5]。孙莉等[6]通过棉花非成像光谱开展了不同生育期高光谱数据与叶面积、叶绿素含量的相关性分析,得到叶面积、叶绿素变化呈抛物线变化趋势,在花期达到最大并且与红边特征显著相关。王克如等[7]使用高光谱数据分析棉花盛蕾期、盛花期、盛铃期和吐絮期冠层反射光谱参数与植株氮含量间的关系,可见光波段为整个生育期氮含量敏感波段。陈兵[8]使用ASD光谱仪对不同胁迫下棉花光谱变化情况进行研究,选取14个光谱特征建立相应的判别式对不同胁迫棉叶进行准确识别。孙中宇等[9]利用无人机遥感对银叶树群健康进行准确快速地诊断。田明璐等[10]采集无人机载成像棉花高光谱影像对不同生育期叶绿素含量的反演取得较好的效果。在识别纯净像元与提取端元波谱曲线方面,提出了多种高光谱影像端元提取算法,如纯净像元指数算法(pixel purity index,PPI)[4],基于最大单形体体积的N-FINDR算法[11],顶点成分分析算法(vertex component analysis ,VCA[12]光谱最小信息熵算法(spectral minimum Shannon entropy,SMSE)[13],正交子空间投影(Orthogonal Subspace Projection, OSP)[14],离散粒子群优化算法(discrete particle swarm optimization,D-PSO)[15]等。其中PPI算法与N-FINDR算法由于其易于实现、提取效果良好得到广泛应用。韩仲志等[16]使用N-FINDR算法分区提取海洋溢油波谱曲线,并计算出溢油覆盖度。唐晓燕等[17]利用N-FINDR算法对飞机与混凝土的光谱进行有效提取。徐君等[18]采用了PPI算法改进的数学形态学端元提取算法进行研究。杨鹏飞等[19]通过将PPI算法与N-FINDR算法相结合,提高了AVIRIS高光谱数据中6种矿石端元提取精度。【本研究切入点】作物冠层光谱与生长期之间有紧密联系,无人机载成像高光谱技术对作物在时序变化方面的遥感监测具有明显优势与显著成效。目前在机载高光谱影像中提取波谱曲线的方法多采用划取感兴趣区域,统计区域内平均光谱,而采用端元提取算法提取农作物的波谱曲线,进行相应波谱分析的研究尚且较少。【拟解决的关键问题】基于无人机载成像高光谱仪的优点,采集棉花从花期到后期的高光谱数据,比较PPI算法和N-FINDR算法在机载高光谱影像上的端元提取效果,以及100 m以内不同航高下提取情况,对比多期提取结果,分析棉花在不同生育期各光谱特征和植被指数变化,为基于无人机多光谱、高光谱遥感的作物监测提供借鉴和参考。
1.1.1 研究区概况
研究区位于新疆塔城地区沙湾县蘑菇湖村(85°51'49.44''E,44°25'26.61''N),地势平坦,平均海拔高度450.8 m。该区域地处天山北麓中段,属典型的温带大陆性气候,冬季长而严寒,夏季短而炎热。1年中的最高气温出现在7月,平均气温25.1~26.1℃。年降水量为125.0~207.7 mm,年日照时数为2 721~2 818 h。研究区棉花于2019年4月11日种植,4月20日第1次灌溉,机采棉种植方式,1 m 15穴,1膜6行,正常水肥管理。
1.1.2 机载成像高光谱数据采集系统
高光谱影像采集系统由飞行平台、增稳云台、高光谱相机组成。飞行平台为大疆无人机经纬M600 PRO,机身质量9.5 kg,最大负载6 kg,最大水平飞行速度65 km/h,可承受最大风速8 m/s,安全飞行时间16 min。增稳云台为大疆RONIN,云台质量4.2 kg,最大工作时间6 h。成像高光谱相机为Rikola高光谱成像仪,横向视场角和纵向视场角为36.5°。图像分辨率1 010×1 010 pixel,光谱范围500~900 nm,焦距9 mm,质量0.72 kg,适宜搭载在无人机飞行平台上作业,在60、80、100 m航高下空间分辨率分别为3.9、5.2、6.5 cm,波段数设为42,光谱分辨率10 nm。
1.1.3 高光谱数据获取
数据采集于2019年6月29日(花期)、7月11日(盛花期)、7月18日(铃期)、8月7日(盛铃期)、8月27日(后期),在研究区共采集5期机载成像高光谱数据。 Rikola高光谱传感器在无人机平台上获取成像高光谱数据时,通过实时获取上行辐射来校正DN值。在航摄范围内放置4块漫反射板(反射率分别为3%、22%、48%、64%)作为辐射校正靶标。测量日天气晴朗、无云无风,为保证有足够的太阳高度角,于12:00~14:00采集。无人机每期飞3个架次,飞行高度分别为60、80、100 m,每架次飞行时间约13 min,多期影像使用同一个任务航线。使用SR-3500便携式光谱仪同步采集棉花冠层光谱数据以作验证,光谱范围350~2 500 nm,光谱分辨率为1 nm。测量时传感器镜头垂直向下,距离冠层顶端垂直高度约为1.5 m[20]。每个测量点采集10条波谱曲线,取平均值作为该区的棉花冠层光谱反射率。
使用端元提取算法提取出非成像高光谱数据,分析棉花在不同生育期的反射光谱特征变化。包括高光谱影像获取与预处理、棉花端元提取、光谱特征分析三部分。
对采集的无人机成像高光谱数据进行暗电流校正、波段配准、影像拼接等,获得正射高光谱影像。通过影像中辐射靶标建立线性拟合模型,将辐亮度值校正为反射率,完成数据预处理。
选取7月11日采集的高光谱影像进行棉花端元提取,对比PPI和N-FINDR算法分别在60、80、100 m航高下端元提取结果,将2种算法提取结果重采样后与SR-3500数据进行波谱角计算,波谱角度越小则数据接近程度越高。为减小误差与对计算机性能的需求,将实验区的影像分为16个17.5 m×19 m的小区域,对每个小区域分别提取端元,求其均值作为该期棉花波谱曲线。棉花影像中含有棉花、土壤2种端元,为使提取的端元有余量,减少图像异常像素带来的误差,每个区域提取3条波谱曲线作为待分析端元。所有程序均在MATLAB2018b中运行。选择N-FINDR算法在每期100 m航高无人机高光谱影像中提取棉花波谱曲线作为光谱分析数据源。
利用提取出的5期非成像高光谱数据,对比分析其植被指数(NDVI、RVI、EVI、DVI、OSAVI)和光谱变量(Rg、Rrg、Vreg、R802)在不同生育期的变化情况,得到在棉花生长过程中其反射光谱特征变化规律。
1.2.1 光谱特征参数选择
选择9个常用的光谱特征参数。表1
表1 光谱参数Table 1 Table of spectral parameters
1.2.2 N-FINDR端元提取算法
N-FINDR算法是根据凸面几何学理论,寻找体积最大的单形体自动提取图像中的所有端元。在高光谱影像的特征空间中,像素点占据了一个由单形体组成的空间,单形体就是能包含所有给定维数的最简单几何体[11],单形体的顶点则是较为纯净的像元,单形体的中心是较为混合的像元。N-FINDR算法需有先验知识获得端元个数n,随机选取n个像元作为一组数据,计算其单形体体积,将每个顶点依次遍历整个影像的像素,若体积增加则用新像素代替原来的像素,最后得到体积最大的几何体的各个顶点即是端元[19]。n个像素e1,e2,…,en形成的凸锥结构体积计算公式[19]如下:
其中ei表示第i个端元向量,V表示由e1,e2,…,en组成的几何单形体体积,n为端元数目个数。
在进行N-FINDR算法端元提取之前,需对原始图像进行MNF降维至n-1维,将原始图像与MNF降维结果放在同一MATLAB程序里面运行。
1.2.3 PPI端元提取算法
PPI算法同样基于凸面几何学理论,凸面单形体的特点决定了位于单形体顶点的样本点在特征空间中任意直线上的投影必在直线的两端[18]。二维空间中投影示意中,圆点是高光谱影像中的所有像素,生成随机向量q、b穿过这个几何单形体,在单形体顶端的纯净像元p投影在这些随机向量上时,总是会落在q、b上投影区间的两端。在特征空间中生成大量的随机直线进行投影,统计每个像元落在线段2个端点的次数作为像元纯净指数,也即是生成大量的随机向量与每个像素代表的向量进行点乘,记录计算结果的最大值与最小值对应的像素与次数,次数作为纯净像元指数来定量表示该像元的纯净度。图2
1.2.4 光谱角匹配(SAM)
光谱角度定义为两地物光谱矢量之间的广义夹角,将像元n个波段的光谱响应(该像元的光谱曲线)视作为n维空间中的矢量,通过计算2条光谱曲线的广义夹角来表征其匹配程度,夹角越接近于0,两者越相似,计算公式[23]为
其中θ表示光谱角度值,t、r表示不同的两条波谱曲线,ti、ri表示2条波谱曲线在第i波段的反射率,n为波段数。
在可见光波段反射率主要取决于冠层叶绿素含量与盖度[24]。
获取的高光谱数据使用Hyperspectral Imager软件进行镜头渐晕校正和暗电流校正,去除感光单元不均匀和边缘光线衰弱及暗电流误差。通过RegMosaic软件进行波段配准,将成像高光谱仪在拍摄影像时多个波段沿飞行方向上的错位进行修正,POS信息与航摄照片输入PhotoScan软件中拼接生成正射影像。将高光谱影像中包含的4块辐射校正靶标辐亮度值通过ENVI用感兴趣区域工具统计,与靶标的实际反射率值进行匹配,计算增益和偏移参数。以辐射靶标为参考,利用最小二乘法建立每个波段的经验线性校正模型,估算每个波段反射率,将高光谱影像由辐亮度值校正为反射率,其计算公式为:
ri=R(i)×bref_i+cref_i.
ri为i波段的地表反射率,R(i)为遥感影像i波段的辐亮度值,bref_i和cref_i为i波段的反射率转换因子,近似理解为高光谱传感器的增益和偏移。
研究表明,L1与L2与土壤反射光谱曲线相似,随波长的增加,反射率呈逐渐升高趋势,但又不完全与土壤相同,土壤与棉花的光谱混合,L3波谱曲线符合绿色植被反射规律,为棉花的波谱曲线,提取出混合像元、棉花的波谱曲线。在PPI算法提取结果中,L1、L2、L3曲线相似,而影像中应该有土壤与棉花2种端元,在7月11日棉花生长较为茂盛,棉花像素个数远大于土壤像素时,PPI算法只能提取出棉花的波谱曲线,未能成功提取出土壤波谱曲线。图3
研究表明,60、80、100 m的航高下提取的3条棉花波谱曲线形态相似,绿峰、红谷、红边等特征相近。3个航高下影像提取结果的光谱角为0.065 8、0.065 9、0.067 7,PPI算法结果为0.070 1、0.072 6、0.071 1。N-FINDR算法对棉花高光谱影像端元提取效果优于PPI算法,2种算法之间波谱角差异在5.02%~10.1%;同一种算法在不同航高下提取结果差异不大,波谱角差异均在2%以内,在100 m高度以内,航高带来的空间分辨率变化对棉花波谱端元提取影响较小。R2均在0.99以上,呈显著相关。图4,图5
研究表明,不同生长期的棉花冠层波谱曲线趋势相同,都具有绿色植被冠层典型的反射光谱特征。在波长550 nm左右呈现反射峰(绿峰),这是由于叶绿素对蓝光、红光强吸收引起的。波谱曲线在650~700 nm波长范围之间呈现吸收谷(红谷),这是由于叶片进行光合作用时对红光强吸收的结果。在780 nm之后的近红外波段呈现强反射现象,这导致在700~780 nm冠层反射率呈现急剧上升趋势。不同生长期波谱曲线在各波长处的反射率具有差别,不同时期的棉花光谱在近红外波段范围(780~850 nm)内差异最为明显,7月18日(铃期)波谱的反射率在此处最高。各期光谱在550 nm左右与680 nm左右处反射率差异也较大,其余波长处反射率差异不明显。一阶微分的峰值分别出现在500~550 nm与730~750 nm,各时期的峰值差异较大,且随时间的增加呈显著的规律性变化,可较好区分棉花不同生长期。图6
研究表明,从6月29日(花期)到7月18日(铃期),随着生长期的推移,棉花的光合作用逐渐增加,对红光、蓝光的吸收逐渐增强。从7月18日到8月27日(后期),光合作用减弱,对红光的吸收逐渐减弱。红谷深度可反映叶片光合作用,红谷越深,代表对红光的吸收越强,光合作用也即越强。随着生育期的推移,红谷值呈抛物线趋势,在7月18号吸收谷最深(谷底反射率0.027 9)且宽度最大。近红外波段反射率均大于50%,由棉花冠层结构与棉叶内部结构共同导致的,从花期到铃期期间,棉花叶面积指数与覆盖度增加,使得近红外波段反射率升高,从铃期到后期期间,棉花养分往棉铃上传输,棉叶内部结构发生变化,使得近红外波段反射率降低,802 nm处反射率呈现出抛物线形,于7月18日(铃期)达到最大值0.625 2,这与红谷的变化一致。红边也呈现先增大,后减小趋势,于7月18日达到峰值0.111。由于叶绿素、叶黄素等生物化学成分的吸收,可见光波段反射率在15%以下,且从花期到后期呈逐渐减小的趋势,绿峰值由0.111 4下降到0.080 5。表2
表2 光谱特征Table 2 Statistical Table of spectral characteristics
研究表明,植被指数均呈现一种先上升,再下降的趋势,在7月11日(盛花期)或7月18日(铃期)达到峰值,棉花在此时进行大量光合作用,对近红外波段的强反射和对蓝光、红光的强吸收现象最强烈,棉花在棉田中分布密度在铃期之后变化不大的特点。表3
表3 植被指数计算Table 3 Calculation results of vegetation index
3.1N-FINDR算法可在棉花无人机高光谱影像中准确地提取出棉花波谱曲线,其在3个航高影像中的结果与SR-3500光谱数据波谱角计算结果分别为0.065 8、0.065 9、0.067 7,R2分别为0.992、0.991、0.991,提取结果可靠。Rikola高光谱仪和SR-3500便携式光谱仪采集的波谱曲线虽然很接近,但依然有细微差异,一方面是由于不同采集仪器的波谱响应函数差异和不同采集方式、不同采集平台的入瞳辐射量差异,辐射校正不能完全消除入瞳辐射量差异引起的波谱曲线变化;另一方面是在机载高光谱影像中,一个像素代表的地面面积较大,棉花像素会受到土地光谱的影像,故2种光谱仪测量的结果肯定有差异。
3.2N-FINDR算法提取结果优于PPI算法结果,两者差异在5%~10%。PPI算法易忽略小样本端元,在实验采集的机载棉花高光谱影像中,棉花像素远远多于土壤像素,使用PPI端元提取算法不能够有效得提取出土壤端元。其次PPI算法易将受到土壤光谱影响的棉花像素视作端元提取出来,这是由于PPI算法原理是统计像素的纯净像元指数的大小,当混合像元大量存在与影像中时,受到少量土壤光谱影响的棉花像素由于其光谱差异性与像素数量两方面因素共同作用,易被视作端元提取出来;N-FINDR算法在每个小区域内遍历所有像素点,找出差异最大的波谱曲线视作端元,不受像素数量因素的影响,有效地将每种地物的端元提取出来。
3.3航高在100 m以下时,其与端元提取结果相关性甚微,2种算法在60、80和100 m航高下提取结果的波谱角差异均在2%以内,表示在100 m以下,由航高变化导致的空间分辨率的变化,对棉花的光谱曲线影响不大,在飞行设置时可优先考虑航高。
3.4对比棉花不同生长期反射光谱特征,可见光波段反射率随生长期推移而降低,红谷、红边等光谱特征和NDVI、RVI等植被指数变化情况呈抛物线形,均在盛花期或铃期达到峰值,表明此阶段光合作用最强烈,对红光的吸收最显著,棉株在大量聚集营养为结铃做准备。铃期之后棉花株高、空间分布已变化不大。这与孙莉等[6]得到的棉花叶面积指数、叶绿素含量变化趋势吻合。反映并验证了棉花不同生长期的反射光谱变化具有规律性,且可作为棉花长势监测、空间分布、理化指标变化的数据参考。
无人机航摄棉花高光谱影像时,100 m以内的航高变化对棉花光谱曲线的影响在2%以内。在棉花无人机高光谱影像的棉花端元提取任务中,N-FINDR算法可实现棉花波谱曲线的提取,其结果与地面测量光谱的R2在0.99以上。不同时期棉花的红谷、红边、NDVI等参数呈抛物线趋势变化,于铃期或盛花期达到极值0.027 9、0.111、0.841 7,棉花的生物化学成分与叶片结构呈规律性变化,使其在7月光合作用最大,对红光的强吸收、近红外波段强反射现象最为明显。