上海海洋大学学报  2017, Vol. 26 Issue (1): 146-153    PDF    
基于最大类间方差法和数学形态学的遥感图像潮沟提取方法
朱言江1, 韩震1,2, 和思海1, 胡旭冉1, 陈佩达1     
1. 上海海洋大学 海洋科学学院, 上海 201306;
2. 远洋渔业协同创新中心, 上海 201306
摘要: 潮沟作为潮滩的主要地貌类型之一,以长江口九段沙为研究对象,利用2015年2月15日获取的Landsat 8分辨率为15 m的全色波段遥感数据为数据源,选取了3条发育不同的潮沟。首先利用顶帽变换来消除直接利用最大类间方差法对图像亮度背景不均匀不能准确分割的问题,然后通过最大类间方差法找到一个最佳的阈值使潮沟和背景之间方差最大,得到二值化图像;接着通过形态膨胀对断裂的潮沟进行连接,并用形态去除方法剔除非目标;最后对潮沟进行骨架化提取和去除短枝处理,得到完整的潮沟信息骨架图。采用视觉分析和定量分析对提取的潮沟信息进行精度评价。结果表明,最大类间方差法和数学形态学的结合对潮沟信息提取有较好的效果,平均准确度达到93.0%,遗漏误差和冗余误差分别为7.0%和0.5%。
关键词潮沟     九段沙     最大类间方差法     数学形态学    

对于线状目标(如河流、道路)的信息提取是遥感图像处理的重要问题之一。随着遥感卫星影像分辨率的提高,在遥感影像中,处理线状目标时背景噪声会增加。同时随着道路、河流等线状目标的扩展,其光谱特性和空间几何特征也会发生相应变化,这就增加了线状目标的提取难度。此外,使用多光谱光学遥感卫星影像时,会受到同谱异物或者同物异谱现象的影响,使提取的目标并不能准确代表其空间分布特征。对于高分辨率遥感影像的“噪声”和“异物同谱”问题使用单一的方法很难得到很好的解决。目前,针对线状目标的提取大多采用多种方法的结合。吴皓等[1]针对桥梁轮廓的几何特征及其在整幅影像中所占面积的大小,利用数学形态学准确地提取出桥梁,该方法简单、有效;郑丽等[2]通过使用阈值分割和数学形态学法实现了道路的完整提取。陈翔等[3]获取长江口九段沙TerraSAR-X微波遥感数据,利用增强LEE滤波和数学形态学完成了九段沙潮沟信息的提取。刘小丹等[4]首先通过Hough变换完成对道路的方向和长度的确定,然后利用形态学对道路进行提取,结果表明两种方法的结合有利于道路信息的提取;郭永飞等[5]通过获取分辨率较高的长江口九段沙潮沟影像,利用灰度形态学完成了对潮沟信息的提取,并对潮沟进行分维研究;SINGH等[6]利用自适应的全局阈值法进行图像的阈值分割,之后通过数学形态学法从遥感图像中完成了完整道路网的提取;SHI等[7]通过使用空间分析的方法和自适应邻域法从遥感图像上准确地提取了道路的中心线。本研究对象是潮沟,它位于海陆交汇的活跃地带,是潮滩与外海进行物质和能量交换的主要通道。潮沟通常也是船只进出航道的避风和停泊场所,通过对潮沟的研究,可以为潮滩资源的合理开发利用提供科学依据。本文在掌握和分析潮沟在遥感图像上的灰度和空间几何特性的基础上,以2015年2月15日的Landsat 8获取的分辨率为15 m的全色波段遥感数据为数据源,提出了基于最大类间方差法和数学形态学法的遥感图像潮沟提取方法,选取长江口九段沙为实验区,使用Landsat 8卫星遥感全色波段数据开展研究。

1 研究区域

九段沙的地理位置位于长江入海口处,坐标为121°46'~122°15'E、31°03'~31°17'N(图 1),其正处在不断的演变和发育中。九段沙上的植物群落主要有3种,分别为互花米草、芦苇和海三棱藨草[8],它已经成为了国家重要的湿地。

图 1 研究区域地理位置 Fig. 1 Geographic position of the study area
2 研究区域数据介绍

本文主要以长江口九段沙为研究对象,利用2015年2月15日Landsat 8全色波段遥感影像为数据源,分辨率为15 m。Landsat 8全色波段遥感影像具有内容丰富、影像清晰等特点,在全色波段遥感影像上,潮沟具有非常明显的形态和线状特征,选取时刻潮情为小潮,这样更有利于潮沟信息的提取。潮沟在影像上的几何特征与道路、桥梁等较为规则的线状目标不一样,由于受到潮汐动力等因素的影响潮沟表现得更加弯曲和复杂,可以看作两条弯曲的曲线所包络的条带状区域,这种条带区域在不同位置表现出的潮沟粗细和弯曲曲率也有显著差异。根据九段沙的地理环境和地貌形态条件,选择了3条不同区域和发育条件不同的潮沟(图 2中A、B、C所示)进行研究。潮沟A位于九段沙中沙和下沙交界处,该区域的潮沟发育较为成熟。潮沟B位于下沙西部,潮沟发育成熟且弯曲曲率较大。潮沟C位于下沙的东北部,发育不是特别成熟。之所以选取这3条不同类型的潮沟,一是为了使数据更具有代表性,二是验证本文算法针对不同类型潮沟的提取效果。

图 2 长江口九段沙Landsat 8全色波段遥感影像 Fig. 2 Landsat 8 panchromatic band remote sensing image in the Jiuduan Shoal in Yangtze River Estuary
3 基于最大类间方差法和数学形态学的潮沟信息提取 3.1 技术路线

首先对读取的图像利用顶帽变换去除斑点噪声,然后利用最大类间方差法计算阈值,使用该阈值对图像进行阈值分割,将其变为二值图像,得到初步的潮沟提取结果。接着利用imdilate函数对二值图像进行形态膨胀,把初步结果中断裂的潮沟连接起来,紧接着再利用bwareopen函数来剔除连通面积小于500像素的孤立目标,最后利用bwmorph函数进行骨架化提取和去除短枝,从而得到完整的潮沟信息。潮沟信息提取技术流程图如图 3所示。

图 3 潮沟提取技术流程 Fig. 3 Technology flow of tidal channel extraction
3.2 顶帽变换

由于潮沟在遥感影像上表现出具有一定弯曲和延伸的特性,相较于周围其他地物其灰度值更低,为了增强潮沟与背景之间差异以及去除图像中的噪声,本文选用顶帽变换对图像进行滤波处理。顶帽变换适合于大范围信息提取,同时可以消除背景亮度不均匀。

srt为要处理的图像,element为所使用的结构元素,dst为结果图像。结构元素的选择决定了顶帽变换滤波的效果[9]。形态学顶帽变换的公式如下:

    (1)

下面以潮沟C为例说明顶帽变换效果。选取半径为10的圆盘结构元素,对包含潮沟C的图像(图 4a)进行顶帽变换,结果如图 4(b)所示。对比图 4(a)图 4(b),可以看出有效地抑制了背景信息,增强了潮沟和非潮沟的对比度,这有利于后续的潮沟信息提取。

图 4 潮沟C的原始影像(a)和顶帽变换影像(b) Fig. 4 Original image of tidal channel C(a)and top-hat transformation image(b)
3.3 最大类间方差法阈值分割

阈值分割是一种较为有效的区分水体和非水体目标的方法。本文利用最大类间方差法,该方法通过找出一个最佳阈值把一幅图像分为两个部分即目标和背景,这个阈值使两个部分之间的方差达到最大[10]

首先假设图像包含L个灰度级[1,2,3.....,L],灰度值为i的像素点数为ni。因此,整幅图像的像素和为,灰度值为i的像素点的概率可以表示为:

    (2)

门限k将整幅图像分为两类,将潮沟目标定为 a1,与潮沟目标灰度差异较大的目标定为a2。因此,两类目标所对应出现的概率分别用w1w2表示:

    (3)

因此,潮沟和非潮沟目标a1,a2所表示的平均灰度μ1,μ2以及整幅图像的平均灰度μT分别为:

    (4)

整幅图像的总方差σT2及类间方差σB2分别为:

    (5)

k*可以表示为:

    (6)

以潮沟C为例,对顶帽变换后的图像利用最大类间方差法计算出阈值为92,使用该阈值进行分割,可以较好地从背景中提取出潮沟目标(图 5)。在裁剪的图像中,一些不属于潮沟C的潮沟也被包含在里面,阈值分割后,这些潮沟也在分割结果中,下一步需要将它们剔除。从图 5我们还可以看出,潮沟在提取出来以后,潮沟C也发生了一些断裂。另外,结果中还包含一些大小不一的噪声。为了能够得到完整的潮沟信息,同时剔除非目标的干扰,本次研究采用数学形态学来进一步处理。

图 5 潮沟C的阈值化图像 Fig. 5 The threshold image of tidal channel C
3.4 数学形态学处理

数学形态学主要是以膨胀、腐蚀为基础,使用这两种简单基础的运算来推导其他较为复杂的运算以实现对图像的处理,具有计算快、算法较为简单的特点。本文简单介绍一下4种基本的算子:

(1) 膨胀和腐蚀,记A 为要处理的图像,B为结构元素:

    (7)
    (8)

(2) 结合膨胀和腐蚀的其它运算主要有:

    (9)
    (10)

图 5可以看出,使用阈值分割能较好地提取出潮沟C,但同时也使得有些潮沟发生断裂,还会出现噪声。考虑到潮沟主要表现为细长状,它们的曲率变化较大,所以选取尺寸较小的结构元素进行膨胀处理,这样一是可以使断裂的潮沟得到连接,二是不会导致潮沟真实信息的丢失。这里针对潮沟C采用正方形结构元素,如图 6所示。

图 6 正方形结构元素 Fig. 6 Square structural elements

利用正方形结构元素对二值化后的潮沟图像进行膨胀,将原先断开的潮沟进行连接,然后利用形态去除bwareaopen函数来剔除图像中孤立的较小的目标,一般用它去除连通面积小于500个像素的孤立目标,500个像素的设定主要根据选取潮沟像元和非潮沟像元数目差异。从图 5中可以看出,提取出来的潮沟C目标的连通面积中的像素远远大于非目标的连通面积,经计算对于潮沟C此处孤立像素的值设定为500(图 7);从图 7中可以看到一些细小断裂的潮沟经过膨胀实现了连接,之后的孤立小目标去除很好地去除了不在研究范围之内的潮沟和一些非潮沟目标。

图 7 潮沟C形态膨胀和去除后信息图 Fig. 7 The information image of tidal channel C morphological dilation and removal

经过膨胀和形态去除后潮沟的宽度发生了变化,比真实的潮沟更宽了。本文通过骨架算法[11]和去除短枝[12]来处理。以下是处理思路和过程,通过提取潮沟的完整骨架信息来研究其形态变化,将经过处理后的潮沟特征图像进行骨架化处理以得到潮沟的中心线。骨架化所用的函数为BW2=bwmorph(BW,‘skel’,n),n=Inf。其主要计算思想是通过提取次数来去除目标边缘上的像素点,同时不至于导致其断裂。提取骨架的次数n同样根据选取潮沟像元和非潮沟像元数目差异进行设定。本次潮沟C的骨架提取n的次数为3。最后使用去除短枝来完成由于骨架化结果所产生的非潮沟像素点。其函数表达式为BW2=bwmorph(BW,'spur',n),n为运行的次数。本文中潮沟C选取n的运行次数为4次。最后得到的潮沟C的提取效果和骨架信息如图 8所示。图 8(a)是经过阈值分割和孤立小目标去除后的图像,并没有改变潮沟的真实宽度和长度信息。图 8(b)是潮沟C的骨架信息图。图 8(b)相对于图 8(a),潮沟更加完整,例如图 8中圆圈标记处。

图 8 潮沟C的提取效果图(a)和骨架信息图(b) Fig. 8 Extracted tidal channel C(a) and skeleton image(b)

利用同样的方法对潮沟A和潮沟B进行提取,得到潮沟A和潮沟B的提取效果和骨架信息如图 9图 10所示。

图 9 潮沟A的提取效果图(a)和骨架信息图(b) Fig. 9 Extracted tidal channel A(a) and skeleton image(b)
图 10 潮沟B的提取效果图(a)和骨架信息图(b) Fig. 10 Extracted tidal channel B(a) and skeleton image(b)
4 实验结果的精确度评价

主要采用两种评价方法对潮沟进行精度评价,一是视觉分析的评价方法;二是定量分析的评价方法。

4.1 视觉分析方法

该方法是通过把提取出来的3条潮沟骨架信息图和原始潮沟图像叠合在一起。通过目视解译的方法来分析本文方法的提取效果。为了提高视觉分析的效果,把提取的潮沟骨架信息图转变为红色。对3条潮沟骨架信息图和原始影像进行叠加的结果如图 11所示。由图 11可以看出,采用本文方法对潮沟信息提取是有效的,并且提取出的潮沟骨架信息和原始影像中的潮沟吻合度很高。

图 11 潮沟A、B、C视觉分析评价 Fig. 11 Visual analysis and evaluation of tidal channels A,B and C
4.2 定量分析法

由于潮沟是线状目标之一,与道路相比,潮沟更加弯曲。这里把潮沟的长度作为其中一个效果评价参数,有3个评价标准分别为准确度(P)、遗漏误差(Q)和冗余误差(R)。其公式表达如下:

    (11)
    (12)
    (13)

式中:线性目标总长度(i)表示所提取出的每条潮沟的总长度。遗漏的线性目标长度(b)指的是把潮沟作为背景噪声没有被提取出来的长度,多余的线性目标长度(c)指的是把背景噪声当作潮沟提取出来的长度,正确提取的潮沟总长度用a表示。利用软件Arcgis对提取出来的潮沟进行长度统计,统计结果如表 1所示。从表 1中可以看出,提取平均准确度可以达到93.0%,遗漏误差和冗余误差分别为7.0%和0.5%。造成3条潮沟准确度不一样且差异较大的原因可能一是由于图像在顶帽变换时对结构元素和元素尺寸的选择,二是由于阈值分割影响。

表 1 潮沟提取精度统计值 Tab.1 Tidal channel extraction accuracy statistics
5 结论

本次研究使用Landsat 8全色波段遥感数据,利用最大类间方差法和数学形态学法对长江口九段沙潮沟信息进行了提取。在使用圆盘结构元素滤波的基础上,考虑潮沟的几何特征,结合多种结构元素对图像进行处理,最终较好地将潮沟目标信息提取出来。将提取的潮沟骨架信息图与原始图像进行视觉对比分析和定量分析,结果表明提取的结果与潮沟实际分布相吻合且准确度达到93.0%。

本文的难点有以下几点:一是如何选择较为理想的滤波器来去除背景噪声的影响;二是如何正确利用不同的结构元素和尺寸对图像进行处理,同时如何使用结构元素来处理由于膨胀后目标发生的形态变化,这几个难点都值得进一步去分析和研究。

参考文献
[1] 吴皓, 刘政凯, 张荣.TM图像中桥梁目标识别方法的研究[J]. 遥感学报, 2003, 7(6): 478–484. WU H, LIU Z K, ZHANG R.A study of bridge recognition from Landsat TM images[J]. Journal of Remote Sensing, 2003, 7(6): 478–484.
[2] 郑丽, 潘建平.基于数学形态学的遥感图像道路提取[J]. 铁道勘察, 2010, 36(1): 12–15. ZHENG L, PAN J P.Extraction of roads from remote sensing images based on mathematical morphology[J]. Railway Investigation and Surveying, 2010, 36(1): 12–15.
[3] 陈翔, 韩震.TerraSAR-X在长江口九段沙潮沟信息提取中的应用[J]. 海洋湖沼通报, 2012: 25–30. CHEN X, HAN Z.Tidal channels information extraction of Jiuduansha in the Yangtze River Estuary by TerraSAR-X data[J]. Transactions of Oceanology and Limnology, 2012: 25–30.
[4] 刘小丹, 刘岩.基于Hough变换和路径形态学的城区道路提取[J]. 计算机工程, 2012, 38(6): 265–268. LIU X D, LIU Y.Urban road extraction based on Hough transform and path morphology[J]. Computer Engineering, 2012, 38(6): 265–268.
[5] 郭永飞, 韩震.基于SPOT遥感影像的九段沙潮沟信息提取及分维研究[J]. 海洋与湖沼, 2013, 44(6): 1436–1441. GUO Y F, HAN Z.Information extraction and fractal dimensions based on SPOT remote sensing image for tidal channels in Jiuduan Shaol[J]. Oceanologia et Limnologia Sinica, 2013, 44(6): 1436–1441.
[6] SINGH P P, GARG R D.Automatic road extraction from high resolution satellite image using adaptive global thresholding and morphological operations[J]. Journal of the Indian Society of Remote Sensing, 2013, 41(3): 631–640. DOI:10.1007/s12524-012-0241-4
[7] SHI W Z, MIAO Z L, DEBAYLE J.An integrated method for urban main-road centerline extraction from optical remotely sensed imagery[J]. IEEE Transactions on Geoscience and Remote Sensing, 2014, 52(6): 3359–3372. DOI:10.1109/TGRS.2013.2272593
[8] 唐承佳, 陆健健.长江口九段沙植物群落研究[J]. 生态学报, 2003, 23(2): 399–403. TANG C J, LU J J.Studies on plant community on the Jiuduansha Shoals at the Yangtze Estuary[J]. Acta Ecologica Sinica, 2003, 23(2): 399–403.
[9] BAI X Z, ZHOU F G. Multi scale top-hat transform based algorithm for image enhancement[C]//Proceedings of the 10th International Conference on Signal Processing (ICSP). Beijing: IEEE, 2010: 797-800.
[10] 钟雪君.一种改进的Otsu双阈值二值化图像分割方法[J]. 电子世界, 2003: 104. ZHONG X J.Improved Otsu dual-threshold binary image segmentation method[J]. Electronics World, 2003: 104.
[11] GAETANO R, ZERUBIA J, SCARPA G, et al. Morphological road segmentation in urban areas from high resolution satellite images[C]//Proceedings of the 17th International Conference on Digital Signal Processing (DSP). Corfu Greece: IEEE, 2011: 1-8.
[12] ZHU C, SHI W, PESARESI M, et al.The recognition of road network from high-resolution satellite remotely sensed data using image morphological characteristics[J]. International Journal of Remote Sensing, 2005, 26(24): 5493–5508. DOI:10.1080/01431160500300354
Remote sensing image extraction of tidal channels based on Otsu and mathematical morphology
ZHU Yanjiang1, HAN Zhen1,2, HE Sihai1, HU Xuran1, CHEN Peida1     
1. College of Marine Sciences, Shanghai Ocean University, Shanghai 201306, China;
2. Collaborative Innovation Center for Distant-water Fisheries, Shanghai 201306, China
Abstract: Tidal channel is one of the major tidal flat landforms. The Jiuduan Shoal in Yangtze River Estuary was taken as the object in this paper. Three different types of tidal channels in panchromatic band remote sensing data received from Landsat 8 satellite on February 15, 2015 were selected. Firstly, the top-hat transformation is to eliminate the problem of the backgrounds of uneven brightness accurately, if the method of maximum between-class variance(Otsu) is used directly; Then, Otsu is used to find an optimal threshold to make the maximum variance between tidal channel and background, obtaining a binary image; Thirdly, the broken tidal channel is connected by morphological dilation, and non-target is removed by morphological removal. Finally the skeleton of tidal channel is extracted and cropped, and the complete information of tidal channel is gained. Visual analysis and quantitative analysis are used to evaluate the accuracy of the information of tidal channel. The result shows that the method combining Otsu with mathematical morphology can be used to extract tidal channel information efficiently. The average accuracy can reach 93.0%, missing error and redundance error are 7.0% and 0.5% respectively.
Key words: tidal channel     Jiuduan Shoal     Otsu     mathematical morphology