如何使用 python 生成规则的地理网格？

Question

我想检索特定地图区域上规则网格的所有 lat/lon 坐标对。我找到了 geopy 库，但根本无法解决问题。

例如，我有一个矩形地理区域，其四个角在 lat/lon 坐标中描述，我试图计算间距为例如覆盖此区域 1 公里。

Answer 1

初步考虑

特定区域的定义方式略有不同。如果它只是一个矩形区域（注意：投影中的矩形在地球表面上不一定是矩形！），您可以使用所需的步长在两个坐标维度上简单地从最小值迭代到最大值。如果您手头有一个任意多边形形状，则需要测试生成的哪些点 与该多边形相交 并且仅 return 满足此条件的坐标对。

计算规则网格

规则网格不等于跨投影的规则网格。你说的是 latitude/longitude 对，这是一个极坐标系，以地球表面形状的近似度为单位测量。在 lat/lon (EPSG:4326) 中，距离的测量单位不是 meters/kilometers/miles，而是度数。

此外，我假设您想计算一个网格，其“水平”步长与赤道（即纬度）平行。对于其他网格（例如旋转的矩形网格，垂直与经度平行等），您需要花费更多的精力来变换形状。

问问自己：您要创建以度为单位还是以米为单位的规则间隔网格？

以度为单位的网格

如果你想要度数，你可以简单地迭代：

stepsize = 0.001
for x in range(lonmin, lonmax, stepsize):
    for y in range(latmin, latmax, stepsize):
        yield (x, y)

但是： 请务必了解以度为单位的台阶的长度（以米为单位）在整个地球表面是不一样的。例如，靠近赤道的 0.001 增量纬度在表面上覆盖的距离（以米为单位）与靠近两极的不同。

以米为单位的网格

如果您想要以米为单位，您需要将输入区域（地图上的特定区域）的lat/lon边界投影到坐标系中支持以米为单位的距离。您可以使用 Haversine formula 作为粗略的近似值来计算 lat/lon 对之间的距离，但这不是您可以使用的最佳方法。

更好的方法是搜索合适的投影，将您感兴趣的区域转换为该投影，通过直接迭代创建网格，获取点，并将它们投影回 lat/lon 对。例如，适合欧洲的投影是 EPSG:3035。 Google 顺便说一下，地图使用 EPSG:900913 作为他们的网络地图服务。

在 python 中，您可以使用库 shapely 和 pyproj 处理地理形状和投影：

import shapely.geometry
import pyproj

# Set up transformers, EPSG:3857 is metric, same as EPSG:900913
to_proxy_transformer = pyproj.Transformer.from_crs('epsg:4326', 'epsg:3857')
to_original_transformer = pyproj.Transformer.from_crs('epsg:4326', 'epsg:3857')

# Create corners of rectangle to be transformed to a grid
sw = shapely.geometry.Point((-5.0, 40.0))
ne = shapely.geometry.Point((-4.0, 41.0))

stepsize = 5000 # 5 km grid step size

# Project corners to target projection
transformed_sw = to_proxy_transformer.transform(sw.x, sw.y) # Transform NW point to 3857
transformed_ne = to_proxy_transformer.transform(ne.x, ne.y) # .. same for SE

# Iterate over 2D area
gridpoints = []
x = transformed_sw[0]
while x < transformed_ne[0]:
    y = transformed_sw[1]
    while y < transformed_ne[1]:
        p = shapely.geometry.Point(to_original_transformer.transform(x, y))
        gridpoints.append(p)
        y += stepsize
    x += stepsize

with open('testout.csv', 'wb') as of:
    of.write('lon;lat\n')
    for p in gridpoints:
        of.write('{:f};{:f}\n'.format(p.x, p.y))

这个例子生成了这个均匀分布的网格：

如何使用 python 生成规则的地理网格？

How can I generate a regular geographic grid using python?

python

google-maps

geopy

初步考虑

计算规则网格