滑動視窗軌跡壓縮Python實現
阿新 • • 發佈:2019-02-19
軌跡壓縮
問題描述:
已知計程車的運動軌跡點檔案,在誤差範圍內,壓縮軌跡。
平臺:
Linux
語言:
Python3.5
結果截圖:
解決思路:
最容易想到的應該是DP演算法,即取初始軌跡的起點A和終點B連線,計算每個點到這條線的距離,距離最大的點C若小於要求誤差則結束;否則將C點加入壓縮後的資料集,對AC和CB重複以上過程直至滿足誤差要求。
DP演算法Java實現:http://www.cnblogs.com/xdlwd086/p/5100425.html
正是參考了這篇文章,作業才得以完成,在此謝過~~~
此篇文章將實現另一種軌跡壓縮演算法一種改進的滑動視窗軌跡資料壓縮演算法
程式碼實現
from math import *
import sys
sys.setrecursionlimit(10000)
def dfTodu(str):
str = str.split('.')
d = (float(str[0][:-2])) + (float(str[0][-2:] + str[1]) / 600000)
return d
#將檔案中的經緯度提取出儲存到陣列
def createDataSet():
rfile = open("2007-10-14-GPS.log" )
time = []
jwd = []
index = 0
for line in rfile:
strlist = line.split(' ', 8)
time.append(strlist[2])
jwd.append([dfTodu(strlist[5]), dfTodu(strlist[3]), index])
index += 1
rfile.close()
return time, jwd
#將壓縮後的軌跡寫入陣列
def Write(enpArrayFilter):
wfile = open('2016-10-28.log', "w")
for i in enpArrayFilter:
wfile.write(str(i) + "\n")
def geoDist(pA, pB):
radLat1 = Rad(pA[0])
radLat2 = Rad(pB[0])
delta_lon = Rad(pB[1] - pA[1])
top_1 = cos(radLat2) * sin(delta_lon)
top_2 = cos(radLat1) * sin(radLat2) - sin(radLat1) * cos(radLat2) * cos(delta_lon)
top = sqrt(top_1 * top_1 + top_2 * top_2)
bottom = sin(radLat1) * sin(radLat2) + cos(radLat1) * cos(radLat2) * cos(delta_lon)
delta_sigma = atan2(top, bottom)
distance = delta_sigma * 6378137.0
return distance
def Rad(d):
return d * pi / 180.0
def distToSegment(pA, pB, pX):
a = abs(geoDist(pA, pB))
b = abs(geoDist(pA, pX))
c = abs(geoDist(pB, pX))
p = (a + b + c) / 2.0
s = sqrt(abs(p * (p - a) * (p - b) * (p - c)))
d = s * 2.0 / a
return d
#滑動視窗演算法
#enpInitenpInit初始軌跡點陣列
#enpArrayFilter過濾陣列
#start視窗內的起始點
#end視窗內的終點
#cur當前點
#m目前誤差最大的點
#DMax最大誤差
#count
def SlideWindow(enpInit, enpArrayFilter, start, end, cur, m, DMax, count):
if (end < count):
d_cur = distToSegment(enpInit[start], enpInit[end], enpInit[cur]) # 當前點到對應線段的距離
d_m = distToSegment(enpInit[start], enpInit[end], enpInit[m]) # 當前點到對應線段的距離
if (d_cur > DMax or d_m > DMax):
enpArrayFilter.append(enpInit[cur]) # 將當前點加入到過濾陣列中
start = cur
cur = start + 1
end = start + 2
m = cur
d_m = 0
SlideWindow(enpInit, enpArrayFilter, start, end, cur, m, DMax, count)
elif ((d_cur <= DMax) and (d_m <= DMax)):
if (d_cur > d_m):
m = cur
cur = end
end = end + 1
SlideWindow(enpInit, enpArrayFilter, start, end, cur, m, DMax, count)
#平均誤差計算函式
def getMeanDistError(enpInit, enpArrayFilter):
sumDist = 0
for i in range(1, len(enpArrayFilter)):
start = enpArrayFilter[i - 1][2]
end = enpArrayFilter[i][2]
for j in range(start + 1, end):
sumDist += distToSegment(enpInit[start], enpInit[end], enpInit[j])
return sumDist / len(enpInit)
time, jwd = createDataSet()
start = 0
end = 2
cur = 1
DMax = 30
enpInit = jwd
m = 1
count = len(enpInit) - 1
enpArrayFilter = []
enpArrayFilter.append(enpInit[0]) # 獲取第一個原始經緯度點座標並新增到過濾後的陣列中
SlideWindow(enpInit, enpArrayFilter, start, end, cur, m, DMax, count)
Write(enpArrayFilter)
c_enpInit = len(enpInit)
c_enpArrayFilter = len(enpArrayFilter)
print("壓縮前的點數=" + str(c_enpInit))
print("壓縮後的點數=" + str(c_enpArrayFilter))
print("壓縮率=" + str(c_enpArrayFilter / c_enpInit * 100) + "%")
print("平均誤差=" + str(getMeanDistError(enpInit, enpArrayFilter)))