我有一个数据可以在美国地图上显示.我不希望在美国以外的土地上有任何颜色(在西南部的墨西哥和东北部的加拿大).如何从contourf中屏蔽这些区域?请注意,是否绘制了这些状态边界并不重要.
代码是:
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap, cm, maskoceans
from scipy.interpolate import griddata
np.random.seed(77)
lllat = 24.396308
lllon = -124.848974
urlat = 49.384358
urlon = -66.885444
m = Basemap(llcrnrlat=lllat,
urcrnrlat=urlat,
llcrnrlon=lllon,
urcrnrlon=urlon,
resolution='i', projection='cyl')
m.drawcountries(linewidth=1.0)
m.drawstates(linewidth=1.0, color='lightgray')
m.drawlsmask(land_color='gray',ocean_color="#b0c4de", lakes=True)
#create 100 random latitudes
lats = np.random.randint(low=lllat-1, high=urlat+1, size=1000) + np.random.ranf(size=1000)
#create 100 random longitudes
lons = np.random.randint(low=lllon-1, high=urlon+1, size=1000) + np.random.ranf(size=1000)
#create 100 random values/probabilities
probabilities = np.random.random(size=1000)
#now use meshgrid and contourf to visualize it
mlon, mlat = m(*(lons, lats))
# grid data
numcols, numrows = 1000, 1000
xi = np.linspace(mlon.min(), mlon.max(), numcols)
yi = np.linspace(mlat.min(), mlat.max(), numrows)
xi, yi = np.meshgrid(xi, yi)
# interpolate
x, y, z = mlon, mlat, probabilities
zi = griddata((mlon, mlat), probabilities, (xi, yi), method='nearest', rescale=False)
data = maskoceans(xi, yi, zi)
con = m.contourf(xi, yi, data, cmap=plt.get_cmap('YlOrRd'))
cbar = m.colorbar(con,location='right',pad="3%")
plt.show()
我尝试从美国不存在的lats和lons中删除点,但仍然在边界外的某些部分使用contourf颜色.
我可以用海洋面具掩盖水域,但不能掩盖墨西哥和加拿大的那些点.
注意:有一个解决方案可以循环遍历所有网格点,并为所有那些不在美国的点设置zi为nan,但由于我的实际数据大小对我来说不是解决方案,因此计算成本太高.
最佳答案 我最终在轮廓上再次绘制墨西哥和加拿大多边形.已下载的国家/地区从
here开始编码文件.代码如下:
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap, cm, maskoceans
from scipy.interpolate import griddata
from matplotlib.patches import Polygon as MplPolygon
import shapefile
import pdb
np.random.seed(77)
lllat = 24.396308
lllon = -124.848974
urlat = 49.384358
urlon = -66.885444
m = Basemap(llcrnrlat=lllat,
urcrnrlat=urlat,
llcrnrlon=lllon,
urcrnrlon=urlon,
resolution='i', projection='cyl')
m.drawcountries(linewidth=1.0)
m.drawstates(linewidth=1, color='lightgray')
m.drawcoastlines()
m.drawlsmask(land_color='gray',ocean_color="#b0c4de", lakes=True)
#create 100 random latitudes
lats = np.random.randint(low=lllat-1, high=urlat+1, size=1000) + np.random.ranf(size=1000)
#create 100 random longitudes
lons = np.random.randint(low=lllon-1, high=urlon+1, size=1000) + np.random.ranf(size=1000)
#create 100 random values/probabilities
probabilities = np.random.random(size=1000)
#now use meshgrid and contourf to visualize it
mlon, mlat = m(*(lons, lats))
# grid data
numcols, numrows = 1000, 1000
xi = np.linspace(mlon.min(), mlon.max(), numcols)
yi = np.linspace(mlat.min(), mlat.max(), numrows)
xi, yi = np.meshgrid(xi, yi)
# interpolate
x, y, z = mlon, mlat, probabilities
zi = griddata((mlon, mlat), probabilities, (xi, yi), method='nearest', rescale=False)
data = maskoceans(xi, yi, zi)
con = m.contourf(xi, yi, data, cmap=plt.get_cmap('YlOrRd'))
cbar = m.colorbar(con,location='right',pad="3%")
world_shp_info = m.readshapefile('./data/CNTR_2014_10M_SH/Data/CNTR_RG_10M_2014','world',drawbounds=False)
ax = plt.gca()
for shapedict,state in zip(m.world_info, m.world):
if shapedict['CNTR_ID'] not in ['CA', 'MX']: continue
poly = MplPolygon(state,facecolor='gray',edgecolor='gray')
ax.add_patch(poly)
plt.show()
我认为应该有一种更聪明的方法.