3D CartoPy 类似于 Matplotlib-Basemap [英] 3D CartoPy similar to Matplotlib-Basemap
问题描述
我是 Python 新手,有一个关于 Cartopy 能够用于 3D 绘图的问题.下面是一个使用 matplotlibBasemap
的例子.
I'm new to Python with a question about Cartopy being able to be used in a 3D plot. Below is an example using matplotlibBasemap
.
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from mpl_toolkits.basemap import Basemap
m = Basemap(projection='merc',
llcrnrlat=52.0,urcrnrlat=58.0,
llcrnrlon=19.0,urcrnrlon=40.0,
rsphere=6371200.,resolution='h',area_thresh=10)
fig = plt.figure()
ax = Axes3D(fig)
ax.add_collection3d(m.drawcoastlines(linewidth=0.25))
ax.add_collection3d(m.drawcountries(linewidth=0.35))
ax.add_collection3d(m.drawrivers(color='blue'))
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Height')
fig.show()
这会在 3D 轴内创建地图,以便您可以在表面上绘制对象.但是 Cartopy 返回一个 matplotlib.axes.GeoAxesSubplot
.不清楚如何使用 matplotlib-basemap
将其添加到 3D 图形/轴,如上所述.
This creates a map within a 3D axis so that you can plot objects over the surface. But with Cartopy returns a matplotlib.axes.GeoAxesSubplot
. Not clear how to take this and add to a 3D figure/axis as above with matplotlib-basemap
.
那么,有人可以就如何使用 Cartopy 进行类似的 3D 绘图提供任何指示吗?
So, can someone give any pointers on how to do a similar 3D plot with Cartopy?
推荐答案
底图 mpl3d 是一个非常巧妙的 hack,但它没有被设计为以描述的方式运行.因此,除了简单的海岸线之外,您目前无法使用相同的技术.例如,填充的大陆在 AFAICT 中不起作用.
The basemap mpl3d is a pretty neat hack, but it hasn't been designed to function in the described way. As a result, you can't currently use the same technique for much other than simple coastlines. For example, filled continents just don't work AFAICT.
也就是说,在使用 cartopy 时可以使用类似的 hack.由于我们可以通用地访问 shapefile 信息,因此该解决方案应该适用于任何折线 shapefile,例如海岸线.
That said, a similar hack is available when using cartopy. Since we can access shapefile information generically, this solution should work for any poly-line shapefile such as coastlines.
第一步是获取 shapefile 和相应的几何图形:
The first step is to get hold of the shapefile, and the respective geometries:
feature = cartopy.feature.NaturalEarthFeature('physical', 'coastline', '110m')
geoms = feature.geometries()
接下来,我们可以将这些转换为所需的投影:
Next, we can convert these to the desired projection:
target_projection = ccrs.PlateCarree()
geoms = [target_projection.project_geometry(geom, feature.crs)
for geom in geoms]
由于这些是匀称的几何图形,因此我们希望将它们转换为 matplotlib 路径:
Since these are shapely geometries, we then want to convert them to matplotlib paths with:
from cartopy.mpl.patch import geos_to_path
import itertools
paths = list(itertools.chain.from_iterable(geos_to_path(geom)
for geom in geoms))
对于路径,我们应该能够在 matplotlib 中创建一个 PathCollection,并将其添加到坐标区,但遗憾的是,Axes3D 似乎无法处理 PathCollection 实例,因此我们需要通过构建 LineCollection 来解决此问题(就像底图一样).遗憾的是 LineCollections 不采用路径,而是采用段,我们可以用它来计算:
With paths, we should be able to just create a PathCollection in matplotlib, and add it to the axes, but sadly, Axes3D doesn't seem to cope with PathCollection instances, so we need to workaround this by constructing a LineCollection (as basemap does). Sadly LineCollections don't take paths, but segments, which we can compute with:
segments = []
for path in paths:
vertices = [vertex for vertex, _ in path.iter_segments()]
vertices = np.asarray(vertices)
segments.append(vertices)
将所有这些放在一起,我们最终得到与您的代码生成的底图相似的结果:
Pulling this all together, we end up with a similar result to the basemap plot which your code produces:
import itertools
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from matplotlib.collections import LineCollection
import numpy as np
import cartopy.feature
from cartopy.mpl.patch import geos_to_path
import cartopy.crs as ccrs
fig = plt.figure()
ax = Axes3D(fig, xlim=[-180, 180], ylim=[-90, 90])
ax.set_zlim(bottom=0)
target_projection = ccrs.PlateCarree()
feature = cartopy.feature.NaturalEarthFeature('physical', 'coastline', '110m')
geoms = feature.geometries()
geoms = [target_projection.project_geometry(geom, feature.crs)
for geom in geoms]
paths = list(itertools.chain.from_iterable(geos_to_path(geom) for geom in geoms))
# At this point, we start working around mpl3d's slightly broken interfaces.
# So we produce a LineCollection rather than a PathCollection.
segments = []
for path in paths:
vertices = [vertex for vertex, _ in path.iter_segments()]
vertices = np.asarray(vertices)
segments.append(vertices)
lc = LineCollection(segments, color='black')
ax.add_collection3d(lc)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Height')
plt.show()
最重要的是,mpl3d 似乎可以很好地处理 PolyCollection,这将是我研究填充几何图形的路线,例如陆地轮廓(与海岸线相反,严格来说是轮廓).
On top of this, mpl3d seems to handle PolyCollection well, which would be the route I would investigate for filled geometries, such as the land outline (as opposed to the coastline, which is strictly an outline).
重要的一步是将路径转换为多边形,并在 PolyCollection 对象中使用它们:
The important step is to convert the paths to polygons, and use these in a PolyCollection object:
concat = lambda iterable: list(itertools.chain.from_iterable(iterable))
polys = concat(path.to_polygons() for path in paths)
lc = PolyCollection(polys, edgecolor='black',
facecolor='green', closed=False)
此案例的完整代码如下所示:
The complete code for this case would look something like:
import itertools
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from matplotlib.collections import LineCollection, PolyCollection
import numpy as np
import cartopy.feature
from cartopy.mpl.patch import geos_to_path
import cartopy.crs as ccrs
fig = plt.figure()
ax = Axes3D(fig, xlim=[-180, 180], ylim=[-90, 90])
ax.set_zlim(bottom=0)
concat = lambda iterable: list(itertools.chain.from_iterable(iterable))
target_projection = ccrs.PlateCarree()
feature = cartopy.feature.NaturalEarthFeature('physical', 'land', '110m')
geoms = feature.geometries()
geoms = [target_projection.project_geometry(geom, feature.crs)
for geom in geoms]
paths = concat(geos_to_path(geom) for geom in geoms)
polys = concat(path.to_polygons() for path in paths)
lc = PolyCollection(polys, edgecolor='black',
facecolor='green', closed=False)
ax.add_collection3d(lc)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Height')
plt.show()
屈服:
这篇关于3D CartoPy 类似于 Matplotlib-Basemap的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!