python - 订购用于绘图的多边形坐标

Question

我有一个由许多单元格组成的模型网格，我想在 matplotlib basemap 上为其绘制阴影多边形。

使用 pyproj，我首先投影点，然后使用 shapely.geometry 的 Polygon 类创建多边形以提取网格的外部坐标从。然后我将它们还原为 WGS84 以传递给我的绘图函数:

grid_x_mesh, grid_y_mesh = pyproj.transform(wgs84, nplaea, grid_lons, grid_lats)
grid_x = grid_x_mesh.ravel()
grid_y = grid_y_mesh.ravel()
grid_poly = Polygon(zip(grid_x, grid_y))
grid_x, grid_y = grid_poly.exterior.coords.xy
grid_plons, grid_plats = pyproj.transform(nplaea, wgs84, grid_x, grid_y)

然后，使用 matplotlib.basemap 方法，我将 WSG84 坐标投影到 map 投影(在本例中为 nplaea)和

grid_poly_x, grid_poly_y = m(grid_plons, grid_plats)
grid_poly_xy = zip(grid_poly_x, grid_poly_y)
grid_poly = Polygon(grid_poly_xy, facecolor='red', alpha=0.4)
plt.gca().add_patch(grid_poly)

尝试这样做时，我得到了一个十字交叉图案，我认为它必须对我提供给多边形函数的坐标进行排序。

我认为这与我提取外部坐标的方式有关，或者与我创建要绘制的最终多边形时坐标列表的顺序有关。

如果这是问题所在，是否有一种聪明的方法可以正确地对这些进行排序？

绘制的多边形

特写

Answer 1

我同意网格坐标有一些排列错误。 grid_lons 是如何创建的？将 Pyproj 与 Shapely 几何一起使用的一种更简洁的方法可能是使用相对较新的函数 shapely.ops.transform。例如:

import pyproj
from shapely.geometry import Polygon, Point
from shapely.ops import transform
from functools import partial

project = partial(
    pyproj.transform,
    pyproj.Proj(init='epsg:4326'),  # WGS84 geographic
    pyproj.Proj(init='epsg:3575'))  # North Pole LAEA Europe

# Example grid cell, in long/lat
poly_g = Polygon(((5, 52), (5, 60), (15, 60), (15, 52), (5, 52)))

# Transform to projected system
poly_p = transform(project, poly_g)

应该通过转换保持坐标的健全性(假设它们一开始是健全的)。