查找geopandas df中所有线串交叉点的每条线上的最近点

Question

我有一个geopandas数据框，其中包含几个从lat、lon点数据创建的线串。对于所有直线交点，我需要在每个直线串中找到距离该交点最近的点。

因此，如果dataframe中的两条线相交，我需要在每个linestring中找到距离该相交点最近的点。我使用itertools找到了所有可能的交叉点，这些交叉点与本文中公认的答案类似：https://gis.stackexchange.com/questions/137909/intersecting-lines-to-get-crossings-using-python-with-qgis

对于geopandas数据帧中的所有线交点，是否有更简单的方法来查找每个线串中距离交点最近的点？

我的dataframe看起来是这样的，每个文件名都包含自己的行字符串：

                                                            geometry
file                                                            
2015_may14_10  LINESTRING (-140.43855 59.80302, -140.44101 59...
2015_may14_11  LINESTRING (-140.84909 59.83433, -140.84758 59...
2015_may14_12  LINESTRING (-140.66859 59.79890, -140.66600 59...
2015_may14_15  LINESTRING (-140.19642 59.86655, -140.19795 59...
2015_may14_16  LINESTRING (-141.08783 59.94741, -141.08610 59...

Answer 1

让我们创建n条随机线：

import geopandas as gpd
from shapely.geometry import LineString, Point, Polygon
from shapely import wkt
import numpy as np
xmin, xmax, ymin, ymax = 0, 10000, 0, 10000
n = 100
xa = (xmax - xmin) * np.random.random(n) + xmin
ya = (ymax - ymin) * np.random.random(n) + ymin
xb = (xmax - xmin) * np.random.random(n) + xmin
yb = (ymax - ymin) * np.random.random(n) + ymin

lines = gpd.GeoDataFrame({'index':range(n),'geometry':[LineString([(a, b), (k, l)]) for a, b, k, l in zip(xa, ya, xb, yb)]})

这给了：

>>> lines
index   geometry
0   0   LINESTRING (4444.630 3081.439, 6132.674 5849.463)
1   1   LINESTRING (7015.940 6378.245, 4568.386 757.205)
2   2   LINESTRING (8766.417 6070.131, 690.359 7511.385)
3   3   LINESTRING (4245.544 4009.196, 8496.307 1557.175)
4   4   LINESTRING (1489.436 9364.784, 2109.740 5923.480)
...     ...     ...
95  95  LINESTRING (4783.454 7840.857, 1935.396 2435.260)
96  96  LINESTRING (1884.455 4982.662, 6257.958 3580.912)
97  97  LINESTRING (7072.811 7843.319, 4811.589 2486.040)
98  98  LINESTRING (6933.272 6427.046, 7528.579 2064.067)
99  99  LINESTRING (3876.400 5183.790, 5360.753 1901.207)

让我们来看看我们的十字路口：

res = []
for i in lines.loc[:, 'geometry']:
    for j in lines.loc[:, 'geometry']:
        inter = i.intersection(j)
        if inter.geom_type != 'LineString':
            res.append(inter)

这里我只是有一点误解，有时inter=I.intersection（j）返回一个LineString对象，我不知道两条不同的线如何作为一个交点输出另一条线（除非它们相同）。我让你来决定。

现在，我们可以创建带有结果点的df：

points = gpd.GeoDataFrame({'geometry':res})
>>>points

    geometry
0   POINT (4811.366 3682.806)
1   POINT (5149.727 4237.644)
2   POINT (4607.312 3348.202)
3   POINT (6026.639 5675.588)
4   POINT (4514.359 3195.779)
...     ...
2215    POINT (4788.793 3166.070)
2216    POINT (4704.895 3351.608)
2217    POINT (4581.390 3624.734)
2218    POINT (4320.392 4201.921)
2219    POINT (4949.041 2811.691)

2220 rows × 1 columns

我们可以看到，我们更多地使用线段而不是纯直线，因为交点（即点）的数量是2220。我不同意认为我们有幸拥有<代码> 7880 >代码>平行线。

然后，我们导入我们最好的朋友进行操作：

from shapely.ops import nearest_points

我们计算期望的输出：

intersection = []
line = []
my_point = []

for i in points.index:
    for j in lines.index:
        intersection.append(points.loc[i, 'geometry'])
        line.append(lines.loc[j, 'geometry'])
        my_point.append([p.wkt for p in nearest_points(points.loc[i, 'geometry'], lines.loc[j, 'geometry'])][1])


result = gpd.GeoDataFrame({'intersection':intersection, 'line':line, 'nearest_point':my_point})

result.geometry = result.loc[:, 'nearest_point'].apply(wkt.loads)
result.drop(columns=['nearest_point'], inplace=True)

>>>result

intersection    line    geometry
0   POINT (4811.365980053641 3682.805619834874)     LINESTRING (4444.630325108094 3081.43918610815...   POINT (4811.366 3682.806)
1   POINT (4811.365980053641 3682.805619834874)     LINESTRING (7015.939846319573 6378.24453843603...   POINT (5677.967 3305.464)
2   POINT (4811.365980053641 3682.805619834874)     LINESTRING (8766.416847858662 6070.13073873083...   POINT (5346.331 6680.480)
3   POINT (4811.365980053641 3682.805619834874)     LINESTRING (4245.544341245415 4009.19558793877...   POINT (4811.366 3682.806)
4   POINT (4811.365980053641 3682.805619834874)     LINESTRING (1489.4355376526 9364.784164867619,...   POINT (2109.740 5923.480)
...     ...     ...     ...
221995  POINT (4949.040525093341 2811.690701237854)     LINESTRING (4783.453909575222 7840.85687296287...   POINT (2745.435 3972.709)
221996  POINT (4949.040525093341 2811.690701237854)     LINESTRING (1884.454611847149 4982.66168904636...   POINT (5294.551 3889.693)
221997  POINT (4949.040525093341 2811.690701237854)     LINESTRING (7072.811488307434 7843.31900543939...   POINT (4949.041 2811.691)
221998  POINT (4949.040525093341 2811.690701237854)     LINESTRING (6933.272054846982 6427.04550331467...   POINT (7381.288 3143.559)
221999  POINT (4949.040525093341 2811.690701237854)     LINESTRING (3876.399925481877 5183.78974899146...   POINT (4949.041 2811.691)

222000 rows × 3 columns

希望，这回答了你的问题，如果你有更好的答案，请告诉我。