scipy 如何在Pandas内插GPS坐标?

ca1c2owp  于 2023-06-06  发布在  其他
关注(0)|答案(1)|浏览(176)

我有一个数据集,包含大约150K个GPS坐标条目,看起来像这样:

log_time    latitude    longitude
0   1.555840e+09    45.429597   11.974981
1   1.555869e+09    45.429597   11.974981
3   1.555869e+09    45.429596   11.974984
4   1.555869e+09    45.429490   11.975089
5   1.555869e+09    45.429092   11.975478
count                     147538
mean      0 days 00:02:27.234798
std       0 days 02:34:54.243149
min              0 days 00:00:00
25%              0 days 00:00:03
50%              0 days 00:00:05
75%              0 days 00:00:08
max      39 days 12:25:39.551000
Name: log_time, dtype: object

我希望这样的Dataframe在不久的将来可以扩展到数百万条记录,因此可扩展性是一个优先事项。
我想对移动进行插值,这样对于更宽的间隙,至少每60秒有一个GPS记录。
标准方法为:

dff = dff.set_index(dff.pop('log_time'))
dff = dff.reindex(np.arange(dff.index.min(), dff.index.max()+1))

生产:

latitude    longitude
log_time        
1.555840e+09    45.429597   11.974981
1.555840e+09    NaN NaN
1.555840e+09    NaN NaN
1.555840e+09    NaN NaN
1.555840e+09    NaN NaN

这将被插入类似dff.interpolate().reset_index()的东西。
但我有一个很大的问题:scipy(以及pandas)提供的插值函数都不适用于GPS距离,即圆弧,而不是直线。从I've seen开始,没有简单的方法来扩展插值函数
我已经有了要使用的distance函数,但我认为没有简单的方法来部署它,而不必求助于嵌套的for循环。

from geographiclib.geodesic import Geodesic
geod = Geodesic.WGS84

def custom_interpolation(starting_value, ending_value, number_of_missing_values):
    filled_array = [starting_value]

    # 1. create a line between starting_value and ending_value 
    # by solving the inverse geodesic problem
    line = geod.InverseLine(starting_value.lat, starting_value.lon, ending_value.lat, ending_value.long)

    # 2. Determine the length of the steps needed to fill 
    # the missing values between the two extremes; 
    # s13 is the total arc length of the line
    step_lenght = line.s13 / number_of_missing_values

    # 3. Add mid values between the two arrays
    for i in range(1, n + 1):
        distance = min(step_lenght * i, line.s13)
        g = line.Position(distance, Geodesic.STANDARD | Geodesic.LONG_UNROLL)
        filled_array.append(g['lat2'], g['lon2'])

    filled_array.append(ending_value)
    return filled_array

所以像[(LAT1, LON1), None, None, None, (LAT2, LON2)]这样的东西可以变成[(LAT1, LON1), (LAT, LON), (LAT, LON), (LAT2, LON2)]

nbnkbykc

nbnkbykc1#

可以使用最近的数据点进行插值;对于每个点x、y,找到四个其他坐标,并使用距离作为权重来平均相关联的数据点、时间、海拔等。下面的代码使用欧几里德,它在短距离上足够好,但它可以很容易地改变为更精确的距离度量,

class QuadTreeInterpolator:
    def __init__(self, x, y, z):
        self.x = x
        self.y = y
        self.z = z
        self.tree = quads.QuadTree((np.mean(x), np.mean(y)), np.std(x)*4, np.std(y)*4)
        for xx,yy,zz in zip(x,y,z):
            self.tree.insert((xx,yy),data=zz)

    def interp_cell(self, x, y, points):
        a = np.array([x,y]).reshape(-1,2)
        b = np.array(points)[:,:2]
        ds = cdist(a,b)
        ds = ds / np.sum(ds)
        ds = 1. - ds
        c = np.array(points)[:,2]
        iz = np.sum(c * ds) / np.sum(ds)
        return iz
            
    def interpolate(self,x,y):
        res = self.tree.nearest_neighbors((x,y), count=4)
        points = np.array([[c.x, c.y, c.data] for c in res])
        return self.interp_cell(x, y, points)

相关问题