如何找到线与网格的交点?

2022-01-18 00:00:00 python numpy algorithm intersection grid

问题描述

我有轨迹数据，其中每个轨迹由一系列坐标(x，y 点)组成，每个轨迹都由唯一的 ID 标识.

I have trajectory data, where each trajectory consists of a sequence of coordinates(x, y points) and each trajectory is identified by a unique ID.

这些轨迹在 x - y 平面上，我想将整个平面划分为大小相等的网格(方形网格).该网格显然是不可见的，但用于将轨迹划分为子段.每当轨迹与网格线相交时，它就会在此处分段，并成为具有 new_id 的新子轨迹.

These trajectories are in x - y plane, and I want to divide the whole plane into equal sized grid (square grid). This grid is obviously invisible but is used to divide trajectories into sub-segments. Whenever a trajectory intersects with a grid line, it is segmented there and becomes a new sub-trajectory with new_id.

我已经包含了一个简单的手工图表，以明确我的期望.

I have included a simple handmade graph to make clear what I am expecting.

可以看出轨迹在网格线的交点处是如何划分的，并且这些线段中的每一个都有新的唯一id.

It can be seen how the trajectory is divided at the intersections of the grid lines, and each of these segments has new unique id.

我正在研究 Python，并寻求一些 Python 实现链接、建议、算法，甚至是相同的伪代码.

I am working on Python, and seek some python implementation links, suggestions, algorithms, or even a pseudocode for the same.

如果有什么不清楚的地方请告诉我.

Please let me know if anything is unclear.

更新

为了将平面划分为网格，单元格索引如下:

In order to divide the plane into grid , cell indexing is done as following:

#finding cell id for each coordinate
#cellid = (coord / cellSize).astype(int)
cellid = (coord / 0.5).astype(int)
cellid
Out[] : array([[1, 1],
              [3, 1],
              [4, 2],
              [4, 4],
              [5, 5],
              [6, 5]])
#Getting x-cell id and y-cell id separately 
x_cellid = cellid[:,0]
y_cellid = cellid[:,1]

#finding total number of cells
xmax = df.xcoord.max()
xmin = df.xcoord.min()
ymax = df.ycoord.max()
ymin = df.ycoord.min()
no_of_xcells = math.floor((xmax-xmin)/ 0.5)
no_of_ycells = math.floor((ymax-ymin)/ 0.5)
total_cells = no_of_xcells * no_of_ycells
total_cells
Out[] : 25

由于平面现在被分为 25 个单元格，每个单元格都有一个 cellid.为了找到交叉点，也许我可以检查轨迹中的下一个坐标，如果 cellid 保持不变，那么轨迹的那段在同一个单元格中并且与网格没有交叉点.比如说，如果 x_cellid[2] 大于 x_cellid[0]，那么段与垂直网格线相交.不过，我仍然不确定如何找到与网格线的交点，并在交点上分割轨迹，为它们提供新的 id.

Since the plane is now divided into 25 cells each with a cellid. In order to find intersections, maybe I could check the next coordinate in the trajectory, if the cellid remains the same, then that segment of the trajectory is in the same cell and has no intersection with grid. Say, if x_cellid[2] is greater than x_cellid[0], then segment intersects vertical grid lines. Though, I am still unsure how to find the intersections with the grid lines and segment the trajectory on intersections giving them new id.

解决方案

这个可以通过shapely来解决:

This can be solved by shapely:

%matplotlib inline
import pylab as pl
from shapely.geometry import MultiLineString, LineString
import numpy as np
from matplotlib.collections import LineCollection

x0, y0, x1, y1 = -10, -10, 10, 10
n = 11

lines = []
for x in np.linspace(x0, x1, n):
    lines.append(((x, y0), (x, y1)))

for y in np.linspace(y0, y1, n):
    lines.append(((x0, y), (x1, y)))

grid = MultiLineString(lines)

x = np.linspace(-9, 9, 200)
y = np.sin(x)*x
line = LineString(np.c_[x, y])

fig, ax = pl.subplots()
for i, segment in enumerate(line.difference(grid)):
    x, y = segment.xy
    pl.plot(x, y)
    pl.text(np.mean(x), np.mean(y), str(i))

lc = LineCollection(lines, color="gray", lw=1, alpha=0.5)
ax.add_collection(lc);

结果:

不使用shapely，自己动手:

To not use shapely, and do it yourself:

import pylab as pl
import numpy as np
from matplotlib.collections import LineCollection

x0, y0, x1, y1 = -10, -10, 10, 10
n = 11
xgrid = np.linspace(x0, x1, n)
ygrid = np.linspace(y0, y1, n)
x = np.linspace(-9, 9, 200)
y = np.sin(x)*x
t = np.arange(len(x))

idx_grid, idx_t = np.where((xgrid[:, None] - x[None, :-1]) * (xgrid[:, None] - x[None, 1:]) <= 0)
tx = idx_t + (xgrid[idx_grid] - x[idx_t]) / (x[idx_t+1] - x[idx_t])

idx_grid, idx_t = np.where((ygrid[:, None] - y[None, :-1]) * (ygrid[:, None] - y[None, 1:]) <= 0)
ty = idx_t + (ygrid[idx_grid] - y[idx_t]) / (y[idx_t+1] - y[idx_t])

t2 = np.sort(np.r_[t, tx, tx, ty, ty])

x2 = np.interp(t2, t, x)
y2 = np.interp(t2, t, y)

loc = np.where(np.diff(t2) == 0)[0] + 1

xlist = np.split(x2, loc)
ylist = np.split(y2, loc)


fig, ax = pl.subplots()
for i, (xp, yp) in enumerate(zip(xlist, ylist)):
    pl.plot(xp, yp)
    pl.text(np.mean(xp), np.mean(yp), str(i))


lines = []
for x in np.linspace(x0, x1, n):
    lines.append(((x, y0), (x, y1)))

for y in np.linspace(y0, y1, n):
    lines.append(((x0, y), (x1, y)))

lc = LineCollection(lines, color="gray", lw=1, alpha=0.5)
ax.add_collection(lc);

相关文章