Geographical Science Research
Vol.05 No.02(2016), Article ID:17510,11 pages
10.12677/GSER.2016.52008

Research on Method of “Abnormal Stripe Problem” of Vegetation Index Product of FY3A MERSI

Chuanshang Wu

School of Resources and Environment, University of Electronic Science and Technology of China, Chengdu Sichuan

Received: Apr. 15th, 2016; accepted: May 6th, 2016; published: May 9th, 2016

Copyright © 2016 by author and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY).

http://creativecommons.org/licenses/by/4.0/

ABSTRACT

There are so many problems in the NVI product that affect the availability of the NVI product. The main problem includes abnormal stripe, which is summarized in the NVI product. Apply manual method and auto method, which could solve the significant problems of the NVI and improve its quality, has been discussed. Auto method avoids wasting a lot of time to search the area of abnormal stripes. Thus the abnormal stripes are removed rapidly by using mathematical morphology to do boundary detection. Finally, the quality of FY3A MERSI has been improved and the data availability has been ensured.

Keywords:Remote Sensing, MERSI, Abnormal Stripe, Mathematical Morphology, North China

FY3A MERSI数据的植被指数产品“异常条带”问题处理方法探究

吴传赏

电子科技大学资源与环境学院,四川 成都

收稿日期:2016年4月15日;录用日期:2016年5月6日;发布日期:2016年5月9日

摘 要

目前,FY3A MERSI数据的植被指数(NVI)产品没有进行数据质量的严格控制。相应的,数据产品出现各种问题,严重影响了NVI数据产品的可用性。通过观察可发现在NVI产品之中,最严重的显著的问题是“异常条带”现象。本文首先研究手动去除“异常条带”的方法,然后根据实际需要基于数学形态学方法探究自动化“异常条带”检测和去除。该自动检测方法通过快速检测异常条带边界,极大地提高了检测和去除“异常条带”的时间。从而FY3A MERSI数据的植被指数(NVI)产品的显著问题,提升FY3A MERSI数据的植被指数(NVI)产品的质量,保证该数据的可用性。

关键词 :遥感,MERSI,异常条带,数学形态学,华北

1. 引言

归一化植被指数(Normalized Difference Vegetation Index即NDVI)是检测植被生长状态、植被覆盖度和消除部分辐射误差等。NDVI能反映出植物冠层的背景影响,如土壤、潮湿地面、雪、枯叶、粗超度等,且与植被覆盖有关。NDVI是目前使用最广泛的植被指数,它于1974年由Rouse等提出(植被指数研究进展,田久庆)。其应用领域包括:土地覆盖变化 [1] ,植被与环境因子变化 [2] - [4] ,叶面积植数(LAI) [5] 和农业产量预报 [6] 等。虽然NDVI应用研究富有成效,但NDVI本身存在着不足,例如NDVI饱和问题;对大气影响纠正不彻底;对低植被覆盖区土壤背景影响未作处理;“最大值合成法”不能保证最佳像元等 [7] 。增强型植被值数(Enhanced Vegetation Index, EVI)是对NDVI的继承和改进,它综合了搞大气植被指数(ARVI)(Kaufman et al., 1992)和土壤调节植被指数(SAVI)(Huete et al. ,1988)的优点,可以提高对高生物量地区的敏感度 [8] ,通过削弱叶冠北景信号和降低大气的影响来改善对植被的监测。增强型植被值数(EVI)依托“中分辨率成像光谱仪(MODIS)”或“(国产)中分辨率成像光谱仪(MERSI)”的丰富信息,以及过去二十年对植被指数的研究,在植被指数和合成算法方面都作了改进 [9] 。归一化植被指数和增强型植被指数的公式(1)和式(2)。当前的归一化植被指数产品如NOAA/AVHRR NDVI以及SPOT VEGETATION NDVI,增强型指被指数产品如MODIS EVI。FY3A-MERSI植被指数(NVI)产品是同时综合NDVI和EVI这两种指数产品。该指数产品对遥感相关科研及生产活动都极为重要。所以对MERSI植被指数产品之中存在的问题进行总结是有非常有必要的。

(1)

(2)

2. 研究数据及研究区介绍

FY3A MERSI数据的植被指数产品投所采用的投影系统为汉麦尔–埃托夫(Hammer-Aitoff)投影,该投影是兰勃特(Lambert)方位等积投影的改良型投影。其中心投影经度为150˚,中心纬度为0˚。MERSI 250米植被指数产品包括全球10˚ × 10˚分幅的250米Hammer投影均匀网格旬、月产品˚,二者区别在于合成的周期分别是旬和月。

卫星传感器进行测时,太阳光照角度和卫星测视角以及云的变化都很大,因此得到的是来自地表的双向反射率信息。要构造植被指数(VI)的季节性的时间曲线,需要把给定时间段内的几幅VI影像合成一张晴空的VI影像,并且要使大气效应和角度效应影响最小。目前为人们所接受的VI合成产品处理方法是最大值合成法(MVC),该方法通过云检测和质量检查后,再逐像元地比较多同VI影像,并选取最大的VI值。一般人们认为MVI倾向于选择最“晴空”的(最小光学路径)、最接近于星下点和太阳天顶角最小的像元。在MERSI植被指数合成算法中采用得主要是最大值合成法,辅以BRDF合成以及约束视角法 [10] 。

多光谱遥感数据经线性和非线性组合构成的对植被状态有一定指示意义的数值,叫植被指数(VI),是对地表植被活动的简单、有效和经验的度量。FY3A MERSI数据的植被指数(NVI)产品包含归一化植被指数(NDVI)、增强植被指数(EVI),二者皆为无量纲量。FY3A MERSI数据的植被指数(NVI)产品分辨率为250米(旬产品和月产品)。主要含有七个数据层:包括了蓝(CH1)、绿(CH2)、红(CH3)、红外(CH4)、热红(CH5)五个波段(详细信息见表1)的观测物理量值及归一化值被指数(NDVI)数据和增强型植被植数(EVI)数据。另外还有太阳天顶角、卫星天顶角、卫星方位角和质量控制信息等。本产品主要用途是用于陆地表面覆盖变化监测、农业遥感和气候变化研究。

研究区为北纬25˚~45˚,东经100˚~125˚。该区域主要包括了华北平原。华北平原是中国东部大平原的重要组成部分,位于北抵燕山,南达大别山,西倚太行山、伏牛山,东临渤海和黄海,跨越京、津、冀、鲁、豫、皖、苏7省市,面积30万平方千米。平原地势平坦,河湖众多,交通便利,经济发达,自古即为中国政治、经济、文化、交通中心。华北平原是全国主要的粮食产地,在针对华北平原进行农业遥感的研究中非常需要获取可靠的NVI数据,所以选作研究区。

3. “异常条带”问题分析及处理方法究

3.1. “异常条带”问题分析

目前,FY3A MERSI数据的植被指数(NVI)产品没有进行数据质量的严格控制,相应的,数据产品出现了显著问题,严重影响NVI数据产品的可用性。所以,有必要针对FY3A_MERSI数据的植被指数(NVI)产品中所存在的“异常条问题”进行研究并探究其处理方式。

具体的,选取时间为2014年5月下旬,分辨率为250米的FY3A MERSI数据的植被指数旬产品进行研究。FY3A MERSI数据的植被指数旬产品 采用设影系统为Hammer投影。为了便于研究我们将30B0 (Hammer投影分块位置代码)、30C0、40B0、40C0共4个区域的数据进行拼接。该区域完全覆盖华北平原。图1即为指数产品其中的红光(CH3)、红外光(CH4)、蓝光(CH1)波段进行RGB假彩色合成的影像。该影像包含NDVI和EVI所包含有的所有波段信息,且影像的对比度较高,色彩饱和度适中,易于反映和说明影像中所存在的相关问题。研究的重点是考虑被用于合成NVI产品的原始FY3A MERSI遥感卫星数据对官方植被指数(NVI)产品的影响。(研究的FY3A MERSI数据的植被指数(NVI)产品为十天合成产品。其中NDVI数据采用的的最大值合成法。该方法为国际公认方法之一,此处不再讨论其正确性。)

图1中不难发现NVI数据产品本身存在的“异常条带”问题。通过手动溯回查找2014年5月下旬FY3A MERSI数据,即2014年5月21日到2014年5月31日的所有涵盖研究区、原始分辨率为250米的MERSI数据中可以找出对应于问题条带的原始FY3A MERSI影像:FY3A_MERSI_GBAL_L1_2014 0530_0235_0250M_MS.HDF。

图2的右上部及图3分别是UTC时间2014年5月30日2点30分MERSI传感器的第1通道(蓝波段)、第3通道(红波段)以及第4通道(红外波段)获取的原始影像及投影转换后的影像。由于第3通道(红波段)以及第4通道(红外波段)出现了异常条带,因而造成了假彩色影像(3,4,1合成)上观察出明显的“异

Table 1. The 250 m resolution channel characteristics of MERSI

表1. MERSI 250米分辨率通道特性

Figure 1. The NVI product of the last ten-day of May 2014 of FY3A MERSI (3, 4, 1 false color composite image)

图1. 2014年5月下旬FY3A MERSI数据的NVI产品(3,4,1假彩色合成影像)

常条带”问题。这种连续大范围的奇异条带现象主要是由遥感影像波段数值异常造成的。遥感卫星在记录扫描行的数据时,可能受到信号干扰抑或是系统错误影响,造成了某行或若干行扫描行数据发生错误输入或者数据丢失。在合成NVI产品之前未对该遥感影像数据进行数据质量检验。于是在最终的NVI合成产品之中,会保留下原始影像中的扫描行数据错误。

3.2. “异常条带”问题处理方法探讨

3.2.1. “异常条带”的手动去除方法

图2图3进行对比可得出如下结论。

1) “异常条带”在未进行预处理的初始影像中规律性很强,“异常条带”的表现形式为行条带。而

Figure 2. The original image of FY3A MERS (red, near-infrared, blue, and composite)

图2. FY3A MERSI原始影像(红、近红、蓝及其合成图)

Figure 3. The image of FY3A MERSI after projection

图3. 投影后的FY3A MERSI影像

经过投影转换和几何校正之后的MERSI影像中,“异常条带”表现形式则较为复杂。所以在初始影像中更容易确定“异常条带”的“空间范围”,从而将该“异常条带”区域的异常值抹除。

2) “异常条带”的在初始影像中呈现方式主要是分成两种形式,如图2右下部所示,通过观察,不难发现每条异常条带几乎都有着显著的边缘水平分隔线。具体的分成两种主要的情况。一是标准模式:“异常条带”以行条带的形式存在,每个“异常条带”覆盖一个或多个扫描观测带(SCAN) (MERSI的一个扫描观测带:40个像素宽度,约两万米)。二是跳跃模式:异常条带由左至右会生成“Z字形跳跃阶梯”的形态。异常条带的边缘区别于标准模式下“异常条带”的整齐状。即条带呈现由左至右按单像素平行下移跃变。一般每个异常条带覆多个扫描观测带。在图2中即同时存在着标准模式和跳跃模式的“异常条带”。

对“异常条带”的去除可以理解为在进行预处理的过程中将MERSI原始数据各波段的原始DN值矩阵中含有“异常条带”的数据置为零值或背景值。而在手动去除“异常条带”的方法中,重点在于将“异常条带”所跨越的行范围找出来,达到删除相应的异常的“行条带”的目的。有效地去除FY3A MERSI数据的植被指数(NVI)产品的条带问题。可以通过对源文件HDF直接进行查看,查到“异常条带”的起始和结束行并读出相应的行号。然后在利用ENVI/IDL对源文件HDF进行几何较正和辐射较正的时候将原始的波段影像矩阵修正为排除掉异常的行条带的有效影像矩阵(把起始行和结束行之间的矩阵元素值定义为零值或者背景值)。在预处理过程中就不会将带有“异常条带”的错误数据带到可用的源数据之中,从而也就避免了“异常条带”问题累积到合成的MERSI植被指数产品之中。图4中云除条带后的MERSI影像。

3.2.2. “异常条带”的自动检测去除方法(数学形态学方法)

区别于手动去除“异常条带”,“异常条带”的自动检测去除方法,则更加实用。因为在做FY3A MERSI数据的植被指数产品过程,不可能一幅一幅去对影像手动进行检查,这样显然不经济,效率低。而且容易出错。

在“异常条带”的自动检测去除方法中,主要基于数学形态学,数学形态学在图像处理方面得到了日益广泛的应用。例如在边缘检测方面。边缘检测是大多数图像处理必不可少的一步,提供了物体形状的重要信息。对于二值图像,边缘检测是求一个集合A的边界,记为

(3)

对于灰度图像,边缘检测是求一幅图像的形态学梯度,记为g:

数学形态学运算用于边缘检测,存在着结构元素单一的问题。它对与结构元素同方向的边缘敏感,而与其不同方向的边缘(或噪声)会被平滑掉,即边缘的方向可以由结构元素的形状确定。但如果采用对称的结构元素,又会减弱对图像边缘的方向敏感性。所以在边缘检测中,可以考虑用多方位的形态结构元素,运用不同的结构元素的逻辑组合检测出不同方向的边缘。

梁勇等人构造了8个方向的多方位形态学结构元素,应用基本形态运算,得到8个方向的边缘检测结果,再把这些结果进行归一化运算、加权求和,得到最终的图像边缘。该算法在保持图像细节特征和平滑边缘等方面,取得了较好的效果。

在本研究中,先做出一条假设:那就是认为异常条带与正常值区域之间有一条明显的边界分隔线,在该边界线上面相领上下两行之间的行像素反射率值差异性很显著,且该显著差异能够在够在水平方向上保持稳定连续。而正常像素的反射率则呈现着截然相反的情况,即无此类型的稳定性连续显著差异。

异常条带的简化模型如图5所示,可将该问题转化为水平行边界分隔线的检测问题。代表影像中的反射率值矩阵,x代表行标,y代表列标,矩阵行数为,列数为。具体算法流程如下:

Figure 4. The image of MERSI of strip removed after projection

图4. 去条带的投影后MERSI影像

Figure 5. Sample model of abnormal strip

图5. 异常行条带简化模型

第一步:定义“相邻行差异”矩阵DIFF。

(4)

该矩阵中每一个元素值的实际含义为:显示影像矩阵里相邻的上下两行间的数值(一般为反射率值)的差值,即“相邻行差异”。

第二步:定义零一标记矩阵FLAG1。矩阵行数为,列数为。并设置一个典型阈值。若DIFF达到该阈值,则,否则,即:

(5)

该零一标记矩阵FLAG1表示含义为每个像元值是否与下一行同一列对应的像元值存在着显著差异。

第三步:定义连续数矩阵Num。矩阵行数为,列数为。同时固定行标

(6)

表示将某行连续有效的像素数累加求和。值开始随行数步长为一递增滑动。从而得到连续数矩阵Num。

第四步:定义零一标记矩阵FLAG2。矩阵行数为,列数为1。并设置一个另一个典型的阈值。若Num达到该阈值,则判定,认定为有效的异常条带分隔线,否则。即

(7)

得到一维矩阵FLAG2,即异常条带的分隔线的标志矩阵。

第五步:得到异常条带的分隔线位置之后,在相应的二维影像矩阵里,将处于分隔线之间的像值判定为异常条带区域像素即可,即将其设置为零值或背景值

上述为单波段的“异常条带”自动识别和去除方法,多波段则需要对每一个波段依次进行识别和去除(图6)。

4. 结果验证与对比分析

参照FY3A_MERSI_GBAL_L1_20140530_0235_0250M_MS.IMG (正常预处理生成的影像)同“异常条带”处理后的预处理影像进行验证对比。如表2,不难得出上述“异常条带”处理方法是正确有效的。在标准模式下,不论是手动还是自动方法,识别的正确判决概率都达到百分之百,漏检概率和虚检概率都为零。在跳跃模式下,手动方法所能检测出的“异常条带”像素数略多,但同时伴随着正确判决概率的稍微回落。主要原因是检测方法只检测结果为“异常条带”的初始行跟结束行,在跳跃模式下,“Z字形跳跃阶梯”状的异常条带,检测出的起始行和结束行为该条带所跨越的最大范围,所以必然会造成虚检问题,此模式下虚检概率在10%到15%。

表3所示,通过比较不难发现,手动去除跟自动去除方法都有很高的精确度,但是手动方法耗时很长,因为需要将参入合成的相关原始MERSI数据一一进行查看,目视判断MERSI数据是否存在“异常条带”,若存在则需要打开HDF文件,数出“异常条带”的开始和结束的行数。显然这个过程非常费时。所以只有自动化的“异常条带”检测去除方法能适应快速业务化的需要。

图7所示技术流程,对上述2014年5月下旬FY3A MERSI数据的NVI产品中存在的各种问题进

Figure 6. Algorithm flow chart of abnormal strip auto removed

图6. 自动去除“异常条带”的算法流程图

Table 2. The comparison of strip detection of manual and auto method

表2. 手动和自动“异常条带”检测方法的比较

Table 3. The difference of strip detection of manual and auto method

表3. 手动和自动“异常条带”检测方法特点对比

Figure 7. Techenique flow chart

图7. 技术流程图

行处理。首先是由FY3A MERSI NVI官方250米数据中发现存在的各种问题,然后,回溯到原始FY3A MERSI的HDF数据,找到五月下旬之中具体出问题的相关若干影像如5月30号2点35分生成的FY3A MERSI HDF影像接着原始的HDF数据经过包括地理几何校正、几何辐射校正、质量控制、拼接等一系列的预处理工作之后,在预处理的同时,加入了去“异常条带”能以确保相关的问题得到解决。再者,遵循最大合成法进行旬数据合成,之后是转换影像的投影方式即由等经纬度投影转换为HAMMER投影。得到最终的FY3A MERSI数据的NVI产品同FY3A MERSI NVI官方250米数据进行论证对比。

图8为进行“异常条带”去除之后,重新进行处理所生成的2014年5月下旬FY3A MERSI数据的NVI产品,对比官方的FY3A MERSI 植被植数产品(如图1所示)可发现显著影响产品质量的“异常条带”不存在了。图9为对应2014年5月下旬FY3A MERSI数据的NVI产品(NDVI)。同样可以发现NVI产品

Figure 8. The processed NVI product of the last ten-day of May 2014 of FY3A MERSI (false color composite image)

图8. 处理后的2014年5月下旬FY3A MERSI数据的NVI产品(假彩色合成图)

Figure 9. The processed NVI product of the last ten-day of May 2014 of FY3A MERSI (NDVI)

图9. 2014年5月下旬FY3A MERSI数据的NVI产品(NDVI)

质量明显提升。所以针对FY3A MERSI影像进行“异常条带”去除的处理是正确且必要的。

5. 结论

1) 为提升FY3A MERSI 植被指数产品质量,针对L1级HDF文件进行“异常条带”的去除是正确且必要的。

2) 为适业务化快自动化需求,进行基于形态学方法的自动“异常条带”检测和去除是必要的。

3) 本研究中的“异常条带检测方法”,在针对跳跃模式的“异常条带的检测”,会出来虚检现象,其准确率还有待近一步提高。

基金项目

国家自然科学基金(41271174)。

文章引用

吴传赏. FY3A MERSI数据的植被指数产品“异常条带”问题处理方法探究
Research on Method of “Abnormal Stripe Problem” of Vegetation Index Product of FY3A MERSI[J]. 地理科学研究, 2016, 05(02): 64-74. http://dx.doi.org/10.12677/GSER.2016.52008

参考文献 (References)

  1. 1. 赵茂盛, 符淙斌, 延晓冬, 等. 应用遥感数据研究中国植被生态系统与气候的关系[J]. 地理学报, 2001, 56(3): 287-296.

  2. 2. 张军, 葛剑平, 国庆喜. 中国东北地区主要植被类型NDV I变化与气候因子的关系[J]. 生态学报, 2001, 21(4): 522-527.

  3. 3. 李本纲, 陶澍. AVHRR NDV I与气候因子的相关分析[J]. 生态学报, 2000, 20(5): 898-902.

  4. 4. 陈云浩, 李小兵, 史培军. 1983~1992年中国陆地NDV I变化的气候因子驱动分析[J]. 植物生态学报, 2001, 25(6): 716-720.

  5. 5. 孙睿, 刘昌明, 朱启疆. 黄河流域植被覆盖度动态变化与降水的关系[J]. 地理学报, 2001, 56(6): 667-672.

  6. 6. 江东, 王乃斌, 杨小唤, 等. NDVI 曲线与农作物长势的时序互动规律[J]. 生态学报, 2002, 22(2): 247-252.

  7. 7. 李晓兵, 史培军. 基于NOAA/AVHRR 数据的中国主要植被类型NDV I变化规律研究[J]. 植物学报, 1999, 41(3): 314-324.

  8. 8. 郭志华, 彭少麟, 王伯荪. GIS和RS支持下广东省植被吸收PAR的估算及其时空分布[J]. 生态学报, 1999, 19(4): 1444-1449.

  9. 9. 王正兴, 刘闯, 等. 植被指数研究进展: 从AVHRR-NDVI到MOD IS-EVI [J]. 生态学报, 2002, 23(5): 979-987.

  10. 10. 杨军, 董超华, 等. 新一代风云极轨气象卫星业务产品及应用[M]. 科学出版社, 2010: 253-257.

期刊菜单