0

我正在尝试创建显示高度的马拉维地图。像这样的东西,但当然是马拉维:

在此处输入图像描述

我从这里下载了一些海拔数据:http ://research.jisao.washington.edu/data_sets/elevation/

这是我创建多维数据集后该数据的打印:

meters, from 5-min data / (unknown) (time: 1; latitude: 360; longitude: 720)
     Dimension coordinates:
          time                           x            -               -
          latitude                       -            x               -
          longitude                      -            -               x
     Attributes:
          history: 
Elevations calculated from the TBASE 5-minute
latitude-longitude resolution...
          invalid_units: meters, from 5-min data

我开始导入我的数据,形成一个立方体,删除额外的变量(时间和历史)并将我的数据限制在马拉维的纬度和经度。

import matplotlib.pyplot as plt
import matplotlib.cm as mpl_cm
import numpy as np
import iris
import cartopy
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
import iris.analysis.cartography


def main():

    #bring in altitude data
    Elev = '/exports/csce/datastore/geos/users/s0899345/Climate_Modelling/Actual_Data/elev.0.5-deg.nc'

    Elev= iris.load_cube(Elev)

    #remove variable for time
    del Elev.attributes['history']
    Elev = Elev.collapsed('time', iris.analysis.MEAN)

    Malawi = iris.Constraint(longitude=lambda v: 32.0 <= v <= 36., latitude=lambda v: -17. <= v <= -8.)      
    Elev = Elev.extract(Malawi)

    print 'Elevation'
    print Elev.data
    print 'latitude'
    print Elev.coord('latitude')
    print 'longitude'
    print Elev.coord('longitude')

这很好用,输出如下:

Elevation
[[  978.  1000.  1408.  1324.  1080.  1370.  1857.  1584.]
 [ 1297.  1193.  1452.  1611.  1354.  1480.  1350.   627.]
 [ 1418.  1490.  1625.  1486.  1977.  1802.  1226.   482.]
 [ 1336.  1326.  1405.   728.  1105.  1559.  1139.   789.]
 [ 1368.  1301.  1463.  1389.   671.   942.   947.   970.]
 [ 1279.  1116.  1323.  1587.   839.  1014.  1071.  1003.]
 [ 1096.   969.  1179.  1246.   855.   979.   927.   638.]
 [  911.   982.  1235.  1324.   681.   813.   814.   707.]
 [  749.   957.  1220.  1198.   613.   688.   832.   858.]
 [  707.  1049.  1037.   907.   624.   771.  1142.  1104.]
 [  836.  1044.  1124.  1120.   682.   711.  1126.   922.]
 [ 1050.  1204.  1199.  1161.   777.   569.   999.   828.]
 [ 1006.   869.  1183.  1230.  1354.   616.   762.   784.]
 [  838.   607.   883.  1181.  1174.   927.   591.   856.]
 [  561.   402.   626.   775.  1053.   726.   828.   733.]
 [  370.   388.   363.   422.   508.   471.   906.  1104.]
 [  504.   326.   298.   208.   246.   160.   458.   682.]
 [  658.   512.   334.   309.   156.   162.   123.   340.]]
latitude
DimCoord(array([ -8.25,  -8.75,  -9.25,  -9.75, -10.25, -10.75, -11.25, -11.75,
       -12.25, -12.75, -13.25, -13.75, -14.25, -14.75, -15.25, -15.75,
       -16.25, -16.75], dtype=float32), standard_name='latitude', units=Unit('degrees'), var_name='lat', attributes={'title': 'Latitude'})
longitude
DimCoord(array([ 32.25,  32.75,  33.25,  33.75,  34.25,  34.75,  35.25,  35.75], dtype=float32), standard_name='longitude', units=Unit('degrees'), var_name='lon', attributes={'title': 'Longitude'})

但是,当我尝试绘制它时,它不起作用......这就是我所做的:

#plot map with physical features 
ax = plt.axes(projection=cartopy.crs.PlateCarree())
ax.add_feature(cartopy.feature.COASTLINE)   
ax.add_feature(cartopy.feature.BORDERS)
ax.add_feature(cartopy.feature.LAKES, alpha=0.5)
ax.add_feature(cartopy.feature.RIVERS)
#plot altitude data
plot=ax.plot(Elev, cmap=mpl_cm.get_cmap('YlGn'), levels=np.arange(0,2000,150), extend='both') 
#add colour bar index and a label
plt.colorbar(plot, label='meters above sea level')
#set map boundary
ax.set_extent([32., 36., -8, -17]) 
#set axis tick marks
ax.set_xticks([33, 34, 35]) 
ax.set_yticks([-10, -12, -14, -16]) 
lon_formatter = LongitudeFormatter(zero_direction_label=True)
lat_formatter = LatitudeFormatter()
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)
#save the image of the graph and include full legend
plt.savefig('Map_data_boundary', bbox_inches='tight')
plt.show() 

我得到的错误是'Attribute Error: Unknown property type cmap'和下面的整个世界地图......

在此处输入图像描述

有任何想法吗?

4

3 回答 3

4

我会像你一样准备数据,除了删除time我将使用的维度iris.util.squeeze,这会删除任何长度为 1 的维度。

import iris

elev = iris.load_cube('elev.0.5-deg.nc')
elev = iris.util.squeeze(elev)
malawi = iris.Constraint(longitude=lambda v: 32.0 <= v <= 36.,
                         latitude=lambda v: -17. <= v <= -8.)      
elev = elev.extract(malawi)

正如@ImportanceOfBeingErnest 所说,您需要等高线图。当不确定要使用什么绘图功能时,我建议浏览 matplotlib以查找与您想要生成的内容相似的内容。单击图像,它会显示代码。

因此,要制作等高线图,您可以使用该函数,但您必须以数组matplotlib.pyplot.contourf的形式从立方体中获取相关数据:numpy

import matplotlib.pyplot as plt
import matplotlib.cm as mpl_cm
import numpy as np
import cartopy

cmap = mpl_cm.get_cmap('YlGn')
levels = np.arange(0,2000,150)
extend = 'max'

ax = plt.axes(projection=cartopy.crs.PlateCarree())
plt.contourf(elev.coord('longitude').points, elev.coord('latitude').points, 
             elev.data, cmap=cmap, levels=levels, extend=extend)

matplotlib.pyplot 版本

但是,iris以. 这会自动设置一个带有正确投影的坐标区实例,并将数据从立方体传递到. 所以最后两行可以简单地变成:maplotlib.pyplotiris.plotmatplotlib.pyplot

import iris.plot as iplt
iplt.contourf(elev, cmap=cmap, levels=levels, extend=extend)

iris.plot 版本

还有iris.quickplot,它与 基本相同iris.plot,只是它会在适当的地方自动添加颜色栏和标签:

import iris.quickplot as qplt
qplt.contourf(elev, cmap=cmap, levels=levels, extend=extend)

iris.quickplot 版本

绘制完成后,您可以获取轴实例并添加其他项目(为此我只是复制了您的代码):

from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter

qplt.contourf(elev, cmap=cmap, levels=levels, extend=extend)
ax = plt.gca()

ax.add_feature(cartopy.feature.COASTLINE)   
ax.add_feature(cartopy.feature.BORDERS)
ax.add_feature(cartopy.feature.LAKES, alpha=0.5)
ax.add_feature(cartopy.feature.RIVERS)

ax.set_xticks([33, 34, 35]) 
ax.set_yticks([-10, -12, -14, -16]) 
lon_formatter = LongitudeFormatter(zero_direction_label=True)
lat_formatter = LatitudeFormatter()
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)

iris.quickplot 与添加

于 2017-12-21T19:34:16.717 回答
1

看来您想要等高线图之类的东西。所以而不是

plot = ax.plot(...)

你可能想使用

plot = ax.contourf(...)

很可能您还想将纬度和经度作为参数提供给contourf,

plot = ax.contourf(longitude, latitude, Elev, ...)
于 2017-12-19T09:43:21.430 回答
0

您可以尝试添加以下内容:

import matplotlib.colors as colors

color = plt.get_cmap('YlGn')  # and change cmap=mpl_cm.get_cmap('YlGn') to  cmap=color 

并尝试更新您的 matplotlib:

pip install --upgrade matplotlib

编辑

color = plt.get_cmap('YlGn')  # and change cmap=mpl_cm.get_cmap('YlGn') to  cmap=color 
于 2017-12-19T09:17:32.927 回答