CN105825510B - 一种兴趣点与道路网的自动配准方法 - Google Patents

Abstract

本发明公开了一种兴趣点与道路网的自动配准方法，该方法充分利用兴趣点和道路网中道路的属性相关性来纠正不同来源兴趣点与道路网间的系统偏差；通过构建道路网骨架图和兴趣点骨架图，挖掘兴趣点与道路网的几何模式对应关系，将兴趣点与道路间的点与线的同名控制点识别问题转换为骨架图的匹配问题，基于骨架图的交叉点匹配关系进行不同来源兴趣点与道路网数据的位置偏差自动纠正。极大的提高了各空间数据的相对位置精度，减少了人工干预，实现了兴趣点与道路网的自动精准匹配。

Description

一种兴趣点与道路网的自动配准方法

技术领域

本发明涉及地理信息技术领域，特别涉及一种兴趣点与道路网的自动配准方法。

背景技术

随着传感器和互联网技术的迅猛发展，用户可通过智能终端和开放地图平台参与地理空间数据的发布和生产，这种由用户生成的空间数据称为众源地理空间数据(或志愿者地理信息)。与专业测绘数据相比，众源地理空间数据在几何表达、精度、尺度、数据质量等方面差异较大，导致不同来源空间数据间存在重复采集而信息共享困难等问题。在导航与移动位置服务领域，兴趣点(Point of Interest，POI)和道路网是主要的数据内容，不同来源的兴趣点与道路网数据直接叠加存在一定位置偏差，直接使用存在偏差的兴趣点和道路网数据导致错误的规划路径和用户行为分析结果。因此，不同采集方式获取的兴趣点与道路数据在投入出行应用前需要进行几何配准，以保证兴趣点与道路网数据间的相对位置和空间关系的准确和一致性。相关文献：Sester,M.,Arsanjani,J.J.,Klammer,R,Burghardt,D.,and Haunert,J.H.,2014.Integrating and generalising volunteeredgeographic information[A].In:Burghardt,D.,Duchêne,C.,and Mackaness,W.(eds.).Abstracting geographic information in a data rich world[M].Switzerland:Springer International Publishing,119-155。传统的地图匹配方法利用GPS轨迹点的位置坐标与道路网的几何特征等将轨迹点匹配到对应路段，但是兴趣点并不落在道路区域，地图匹配方法无法建立兴趣点与道路网的匹配关系。基于路网匹配的配准方法利用道路网的匹配关系，将其中一份道路网的附属兴趣点配准到另一份道路网中，但该方法需要参考路网作为中介，而实际应用中很难获取兴趣点的参考路网数据。相关文献：Chen,B.,Yuan,H.,Li,Q.,Lam,W.H.K.,Shaw,S.,and Yan,K.,2014.Map-matching algorithm for large-scale low-frequency floating car data[J].International Journal ofGeographical Information Science,28(1):22-38；Zhang,M.,2009.Methods andimplementations of road-network matching[D]:Germany:Technical University ofMunich.现有方法无法解决不同来源兴趣点与道路网间同名控制点的自动识别问题。

发明内容

本发明利用兴趣点与道路网的几何模式相关性，提供了一种兴趣点与道路网的自动配准方法，有效解决不同来源兴趣点与道路网的位置偏差纠正问题。该方法综合考虑兴趣点与道路网数据的属性内容和几何位置，识别和提取道路链模式和兴趣点聚类模式，自动识别两者的同名控制点对，进行兴趣点与道路网的几何配准。

尽管兴趣点与道路网之间存在不可预测的非刚性偏差，但不同来源的兴趣点与关联的道路之间具有一定的几何模式相关性，这种几何模式相关性对兴趣点与道路网的自动配准具有重要的借鉴意义。

一种兴趣点与道路网的自动配准方法，该方法包括以下步骤，

步骤1：提取兴趣点所在的道路信息，与道路网进行第一次配准；

提取兴趣点地址属性中包含的道路信息，与道路网中道路信息进行判断，确认道路两边关联的兴趣点，消除兴趣点与道路网的系统误差；

比较兴趣点地址信息与道路名称属性的公共字符串，将兴趣点与道路划分至两两对应的集合；

统计所有道路称，遍历每个道路名rn，将“地址属性包含rn”的兴趣点p_i与“名称属性等于rn”的道路r_j分别划分至相互对应的兴趣点集合set_poi和道路集合set_road。通过道路名rn，set_poi与set_road一一对应；

根据最小二乘算法将各个兴趣点集合set_poi拟合为直线段line_poi。由于兴趣点集合set_poi与拟合直线段line_poi一一对应，而根据步骤1.1，兴趣点集合set_poi又与道路集合set_road一一对应，则可推出line_poi与set_road也是一一对应，即拟合直线段line_poi与道路集合set_road也通过道路名称相互关联起来；

识别相互关联的线性特征line_poi与道路集合set_road的两两交叉点，以此作为控制点集，建立仿射变换模型，求解变换参数，将兴趣点粗配准至目标道路网，消除系统偏差；

步骤2：识别道路网中道路链，构建道路骨架图；

道路网中方向连续道路相互连接构成道路链(stroke)，按道路连接夹角之和最大生成单线道路链，通过构建缓冲区识别和简化平行道路链，最后以单线道路链和平行道路链简化结果为边，构建道路骨架图；

道路链是由方向连续的道路连接后构成；

步骤3：将兴趣点与道路骨架图中的道路链进行匹配；

利用基于划分的空间聚类方法，根据兴趣点到邻近道路链的投影矢量，将兴趣点划分至不同的道路链或空道路链，生成与道路链对应的兴趣点聚类；

步骤4：拟合兴趣点聚类的主曲线，构建兴趣点骨架图；

对步骤3得到的任一兴趣点聚类，以最小二乘拟合直线段为初始主曲线，通过迭代添加顶点和优化顶点位置两个过程确定各兴趣点聚类的拟合主曲线，最后以拟合主曲线为边构建兴趣点骨架图；

步骤5：利用兴趣点与道路链的匹配关系，寻找兴趣点骨架图与道路网骨架图的匹配控制点，进行兴趣点与道路的精配准。

所述步骤2的具体过程如下：

步骤2.1，遍历所有道路结点，选择与道路结点连接的道路夹角之和最大的两两道路连接组合，对连接的道路进行道路链编号；

步骤2.2，将道路链编号相同的道路进行首尾连接，得到初始单线道路链；

步骤2.3，构建单线道路链的40米缓冲区，计算与缓冲区内邻近道路链的角度差异和长度比值，将角度差异和长度比值在给定阈值范围内的道路链识别为平行道路链，对各组平行道路链构建40米缓冲多边形，提取缓冲多边形的中心线，简化原有平行道路链；

步骤2.4，在简化后的中心线道路链与其他邻近道路链在相交处打断，相交处的结点添加到道路结点集合；遍历与简化前平行道路链相连接的道路链，判断与简化后中心线道路链是否连接，若连接，则删除悬挂边，若不连接，则延长中心线道路链使之连接。

所述步骤3的具体过程如下：

步骤3.1，计算兴趣点与道路链的投影矢量x_ij，x_ij为兴趣点p_i到道路链S_j的投影点p_j’与p_j构成的矢量，初始化k(k为道路链条数)个聚类及其中心c_j，聚类参数λ₀以及分配矩阵u_ij；

k个聚类与步骤2提取的k条道路链对应；

步骤3.2，将每个兴趣点划分至最邻近道路链，兴趣点与道路链的邻近程度由两者的投影矢量与道路链对应聚类中心的欧氏距离计算得到：|x_ij-c_j|，若该值最小，则u_ij＝1；

步骤3.3，按照下式重新计算各道路链对应的聚类中心c_j和聚类参数λ：

其中，u_ij为分配函数，u_ij＝1表示兴趣点p_i分配到道路链S_j，u_ij＝0表示p_i不分配到S_j，m为兴趣点个数；

步骤3.4，重新将各兴趣点划分至最邻近的道路链或空道路链；

对兴趣点p_i，若其与最邻近道路链S_j的聚类中心的欧氏距离|x_ij-c_j|小于聚类参数λ，则该兴趣点将被划分至道路链S_j，即u_ij＝1，否则被划分至空道路链，即u_ij＝0，j＝1…k；

步骤3.5，重复步骤3.3和3.4，重新计算聚类中心和重新分配兴趣点，直到两次迭代间聚类中心的变化量小于给定阈值。迭代后，根据分配矩阵u_ij将所有兴趣点划分至k条道路链或空道路链。

所述步骤4的具体过程如下：

步骤4.1：通过步骤3的分配函数u_ij得到与k条道路链对应的k个兴趣点聚类，逐一遍历各兴趣点聚类PS＝{p₁…p_q}，q为该聚类中兴趣点个数，通过最小二乘拟合将聚类中兴趣点拟合成直线段，以拟合直线段为该聚类的主曲线的初始值f，f＝{v₁…v_j…v_vn}，v_j为f的顶点，vn为顶点个数，此时vn＝2；

步骤4.2：将聚类PS中所有兴趣点投影到f的顶点v_j或边v_jv_j+1上，选择投影点个数最多的边，计算该边的中点，并将计算得到的中点添加到主曲线f中，此时vn＝vn+1；

步骤4.3：利用梯度下降法，以主曲线中每个顶点v_j(j＝1...vn)的局部惩罚误差函数最小为目标逐一优化主曲线f的各顶点位置；

步骤4.4：重复步骤4.2和4.3，迭代添加新的顶点和优化主曲线顶点位置，直到顶点个数大于给定阈值，返回4.1拟合下一兴趣点聚类，直到拟合完所有的兴趣点聚类，获得兴趣点骨架图skeleton_poi＝{f_t|j＝1…k}，f_t为第t个兴趣点聚类拟合得到的主曲线；

所述步骤4.3中主曲线的顶点v_j(j＝1...vn)的局部惩罚误差函数如下：

①Δ_q(v_j)为投影到顶点v_j或其邻边v_j-1v_j,v_jv_j+1的兴趣点到主曲线f的投影距离平方和，计算公式为：

其中，V_j和S_j-1、S_j-1分别为投影到顶点v_j和邻边v_j-1v_j，v_jv_j+1上的兴趣点集合，d(p_i,f)为兴趣点p_i到f的投影距离；

②Penalty(v_j)为主曲线f本身的形状惩罚因子，计算公式为：

其中，u(v_j)＝[d(v_j,v_j-1)]²,π(v_j)＝r²(1+cosγ_j)，γ_j为主曲线f在顶点v_j处的夹角，d(v_j,v_j-1)为顶点v_j，v_j-1的距离平方，r为聚类PS中所有兴趣点到中心点的最大距离；

③η为惩罚系数，计算公式为：

其中，Δ_q(f)为聚类PS中所有兴趣点到主曲线f的投影距离平方的均值，η₀为初始惩罚系数，设为0.4。

所述步骤5的具体过程如下：

步骤5.1，根据步骤3得到的兴趣点与道路链对应关系，确定道路和兴趣点骨架图的对应交叉点；

逐一遍历道路骨架图中所有道路链，寻找与之相交的道路链及其交点p＝(x_p,y_p)，确定与该道路链及相交道路链对应的兴趣点骨架图中的主曲线，计算这些主曲线的交点p’＝(x_p’,y_p’)，构成一对对应交叉点，其中p’为源控制点，p为目标控制点，相应位置偏差为(dx_p,dy_p),dx_p＝x_p-x_p’,dy_p＝y_p-y_p’，继续寻找下一条相交道路链，直到所有道路链均被遍历；

步骤5.2，以道路骨架图与兴趣点骨架图的对应交叉点作为同名控制点对，对同名控制点对中的源控制点构建三角网，每个兴趣点将落在三角网的某个三角形中，则该兴趣点的位置偏差可通过三角形三个顶点所对应的位置偏差加权得到，权重由该兴趣点到三角形三个顶点的距离确定，根据各兴趣点的位置偏差将兴趣点纠正到原始目标道路网中，完成精配准。

有益效果

本发明提供了一种兴趣点与道路网的自动配准方法，该方法充分利用兴趣点和道路网中道路的属性相关性来纠正不同来源兴趣点与道路网间的系统偏差；通过构建道路网骨架图和兴趣点骨架图，挖掘兴趣点与道路网的几何模式对应关系，将兴趣点与道路间的点与线的同名控制点识别问题转换为骨架图的匹配问题，基于骨架图的交叉点匹配关系进行不同来源兴趣点与道路网数据的位置偏差自动纠正。极大的提高了各空间数据的相对位置精度，减少了人工干预，实现了兴趣点与道路网的自动精准匹配。

附图说明

图1为本发明实施例的流程示意图；

图2为本发明实施例的粗配准过程示意图；

图3为本发明实施例的单线道路链连接规则示意图；

图4为本发明实施例的平行道路链识别示意图；

图5为本发明实施例的平行道路链简化示意图；

图6为本发明实施例的主曲线拟合示意图；

图7为本发明实施例的匹配控制点确定示意图；

图8为本发明实施例的兴趣点位置偏差计算示意图。

具体实施方式

下面将结合附图和实施例对本发明做进一步的说明。

图1为本发明实施例的兴趣点与道路网自动配准流程，该发明实施例首先通过比较两份兴趣点与道路网数据的属性信息，将兴趣点粗配准至目标道路网；然后提取道路网中的道路链模式，构建道路骨架图(skeleton_road)；之后通过基于划分的空间聚类方法建立兴趣点与道路链的匹配关系，生成与道路链对应的兴趣聚类；在此基础上，提取兴趣点聚类的主曲线特征，构建兴趣点骨架图(skeleton_poi)；最后基于兴趣点聚类与道路链的匹配关系，确定兴趣点骨架图与道路骨架图的匹配控制点，建立几何变换模型将兴趣点准确配准至目标路网中。

以下分步骤详细说明实现过程。假设兴趣点集合表达为P＝{p_i|p_i＝(v_i,addr_i,a_i1…a_is),i＝1…m}，道路集合表达为R＝{r_j|r_j＝(pline_j,name_j,b_j1…b_jt),j＝1…n}，其中点v_i和折线pline_j为兴趣点与道路的几何特征，道路间通过结点集N＝{n_h|r_i∩r_j＝n_h,r_i∈R,r_j∈R,h＝1…g}连通，addr_i为兴趣点p_i的地址属性，name_j为道路r_j的名称属性，a_i1与b_j1等为除地址、道路名称外的其他相关属性，如兴趣点类别，道路宽度等。

本发明可以通过计算机软件编程实现自动配准流程，本发明的具体配准过程如下：

步骤1，提取兴趣点所在的道路信息，与道路网进行第一次配准；

步骤1.1，统计道路数据的所有道路名称，遍历每一个道路名rn，按包含的道路名称rn将兴趣点与道路划分为一一对应的不同集合set_poi＝{p_i|i＝1…m,addr_i∩rn＝rn}和set_road＝{r_j|j＝1…n,name_j＝rn}，其中set_poi与set_road通过道路名称rn一一对应，即rn:如图2中，黑色和白色的圆点分别为两个兴趣点集合set_poi1和set_poi2。通过道路名称rn1和rn2，它们分别与道路集合set_road1和set_road2一一对应。

步骤1.2，通过最小二乘线性拟合将集合set_poi中的兴趣点拟合为直线段line_poi，则直线段line_poi与兴趣点集合set_poi对应的道路集合set_road也一一对应，即步骤1.1得到①：且步骤1.2直线段拟合得到②：则可推出③： 如图2中，两条黑色和灰色的虚线分别为兴趣点集合set_poi1和set_poi2拟合的线段line_poi1和line_poi2，它们与兴趣点集合对应的道路集合set_road1和set_road2也一一对应，且具有相似的几何形状。

步骤1.3，分别计算拟合的直线段及所对应道路的对应交叉点对，如图2中，计算两条相交的直线段line_poi1和line_poi2的交点p，同时计算其对应的道路集合set_road1和set_road2的交叉点，若存在多个道路交叉点，取其中心点，如图2中点p’，则对应交叉点对(p,p’)构成一对控制点，同理寻找其他控制点对，建立仿射变换模型，求解变换参数，消除兴趣点与道路数据的系统误差。

步骤2，提取道路链模式，根据连接道路间夹角提取单线道路链模式，通过中心线生成算法将平行道路链简化为单线道路链，得到k条道路链组成的道路骨架图skeleton_road；

步骤2.1，按与结点连接的道路条数从大到小逐一遍历道路结点，计算与之连接的两两道路夹角，选择夹角大于给定阈值且两两连接夹角之和最大的连接组合，按该组合进行道路链编号标识。

步骤2.1.1，所有道路的道路链编号属性stroke id初始化为-1，设current_id＝1。

步骤2.1.2，按与结点连接的道路条数从大到小遍历道路结点，计算与之连接的两两道路夹角，选择道路夹角大于120°且两两夹角总和最大的连接组合。如遍历图3中道路结点n_q，计算两两道路连接夹角，如r₁与r₄的夹角α，选择两两夹角之和最大的连接组合，即(r₁,r₄),(r₂,r₅),(r₃,r₆)。

步骤2.1.3，对最大连接组合中任一对连接道路(a,b),按以下规则进行道路链编号标识：

若a与b的stroke id皆为-1，则a与b的stroke id皆赋为current_id，current_id自加1；

若a与b的stroke id皆不为-1，分别为sa和sb，不妨设s_max＝max(sa,sb)，s_min＝min(sa,sb)，则将所有stroke id为s_min的道路修改为s_max；

若a与b的stroke id有且仅有一个为-1，不妨设a的stroke id为-1，b的stroke id为sb，则将a的stroke id赋为sb。

步骤2.2，将步骤2.1生成的具有相同编号的道路进行首尾连接，生成初始单线道路链。

步骤2.2.1，逐一遍历所有道路，如r，判断其是否已经生成道路链，若生成则继续遍历下一道路，否则，初始化队列Q，添加r到队列Q，以r为初始道路链S，转步骤2.2.2。

步骤2.2.2，弹出队列Q的当前元素q，寻找与q连接的邻接道路中具有相同编号的未连接道路next_r，连接S与next_r得到新的道路链S，添加next_r到队列Q中。如此不断遍历队列Q中的元素，直到Q中元素为空，返回步骤2.2.1遍历下一条道路。

按步骤2.1和2.2提取的道路链中存在平行多车道形式道路链模式，执行步骤2.3对平行道路链进行识别和简化。

步骤2.3，构建单线道路链的40米缓冲区，计算与缓冲区内邻近道路链的角度差异和长度比值，将角度差异和长度比值在给定阈值范围内的道路链识别为平行道路链，对各组平行道路链构建40米缓冲多边形，提取缓冲多边形的中心线，简化原有平行道路链；

步骤2.3.1，对步骤2.1和2.2生成的单线道路链，按道路链长度从大到小逐一遍历，如S_i。以40米为半径建立缓冲区，搜索S_i邻近的道路链，如P_j。如图4所示，若P_j落在S_i缓冲区内的长度l与S_i，P_j本身长度l₁，l₂满足①l/min{l₁,l₂}>0.8；②l/max{l₁,l₂}>0.2；③S_i与P_j的角度差异θ<30°，则P_j为S_i的平行道路链。按此规则选取与S_i平行的道路链。

步骤2.3.2，构建与S_i平行的道路链的40米缓冲区多边形，对该缓冲多边形的所有顶点构建Delaunay三角网，根据邻近三角形个数，将所有三角形分为三种类型。如图5所示，I型三角形仅与一个三角形相邻，具有1条公共边；II型三角形与两个三角形相邻，具有2条公共边；III型三角形与三个三角形相邻，具有三条公共边。

步骤2.3.3，对步骤2.3.2得到的三种类型的三角形按以下规则进行连接，得到该缓冲区多边形的中心线。如图5所示，对I型三角形，连接唯一的公共边中点与公共边所对的三角形顶点；对II型三角形，连接两条公共边的中点；对III型三角形，分别连接三角形外心与三条公共边的中点。按以上规则生成的多边形中心线可能包含不止一条路径，本发明选取路径最长的作为平行道路链的中心线，如图5中粗黑线。

步骤2.4，在简化后的中心线道路链与邻近道路链在相交处打断，相交处的结点添加到道路结点集合；遍历与简化前平行道路链连接的道路链，判断与简化后中心线道路链是否连接，若连接，则删除悬挂边，若不连接，则延长中心线道路链使之连接。

步骤2.4.1，对简化后的中心线道路链与邻近道路链在相交处进行一一打断，相交处的交叉点添加到道路结点集合。

步骤2.4.2，遍历与简化前平行道路链连接的道路链，若该道路链与简化后的中心线道路链连接，则去除悬挂边；若与简化后的中心线道路链不连接，则延长末节点至中心线道路链。

步骤3，兴趣点与道路链匹配，利用基于划分的空间聚类将兴趣点划分至不同聚类，生成k个兴趣点聚类与k条道路链对应，剩余兴趣点与空道路链对应；

步骤3.1，计算兴趣点与道路链的投影矢量，初始化k个聚类、聚类中心和聚类参数。首先计算兴趣点与道路链的投影矢量x_ij，设p_i’为兴趣点p_i到100米范围内的邻近道路链S_j的投影点，则投影矢量x_ij即为投影点p_i’与兴趣点p_i连成的矢量，然后初始化k个聚类中心为c_j＝(0,0)，j＝1…k，设定聚类参数λ₀＝1.5，初始化分配矩阵u_ij＝0,i＝1…m,j＝1…k，m为兴趣点个数，k为步骤2生成的道路链条数。

步骤3.2，按公式(1)将所有兴趣点划分至最邻近的道路链，使投影矢量与聚类中心的欧氏距离最小，分配矩阵中u_ij＝1表示将第i个兴趣点划分至第j条道路链，否则u_ij＝0。

步骤3.3，按公式(2)重新计算聚类中心和聚类参数。

步骤3.4，将所有兴趣点重新分配至最邻近的道路链或空道路链。对某一兴趣点p_i与最邻近道路链S_j，若满足|x_ij-c_j|小于聚类参数λ，则该兴趣点将被重新划分至道路链S_j，否则，该兴趣点将重新匹配至空道路链。

步骤3.5，如此重复执行步骤3.3和3.4，直到两次迭代期间所有聚类中心的变化量|c^(r+1)-c^(r)|皆小于给定阈值。迭代结束后，根据分配矩阵u_ij可确定每条道路链对应的兴趣点聚类为{p_i|u_ij＝1,i＝1…m}，j＝1…k。

步骤4，主曲线拟合，将步骤3得到的与道路链对应的兴趣点聚类中的兴趣点分别拟合为主曲线，这样每个兴趣点聚类拟合的主曲线也与道路链一一对应，并构成兴趣点骨架图skeleton_poi；

步骤4.1，初始化主曲线，逐一遍历步骤3得到的k个兴趣点聚类，不妨设遍历某一兴趣点聚类为PS＝{p₁…p_q}，该聚类包含q个兴趣点，利用最小二乘算法将PS中兴趣点拟合为直线段，以拟合直线段为兴趣点聚类PS的主曲线的初始值，f＝{v₁…v_j…v_vn}，v_j为f的顶点，vn为顶点个数，此时vn＝2。

步骤4.2，添加新的顶点到主曲线f，将所有兴趣点投影到主曲线f的顶点或边上，选择投影点个数最多的线段，计算其中点，添加到f中，此时vn＝vn+1。

步骤4.2.1，按公式(3)将PS中的所有兴趣点投影至f的各边s_j＝v_jv_j+1和顶点v_j，其中d(p_i,v_j)表示兴趣点p_i到顶点v_j的欧氏距离,d(p_i,s_j)或d(p_i,f)分别表示该点到边s_j或主曲线f的投影距离，V_j表示投影到主曲线f的顶点v_j上的兴趣点集合，S_j表示投影到主曲线f的边s_j上的兴趣点集合。

若兴趣点p_i到f的投影距离等于该点到主曲线f的顶点v_j的欧氏距离，则兴趣点p_i添加到集合V_j，如图6中兴趣点p_b。

若兴趣点p_i到f的投影距离等于该点到主曲线f的边s_j的投影距离，则兴趣点p_i添加到集合S_j，如图6中兴趣点p_a。

步骤4.2.2，对主曲线f中的vn-1条边，选择投影的兴趣点个数最多(即|S_j|最大)的边s_j＝v_jv_j+1，计算s_j的中点v，将v添加到f中，此时vn＝vn+1。

步骤4.3，优化主曲线f的各顶点v_j(j＝1…vn)，使局部惩罚误差最小。

步骤4.3.1，逐一遍历主曲线f的各顶点，如v_j。按公式(4)计算顶点v_j处的局部惩罚误差目标函数。

其中，①d(p_i,f)为兴趣点p_i到主曲线f的投影距离，S_j-1,S_j和V_j分别为投影到主曲线f的边v_j-1v_j,v_jv_j+1和顶点v_j的兴趣点集合(见步骤4.2.1)。

②Penalty(v_j)为f本身的形状惩罚因子，按公式(5)计算得到，其中，u(v_j)＝[d(v_j,v_j-1)]²,π(v_j)＝r²(1+cosγ_j),γ_j为f在顶点v_j处的夹角，d(v_i,v_j-1)为主曲线f中顶点v_i和v_j-1间距离。

③η为惩罚系数，Δ_q(f)为PS中各点到f的投影距离平方的均值，η₀设为0.4。

步骤4.3.2，采用梯度下降法迭代寻求各顶点的局部惩罚误差函数G_q(v_j)最小来优化该顶点v_j的位置。返回步骤4.3.1继续遍历和优化下一顶点v_j+1，直到所有顶点皆已优化。

步骤4.4，重复步骤4.2和4.3不断添加新的顶点和优化各顶点位置，直到主曲线顶点个数大于给定阈值(本发明设为30)。返回步骤4.1遍历下一兴趣点聚类，直到k个兴趣聚类皆被拟合为主曲线，构成兴趣点骨架图sekleton_poi＝{f_t|j＝1…k}，f_t为第t个兴趣点聚类拟合得到的主曲线。

步骤5，根据步骤3建立的兴趣点与道路链的匹配关系，提取道路骨架图skeleton_road与兴趣点骨架图sekleton_poi的对应交叉点，进行兴趣点与道路网的精配准。

步骤5.1，根据步骤3得到的兴趣点与道路链对应关系，确定道路和兴趣点骨架图的对应交叉点。具体地，逐一遍历所有道路链，寻找与之相交的道路链及其交点p＝(x_p,y_p)，确定与该道路链及相交道路链对应的兴趣点主曲线，计算这些主曲线的交点p’＝(x_p’,y_p’)，构成一对对应交叉点，计算该对应交叉点的位置偏差(dx_p,dy_p),dx_p＝x_p-x_p’,dy_p＝y_p-y_p’。如图7所示，道路链S₁与对应的兴趣点聚类所拟合的主曲线f₁也对应，同样地，S₂与f₂，S₃与f₃也一一对应。我们寻找S₁与S₃的交点p，f₁与f₃的交点p’，两者构成一对对应交叉点(p,p’)，同样可以确定S₂与S₃的交点q，f₂与f₃的交点q’，构成另一对对应交叉点(q,q’)。按以上规则依次继续寻找兴趣点骨架图与道路骨架图的对应交叉点。

步骤5.2，以道路骨架图与兴趣点骨架图的对应交叉点作为同名控制点对(如(p,p’)中p为目标控制点，p’为源控制点)，对所有源控制点构建控制三角网，每个兴趣点将落在三角网的某个三角形中，则该兴趣点的位置偏差可通过三角形三个顶点所对应的位置偏差加权得到，根据各兴趣点的位置偏差将兴趣点纠正到原始目标道路网中。

为便于理解本发明实施过程，下面举例说明几何精配准过程。如图8中，兴趣点p落在控制三角网的某个三角形，其位置偏差(dx,dy)通过三角形三个顶点对应的目标控制点与源控制点的位置偏差加权得到，本发明可采用面积加权和反距离加权两种方法，以此位置偏差进一步将兴趣点纠正到相应位置。

以上所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员在不脱离本发明实施例所述原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (5)

1.一种兴趣点与道路网的自动配准方法，其特征在于：该方法包括以下步骤，

步骤1：提取兴趣点所在的道路信息，与道路网进行第一次配准；

提取兴趣点地址属性中包含的道路信息，与道路网中道路信息进行判断，确认道路两边关联的兴趣点，消除兴趣点与道路网的系统误差；

步骤2：识别道路网中道路链，构建道路骨架图；

道路网中方向连续的道路相互连接构成道路链，按道路连接夹角之和最大生成单线道路链，通过构建缓冲区识别和简化平行道路链，最后以单线道路链和平行道路链简化结果为边，构建道路骨架图；

步骤3：将兴趣点与道路骨架图中的道路链进行匹配；

利用基于划分的空间聚类方法，根据兴趣点到邻近道路链的投影矢量，将兴趣点划分至不同的道路链或空道路链，建立兴趣点与道路链对应关系；

步骤4：拟合兴趣点聚类的主曲线，构建兴趣点骨架图；

对任一兴趣点聚类，以最小二乘拟合直线段为初始主曲线，通过迭代添加顶点和优化顶点位置两个过程确定任一兴趣点聚类的拟合主曲线，最后以拟合主曲线为边，构建兴趣点骨架图；

步骤5：利用步骤3确定的兴趣点与道路链对应关系，寻找兴趣点骨架图与道路网骨架图的匹配控制点，进行兴趣点与道路的精配准；

所述步骤2的具体过程如下：

步骤2.1，遍历所有道路结点，选择与道路结点连接的道路夹角之和最大的两两道路连接组合，对连接的道路进行道路链编号；

步骤2.2，将道路链编号相同的道路进行首尾连接，得到初始单线道路链；

步骤2.3，构建单线道路链的40米缓冲区，计算与缓冲区内邻近道路链的角度差异和长度比值，将角度差异和长度比值在给定阈值范围内的道路链识别为平行道路链，对各组平行道路链构建缓冲区多边形，提取缓冲区多边形的中心线，简化原有平行道路链；

步骤2.4，在简化后的中心线道路链与邻近道路链在相交处打断，相交处的结点添加到道路结点集合；遍历简化前道路链连接的道路链，判断与简化后中心线道路链是否连接，若连接，则删除悬挂边，若不连接，则延长中心线道路链使之连接。

2.根据权利要求1所述的方法，其特征在于，所述步骤3的具体过程如下：

步骤3.1，计算兴趣点与道路链的投影矢量x_ij，x_ij为兴趣点p_i到道路链S_j的投影点p_j’与p_j构成的矢量，初始化k个聚类及其中心c_j＝(0,0)，聚类参数λ₀＝1.5以及分配矩阵u_ij＝0；

k为道路链条数，k个聚类与步骤2提取的k条道路链对应；

步骤3.2，将每个兴趣点划分至最邻近道路链，兴趣点与道路链的邻近程度由两者的投影矢量与道路链对应聚类中心的欧氏距离计算得到：|x_ij-c_j|，若该值最小，则u_ij＝1；

步骤3.3，按照下式重新计算各道路链对应的聚类中心c_j和聚类参数λ：

$c_j=\frac{\underset{i=1}{\overset{m}{\&Sigma;}}{u}_{ij}*x_{ij}}{\underset{i=1}{\overset{m}{\&Sigma;}}{u}_{ij}},\mathrm{\&lambda;}={\mathrm{\&lambda;}}_0\frac{\underset{i=1}{\overset{m}{\&Sigma;}}\underset{j=1}{\overset{k}{\&Sigma;}}{u}_{ij}*|x_{ij}-c_j|}{\underset{i=1}{\overset{m}{\&Sigma;}}\underset{j=1}{\overset{k}{\&Sigma;}}{u}_{ij}}$

其中，u_ij为分配函数，u_ij＝1表示兴趣点p_i分配到道路链S_j，u_ij＝0表示p_i不分配到S_j，m为兴趣点个数；

步骤3.4，重新将各兴趣点划分至最邻近的道路链或空道路链；

对兴趣点p_i，若其与最邻近道路链S_j的聚类中心的欧氏距离|x_ij-c_j|小于聚类参数λ，则该兴趣点将被划分至道路链S_j，即u_ij＝1，否则被划分至空道路链，即u_ij＝0，j＝1…k；

步骤3.5，重复步骤3.3和3.4，重新计算聚类中心和重新分配兴趣点，每次迭代后，根据分配矩阵u_ij将所有兴趣点划分至k条道路链或空道路链，直到两次迭代间聚类中心的变化量小于给定阈值。

3.根据权利要求2所述的方法，其特征在于，所述步骤4的具体过程如下：

步骤4.1：通过步骤3的分配函数u_ij得到与道路链对应的k个兴趣点聚类，逐一遍历各兴趣点聚类PS＝{p₁…p_q}，q为该聚类中兴趣点个数，通过最小二乘拟合将聚类中兴趣点拟合成直线段，以拟合直线段为该兴趣点聚类的主曲线f，f＝{v₁…v_j…v_vn}，v_j为f的顶点，vn为顶点个数，此时vn＝2；

步骤4.2：将所有兴趣点投影到f的顶点v_j或边v_jv_j+1上，选择投影点个数最多的边，计算该边的中点，并将计算得到的中点添加到主曲线f中，此时vn＝vn+1；

步骤4.3：利用梯度下降法，以主曲线中每个顶点v_j的局部惩罚误差函数最小为目标逐一优化主曲线f的各顶点位置；

步骤4.4：重复步骤4.2和4.3，迭代添加新的顶点和优化主曲线顶点位置，直到顶点个数大于给定阈值，返回4.1拟合下一兴趣点聚类，直到拟合完所有的兴趣点聚类，获得兴趣点骨架图skeleton_poi＝{f_t|j＝1…k}，f_t为第t个兴趣点聚类拟合得到的主曲线；

4.根据权利要求3所述的方法，其特征在于，所述步骤4.3中主曲线f中顶点v_j的局部惩罚误差函数如下：

$G_q\left(v_j\right)=\frac{1}{q}{\mathrm{\&Delta;}}_q\left(v_j\right)+\mathrm{\&eta;}\&CenterDot;\frac{1}{k}Penalty\left(v_j\right)$

①q为兴趣点聚类PS中包含的兴趣点个数，Δ_q(v_j)为投影到顶点v_j或其邻边v_j-1v_j,v_jv_j+1的兴趣点到主曲线f的投影距离平方和，计算公式为：

${\mathrm{\&Delta;}}_q\left(v_j\right)={\mathrm{\&Sigma;}}_{p_i\&Element;S_{j-1}\&cup;S_j\&cup;V_j}{\&lsqb;d(p_i,f)\&rsqb;}^2$

其中，V_j和S_j-1、S_j-1分别为投影到顶点v_j和邻边v_j-1v_j，v_jv_j+1上的兴趣点集合，d(p_i,f)为兴趣点p_i到f的投影距离；

②Penalty(v_j)为f本身的形状惩罚因子，计算公式为：

$Penalty\left(v_j\right)=\left\{\begin{array}{cc}u\left(v_{j+1}\right)+\mathrm{\&pi;}\left(v_{j+1}\right),& j=1\\ u\left(v_j\right)+\mathrm{\&pi;}\left(v_j\right)+\mathrm{\&pi;}\left(v_{j+1}\right),& j=2\\ \mathrm{\&pi;}\left(v_{j-1}\right)+\mathrm{\&pi;}\left(v_j\right)+\mathrm{\&pi;}\left(v_{j+1}\right),& 2<j<vn-1\\ \mathrm{\&pi;}\left(v_{j-1}\right)+\mathrm{\&pi;}\left(v_j\right)+u\left(v_{j+1}\right),& j=vn-1\\ \mathrm{\&pi;}\left(v_{j-1}\right)+u\left(v_j\right),& j=vn\end{array}\right.$

其中，u(v_j)＝[d(v_j,v_j-1)]²,π(v_j)＝r²(1+cosγ_j)，γ_j为主曲线f在顶点v_j处的夹角，d(v_j,v_j-1)为顶点v_j，v_j-1的距离平方，r为聚类PS中所有兴趣点到中心点的最大距离；

③η为惩罚系数，计算公式为：

$\mathrm{\&eta;}={\mathrm{\&eta;}}_0*vn*q^{-1/3}*\sqrt{{\mathrm{\&Delta;}}_q\left(f\right)}*r^{-1}$

其中，Δ_q(f)为聚类PS中所有兴趣点到f的投影距离平方的均值，η₀为初始惩罚系数，设为0.4。

5.根据权利要求1所述的方法，其特征在于，所述步骤5的具体过程如下：

步骤5.1，根据步骤3得到的兴趣点与道路链对应关系，确定道路和兴趣点骨架图的对应交叉点；

逐一遍历所有道路链，寻找与之相交的道路链及其交点p＝(x_p,y_p)，确定与该道路链及相交道路链对应的兴趣点骨架图中的主曲线，计算主曲线的交点p’＝(x_p’,y_p’)，构成一对对应交叉点，其中p’为源控制点，p为目标控制点，相应位置偏差为(dx_p,dy_p),dx_p＝x_p-x_p’,dy_p＝y_p-y_p’，继续寻找下一条相交道路链，直到所有道路链均被遍历；

步骤5.2，以道路骨架图与兴趣点骨架图的对应交叉点作为同名控制点对，对同名控制点对中源控制点构建三角网，每个兴趣点将落在三角网的某个三角形中，则该兴趣点的位置偏差通过三角形三个顶点所对应的位置偏差加权得到，权重由该兴趣点到三角形三个顶点的距离确定，根据各兴趣点的位置偏差将兴趣点纠正到原始目标道路网中，完成精配准。