我想使用 DEM 文件使用 matplotlib 生成模拟地形表面。但我不知道如何将栅格坐标地理配准到给定的 CRS。我也不知道如何以适合在 3D matplotlib 绘图中使用的格式来表达地理引用栅格,例如作为 numpy 数组。
到目前为止,这是我的 python 代码:
import osgeo.gdal
dataset = osgeo.gdal.Open("MergedDEM")
gt = dataset.GetGeoTransform()
最佳答案
您可以使用 matplotlib 中的普通 plot_surface
方法。因为它需要一个 X 和 Y 数组,所以它已经用正确的坐标绘制。我总是发现很难做出好看的 3D 图,所以视觉方面当然可以改进。 :)
import gdal
from mpl_toolkits.mplot3d import Axes3D
dem = gdal.Open('gmted_small.tif')
gt = dem.GetGeoTransform()
dem = dem.ReadAsArray()
fig, ax = plt.subplots(figsize=(16,8), subplot_kw={'projection': '3d'})
xres = gt[1]
yres = gt[5]
X = np.arange(gt[0], gt[0] + dem.shape[1]*xres, xres)
Y = np.arange(gt[3], gt[3] + dem.shape[0]*yres, yres)
X, Y = np.meshgrid(X, Y)
surf = ax.plot_surface(X,Y,dem, rstride=1, cstride=1, cmap=plt.cm.RdYlBu_r, vmin=0, vmax=4000, linewidth=0, antialiased=True)
ax.set_zlim(0, 60000) # to make it stand out less
ax.view_init(60,-105)
fig.colorbar(surf, shrink=0.4, aspect=20)
关于python - 在 3D matplotlib 中显示地理引用 DEM 表面,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/17713050/