便携设备可以方便的获取经纬度信息,如果按照一定的时间间隔就能获取到具体的行动轨迹。比如:
- 网约车平台通过获取实时的司机位置判定司机是否偏航
- 运动记录平台用户的运动轨迹,比如跑步路径等。
如果记录时间间隔较大,则会丢失很多数据。如果记录时间间隔较小,则会记录很多距离较近的点,比如网约车堵车。所谓的经纬度轨迹的抽稀,其实就是删除哪些不重要的点。以下是整理的一些方案。
为了更好的讲解数据,我们这里采用高德的线路规划接口的数据进行模拟。(苏州站→苏州园区站)
# 获取线路规划 import requests def get_route(start_, end_coordinate): api_url = 'https://restapi.amap.com/v5/direction/driving' payloads = { "key": "", "origin": start_coordinate, # 苏州站 "destination": end_coordinate, # 上海虹桥机场 "show_fields": "polyline" } r = requests.get(api_url, payloads) data_json = r.json() steps = data_json['route']['paths'][0]['steps'] return steps steps = get_route(start_coordinate="120.610868,31.329679", end_coordinate="120.710574,31.341096") print(steps)
线路规划返回的数据:
'polyline':'120.672036,31.331114;120.672201,31.331152;120.672366,31.33119;120.672531,31.331228;120.672696,31.331266;120.672861,31.331304;120.673026,31.331342;120.673191,31.33138;120.673356,31.331418;120.673521,31.331456;120.673686,31.331494;120.673851,31.331532;120.674016,31.33157;120.674181,31.331608;120.674346,31.331646;120.674511,31.331684;120.674676,31.331722;120.674841,31.33176;120.675006,31.331798;120.675171,31.331836;120.675336,31.331874;120.675501,31.331912;120.675666,31.33195;120.675831,31.331988;120.675996,31.332026;120.676161,31.332064;120.676326,31.332102;120.676491,31.33214;120.676656,31.332178;120.676821,31.332216;120.676986,31.332254;120.677151,31.332292;120.677316,31.33233;120.677481,31.332368;120.677646,31.332406;120.677811,31.332444;120.677976,31.332482;120.678141,31.33252;120.678306,31.332558;120.678471,31.332596;120.678636,31.332634;120.678801,31.332672;120.678966,31.33271;120.679131,31.332748;120.679296,31.332786;120.679461,31.332824;120.679626,31.332862;120.679791,31.3329;120.679956,31.332938;120.680121,31.332976;120.680286,31.333014;120.680451,31.333052;120.680616,31.33309;120.680781,31.333128;120.680946,31.333166;120.681111,31.333204;120.681276,31.333242;120.681441,31.33328;120.681606,31.333318;120.681771,31.333356;120.681936,31.333394;120.682101,31.333432;120.682266,31.33347;120.682431,31.333508;120.682596,31.333546;120.682761,31.333584;120.682926,31.333622;120.683091,31.33366;120.683256,31.333698;120.683421,31.333736;120.683586,31.333774;120.683751,31.333812;120.683916,31.33385;120.684081,31.333888;120.684246,31.333926;120.684411,31.333964;120.684576,31.334002;120.684741,31.33404;120.684906,31.334078;120.685071,31.334116;120.685236,31.334154;120.685401,31.334192;120.685566,31.33423;120.685731,31.334268;120.685896,31.334306;120.686061,31.334344;120.686226,31.334382;120.686391,31.33442;120.686556,31.334458;120.686721,31.334496;120.686886,31.334534;120.687051,31.334572;120.687216,31.33461;120.687381,31.334648;120.687546,31.334686;120.687711,31.334724;120.687876,31.334762;120.688041,31.3348;120.688206,31.334838;120.688371,31.334876;120.688536,31.334914;120.688701,31.334952;120.688866,31.33499;120.689031,31.335028;120.689196,31.335066;120.689361,31.335104;120.689526,31.335142;120.689691,31.33518;120.689856,31.335218;120.690021,31.335256;120.690186,31.335294;120.690351,31.335332;120.690516,31.33537;120.690681,31.335408;120.690846,31.335446;120.691011,31.335484;120.691176,31.335522;120.691341,31.33556;120.691506,31.335598;120.691671,31.335636;120.691836,31.335674;120.692001,31.335712;120.692166,31.33575;120.692331,31.335788;120.692496,31.335826;120.692661,31.335864;120.692826,31.335902;120.692991,31.33594;120.693156,31.335978;120.693321,31.336016;120.693486,31.336054;120.693651,31.336092;120.693816,31.33613;120.693981,31.336168;120.694146,31.336206;120.694311,31.336244;120.694476,31.336282;120.694641,31.33632;120.694806,31.336358;120.694971,31.336396;120.695136,31.336434;120.695301,31.336472;120.695466,31.33651;120.695631,31.336548;120.695796,31.336586;120.695961,31.336624;120.696126,31.336662;120.696291,31.3367;120.696456,31.336738;120.696621,31.336776;120.696786,31.336814;120.696951,31.336852;120.697116,31.33689;120.697281,31.336928;120.697446,31.336966;120.697611,31.337004;120.697776,31.337042;120.697941,31.33708;120.698106,31.337118;120.698271,31.337156;120.698436,31.337194;120.698601,31.337232;120.698766,31.33727;120.698931,31.337308;120.699096,31.337346;120.699261,31.337384;120.699426,31.337422;120.699591,31.33746;120.699756,31.337498;120.699921,31.337536;120.700086,31.337574;120.700251,31.337612;120.700416,31.33765;120.700581,31.337688;120.700746,31.337726;120.700911,31.337764;120.701076,31.337802;120.701241,31.33784;120.701406,31.337878;120.701571,31.337916;120.701736,31.337954;120.701901,31.337992;120.702066,31.33803;120.702231,31.338068;120.702396,31.338106;120.702561,31.338144;120.702726,31.338182;120.702891,31.33822;120.703056,31.338258;120.703221,31.338296;120.703386,31.338334;120.703551,31.338372;120.703716,31.33841;120.703881,31.338448;120.704046,31.338486;120.704211,31.338524;120.704376,31.338562;120.704541,31.3386;120.704706,31.338638;120.704871,31.338676;120.705036,31.338714;120.705201,31.338752;120.705366,31.33879;120.705531,31.338828;120.705696,31.338866;120.705861,31.338904;120.706026,31.338942;120.706191,31.33898;120.706356,31.339018;120.706521,31.339056;120.706686,31.339094;120.706851,31.339132;120.707016,31.33917;120.707181,31.339208;120.707346,31.339246;120.707511,31.339284;120.707676,31.339322;120.707841,31.33936;120.708006,31.339398;120.708171,31.339436;120.708336,31.339474;120.708501,31.339512;120.708666,31.33955;120.708831,31.339588;120.708996,31.339626;120.709161,31.339664;120.709326,31.339702;120.709491,31.33974;120.709656,31.339778;120.709821,31.339816;120.709986,31.339854;120.710151,31.339892;120.710316,31.33993;120.710481,31.339968;120.710646,31.340006;120.710811,31.340044;120.710976,31.340082;120.711141,31.34012;120.711306,31.340158;120.711471,31.340196;120.711636,31.340234;120.711801,31.340272;120.711966,31.34031;120.712131,31.340348;120.712296,31.340386;120.712461,31.340424;120.712626,31.340462;120.712791,31.3405;120.712956,31.340538;120.713121,31.340576;120.713286,31.340614;120.713451,31.340652;120.713616,31.34069;120.713781,31.340728;120.713946,31.340766;120.714111,31.340804;120.714276,31.340842;120.714441,31.34088;120.714606,31.340918;120.714771,31.340956;120.714936,31.340994;120.715101,31.341032;120.715266,31.34107;120.715431,31.341108;120.715596,31.341146;120.715761,31.341184;120.715926,31.341222;120.716091,31.34126;120.716256,31.341298;120.716421,31.341336;120.716586,31.341374;120.716751,31.341412;120.716916,31.34145;120.717081,31.341488;120.717246,31.341526;120.717411,31.341564;120.717576,31.341602;120.717741,31.34164;120.717906,31.341678;120.718071,31.341716;120.718236,31.341754;120.718401,31.341792;120.718566,31.34183;120.718731,31.341868;120.718896,31.341906;120.719061,31.341944;120.719226,31.341982;120.719391,31.34202;120.719556,31.342058;120.719721,31.342096;120.719886,31.342134;120.720051,31.342172;120.720216,31.34221;120.720381,31.342248;120.720546,31.342286;120.720711,31.342324;120.720876,31.342362;120.721041,31.3424;120.721206,31.342438;120.721371,31.342476;120.721536,31.342514;120.721701,31.342552;120.721866,31.34259;120.722031,31.342628;120.722196,31.342666;120.722361,31.342704;120.722526,31.342742;120.722691,31.34278;120.722856,31.342818;120.723021,31.342856;120.723186,31.342894;120.723351,31.342932;120.723516,31.34297;120.723681,31.343008;120.723846,31.343046;120.724011,31.343084;120.724176,31.343122;120.724341,31.34316;120.724506,31.343198;120.724671,31.343236;120.724836,31.343274;120.725001,31.343312;120.725166,31.34335;120.725331,31.343388;120.725496,31.343426;120.725661,31.343464;120.725826,31.343502;120.725991,31.34354;120.726156,31.343578;120.726321,31.343616;120.726486,31.343654;120.726651,31.343692;120.726816,31.34373;120.726981,31.343768;120.727146,31.343806;120.727311,31.343844;120.727476,31.343882;120.727641,31.34392;120.727806,31.343958;120.727971,31.343996;120.728136,31.344034;120.728301,31.344072;120.728466,31.34411;120.728631,31.344148;120.728796,31.344186;120.728961,31.344224;120.729126,31.344262;120.729291,31.3443;120.729456,31.344338;120.729621,31.344376;120.729786,31.344414;120.729951,31.344452;120.730116,31.34449;120.730281,31.344528;120.730446,31.344566;120.730611,31.344604;120.730776,31.344642;120.730941,31.34468;120.731106,31.344718;120.731271,31.344756;120.731436,31.344794;120.731601,31.344832;120.731766,31.34487;120.731931,31.344908;120.732096,31.344946;120.732261,31.344984;120.732426,31.345022;120.732591,31.34506;120.732756,31.345098;120.732921,31.345136;120.733086,31.345174;120.733251,31.345212;120.733416,31.34525;120.733581,31.345288;120.733746,31.345326;120.733911,31.345364;120.734076,31.345402;120.734241,31.34544;120.734406,31.345478;120.734571,31.345516;120.734736,31.345554;120.734901,31.345592;120.735066,31.34563;120.735231,31.345668;120.735396,31.345706;120.735561,31.345744;120.735726,31.345782;120.735891,31.34582;120.736056,31.345858;120.736221,31.345896;120.736386,31.345934;120.736551,31.345972;120.736716,31.34601;120.736881,31.346048;120.737046,31.346086;120.737211,31.346124;120.737376,31.346162;120.737541,31.3462;120.737706,31.346238;120.737871,31.346276;120.738036,31.346314;120.738201,31.346352;120.738366,31.34639;120.738531,31.346428;120.738696,31.346466;120.738861,31.346504;120.739026,31.346542;120.739191,31.34658;120.739356,31.346618;120.'polyline': '120.672036,31.331114;120.672299,31.331173;120.672567,31.331231;120.673034,31.331317;120.673828,31.331473;120.673844,31.331478;120.677111,31.332133;120.679583,31.332626;120.681332,31.332991;120.681402,31.333007;120.682995,31.333372;120.684996,31.333817;120.68534,31.333892;120.688006,31.334498;120.689733,31.334858;120.691798,31.335239;120.695425,31.335904' },{ 'instruction': '沿玲珑街立交向东行驶201米靠左', 'orientation': '东', 'road_name': '玲珑街立交', 'step_distance': '201', 'polyline': '120.695425,31.335904;120.69572,31.335893;120.6958,31.335904;120.696562,31.336022;120.697157,31.336108;120.697511,31.336183' },{ 'instruction': '沿玲珑街立交途径玲珑街向东行驶1.5千米到达目的地', 'orientation': '东', 'road_name': '玲珑街立交', 'step_distance': '1548', 'polyline': '120.697511,31.336183;120.697962,31.336258;120.69896,31.33644;120.700961,31.336821;120.701138,31.336864;120.701288,31.336923;120.701411,31.336993;120.701508,31.337057;120.701588,31.337138;120.701658,31.337229;120.701728,31.337347;120.701771,31.337454;120.701798,31.337578;120.701808,31.337712;120.701792,31.337846;120.701744,31.33806;120.701637,31.338463;120.701588,31.338838;120.701605,31.339021;120.701621,31.339128;120.701674,31.339294;120.701696,31.339348;120.701808,31.339541;120.70205,31.339783;120.70228,31.339938;120.702586,31.340078;120.703418,31.340308;120.703509,31.34034;120.703579,31.340383;120.70499,31.340705;120.706406,31.341027;120.70743,31.341247;120.70772,31.341279;120.70794,31.34129;120.708766,31.34122;120.709034,31.341231;120.709576,31.341311;120.710016,31.341403;120.710381,31.341472;120.710465,31.341489' }]
线路规划返回的经纬度坐标点拼接后在地图上呈现如下:
points = [] for step in steps: points.extend(step['polyline'].split(';')) point_list = [(float(point.split(',')[1]), float(point.split(',')[0])) for point in points] m = folium.Map(location=[31.363364, 120.638043], zoom_start=12, width=800, height=400, zoom_control='False', tiles='https://map.geoq.cn/arcgis/rest/services/ChinaOnlineStreetGray/MapServer/tile/{z}/{y}/{x}', attr='GeoQ') folium.PolyLine(locations=point_list, color='green', weight=6, opacity=0.5).add_to(m) for point in point_list: folium.Circle(location=point, weight=1, color='green', fill=True, fill_opacity=0.8, opacity=1).add_to(m)
可以看到点之间非常的密集,有的甚至连成了直线。接下来要做的就是对这些点进行抽稀。
基于道路的抽稀
从上面的数据看,高德返回的数据是按照道路进行返回的,为了能得到抽稀的效果。可以采用选取每条道路起始点+最后道路的终点。代码示例:
points = [] for i, step in enumerate(steps): if i == len(steps) - 1: points.append(step['polyline'].split(';')[0]) points.append(step['polyline'].split(';')[-1]) else: points.append(step['polyline'].split(';')[0]) m = folium.Map(location=[31.363364, 120.638043], zoom_start=12, width=800, height=400, zoom_control='False', tiles='https://map.geoq.cn/arcgis/rest/services/ChinaOnlineStreetGray/MapServer/tile/{z}/{y}/{x}', attr='GeoQ') point_list = [(float(point.split(',')[1]), float(point.split(',')[0])) for point in points] folium.PolyLine(locations=point_list, color='green', weight=6, opacity=0.5).add_to(m) for point in point_list: folium.Circle(location=point, weight=1, color='green', fill=True, fill_opacity=0.8, opacity=1).add_to(m)
由于按道路抽稀,线段的长短与行驶的道路长短相关,市内相较于其他区域更加密集,若涉及到高速可能会更长。
基于距离的抽稀
基于距离的抽稀比较好理解,比如每个500米选取一个点。大致的逻辑如下:
points = [] for step in steps: points.extend(step['polyline'].split(';')) point_list = [(float(point.split(',')[1]), float(point.split(',')[0])) for point in points] import math def haversine(latlon1, latlon2): lat1, lon1 = latlon1 lat2, lon2 = latlon2 lon1, lat1, lon2, lat2 = map(math.radians, [lon1, lat1, lon2, lat2]) dlon = lon2 - lon1 dlat = lat2 - lat1 a = math.sin(dlat/2)**2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon/2)**2 c = 2 * math.asin(math.sqrt(a)) r = 6371 #地球平均半径,单位为公里 return c * r def distance_simplify(coords, epsilon): result_list = [coords[0]] for i in range(1, len(coords)): a = result_list[-1] b = coords[i] if haversine(a, b) < epsilon: pass else: result_list.append(coords[i]) return result_list points_simplify = distance_simplify(point_list, 0.5) m = folium.Map(location=[31.363364, 120.638043], zoom_start=12, width=800, height=400, zoom_control='False', tiles='https://map.geoq.cn/arcgis/rest/services/ChinaOnlineStreetGray/MapServer/tile/{z}/{y}/{x}', attr='GeoQ') folium.PolyLine(locations=points_simplify, color='green', weight=6, opacity=0.5).add_to(m) for point in points_simplify: folium.Circle(location=point, weight=1, color='green', fill=True, fill_opacity=0.8, opacity=1).add_to(m)
抽稀后整体点与点之间的距离差不多。
基于Douglas-Peucker算法的抽稀
道格拉斯-普克算法(Douglas–Peucker algorithm,亦称为拉默-道格拉斯-普克算法、迭代适应点算法、分裂与合并算法)是将曲线近似表示为一系列点,并减少点的数量的一种算法。该算法的原始类型分别由乌尔斯·拉默(Urs Ramer)于1972年以及大卫·道格拉斯(David Douglas)和托马斯·普克(Thomas Peucker)于1973年提出,并在之后的数十年中由其他学者予以完善。
算法逻辑:
- 在轨迹曲线在曲线首尾两点A,B之间连接一条直线AB,该直线为曲线的弦;
- 遍历曲线上其他所有点,求每个点到直线AB的距离,找到最大距离的点C,最大距离记为maxDistance
- 比较该距离maxDistance与预先定义的阈值epsilon大小,如果maxDistance<epsilon,则将该直线AB作为曲线段的近似,舍去AB之间的所有点,曲线段处理完毕;
- 若maxDistance>=epsilon,则使C点将曲线AB分为AC和CB两段,并分别对这两段进行(1)~(3)步处理;
- 当所有曲线都处理完毕时,依次连接各个分割点形成的折线,即为原始曲线的路径。
使用Python实现道格拉斯-普克算法:
import numpy as np #定义一个函数计算点到线的距离 def perpendicular_distance(point, line_start, line_end): if np.all(np.equal(line_start, line_end)): return np.linalg.norm(point-line_start) return np.divide( np.abs(np.linalg.norm(np.cross(line_end-line_start, line_start-point))), np.linalg.norm(line_end-line_start) ) #实现道格拉斯-普克算法 def douglas_peucker(points, tolerance): if len(points) < 3: return points #初始化最远点的距离和索引 max_distance = 0 index_of_farthest = 0 #查找最远点 for i in range(1, len(points)-1): distance = perpendicular_distance(np.array(points[i]), np.array(points[0]), np.array(points[-1])) if distance > max_distance: max_distance = distance index_of_farthest = i #如果最大距离超过容忍度,递归地抽稀 if max_distance > tolerance: #递归处理线段两边的点 left_side = douglas_peucker(points[:index_of_farthest+1], tolerance) right_side = douglas_peucker(points[index_of_farthest:], tolerance) #合并结果 return left_side[:-1] + right_side else: return [points[0], points[-1]] #设置容忍度 tolerance = 0.0001 #调用道格拉斯-普克算法 reduced_points = douglas_peucker(point_list, tolerance) m = folium.Map(location=[31.363364, 120.638043], zoom_start=12, width=800, height=400, zoom_control='False', tiles='https://map.geoq.cn/arcgis/rest/services/ChinaOnlineStreetGray/MapServer/tile/{z}/{y}/{x}', attr='GeoQ') folium.PolyLine(locations=reduced_points, color='green', weight=6, opacity=0.5).add_to(m) for point in reduced_points: folium.Circle(location=point, weight=1, color='green', fill=True, fill_opacity=0.8, opacity=1).add_to(m)
可以看到在拐弯处抽稀后比较密集,直线处比较稀疏。
具体在经纬度抽稀时使用何种方法,还需要看具体的使用场景,如果是界面上线程线路,建议使用格拉斯-普克算法。如果是做轨迹的相似度检验,建议使用按距离抽稀。