2017-04-12 65 views
0

我已經在底圖中繪製了帶有形狀文件的地圖,在此基礎上,我想在DEM中添加DEM數據(.tif)地圖。我已經在SRTM網站上下載了一塊.tif數據(經度和緯度真的在shapefile範圍內),但是最後當我運行該程序時,它沒有顯示特定地圖中的柵格數據如何在gdal中用python添加DEM數據(.tif)

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() 

回答

0

如果您在gis.stackexchange上提問,您可能有更簡潔的答案。

無論如何,這對我來說似乎是一個投影錯誤。我希望糾正,但STRM tiffs是GeoTiffs,所以他們在標題中有投影信息。

看起來您正在使用橫​​軸墨卡託投影機,而SRTM數據使用的是地理座標系。

嘗試重新投影的TIFF,我覺得gdal.Warp應該這樣做:

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