2017-04-12 1 views
0

Ich habe eine Karte mit Formdatei in Grundkarte geplottet, auf dieser Grundlage möchte ich die DEM - Daten (.tif) hinzufügen Karte. Ich habe einen Teil der .tif-Daten (Länge und Breite war wirklich innerhalb des Bereichs der Shapefile) in der SRTM-Website herunterladen, aber schließlich, wenn ich das Programm ausführen, zeigte es nicht die Raster-Daten in der spezifischen Position der KarteWie man DEM - Daten (.tif) in einer Formdatei in gdal mit python hinzufügt

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

Antwort

0

Sie können eine knappe Antwort haben, wenn Sie dies bei gis.stackexchange nachgefragt haben.

Wie auch immer, das scheint mir ein Projektionsfehler zu sein. Ich stehe zu korrigieren, aber die STRM-Tiffs sind GeoTiffs, so dass sie die Projektionsinformationen in der Kopfzeile haben.

Es sieht so aus, als ob Sie Transverse Mercator verwenden, während die SRTM-Daten ein geographisches Koordinatensystem verwenden.

Versuchen Sie TIFF Reprojizierens, denke ich gdal.Warp es tun sollten:

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