CN107067384B - 一种基于射线法的黄土高原梯田提取方法 - Google Patents
一种基于射线法的黄土高原梯田提取方法 Download PDFInfo
- Publication number
- CN107067384B CN107067384B CN201710234686.4A CN201710234686A CN107067384B CN 107067384 B CN107067384 B CN 107067384B CN 201710234686 A CN201710234686 A CN 201710234686A CN 107067384 B CN107067384 B CN 107067384B
- Authority
- CN
- China
- Prior art keywords
- terraced fields
- fields
- ray
- terraced
- point
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 52
- 238000004458 analytical method Methods 0.000 claims abstract description 8
- 239000000284 extract Substances 0.000 claims description 8
- 238000001914 filtration Methods 0.000 claims description 2
- 238000004364 calculation method Methods 0.000 abstract description 2
- 238000003672 processing method Methods 0.000 abstract 1
- 238000003708 edge detection Methods 0.000 description 6
- 238000000605 extraction Methods 0.000 description 6
- 239000002689 soil Substances 0.000 description 6
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 6
- 238000002474 experimental method Methods 0.000 description 5
- VMXUWOKSQNHOCA-UKTHLTGXSA-N ranitidine Chemical compound [O-][N+](=O)\C=C(/NC)NCCSCC1=CC=C(CN(C)C)O1 VMXUWOKSQNHOCA-UKTHLTGXSA-N 0.000 description 4
- 238000012545 processing Methods 0.000 description 3
- 238000013139 quantization Methods 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 238000004422 calculation algorithm Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000005855 radiation Effects 0.000 description 2
- 238000004162 soil erosion Methods 0.000 description 2
- 238000001228 spectrum Methods 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- 238000011426 transformation method Methods 0.000 description 2
- 241001465754 Metazoa Species 0.000 description 1
- 230000006978 adaptation Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 230000003628 erosive effect Effects 0.000 description 1
- 239000003337 fertilizer Substances 0.000 description 1
- 238000003709 image segmentation Methods 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 230000001788 irregular Effects 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 238000000638 solvent extraction Methods 0.000 description 1
- 238000007619 statistical method Methods 0.000 description 1
- 238000003971 tillage Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/20—Image enhancement or restoration using local operators
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/73—Deblurring; Sharpening
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30181—Earth observation
- G06T2207/30188—Vegetation; Agriculture
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Image Processing (AREA)
- Image Analysis (AREA)
Abstract
本发明提供了一种基于射线法的黄土高原梯田提取方法,首先分析梯田田坎的排列方式以及延伸方向;然后根据梯田田坎排布特征,以山顶点为初始位置,分析计算以山顶点为顶点的各条射线经过处像元灰度值的变化规律;最后,基于梯田的不连续性、梯田的延伸方向、田坎宽度特征,利用信号处理方法设定相应的规则,去掉不存在梯田或者过于破碎的梯田像素,保留清晰连续的梯田像素并得到梯田的分布图。本发明基于图像处理和信号分析理论提出的射线法提取方法,能够较为有效的从高分辨率遥感影像上提取梯田范围,根据梯田田坎排列的具体方向、宽度、数量等,对梯田范围进行提取,具有普遍的适应性和创新性。
Description
技术领域
本发明涉及一种梯田提取方法,具体涉及一种基于射线法的梯田提取方法。
背景技术
黄土高原水土流失严重,人均耕地面积少。梯田是为了农业耕作或者进行水土保持而修建的人工地貌,把坡耕地修成水平梯田,可以有效的保水、保土、保肥,增强抗旱能力。梯田对自然地表的改造造成了土地利用类型和土壤侵蚀方式的改变,因此在地貌学、土地利用调查和水土保持等方面,梯田都具有重要研究意义。近年来关于梯田的水土保持作用和土地利用方式的研究在不断深入,但是在梯田提取方面的研究却略显不足。
现有技术中存在采用傅里叶变换和纹理分析的方法来提取梯田。傅里叶变换的方法首先利用傅里叶变换得到频率谱,继而通过对频率谱的方向和数值大小进行分析得到梯田范围。该方法很好地利用了梯田纹理重复变化的特征,但是这种单纯依据方向变量来提取梯田的方法会受到同方向其他地物的干扰。纹理分析方法原理简单并且具有普遍的适应性,但是这种方法得到的梯田位置不够精确,对梯田边缘的识别精度较低。另外,还曾出现了使用面向对象的影像分割提取技术对梯田等地物进行了分割分类,其分割参数的确定主要是通过统计分析得到的。该方法操作简单且具有创新性,但是其分类结果不具有普适性,容易形成破碎斑块并且缺乏对梯田整体形态特征的诠释。
发明内容
发明目的:本发明的目的在于针对现有技的问题,提出了一种基于射线法的黄土高原梯田提取方法。
技术方案:本发明提供了一种基于射线法的黄土高原梯田提取方法,包括以下步骤:
(1)山顶点的选取:若有相应区域的高分辨率DEM数据,采用流域边界线的邻域分析方法提取山顶点,实现山顶点的自动提取;若没有DEM数据,进行人工点选山顶点,相比于使用自动提取山顶点,人工选点具有更大的灵活性,对于位于坡面上,又离山顶点较远的梯田区域可以人为的补充点,从而得到更为完整的梯田。
(2)梯田的识别:如图1所示为梯田的遥感影像图,为了准确地识别梯田,首先总结了梯田的图像特征:1、遥感影像中的梯田在田坎处存在明显的灰度变化,这是梯田识别的关键;2、个别梯田由于人畜破坏或者年久失修,已经难于分辨或者存在缺口如图1(a);3、位于不同坡面上的梯田,其延伸方向不一致,如图1(b),因此在梯田识别过程中方向问题必须解决;4、每层梯田宽度不同,绝大多数小于10m,如图1(c),所以不能用单一的宽度来判断梯田是否存在。
使用射线法对每一条射线来说都是相同的,所以只介绍该方法对于一条射线的处理过程。本发明将山顶点对应到遥感影像中,从山顶点开始引出射线,并对梯田的规则进行提取,将射线经过处像元的灰度值按从山顶到山脚的顺序排列成一个数组,提取这列数组中的极大值和极小值,记录这些点在数组中的位置并根据梯田特征设定规则。
梯田通常位于沟沿线之上的正地形区域,围绕山顶或者在山坡上分布。一般来说梯田上的点与山顶点的连线方向和梯田田坎延伸方向基本垂直或者夹角较大,如果从山顶点出发向周围做若干条射线,根据射线覆盖区域像素点灰度值的变化特征,就能大致确定梯田位置。
(3)梯田范围可视化:在影像上标注梯田位置。
进一步,步骤(2)中的遥感影像进行均值滤波、中值滤波、高斯模糊或表面模糊进行预处理以突出梯田信息。
进一步,步骤(2)设定以下规则,如图4、5所示:
①忽略田坎数小于三阶的破碎梯田,同时也是为了避免与梯田延伸方向相同的其他地物的干扰。三阶梯田的侧视图如图2所示,遥感影像为如图3的俯视图,图中竖直的直线表示所画射线,俯视图中用圆点标出的是灰度变化较大处,这7个点处的灰度值也就是图像灰度值变化的极值点。然后判断极值点数目是否小于7,如果小于7,不再对这条线进行后面的操作,否则继续后面的处理。
②获得相邻极值点下标的差值,如果差值小于1个像元,那么就忽略该差值对应的一个极值点,该步骤可解决边缘检测形成双线的问题。
③将射线经过之处的灰度极值依次连接,形成灰度变化曲线。一般来说田坎间距不大于10m,实验图像分辨率为0.6m,所以将阈值设置为16,如果曲线上相邻极值点的间距大于16个像元,则将曲线在此处断开,以便进行分块处理。该方法主要用于检测不连续梯田。
④一块区域内,每个极值点与其相邻两个极值点灰度值之差不能小于阈值,阈值根据影响中梯田的清晰度确定。在提取梯田之前要对图像灰度值进行归一化处理,所以阈值的取值范围为0~1。
⑤梯田灰度变化的周期性是由于阳光对不同倾斜角度的地面进行照射造成的,在梯田的灰度变化曲线中,所有极大值点的灰度值的标准差不大于0.2,极小值点同理。
进一步,步骤(3)使用扇形面填充存在梯田的区域,图6(a)中的黑色实线表示从山顶点发出的射线,其上的粗线表示经过上述实验得到的存在梯田的地方。假设每两条相邻射线间的角度为a,那么一共存在360/a条射线,如果在第i条射线的第50到100个像素处有梯田,那么就将角度在 [360/a*i-a/2,360/a*i+a/2] 区间内,半径在50到100个像素之间的像元标记为梯田,如图6(b)中的灰色区域。
有益效果:本发明首先分析梯田田坎的排列方式以及延伸方向;然后根据梯田田坎排布特征,以山顶点为初始位置,分析计算以山顶点为顶点的各条射线经过处像元灰度值的变化规律;最后,基于梯田的不连续性、梯田的延伸方向、田坎宽度特征,利用信号处理方法设定相应的规则,去掉不存在梯田或者过于破碎的梯田像素,保留清晰连续的梯田像素并得到梯田的分布图。
本发明基于图像处理和信号分析理论提出的射线法提取方法,能够较为有效的从高分辨率遥感影像上提取梯田范围,该方法解决了傅里叶变换方法容易受到同方向地物干扰,以及面向对象分析方法不具有普适性和地学意义难于解释的问题,能够根据梯田田坎排列的具体方向、宽度、数量等,对梯田范围进行提取,具有普遍的适应性和创新性。
附图说明
图1为梯田遥感影像图;
图2为梯田侧视图;
图3为梯田俯视图;
图4为梯田特征规则流程图;
图5为梯田特征规则流程图;
图6为梯田可视化示意图;
图7为实施例的实验样区;
图8为边缘检测后图像;
图9为山顶点分布图;
图10为射线分布图;
图11为梯田分布示意图;
图12为人工矢量化梯田范围图。
具体实施方式
下面对本发明技术方案进行详细说明,但是本发明的保护范围不局限于所述实施例。
实施例:一种基于射线法的黄土高原梯田提取方法,本实施例的实验样区取位于陕西绥德韭园沟,经度:110°20′、纬度:37°35′附近,属黄土丘陵沟壑区第一副区,是黄河水利委员会水土保持重点实验区,该地区沟谷溯源侵蚀强烈,梯田分布广泛且保留较为完好。实验选取了一段连续且形状不规则的梯田作为实验数据,实验数据为0.6 m分辨率全色波段遥感影像,如图7所示,本次试验是提取图中图框内的梯田。
图像预处理是为了减少遥感影像中的噪声,突出梯田信息。首先对图像中的噪声(比如植被)进行去除,本实施例使用的方法是表面模糊,表面模糊原理与高斯模糊相似,但是增加了阈值变量,杂点像素灰度值与周围像素灰度值之间的差值就是阈值,该方法可以在保留边界的同时模糊图像。之后对执行过表面模糊的图像进行边缘检测,这一步操作是为了使梯田信息更加突出,以便于后续的梯田识别工作,边缘检测结果如图8。这里需要说明一点,由于边缘检测算法本身的原因,执行过边缘检测的图像可能在边缘处形成双线,这将对梯田田坎的识别造成干扰,所以在梯田识别过程中,步骤(2)规则设定了相应的算法来避免这一弊端。另外为了统一图像灰度范围,将处理后的图像灰度值做了归一化处理。
进行设置的参数有:1、从每个山顶点发出的射线数目。2、射线的长度。3、相邻极值点最小灰度差值的阈值。前两个参数主要根据梯田的面积大小和对梯田提取精度来选择。如果梯田面积较大,则需要更长的射线或者增加山顶点数目来覆盖所有梯田区域。射线数目越多,得到的梯田边界越精确,但是也需要考虑程序在执行过程中的时间成本和硬件条件。本实施例规定射线数目为120条,长度为150像素。相邻极值点最小灰度差值的阈值为0.25。
本次实验主要为了提取图7黑框内的梯田,因此选取这部分范围内山顶点,通过观察影像图可知,左边第一个山顶点的左下方有一块离山顶点超过150个像素的梯田,所以手动在这块区域内添加一个顶点,如图9所示。继而可以得到射线分布图,如图10。
通过按照本发明方法制定的规则提取和分析图10中每根射线覆盖处灰度值的变化特征,可以知道射线覆盖处是否存在梯田以及梯田位置。以山顶点为圆心,先绘制出山顶点到梯田终点距离为半径的扇形面,然后擦除以山顶点到梯田起点距离为半径做的扇形面,其中扇形的圆心角等于相邻两条射线间的角度。这样就得到了许多不连续的扇形面,不同顶点四周的扇形面相互组合、重叠,从而形成了图11中黑色区域所示的梯田范围。
通过观察图11中黑色区域的覆盖范围可知,射线法能够提取大部分的真实梯田区域,是较为有效的梯田提取方法。为了验证该区域的梯田提取精度,将人工勾绘的梯田范围与射线法提取的梯田进行比较。首先对原图中的梯田进行人工矢量化,如图12所示。矢量化结果与射线法获得的梯田位置较为吻合,人工提取的梯田像素数目为59120,射线法提取的梯田像素数目为49540,两者重合的像素数目为45001,重合像素数占人工提取梯田像素数目的比例为76.12%。
Claims (3)
1.一种基于射线法的黄土高原梯田提取方法,其特征在于:包括以下步骤:
(1)山顶点的选取:若有相应区域的高分辨率DEM数据,采用流域边界线的邻域分析方法提取山顶点;若没有DEM数据,进行人工点选山顶点;
(2)梯田的识别:将山顶点对应到遥感影像中,从山顶点开始引出射线,山顶点的个数及射线的数目、长度要求覆盖所有梯田区域;对梯田的规则进行提取,将射线经过处像元的灰度值按从山顶到山脚的顺序排列成一个数组,提取这列数组中的极大值和极小值,记录这些点在数组中的位置并根据梯田特征设定以下规则:
①忽略田坎数小于三阶的破碎梯田,判断各射线上的灰度变化极值点数目是否小于7,忽略小于7的射线,其他射线进行后面的操作;
②如果数组中极值点位置间隔小于1个像元,那么忽略其中一个极值点;
③将射线经过之处的灰度极值依次连接,形成灰度变化曲线,如果曲线上相邻极值点的间距大于16个像元,则将曲线在此处断开;
④一块区域内,极值点与相邻两个极值点灰度值之差不能小于阈值,阈值根据影像中梯田的清晰度确定;
⑤在梯田的灰度变化曲线中,所有极大值点的灰度值的标准差不大于0.2,极小值点同理;
(3)梯田范围可视化:在影像上标注梯田位置。
2.根据权利要求1所述的基于射线法的黄土高原梯田提取方法,其特征在于:步骤(2)中的遥感影像进行均值滤波、中值滤波、高斯模糊或表面模糊进行预处理以突出梯田信息。
3.根据权利要求1所述的基于射线法的黄土高原梯田提取方法,其特征在于:步骤(3)使用扇形面填充存在梯田的区域:假设每两条相邻射线间的角度为a,那么一共存在360/a条射线,如果在第i条射线的第50到100个像素处有梯田,那么就将角度在 [360/a*i-a/2,360/a*i+a/2] 区间内,半径在50到100个像素之间的像元标记为梯田。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710234686.4A CN107067384B (zh) | 2017-04-12 | 2017-04-12 | 一种基于射线法的黄土高原梯田提取方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710234686.4A CN107067384B (zh) | 2017-04-12 | 2017-04-12 | 一种基于射线法的黄土高原梯田提取方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN107067384A CN107067384A (zh) | 2017-08-18 |
CN107067384B true CN107067384B (zh) | 2019-08-16 |
Family
ID=59601993
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201710234686.4A Active CN107067384B (zh) | 2017-04-12 | 2017-04-12 | 一种基于射线法的黄土高原梯田提取方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN107067384B (zh) |
Families Citing this family (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107909625B (zh) * | 2017-11-15 | 2020-11-20 | 南京师范大学 | 一种基于等高线的山顶点提取方法 |
CN107958485B (zh) * | 2017-12-08 | 2020-12-08 | 南京师范大学 | 一种伪山顶点剔除方法 |
CN110211128B (zh) * | 2019-05-06 | 2023-02-17 | 南京泛在地理信息产业研究院有限公司 | 基于遥感影像和dem的黄土高原梯田提取方法 |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104452654A (zh) * | 2014-11-17 | 2015-03-25 | 中国水利水电科学研究院 | 一种生态理念下的灌溉梯田设计 |
CN106033611A (zh) * | 2016-05-23 | 2016-10-19 | 西安建筑科技大学 | 一种dem数据中的山顶点提取方法 |
CN106447510A (zh) * | 2016-08-30 | 2017-02-22 | 水利部水土保持监测中心 | 一种梯田信息提取方法和系统 |
-
2017
- 2017-04-12 CN CN201710234686.4A patent/CN107067384B/zh active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104452654A (zh) * | 2014-11-17 | 2015-03-25 | 中国水利水电科学研究院 | 一种生态理念下的灌溉梯田设计 |
CN106033611A (zh) * | 2016-05-23 | 2016-10-19 | 西安建筑科技大学 | 一种dem数据中的山顶点提取方法 |
CN106447510A (zh) * | 2016-08-30 | 2017-02-22 | 水利部水土保持监测中心 | 一种梯田信息提取方法和系统 |
Non-Patent Citations (2)
Title |
---|
基于DEM与遥感影像数据的梯田信息强化表达研究;程德强 等;《干旱区资源与环境》;20161115;第30卷(第11期);第124-130页 |
基于面向对象的遥感影像梯田信息提取研究;张雨果;《中国优秀硕士学位论文全文数据库-农业科技辑》;20161115(第11期);第D043-2页 |
Also Published As
Publication number | Publication date |
---|---|
CN107067384A (zh) | 2017-08-18 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
WO2020258899A1 (zh) | 一种基于高分遥感影像的山地景观建筑提取方法和设备 | |
Yun et al. | Individual tree crown segmentation from airborne LiDAR data using a novel Gaussian filter and energy function minimization-based approach | |
Peng et al. | Identifying the urban-rural fringe using wavelet transform and kernel density estimation: A case study in Beijing City, China | |
CN106683112B (zh) | 一种基于高分辨率图像的道路路域建筑物变化提取方法 | |
CN103824077B (zh) | 一种基于多源遥感数据的城市不透水层率息提取方法 | |
CN107247927B (zh) | 一种基于缨帽变换的遥感图像海岸线信息提取方法及系统 | |
CN107609526A (zh) | 基于规则的精细尺度城市不透水面快速提取方法 | |
CN112164062A (zh) | 一种基于遥感时序分析的抛荒地信息提取方法及装置 | |
CN103226832B (zh) | 基于光谱反射率变化分析的多光谱遥感影像变化检测方法 | |
CN107067384B (zh) | 一种基于射线法的黄土高原梯田提取方法 | |
CN105550685A (zh) | 基于视觉注意机制的大幅面遥感影像感兴趣区域提取方法 | |
CN104361582A (zh) | 一种对象级高分辨率sar影像洪水灾害变化检测方法 | |
CN113408460B (zh) | 一种基于遥感大数据和云平台探测互花米草分布的方法 | |
CN106156758B (zh) | 一种sar海岸图像中海岸线提取方法 | |
CN105139396B (zh) | 一种全自动遥感影像云雾检测方法 | |
CN104850853B (zh) | 城市快速提取方法及装置 | |
CN111091079A (zh) | 基于tls的高寒脆弱区植被优势单株结构参数测定方法 | |
CN110889840A (zh) | 面向地物目标的高分6号遥感卫星数据的有效性检测方法 | |
CN108764132A (zh) | 一种湖泊湿地遥感图像误差检测方法 | |
Zhang et al. | Methods for automatic identification and extraction of terraces from high spatial resolution satellite data (China-GF-1) | |
CN101788664A (zh) | 面向土地退化响应单元的高光谱土地退化信息提取方法 | |
CN102567726A (zh) | 极地冰盖边缘区域浮冰自动提取技术 | |
Uchida | Applicability of satellite remote sensing for mapping hazardous state of land degradation by soil erosion on agricultural areas | |
CN103279954A (zh) | 一种基于土地利用数据库的遥感影像变化检测方法 | |
Meng et al. | Monitoring human-induced surface water disturbance around Taihu Lake since 1984 by time series Landsat images |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |