本发明涉及地理信息技术应用领域，具体涉及一种基于矢量数据的树状水系自动识别方法。
背景技术
自然界中相同类型的水系，其所在区域往往有着相似的地质构造和自然环境。通过对水系的排列形式、平面形态等进行分析，可以大致推断出水系所在区域的地质构造和地壳运动的大致情况。从而，水系类型的快速、准确解译具有重要的研究意义。
常见的水系类型有树状水系、格状水系、羽状水系、平行状水系、放射状水系、向心状水系以及辫状水系等。目前，水系类型的判别方式，主要是人工判别，此种方法虽然简单易行，但效率低下，不利于大范围的判别。为提高处理效率，已有学者开展了水系类型的自动识别研究。如王一川使用二值化水系图像，基于最优统计图象分类器，形成了以水系纹理单元为单位进行水系分类的方法(参见王一川.水系自动分类研究[d].西南交通大学,2006.)；车国泉利用句法模式识别技术，对辫状水系的图像进行了自动识别(参见车国泉.辫状水系的自动识别[d].西南交通大学,2007.)。相关水系类型的自动识别方法，虽能够进行水系类型的一定识别，但由于是利用二值化水系图像，而不是直接利用水系矢量数据，导致算法的复杂度较高，识别效果较差。
一般发育于地形比较平坦的地区，地壳较稳定和水平岩层地区的树状水系，主要具有支流多而不规则、主流与支流之间呈锐角、排列呈树枝状等典型特征。为此，针对树状水系的上述特征，基于矢量水系数据，本文提出了一种基于形态特征的树状水系自动识别方法。
技术实现要素：
发明目的：本发明针对现有技术存在的问题，提供一种树状水系的识别方法，该方法基于矢量数据识别树状水系，算法复杂度较低，自动化程度高，识别准确度较好。
技术方案：本发明所述的树状水系的识别方法包括：
(1)获取shp格式待识别水系的线图层数据，并读入到数据集l＝{li|i＝0,1,2,…,m}中，其中，li表示第i个河流要素，其属性包含河流长度leni和标志属性flagi＝0，m为河流要素的个数；
(2)基于河流要素的端点坐标，计算各河流要素的方位角；
(3)根据各河流要素的方位角计算两两相交的河流要素间的夹角；
(4)若两两相交的河流要素间的夹角中锐角占比大于预设阈值，则判定待识别水系为树状水系，否则判别为非树状水系。
可选的，所述识别方法还可以包括：
(5)根据各河流要素的方位角计算河流主要走向频率，从而对树状水系类型细分。
进一步的，步骤(2)具体包括：
(2-1)获取各河流要素的首尾点坐标，写入坐标集合coord；
(2-2)基于各河流要素的首尾点坐标，根据公式(1)计算对应河流要素的方位角azimuthi:
式中，i∈[0,m-1]，(xis,yis),(xie,yie)分别为第i条河流要素的起点坐标和终点坐标，azimuthi的阈值范围为[0,360)。
进一步的，步骤(3)具体包括：
(3-1)从河流要素中提取出河流长度最长的作为主干流，并标记为n级河流，其中，n为正整数，且n>＝3；
(3-2)设置j＝n；
(3-3)检索与j级河流相交的所有河流要素，作为j-1级河流，并基于公式(2)，依次计算每一j-1级河流与j级河流之间的夹角，并存入数组angle[m]，同时设置该j-1级河流的flag＝1；
anglei＝|azimuthj-1-azimuthj|式(2)
式中，anglei表示第i个夹角，i∈[0,m-1]，azimuthj-1、azimuthj表示相交的j-1级河流与j级河流的方位角；
(3-4)将j＝j-1，循环执行步骤(3-3)，直至完成所有夹角的计算。
进一步的，步骤(4)具体包括：
计算夹角数组angle[m]中值属于(0,90)的元素个数number，如果number/m大于预设阈值0.9，则判定待识别水系为树状水系；否则，为非树状水系。
进一步的，步骤(5)具体包括：
(5-1)以10度为间隔划分方位角区间，并根据每一河流要素的方位角azimuthi，统计各方位角区间内的河流要素数量，统计结果存入列表dl＝{dlk|k＝0,1,..,35}，其中，k为方位角区间序号，其区间范围为[k*10,k*10+10)，dlk为相应区间中各河流要素的数量；
(5-2)获取列表dl中的元素的最大值dlmax＝max{dlk}，以及与dlmax邻近的四个元素，分别为dlmax-2、dlmax-1和dlmax+1、dlmax+2；
(5-3)根据以下公式3，计算河流要素主要走向占比mrate；
(5-4)若mrate≥0.5，则判定待识别水系为平行树状水系，否则，为一般树状水系。
有益效果：本发明与现有技术相比，其显著优点是：本发明基于矢量数据识别树状水系，算法复杂度较低，自动化程度高，识别准确度较好。
附图说明
图1为本发明方法的流程图；
图2为实施例中实验数据1；
图3为实施例中实验数据2。
具体实施方式
本实施例提供了一种树状水系的识别方法，如图1所示，包括以下步骤：
(1)获取shp格式待识别水系的线图层数据，并读入到数据集l＝{li|i＝0,1,2,…,m}中，其中，li表示第i个河流要素，其属性包含河流长度leni和标志属性flagi＝0，m为河流要素的个数。
例如，实验数据选取图2中shp格式的待识别水系线图层数据，可以得到数据集l＝{li|i＝0,1,2,…,32}。
(2)基于河流要素的端点坐标，计算各河流要素的方位角。
步骤(2)具体包括：
(2-1)获取各河流要素的首尾点坐标，写入坐标集合coord。
(2-2)基于各河流要素的首尾点坐标，根据公式(1)计算对应河流要素的方位角azimuthi:
式中，i∈[0,m-1]，(xis,yis),(xie,yie)分别为第i条河流要素的起点坐标和终点坐标，azimuthi的阈值范围为[0,360)。
对图2中进行计算后得到的河流首尾两个端点坐标coordi和河流方位角azimuthi具体如表1所示。
表1
(3)根据各河流要素的方位角计算两两相交的河流要素间的夹角。
步骤(3)具体包括：
(3-1)从河流要素中提取出河流长度最长的作为主干流，并标记为n级河流，其中，n为正整数，且n>＝3；本实施例中设置n＝3。
(3-2)设置j＝n；
(3-3)检索与j级河流相交的所有河流要素，作为j-1级河流，并基于公式(2)，依次计算每一j-1级河流与j级河流之间的夹角，并存入数组angle[m]，同时设置该j-1级河流的flag＝1；
anglei＝|azimuthj-1-azimuthj|式(2)
式中，anglei表示第i个夹角，i∈[0,m-1]，azimuthj-1、azimuthj表示相交的j-1级河流与j级河流的方位角；
(3-4)将j＝j-1，循环执行步骤(3-3)，直至完成所有夹角的计算。
接上例，计算得到的所有夹角如表2所示。
表2
(4)若两两相交的河流要素间的夹角中锐角占比大于预设阈值，则判定待识别水系为树状水系，否则判别为非树状水系。
具体为，计算夹角数组angle[m]中值属于(0,90)的元素个数number，如果number/m大于预设阈值0.9，则判定待识别水系为树状水系；否则，为非树状水系。上例中，总的夹角共有32个，其中锐角有30个，锐角占比为0.9375。number/m大于0.9，故为树状水系。
(5)根据各河流要素的方位角计算河流主要走向频率，从而对树状水系类型细分。
步骤(5)具体包括：
(5-1)以10度为间隔划分方位角区间，并根据每一河流要素的方位角azimuthi，统计各方位角区间内的河流要素数量，统计结果存入列表dl＝{dlk|k＝0,1,..,35}，其中，k为方位角区间序号，其区间范围为[k*10,k*10+10)，dlk为相应区间中各河流要素的数量。
接上例，其统计到的河流走向区间的值valuei和个数counti具体如表3：
表3
(5-2)获取列表dl中的元素的最大值dlmax＝max{dlk}，以及与dlmax邻近的四个元素，分别为dlmax-2、dlmax-1和dlmax+1、dlmax+2；
(5-3)根据以下公式3，计算河流要素主要走向占比mrate；
(5-4)若mrate≥0.5，则判定待识别水系为平行树状水系，否则，为一般树状水系。
接上例，可以计算得到mrate＝0.2727，故判定为一般树状水系。
图3是另外一组shp格式待识别水系的线图层数据，对其按照本实施例方法进行同样的处理，可以得到锐角占比为1，主要河流走向为mrate＝0.6786，故而其判别结果为平行树状水系。
以上所揭露的仅为本发明一种较佳实施例而已，不能以此来限定本发明之权利范围，因此依本发明权利要求所作的等同变化，仍属本发明所涵盖的范围。