又是地图,用射线法(点在多边形内部)计算地图点阵,当边界点多起来的时候,计算就很难了。
高德地图给出的基础中国经纬度边界有8万多个点,运行起来基本瘫痪。所以希望能够在尽量保持形状的基础上减少点的个数。下面是最常用的 道格拉斯-普克算法
道格拉斯-普克(Douglas-Peuker)算法
Douglas-Peuker算法(DP算法)过程如下:
1、连接曲线首尾两点A、B;
2、依次计算曲线上所有点到A、B两点所在曲线的距离;
3、计算最大距离D,如果D小于阈值threshold,则去掉曲线上出A、B外的所有点;如果D大于阈值threshold,则把曲线以最大距离分割成两段;
4、对所有曲线分段重复1-3步骤,知道所有D均小于阈值。即完成抽稀。
这种算法的抽稀精度与阈值有很大关系,阈值越大,简化程度越大,点减少的越多;反之简化程度越低,点保留的越多,形状也越趋于原曲线。下面Python代码是抄来的:
# -*- coding: utf-8 -*-"""-------------------------------------------------File Name DouglasPeukerDescription : -Author : J_haodate 2017/8/16-------------------------------------------------Change Activity:2017/8/16: --------------------------------------------------"""from __future__ import divisionfrom math import sqrt, pow__author__ = 'J_hao'THRESHOLD = 0.0001 #def point2LineDistance(point_a, point_b, point_c):"""ab c:param point_a::param point_b::param point_c::return:"""# b cif point_b[0] == point_c[0]:return 9999999slope = (point_b[1] - point_c[1]) / (point_b[0] - point_c[0])intercept = point_b[1] - slope * point_b[0]# ab cdistance = abs(slope * point_a[0] - point_a[1] + intercept) / sqrt(1 + pow(slope, 2))return distanceclass DouglasPeuker(object):def __init__(self):self.threshold = THRESHOLDself.qualify_list = list()self.disqualify_list = list()def diluting(self, point_list):""":param point_list::return:"""if len(point_list) < 3:self.qualify_list.extend(point_list[::-1])else:#max_distance_index, max_distance = 0, 0for index, point in enumerate(point_list):if index in [0, len(point_list) - 1]:continuedistance = point2LineDistance(point, point_list[0], point_list[-1])if distance > max_distance:max_distance_index = indexmax_distance = distance#if max_distance < self.threshold:self.qualify_list.append(point_list[-1])self.qualify_list.append(point_list[0])else:#sequence_a = point_list[:max_distance_index]sequence_b = point_list[max_distance_index:]for sequence in [sequence_a, sequence_b]:if len(sequence) < 3 and sequence == sequence_b:self.qualify_list.extend(sequence[::-1])else:self.disqualify_list.append(sequence)def main(self, point_list):self.diluting(point_list)while len(self.disqualify_list) > 0:self.diluting(self.disqualify_list.pop())print self.qualify_listprint len(self.qualify_list)if __name__ == '__main__':d = DouglasPeuker()d.main([[104.066228, 30.644527], [104.066279, 30.643528]])
垂距限值法
垂距限值法其实和DP算法原理一样,但是垂距限值不是从整体角度考虑,而是依次扫描每一个点,检查是否符合要求。算法过程如下:
1、以第二个点开始,计算第二个点到前一个点和后一个点所在直线的距离d;
2、如果d大于阈值,则保留第二个点,计算第三个点到第二个点和第四个点所在直线的距离d;若d小于阈值则舍弃第二个点,计算第三个点到第一个点和第四个点所在直线的距离d;
3、依次类推,直线曲线上倒数第二个点。
# -*- coding: utf-8 -*-"""-------------------------------------------------File Name LimitVerticalDistanceDescription :Author : J_haodate 2017/8/17-------------------------------------------------Change Activity:2017/8/17:-------------------------------------------------"""from __future__ import divisionfrom math import sqrt, pow__author__ = 'J_hao'THRESHOLD = 0.0001 #def point2LineDistance(point_a, point_b, point_c):"""ab c:param point_a::param point_b::param point_c::return:"""# b cif point_b[0] == point_c[0]:return 9999999slope = (point_b[1] - point_c[1]) / (point_b[0] - point_c[0])intercept = point_b[1] - slope * point_b[0]# ab cdistance = abs(slope * point_a[0] - point_a[1] + intercept) / sqrt(1 + pow(slope, 2))return distanceclass LimitVerticalDistance(object):def __init__(self):self.threshold = THRESHOLDself.qualify_list = list()def diluting(self, point_list):""":param point_list::return:"""self.qualify_list.append(point_list[0])check_index = 1while check_index < len(point_list) - 1:distance = point2LineDistance(point_list[check_index], self.qualify_list[-1], point_list[check_index + 1])if distance < self.threshold:check_index += 1else:self.qualify_list.append(point_list[check_index])check_index += 1return self.qualify_listif __name__ == '__main__':l = LimitVerticalDistance()diluting = l.diluting([[104.066228, 30.644527], [104.066279,30.643528]])
其实DP算法和垂距限值法原理一样,DP算法是从整体上考虑一条完整的曲线,实现时较垂距限值法复杂,但垂距限值法可能会在某些情况下导致局部最优。另外在实际使用中发现采用点到另外两点所在直线距离的方法来判断偏离,在曲线弧度比较大的情况下比较准确。如果在曲线弧度比较小,弯曲程度不明显时,这种方法抽稀效果不是很理想,建议使用三点所围成的三角形面积作为判断标准。
