如何使用python在gdal的形状文件中添加DEM数据(.tif)

Question

我在底图中绘制了一张带有形状文件的地图，在此基础上，我想在地图中添加 DEM 数据(.tif(。 我在SRTM网站上下载了一段.tif数据(经度和纬度确实在shapefile的范围内(，但是最后当我运行该程序时，它未在地图的特定位置显示栅格数据

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.patches import Polygon
from osgeo import gdal
from numpy import linspace
from numpy import meshgrid
fig = plt.figure(figsize=(8,6), dpi=80)
ax1 = fig.add_axes([0.1,0.1,1.0,1.0])
map = Basemap(llcrnrlon=114.7, 
          llcrnrlat=29.3, 
          urcrnrlon=120, 
          urcrnrlat=34.6,
         resolution='h', projection='tmerc', lat_0 =31.5,lon_0=116.5,ax=ax1)
shp_info =  map.readshapefile("bou2_4l",'state',color='k',linewidth='1',drawbounds=True)
ds = gdal.Open("srtm_60_06.tif")
data = ds.ReadAsArray()
loc3=[30,35]
lat3=[115,120]
x1,y1 = map(loc3,lat3)
x = linspace(x1[0],x1[1], data.shape[1])
y = linspace(y1[0],y1[1], data.shape[0])
xx, yy = meshgrid(x, y)
map.pcolormesh(xx, yy, data)
plt.show()

Answer 1

如果你在gis.stackexchange上问这个问题，你可能会有一个更简洁的答案。

无论如何，这对我来说似乎是一个投影错误。我愿意纠正，但 STRM tiffs 是 GeoTiffs，因此它们在标题中有投影信息。

看起来您使用的是横轴墨卡托，而 SRTM 数据使用的是地理坐标系。

尝试重新投影你的TIFF，我想是gdal。Warp应该这样做：

gdal.Warp(output_raster,input_raster,dstSRS='EPSG:your_basemap_projection')