坡体因子的融合方法及系统

1.本发明涉及数字地形分析技术领域，尤其涉及一种坡体因子的融合方法及系统。


背景技术：

2.基于栅格数字高程模型(digital elevation model，以下简称dem)地形因子的自动提取与计算是数字地形分析(digital terrainanalysis，以下简称dta)的重要研究内容，在水土保持、城乡建设、土木工程以及军事等多领域有着广阔的应用前景。地形因子间的相互组合及其结构在相当程度上映射这该样区地貌的发育机制。因此，挖掘不同地貌形态的地形因子组合特征和表达方式以及地形因子之间的关联与分布特征对于探索深层次的地貌特征与过程，以及揭示地貌发育演化规律与驱动机制具有重要意义。仅依靠单一或简单地形因子无法有效刻画地表复杂的变化，多因子综合评价研究对不用的研究对象从多角度出发，遴选初步的目标地形因子，基于多因子融合理论与分析方法更好的表达地形信息。
3.坡体是自然地形实体的基本组成单元，是构成地面的基本要素，地貌的变化实际上是坡体特征及其组合关系的变化，一个坡体完整的表达应涵盖坡度、坡长、坡宽三种坡体地形因子，坡体地形因子可以反映坡体形态及其氧化过程，在黄土地貌、土壤侵蚀、生态建设等方面都发挥着重要作用。常见单一坡体因子能够通过栅格dem直接或者间接计算得到，是进行地形分析及模拟的依据，当前坡体研究通常选用s因子、l因子、ls因子、坡形指数以及曲率等地形因子单独或组合表达来表示坡体形态。地貌学研究、水文侵蚀学以及地址灾害学等领域均把坡宽作为描述坡体形态的关键要素，目前，仅通过坡体综合表达因子ls作为地形计算的重要表达参数，广泛应用于描述坡面形态并进行坡面侵蚀预测，尚未将坡宽纳入坡体综合表达体系。


技术实现要素：

4.本发明的主要目的是提供一种坡体因子的融合方法及系统，旨在解决现有的坡体综合表达因子ls未将坡宽因子纳入综合表达体系、应用的局限性以及在地形隐含信息分析和挖掘上不够全面和精细的问题。
5.为实现上述目的，本发明提供一种坡体因子的融合方法及系统，其中，所述坡体因子的融合方法包括以下步骤：
6.s1、根据栅格dem数据获取坡体因子，所述坡体因子包括：l因子、s因子、ls因子以及w因子；
7.s2、将所述坡体因子极值归一化处理，得到归一化后的地形因子；
8.s3、将所述归一化后的地形因子进行多因子融合，得到多个坡体融合因子；
9.s4、对多个所述坡体融合因子进行综合分析、评定，筛选出目标坡体融合因子。
10.优选地，所述l因子为坡长因子；所述s因子为坡度因子；所述ls因子为坡度坡长因子；所述w因子为坡宽因子。
11.优选地，所述步骤s1中w因子获取的具体步骤为：
12.s111、对栅格dem数据进行填洼处理，得到无洼地栅格dem数据；
13.s112、根据无洼地栅格dem数据获取汇流累积量；
14.s113、采用河网密度二阶导数分析法确定所述汇流累积量的目标阈值；
15.s114、提取所述汇流累积量目标阈值下的栅格河网，利用斯特拉勒(strahler)方法分级河网并进行子流域的划分；
16.s115、将划分好的流域数据利用水蚀预报模型(water erosion prediction project，以下简称为wepp)/基于地理信息系统的wepp模型(geo-spatial interface for wepp，以下简称为geowepp)进行坡面概化，最终确定w因子。
17.优选地，所述目标阈值为1800。
18.优选地，所述步骤s1中l因子、s因子和ls因子获取的具体步骤为：
19.s121、根据无洼地栅格dem数据采用水文分析工具生成多流向栅格图层；
20.s122、对所述栅格图层进行循环迭代处理，获取坡长λ以及l因子；
21.s123、根据所述无洼地栅格dem数据提取坡度θ，并根据坡度分级函数计算出s因子；
22.s124、通过提取所述l因子和s因子的地形指标，得到ls因子。
23.优选地，所述对栅格dem数据进行填洼处理之前，还包括：
24.判断所述栅格dem数据是否存在洼地，若有，则对所述栅格dem数据进行洼地填充，得到无洼地栅格dem数据。
25.优选地，所述步骤s122具体为：
26.根据所述多流向栅格图层，定义局部高点作为坡长累计计算起点，连续提取径流结束点，进行循环迭代处理，通过坡长λ计算公式以及l因子计算公式获取坡长λ以及l因子。
27.优选地，所述步骤s2中坡体因子极值归一化处理的具体步骤为：
28.将步骤s1中获取的坡体因子采用极值归一法进行计算，得到归一化后地形因子的值；
29.将所述归一化后地形因子的值映射到0-1。
30.优选地，所述坡体因子和目标坡体融合因子均为栅格dem数值矩阵。
31.一种坡体因子的融合系统，包括处理器、存储器以及存储在所述存储器上并可在所述处理器上运行的坡体因子融合的应用程序，所述坡体因子融合的应用程序执行时实现所述的坡体因子融合方法的步骤。
32.本发明的上述技术方案中，该坡体因子的融合方法根据栅格dem数据获取坡体因子所述坡体因子包括：l因子、s因子、ls因子以及w因子；将所述坡体因子极值归一化处理，得到归一化后的地形因子；将所述归一化后的地形因子进行多因子融合，得到多个坡体融合因子；对多个所述坡体融合因子进行综合分析、评定，筛选出目标坡体融合因子。本发明在原坡体综合表达因子中融入坡宽，形成一种新的坡体融合因子，解决现有的坡体综合表达因子ls应用的局限性以及不够精细的问题，通过坡体融合因子可以实现地形隐含信息的分析和挖掘，有助于区域地形地貌的分析和规划。
33.在本发明中，采用河网密度二阶导数分析法确定汇流累积量的目标阈值，实现生成的河流网络与实际流域的河流网络相吻合，以及子流域高效精准的划分，减少误差。
附图说明
34.为了更清楚地说明本发明实施方式或现有技术中的技术方案，下面将对实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图示出的结构获得其他的附图。
35.图1为本发明实施例坡体因子的融合方法的流程示意图；
36.图2为本发明实施例坡体因子的融合系统的结构示意图；
37.图3为本发明实施例坡体因子的融合方法中汇流累积量阈值与河网密度二阶导数关系的示意图；
38.图4为本发明实施例坡体因子的融合方法中不同阈值下提取的河流网络的结构示意图；
39.图5为本发明实施例坡体因子的融合方法中w因子获取的流程示意图；
40.图6为本发明实施例坡体因子的融合方法中坡体融合因子与归一化后的地形因子关系的示意图；
41.图7为本发明实施例坡体因子的融合方法中实验区dem栅格图；
42.图8为本发明实施例坡体因子的融合方法中实验区dem渲染图；
43.图9-a、b、c、d分别本发明实施例坡体因子的融合方法中实验区坡体l、s、ls、w因子结构示意图；
44.图10-a、b、c、d均为本发明实施例坡体因子的融合方法中融合后的lsw因子示意图。
45.本发明目的的实现、功能特点及优点将结合实施方式，参照附图做进一步说明。
具体实施方式
46.下面将结合本发明实施方式中的附图，对本发明实施方式中的技术方案进行清楚、完整地描述，显然，所描述的实施方式仅仅是本发明的一部分实施方式，而不是全部的实施方式。基于本发明中的实施方式，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施方式，都属于本发明保护的范围。
47.本发明各个实施方式之间的技术方案可以相互结合，但是必须是以本领域普通技术人员能够实现为基础，当技术方案的结合出现相互矛盾或无法实现时应当认为这种技术方案的结合不存在，也不在本发明要求的保护范围之内。
48.实施例1：
49.参见图1，根据本发明的一方面，本发明提供一种坡体因子的融合方法及系统，其中，所述坡体因子的融合方法包括以下步骤：
50.s1、根据栅格dem数据获取坡体因子，所述坡体因子包括：l因子、s因子、ls因子以及w因子；
51.s2、将所述坡体因子极值归一化处理，得到归一后的化地形因子；
52.s3、将所述归一化后的地形因子进行多因子融合，得到多个坡体融合因子；
53.s4、对多个所述坡体融合因子进行综合分析、评定，筛选出目标坡体融合因子。
54.本实施例中根据栅格dem数据获取坡体单向地形因子，对坡体因子进行融合、分析
和处理，进而得到目标坡体融合因子，该方法解决了现有的坡体综合表达因子ls应用的局限性以及不够精细的问题，通过坡体融合因子可以实现地形隐含信息的分析和挖掘，有助于区域地形地貌的分析和规划。
55.具体地，在本实施例中，所述坡体因子包括：l因子、s因子、ls因子以及w因子，所述l因子为坡长因子；所述s因子为坡度因子；所述ls因子为坡度坡长因子；所述w因子为坡宽因子。
56.具体地，在本实施例中，所述步骤s1中l因子、s因子和ls因子获取的具体步骤为：
57.s121、将无洼地栅格dem数据采用水文分析工具生成多流向栅格图层；
58.s122、对所述栅格图层进行循环迭代处理，获取坡长λ以及l因子；
59.s123、根据所述无洼地栅格dem数据提取坡度θ，并根据坡度分级函数计算出s因子；
60.s124、通过提取所述l因子和s因子的地形指标，得到ls因子。
61.以高分辨率的栅格dem为数据源，对栅格dem数据进行填洼处理，生成无洼地的栅格dem数据，将所述无洼地栅格dem数据采用水文分析工具生成多流向水流方向的栅格图层；根据水流方向的栅格图层，定义局部高点作为坡长累计计算起点，连续提取径流结束点，进行循环迭代处理，利用坡长λ计算公式(1)以及l因子计算公式(2)获取坡长λ和l因子；根据生成的无洼地栅格dem数据，利用planar算法提取坡度，并根据s因子计算公式(3)利用坡度分级函数进行计算，获取s因子；通过提取l因子和s因子的地形指标，经过公式(4)计算得到一个分布式、栅格式图形格式的复合地形指标ls因子。
[0062][0063][0064][0065]
ls＝l*s
ꢀꢀꢀ
(4)
[0066]
具体地，在本实施例中，所述对栅格dem数据进行填洼处理之前，还包括：
[0067]
s110、判断所述栅格dem数据是否存在洼地，若有，则对所述栅格dem数据进行洼地填充，得到无洼地栅格dem数据。
[0068]
具体地，在本实施例中，如图6所示，所述步骤s1中w因子获取的具体步骤为：
[0069]
s111、对栅格dem数据进行填洼处理，得到无洼地栅格dem数据；
[0070]
s112、根据无洼地栅格dem数据获取汇流累积量；
[0071]
s113、采用河网密度二阶导数分析法确定所述汇流累积量的目标阈值；
[0072]
s114、提取所述汇流累积量目标阈值下的栅格河网，利用strahler方法分级河网并进行子流域的划分；
[0073]
s115、将划分好的流域数据利用wepp/geowepp模型进行坡面概化，最终确定w因子。
[0074]
具体地，在本实施例中，首先基于ls因子阶段生成的无洼地栅格dem数据获取汇流累积量，并采用河网密度二阶导数分析法确定汇流累积量的目标阈值，利用gis水文分析工具，分别设定600、800、1000、1200、1400、1600、1800、2000、2200、2400共10个汇流累积量阈值生成栅格河网，阈值的起始设置原则遵循相似大小研究区域目标阈值成果并辅以较大研究区等比例相似性设定，从汇流累积量图层中提取不同阈值下试验区流域河网，通过属性表计算其河流密度、河流密度下降率、河流长度、河流长度下降率、河源数量、河源数量下降率等水文信息，如表1所示：
[0075][0076]
表1流域河网特征随汇流累积量阈值的变化关系
[0077]
随着阈值的增大，河流长度、河网密度以及河源数量都呈下降趋势，但降低的趋势减缓，最终合为主干河流，趋于稳定。当阈值由600升至800时，河网密度、河流长度以及河源数量分别减少了18.29％、18.29％以及40.63％；当阈值由600升至1000时，河网密度、河流长度以及河源数量分别减少了7.53％、7.53％以及17.11％，随着阈值的增大，河网特征均处于缓慢下降，最终趋于平稳；当阈值在1600-1800时，下降速率突然增大后平稳下降；当阈值增大至2400时，众多支流消失，仅剩下主干河道，与实际情况不相符合。
[0078]
具体地，如图3所示，拟合函数为y＝4.7048x-0.285
，拟合度r2为0.991，y为河网密度(km/km2),x为阈值，对所述河网密度的幂函数关系求二阶导数得到y”＝1.723x-2.285
，拟合度为0.997，将阈值代入所述y”＝1.723x-2.285
中，得到汇流累积量阈值与与河网密度二阶导数关系，如图3所示，河网密度从阈值600处降低，最后趋于平缓，当阈值为1800时出现拐点，河网密度与河源数量变化趋于稳定，当阈值为1800时，河网密度区域平缓，其生成的河流网络与实际流域的河流网络相吻合，因此，选取阈值1800为目标阈值，在所述目标阈值下提取栅格河网，利用strahler方法分级河网并划分子流域，在流域划分的基础上，将划分好的流域数据导入wepp/geowepp模型，确定关键源区域csa和最小源通道长度mscl两个参数，模拟坡面流水过程概化坡面，通过手工输入属性统计坡宽因子，并通过gis软件实现坡宽可视化。
[0079]
具体地，在本实施例中，所述步骤s2中坡体因子极值归一化处理的步骤具体包括：
[0080]
将步骤s1中获取的坡体因子l、s、ls以及w因子值进行归一化处理，也即将各坡体因子采用极值归一法计算得到归一化后的地形因子的值，再将各地性因子映射到0-1；计算式为g＝(x-x
min
)/(x
max-x
min
)，最终得到s因子归一化后的地形因子的值gs、l因子归一化后
的地形因子的值gs、ls因子归一化后的地形因子的值g
ls
以及w因子归一化后的地形因子的值gw。
[0081]
具体地，在本实施例中，所述坡体因子为栅格dem数值矩阵，相互不存在特殊的逻辑关联，权重融合方法需考虑具体应用的关联性，但采用几何代数运算融合方法不仅可以完成两个及以上的栅格单元运算，同时，融合后的综合评价因子满足原始坡体形态信息量的条件，因此，对极值归一化后的地形因子采用l*s*w、ls*w、l+s+w以及ls+w四种方式进行坡体地形因子的融合，计算式分别如下：
[0082]
lsw＝gs*g
l
*gwꢀꢀ
(5)
[0083]
lsw＝g
ls
*gwꢀꢀ
(6)
[0084]
lsw＝gs+g
l
+gwꢀꢀ
(7)
[0085]
lsw＝g
ls
+gwꢀꢀ
(8)
[0086]
式中，gs、g
l
、g
ls
和gw分别为s因子、l因子、ls因子、w因子的归一化后的地形因子的值，lsw为坡体融合因子。
[0087]
具体地，在本实施例中，步骤s4中对多个所述坡体融合因子进行综合分析、评定，筛选出目标坡体融合因子，具体为，融合后的坡体因子仍然是栅格dem数值矩阵，所述目标坡体融合因子为栅格dem数值矩阵，设置足够大的统计样本，所有的统计数据都呈正态分布，采用方差矩阵值和标准差矩阵均值作为数理统计分析指标，表示原始因子与融合因子的偏离程度，所述原始因子为坡体因子，所述融合因子为坡体融合因子，其均值越大表示所述融合因子与原始因子存在较大偏差；采用信息熵表示每个栅格所提供信息的总和，可有效评价所述融合因子所携带信息量的大小，信息熵计算公式为h＝-∑pi(x)lnpi(x)，其中h表示信息熵，x为随机变量，p为随机变量出现的概率。
[0088]
具体地，假设所述l因子、s因子、ls因子以及w因子贡献度一致，归一化后各坡体因子数值范围为0-1，通过分析软件获得坡体融合因子与坡体因子关系如图6所示，i表示lsw＝gs+g
l
+gw和lsw＝g
ls
+gw所获得的融合因子的与坡体因子之间的关系，ii表示lsw＝g
ls
*gw所获得的融合因子的与坡体因子之间的关系，iii表示lsw＝gs*g
l
*gw所获得的融合因子的与坡体因子之间的关系，从图6中可以知道，通过lsw＝g
ls
*gw和lsw＝gs*g
l
*gw方式获得的坡体因子存在不同程度的压缩，不利于平衡各坡体因子对片坡体融合因子的贡献，而通过lsw＝gs+g
l
+gw和lsw＝g
ls
+gw方式获得的坡体融合因子不会出现压缩的情况；结合对采用4种方式融合提取的lsw因子进行采用分析软件进行统计分析，可以得到通过lsw＝gs+g
l
+gw方式获得的坡体融合因子所携带的信息量最大，即通过信息熵的计算与比较、方差和标准差矩阵均值以及分析坡体因子与坡体融合因子之间的关系以及其他基本统计量进行综合分析、评定，通过所述l+s+w方式可筛选出目标坡体融合因子。
[0089]
具体地，在本实施例中，以5m栅格dem为基础性数据源，栅格行列为400
×
400，选取黄土高原某小流域为典型进行研究分析，实验区dem及晕渲图如图7、8所示，实验区高程最大值为1062.6m，最小值为851.4m，高程差为211.3m，栅格dem分辨率为5m，共计160000个像元，实验区面积共4km2。根据栅格dem数据获取坡体因子如图9所述，其中a表示实验区l因子，b表示试验区s因子，c表示实验区ls因子，d表示试验区w因子，再对获得的所述坡体因子进行归一化处理以及多因子融合，根据l*s*w、ls*w、l+s+w以及ls+w四种融合方式获取坡体融合因子，如图10所述，其中a表示lsw＝gs+g
l
+gw融合方式获得的lsw的结果，b表示lsw＝g
ls
+gw融合方式获得的lsw的结果，c表示lsw＝gs*g
l
*gw融合方式获得的lsw的结果，d表示lsw＝g
ls
*gw融合方式获得的lsw的结果。对所述融合结果进行综合分析、评价，通过信息熵表示每个栅格所提供信息的总和，对采用4种融合方式提取的lsw因子进行统计分析，获得每个坡体栅格的像元统计值和坡体融合因子lsw的统计规律，具体数值如表2所示：
[0090][0091]
表2采用四种融合方式提取的lsw因子的结果统计
[0092]
通过lsw＝gs+g
l
+gw和lsw＝g
ls
+gw获得的坡体融合因子其极值均差比lsw＝gs*g
l
*gw、lsw＝g
ls
*gw方式获得的坡体融合因子的极值均差大，也即通过lsw＝gs+g
l
+gw和lsw＝g
ls
+gw获得的坡体融合因子数据分布更均匀，对于lsw＝gs*g
l
*gw和lsw＝g
ls
*gw方式均会在总体上造成坡体融合因子值不同程度的压缩，融合坡体因子统计值不利于平衡各因子对坡体融合因子的贡献；通过图10表示的栅格数据分布和融合结果图效果可知，采用lsw＝g
ls
+gw的方式，其中坡宽因素占据主导，导致融合结果与坡宽因子相差无异，其栅格数据的标准差相比其他方式获得的融合因子偏大，也即通过所述方式获得的坡体融合因子其数据波动性较大，不具有代表性。因此，通过lsw＝gs+g
l
+gw融合方式获得栅格数据符合正态分布，标准差较小，较为稳定，代表性较强，其方差矩阵均值越小说明融合后得到的坡体融合因子与原始地形因子的相关性越高，更能有效且综合的表征多项原始地形因子，反之，均值越大表明该融合因子与原始因子的多项特征有较大的偏差，无法均衡地表征多项原始地形因子，通过lsw＝gs+g
l
+gw融合方式获得的坡体融合因子其h值最大，为1.0003，表明所包含的信息量最为丰富。综上，通过lsw＝gs+g
l
+gw融合方式可获得目标坡体融合因子，为均衡融合坡体因子量级，坡体融合因子值等于坡体因子归一化和平均值lsw＝1/3(gs+g
l
+gw)。
[0093]
实施例2：
[0094]
如图2所示，图2是本发明实施例方案涉及的一种坡体因子的融合系统的结构示意图。
[0095]
如图2所示，该系统可以包括：处理器10，例如cpu，通信总线12，用户接口13，网络接口14，存储器15。其中，通信总线12用于实现这些组件之间的连接通信。用户接口13可以红外接收模块，用于接收用户通过遥控器触发的控制指令，可选的用户接口13还可以包括标准的有线接口、无线接口。网络接口14可选的可以包括标准的有线接口、无线接口(如wi-fi接口)。存储器15可以是高速ram存储器，也可以是稳定的存储器(non-volatile memory)，例如磁盘存储器。存储器15可选的还可以是独立于前述处理器11的存储装置。
[0096]
本领域技术人员可以理解，图2中示出的坡体因子的融合系统的结构并不构成对坡体因子的融合系统的限定，可以包括比图示更多或更少的部件，或者组合某些部件，或者不同的部件布置。
[0097]
本发明坡体因子的融合系统的具体实施例与上述坡体因子的融合方法各实施例基本相同，在此不作赘述。
[0098]
以上仅为本发明的优选实施例，并非因此限制本发明的专利范围，凡是在本发明的发明构思下，利用本发明说明书及附图内容所作的等效结构变换，或直接/间接运用在其
他相关的技术领域均包括在本发明的专利保护范围。