Cartopy教学
时间:2020-05-28
本文章向大家介绍Cartopy教学,主要包括Cartopy教学使用实例、应用技巧、基本知识点总结和需要注意事项,具有一定的参考价值,需要的朋友可以参考一下。
Cartopy 教学
本章全部转载于‘微信公众号”气象科研猫“。
主要介绍了Cartopy的13个示例,如下。
1
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import cartopy.crs as ccrs
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
dlon, dlat = 60, 30 #设置步长
xticks = np.arange(0, 360.1, dlon) #设置绘图范围
yticks = np.arange(-90, 90.1, dlat)
fig = plt.figure(figsize=(6,5)) #设置画布大小
ax = fig.add_subplot(1,1,1, projection=ccrs.PlateCarree(central_longitude=180)) #添加图幅本体
ax.coastlines() #海岸线
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=False,linewidth=1, linestyle=':', color='k', alpha=0.8) #格网
gl.xlocator = mticker.FixedLocator(xticks) #格网设置x轴刻度
gl.ylocator = mticker.FixedLocator(yticks)
ax.set_xticks(xticks, crs=ccrs.PlateCarree()) #图幅设置坐标轴刻度
ax.set_yticks(yticks, crs=ccrs.PlateCarree())
ax.xaxis.set_major_formatter(LongitudeFormatter(zero_direction_label=True)) #设置坐标轴刻度标签格式
ax.yaxis.set_major_formatter(LatitudeFormatter())
fig_fname = "全球地图.png"
#plt.savefig(fig_fname, dpi=500, bbox_inches='tight')
plt.savefig('F:/Rpython/lp11/plot2.png',dpi=1200)
plt.show()
2
import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
plt.clf()
fig = plt.figure(dpi=200)
# set the projection to PlateCarree
proj = ccrs.PlateCarree()
ax = fig.add_subplot(1, 1, 1, projection = proj)
ax.set_global()
# set the gridlines to dashed line and the transparency to 0.7
gl = ax.gridlines(ylocs=np.arange(-90,90+30,30),xlocs=np.arange(-180,180+60,60),draw_labels=True,linestyle='--',alpha=0.7)
gl.xlabels_top = False
gl.ylabels_right = False
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
# set background image to 50-natural-earth-1-downsampled.png
ax.stock_img()
ax.add_feature(cfeature.COASTLINE)
ax.add_feature(cfeature.RIVERS)
ax.add_feature(cfeature.LAKES)
plt.savefig('F:/Rpython/lp11/plot11.png',dpi=600)
plt.show()
3
import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
# Define the interested area
region = [70,140,15,60] # [west,east,south,north]
plt.clf()
fig = plt.figure(dpi=200)
# set the projection to PlateCarree
proj = ccrs.PlateCarree()
ax = fig.add_subplot(1, 1, 1,projection = proj)
# set the range of the interested area
ax.set_extent(region,crs = ccrs.PlateCarree())
# set the background image to the 1cm:500km raster
fname = 'F:/Rpython/lp3/HYP_50M_SR_W/HYP_50M_SR_W.tif'
ax.imshow(plt.imread(fname), origin='upper',transform=ccrs.PlateCarree(),extent=[-180, 180, -90, 90])
# set the gridlines to dashed line and the transparency to 0.7
gl = ax.gridlines(ylocs=np.arange(region[2],region[3]+5,5),xlocs=np.arange(region[0],region[1]+10,10),draw_labels=True,linestyle='--',alpha=0.7)
gl.xlabels_top = False
gl.ylabels_right = False
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
railroads = cfeature.NaturalEarthFeature('cultural','railroads','10m',facecolor='None')
urban_areas = cfeature.NaturalEarthFeature('cultural','urban_areas','10m')
# add borders, coastline, rivers, lakes, railroads,and urban areas
ax.add_feature(cfeature.BORDERS.with_scale('50m')) # '50m' for mediate resolution and '10m' for high resolution
ax.add_feature(cfeature.COASTLINE.with_scale('50m'))
ax.add_feature(cfeature.RIVERS) # low resolution
ax.add_feature(cfeature.LAKES) # low resolution
ax.add_feature(railroads,edgecolor='gray')
ax.add_feature(urban_areas, edgecolor='m')
plt.savefig('F:/Rpython/lp11/plot12.png',dpi=600)
plt.show()
4
import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
# Read VTEC data from a file
with open('F:/Rpython/lp5/1st-tecmap.dat') as src:
data = [line for line in src if not line.endswith('LAT/LON1/LON2/DLON/H\n')]
tec = np.fromstring(''.join(data), dtype=float, sep=' ')
# Generate a geodetic grid
nlats, nlons = 71, 73
lats = np.linspace(-87.5, 87.5, nlats)
lons = np.linspace(-180, 180, nlons)
lons, lats = np.meshgrid(lons, lats)
tec.shape = nlats, nlons
fig = plt.figure(figsize=(8, 3))
# Set projection
ax = plt.axes(projection=ccrs.PlateCarree())
# Plot contour
cp = plt.contourf(lons, lats, tec, transform=ccrs.PlateCarree(), cmap='jet')
ax.coastlines()
# Add a color bar
plt.colorbar(cp)
# Set figure extent & ticks
ax.set_extent([-180, 180, -87.5, 87.5])
ax.set_xticks([-180, -150, -120, -90, -60, -30, 0, 30, 60, 90, 120, 150, 180])
ax.set_yticks([-80, -60, -40, -20, 0, 20, 40, 60, 80])
plt.savefig('F:/Rpython/lp11/plot15.png',dpi=600)
plt.show()
5
import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeat
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
from cartopy.io.shapereader import Reader
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import warnings
import os
from cartopy import config
from pylab import imread
#from cv2 import imread, imshow
def create_map():
extent = [70, 140, 0, 60]
shp_path = 'F:/Rpython/lp11/'
# --创建画图空间
proj = ccrs.PlateCarree() # 创建坐标系
fig = plt.figure(figsize=(6, 8), dpi=350) # 创建页面
ax = fig.subplots(1, 1, subplot_kw={'projection': proj}) # 创建子图
# --设置地图属性
reader = Reader(shp_path + 'shp/Province_9.shp')
provinces = cfeat.ShapelyFeature(reader.geometries(), proj, edgecolor='k', facecolor='none')
ax.add_feature(provinces, linewidth=1)
ax.set_extent(extent, crs=proj)
# --增加低分辨率地形图(cartopy自带)
#ax.stock_img() # 增加地形图
# --增加高分辨率地形图(需自行下载)
#https://www.naturalearthdata.com/downloads/10m-raster-data/10m-natural-earth-1/
#fname = os.path.join(config["repo_data_dir"], 'raster', 'natural_earth', 'NE1_HR_LC.tif')
fname='F:/Rpython/lp11/shp/NE1_LR_LC.tif'
ax.imshow(imread(fname), origin='upper', transform=proj, extent=[-180, 180, -90, 90])
# --设置网格点属性
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True, linewidth=1.2, color='k', alpha=0.5, linestyle='--')
gl.xlabels_top = False # 关闭顶端的经纬度标签
gl.ylabels_right = False # 关闭右侧的经纬度标签
gl.xformatter = LONGITUDE_FORMATTER # x轴设为经度的格式
gl.yformatter = LATITUDE_FORMATTER # y轴设为纬度的格式
gl.xlocator = mticker.FixedLocator(np.arange(extent[0], extent[1]+10, 10))
gl.ylocator = mticker.FixedLocator(np.arange(extent[2], extent[3]+10, 10))
return ax
if __name__ == '__main__':
warnings.filterwarnings('ignore')
ax = create_map()
plt.savefig('F:/Rpython/lp11/plot17.png',dpi=600)
plt.show()
6
import warnings
warnings.filterwarnings('ignore')
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np
fig = plt.figure(figsize=(6, 4))
plt.subplots_adjust(left=0.05, right=0.95, top=0.90, bottom=0.05, wspace=0.15, hspace=0.05)
ax = plt.subplot(111)
m = Basemap(resolution='i', projection='merc', llcrnrlat=10.0, urcrnrlat=55.0, llcrnrlon=60., urcrnrlon=140.0)
m.drawcountries(linewidth=0.5)
m.drawcoastlines(linewidth=0.5)
m.drawparallels(np.arange(10., 55., 10.), labels=[1, 0, 0, 0], color='black',linewidth=0.2,dashes=(None, None)) # draw parallels,dashes=[1,0],
m.drawmeridians(np.arange(60., 140., 10.), labels=[0, 0, 0, 1], color='black',linewidth=0.2,dashes=(None, None)) # draw meridians ,dashes=[1,0]
x, y = m(116.4204, 40.21244) # Bejing
x2, y2 = m(125.27538, 43.83453) # Changchun
plt.annotate('Beijing', xy=(x, y), xycoords='data',
# xytext=(x2, y2), textcoords='offset points',
xytext=(x2, y2), textcoords='data',
color='r',arrowprops=dict(arrowstyle="fancy", color='g'))
plt.savefig('F:/Rpython/lp11/plot21.png',dpi=600)
plt.show()
7
from netCDF4 import Dataset
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import cmaps
meteo_file = r'F:/Rpython/lp3/hls/ECMWF.nc'
ds = Dataset(meteo_file, mode='r')
# 获取每个变量的值
lons = ds.variables['longitude'][:]
lats = ds.variables['latitude'][:]
# surface_air_pressure
sp = ds.variables['sp'][:]
sp_units = ds.variables['sp'].units
scale_factor = ds.variables['sp'].scale_factor
add_offset = ds.variables['sp'].add_offset
sp = scale_factor * sp + add_offset
# 2 metre temperature
t2m = ds.variables['t2m'][:]
t2m_units = ds.variables['t2m'].units
scale_factor = ds.variables['t2m'].scale_factor
add_offset = ds.variables['t2m'].add_offset
t2m = scale_factor * t2m + add_offset
# Total column ozone
tco3 = ds.variables['tco3'][:]
tco3_units = ds.variables['tco3'].units
scale_factor = ds.variables['tco3'].scale_factor
add_offset = ds.variables['tco3'].add_offset
tco3 = scale_factor * tco3 + add_offset
# 经纬度平均值
lon_0 = lons.mean()
lat_0 = lats.mean()
# 画图大小设置
fig = plt.figure(figsize=(16, 9))
plt.rc('font', size=15, weight='bold')
ax = fig.add_subplot(111)
m = Basemap(lat_0=lat_0, lon_0=lon_0)
lon, lat = np.meshgrid(lons, lats)
xi, yi = m(lon, lat)
# 这里数据时间是UTC 00:00,2018年1月的日平均数据,只展示1月1号的数据
sp_01 = sp[0:1, :, :]
t2m_01 = t2m[0:1, :, :]
tco3_01 = tco3[0:1, :, :]
levels = m.pcolor(xi, yi, np.squeeze(tco3_01), cmap=cmaps.GMT_panoply)
# 添加格网与绘制经纬线
m.drawparallels(np.arange(-90., 91., 20.), labels=[1, 0, 0, 0], fontsize=15)
m.drawmeridians(np.arange(-180., 181., 40.), labels=[0, 0, 0, 1], fontsize=15)
# 添加海岸线,省州边界以及国家行政边界
m.drawcoastlines()
m.drawstates()
m.drawcountries()
# 添加colorbar
cbar = m.colorbar(levels, location='bottom', pad="10%")
cbar.set_label(tco3_units, fontsize=15, weight='bold')
# 添加图的标题
plt.title('Total column ozone')
#----------------------------------------------------------------------------------
plt.savefig('F:/Rpython/lp7/plot22.png',dpi=1200)
plt.show()
ds.close()
8
import numpy as np
import matplotlib.pyplot as plt
import maskout1
from mpl_toolkits.basemap import Basemap
from scipy.interpolate import griddata
from matplotlib.font_manager import FontProperties
font = FontProperties(fname=r'F:/Rpython/lp3/guangxi/SimHei/SimHei.ttf')
#创建三个列表用于存放点的数据
x=[]
y=[]
h=[]
#with open(r'F:/Rpython/lp3/guangxi/广西站点修正.csv','r') as f:
with open(r'F:/Rpython/lp3/guangxi/广西海拔高度.csv','r') as f:
skip1=f.readline()
datas=f.readlines()
for line in datas:
p=line[:-1] #去掉每行最后的换行符
a=p.split(',') #因为每行每个数据之间用","隔开,去掉逗号。
x.append(float(a[2]))
y.append(float(a[3]))
h.append(float(a[4]))
points=[]
for i in range(len(x)):
point=[]
point.append(x[i])
point.append(y[i])
points.append(point)
points=np.array(points)
xi=np.linspace(min(x),max(x),200)
yi=np.linspace(min(y),max(y),200)
xi,yi=np.meshgrid(xi,yi)#网格化
#scipy.interpolate.griddata(points, values, xi, method=‘linear’, fill_value=nan, rescale=Fale)
#zi=griddata(points,h,(xi,yi),method='cubic')
zi=griddata(points,h,(xi,yi),method='linear')
fig = plt.figure(figsize=(8,6))
ax = fig.add_subplot(111)
# [west,east,south,north]
map = Basemap(
# llcrnrlon=72,
# llcrnrlat=0,
# urcrnrlon=136,
# urcrnrlat=55,
llcrnrlon=104,
llcrnrlat=21,
urcrnrlon=112.5,
urcrnrlat=27,
resolution = None,
projection = 'cyl')
#map.readshapefile('e:\python\map\world0','whatevername',color='black',linewidth=1.2)
#map.readshapefile('e:\python\map\china1','whatevername',color='black',linewidth=0.3)
map.readshapefile(r'F:/Rpython/lp3/guangxi/map/guangxi','whatevername',color='black',linewidth=1.0)
minval,maxval=int(np.amin(h)),int(np.amax(h))
print(minval,maxval)
#c2 = plt.contourf(xi,yi,zi,range(minval,maxval,1),cmap= 'jet')
#c2 = plt.contourf(xi,yi,zi,18,cmap= 'jet')
c2 = plt.contourf(xi,yi,zi,range(0,2000,10),cmap= 'jet')
clip=maskout1.shp2clip(c2,ax,r'F:/Rpython/lp3/guangxi/map/world0',r'F:/Rpython/lp3/guangxi/map/guangxi')
#plt.rcParams['font.family'] = 'Times New Roman' #设置colorbar的label字体
plt.rcParams['font.family'] = 'DejaVu Sans' #设置colorbar的label字体
plt.rcParams['font.size'] = 14 #设置colorbar的label字体大小
bar=map.colorbar(c2,label='height / m')
bar.ax.tick_params(labelsize=14) #设置colorbar刻度字体大小。
# 以下标经纬度
map.drawparallels(np.arange(-90.,90.,2),labels=[1,0,0,0],fontsize=14,linewidth=0.1)
map.drawmeridians(np.arange(-180.,180.,2),labels=[0,0,0,1],fontsize=14,linewidth=0.1)
#u以下写汉字标题
plt.title(u'广西海拔高度',color='red',fontproperties=font,fontsize=14)
#plt.xlabel(u'经度',color='blue',fontproperties=font,fontsize=14)
#plt.ylabel(u'纬度',color='blue',fontproperties=font,fontsize=14)
plt.savefig(r'F:/Rpython/lp3/guangxi/gxfig.png',format='png',dpi=1200)
plt.show()
9
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import cartopy.crs as ccrs
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
from cartopy.util import add_cyclic_point
#from palettable.colorbrewer.diverging import RdBu_11_r
import netCDF4 as nc
sns.set_context('talk', font_scale=1.2) # 设置图形属性
# read NetCDF file
fn = 'F:/Rpython/lp11/shp/air.sig995.2012.nc'
data = nc.Dataset(fn, 'r') # 默认为读文件,此处 'r' 可省略
# 读取相关变量
lat = data.variables['lat'][:].data
lon = data.variables['lon'][:].data
time = data.variables['time'][:].data
air = data.variables['air'][:].data
# 添加数据循环,以防止在0和360时出现白色条状
# 添加数据循环和不添加数据循环的效果见后文两张图
cycle_air, cycle_lon = add_cyclic_point(air, coord=lon)
cycle_LON, cycle_LAT = np.meshgrid(cycle_lon, lat)
projection = ccrs.PlateCarree() # 设置投影
fig, ax = plt.subplots(figsize=(16, 9), subplot_kw=dict(projection=projection))
con = ax.contourf(cycle_LON, cycle_LAT, cycle_air[0, ...],np.arange(220, 321), cmap='Spectral_r')
ax.coastlines(linewidth=1.5) # 添加海岸线
ax.set_xticks(np.arange(-180, 181, 60), crs=projection)
ax.set_yticks(np.arange(-90, 91, 30), crs=projection)
# 设置 ticklabels 格式
lon_formatter = LongitudeFormatter(number_format='.0f',degree_symbol='',dateline_direction_label=True)
lat_formatter = LatitudeFormatter(number_format='.0f',degree_symbol='')
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)
# shrink 控制 colorbar 长度,pad 控制colorbar和图的距离
cb = fig.colorbar(con, shrink=0.73, pad=0.02)
# 调整 ticklabels
cb.set_ticks(np.arange(220, 321, 20))
cb.set_ticklabels(np.arange(220, 321, 20))
cb.ax.tick_params(direction='in', length=5) # 控制 colorbar tick
# 保存文件, dpi用于设置图形分辨率, bbox_inches 尽量减小图形的白色区域
#fig.savefig(F:/Rpython/lp11/plot24.png', dpi=600, bbox_inches='tight')
plt.savefig('F:/Rpython/lp11/plot24.png',dpi=600)
plt.show()
10
import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeat
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
from cartopy.io.shapereader import Reader
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
extent = [70, 140, 0, 60]
Chinese_land_territory_path = 'F:/Rpython/lp11/shp/China_land_territory'
Chinese_10dash_line_path = 'F:/Rpython/lp11/shp/China_10-dash_line'
world_countries_path = 'F:/Rpython/lp11/shp/world_countries'
# 创建坐标系
prj = ccrs.PlateCarree()
fig = plt.figure(figsize=(10, 8), dpi=350)
ax = fig.subplots(1, 1, subplot_kw={'projection': prj})
ax.set_extent(extent, crs=prj)
# 绘制中国陆地领土边界
# reader = Reader(shp_path + 'shp/Province_9.shp')
# provinces = cfeat.ShapelyFeature(reader.geometries(), proj, edgecolor='k', facecolor='none')
# ax.add_feature(provinces, linewidth=1)
Chinese_land_territory2 = Reader('F:/Rpython/lp11/shp/China_land_territory/China_land_territory.shp')
Chinese_land_territory = cfeat.ShapelyFeature(Chinese_land_territory2.geometries(),prj, edgecolor='k',facecolor='none')
ax.add_feature(Chinese_land_territory, linewidth=1)
# 绘制中国十段线
Chinese_10dash_line2 = Reader('F:/Rpython/lp11/shp/China_10-dash_line/China_10-dash_line.shp').geometries()
Chinese_10dash_line = cfeat.ShapelyFeature(Chinese_10dash_line2,prj, edgecolor='r')
ax.add_feature(Chinese_10dash_line, linewidth=2)
# 绘制世界各国领土边界
world_countries2 = Reader('F:/Rpython/lp11/shp/world_countries/world_countries.shp').geometries()
world_countries = cfeat.ShapelyFeature(world_countries2,prj, edgecolor='k', alpha=0.3,facecolor='none')
ax.add_feature(world_countries, linewidth=0.5)
# 添加自带地形数据
ax.stock_img()
# 绘制网格点
gl = ax.gridlines(crs=prj, draw_labels=True, linewidth=1.2, color='k',alpha=0.5, linestyle='--')
gl.xlabels_top = False # 关闭顶端的经纬度标签
gl.ylabels_right = False # 关闭右侧的经纬度标签
gl.xformatter = LONGITUDE_FORMATTER # x轴设为经度的格式
gl.yformatter = LATITUDE_FORMATTER # y轴设为纬度的格式
gl.xlocator = mticker.FixedLocator(np.arange(extent[0], extent[1]+10, 10))
gl.ylocator = mticker.FixedLocator(np.arange(extent[2], extent[3]+10, 10))
plt.savefig('F:/Rpython/lp11/plot25.png',dpi=600)
plt.show()
11
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
data = pd.read_csv('F:/Rpython/lp11/shp/1950-2018_all_tornadoes.csv')
fig, ax = plt.subplots(figsize=(12,8))
m = Basemap(projection='mill',llcrnrlat=25,urcrnrlat=50,llcrnrlon=-125,urcrnrlon=-65,rsphere=6371200.,resolution='l',area_thresh=10000)
m.drawcoastlines()
m.drawstates()
parallels = np.arange(25,51,5.)
m.drawparallels(parallels,labels=[1,0,0,0],fontsize=10)
meridians = np.arange(-125.,-64.,10.)
m.drawmeridians(meridians,labels=[0,0,0,1],fontsize=10)
# 这里要转换坐标,否则可能无法显示图形
# 由于 data 为 DataFrame,所以要直接提取值,因为 basemap 在转换
# 坐标时只能是列表,元组或标量
x, y = m(data.slon.values, data.slat.values)
hex = m.hexbin(x, y, bins = 150, cmap = plt.get_cmap('Reds'))
fig.colorbar(hex, ax = ax, shrink = 0.6)
plt.savefig('F:/Rpython/lp11/plot28.png',dpi=600)
plt.show()
12
import numpy as np
import xarray as xr
import cartopy.crs as ccrs
import cartopy.feature as cfeat
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
import matplotlib.pyplot as plt
# 数据读取及时间平均处理
ds = xr.open_dataset('F:/Rpython/lp11/shp/EC-Interim_monthly_2018.nc')
temp = (ds['t2m'] - 273.15).mean(dim='time') #把温度转换为℃并对其时间纬求平均
temp.attrs['units'] = 'deg C' #温度单位转换为℃
# 创建画图空间
proj = ccrs.PlateCarree() #创建投影
fig = plt.figure(figsize=(9,6)) #创建页面
ax = fig.subplots(1, 1, subplot_kw={'projection': proj}) #子图
# 设置地图属性:加载国界、海岸线、河流、湖泊
ax.add_feature(cfeat.BORDERS.with_scale('50m'), linewidth=0.8, zorder=1)
ax.add_feature(cfeat.COASTLINE.with_scale('50m'), linewidth=0.6, zorder=1)
ax.add_feature(cfeat.RIVERS.with_scale('50m'), zorder=1)
ax.add_feature(cfeat.LAKES.with_scale('50m'), zorder=1)
# 设置网格点属性
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True,linewidth=1.2, color='k', alpha=0.5, linestyle='--')
gl.xlabels_top = False #关闭顶端标签
gl.ylabels_right = False #关闭右侧标签
gl.xformatter = LONGITUDE_FORMATTER #x轴设为经度格式
gl.yformatter = LATITUDE_FORMATTER #y轴设为纬度格式
# 设置colorbar
cbar_kwargs = {'orientation': 'horizontal',
'label': '2m temperature (℃)',
'shrink': 0.8,
'ticks': np.arange(-30,30+5,5)}
# 画图
levels = np.arange(-30,30+1,1)
temp.plot.contourf(ax=ax, levels=levels, cmap='Spectral_r',cbar_kwargs=cbar_kwargs, transform=ccrs.PlateCarree())
plt.savefig('F:/Rpython/lp11/plot26.png',dpi=600)
plt.show()
13
import os
import maskout
import numpy as np
import xarray as xr
from scipy.interpolate import Rbf
import matplotlib as mpl
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.io.shapereader import Reader
from matplotlib.font_manager import FontProperties
SHP = 'F:/Rpython/lp1/wchina/china_shp'
# 数据处理
ds = xr.open_dataset('F:/Rpython/lp1/wchina/r160.nc')
da = ds.PRE.sel(month=slice(6, 8)).mean(dim='month').mean(dim='year')
# 插值
olon = np.linspace(72, 137, 260)
olat = np.linspace(15, 55, 120)
olon, olat = np.meshgrid(olon, olat)
func = Rbf(ds.station_lon, ds.station_lat, da, function='linear')
pre_grid = func(olon, olat)
# 字体
#mpl.rcParams['font.sans-serif'] = ['Times New Roman']
#mpl.rc('font', size=15, weight='normal')
#font = {'family': 'NSimSun','weight': 'normal','size': 15,}
font = FontProperties(fname=r"F:/Rpython/fonts/simsun.ttf")
# 颜色
levels = [0, 0.1, 10, 25, 50, 100, 250, 500]
cmaps = mpl.colors.ListedColormap([ '#FFFFFF', '#A9F18D', '#3DB93D', '#60B8FF', '#0001E2', '#F900FA', '#810045'])
norm = mpl.colors.BoundaryNorm(levels, 7)
# 画图
proj = ccrs.LambertConformal(central_latitude=90, central_longitude=105)
fig, ax = plt.subplots(figsize=(15, 15), subplot_kw=dict(projection=proj))
ax.set_extent([80, 130, 13, 55], ccrs.PlateCarree())
ax.gridlines(linestyle='--')
ax.add_feature(cfeature.OCEAN.with_scale('50m'))
ax.add_feature(cfeature.LAND.with_scale('50m'))
ax.add_feature(cfeature.RIVERS.with_scale('50m'))
ax.add_feature(cfeature.LAKES.with_scale('50m'))
ax.add_geometries(Reader(os.path.join(SHP, 'cnhimap.shp')).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k', linewidth=1)
cf = ax.contourf(olon,olat,pre_grid,levels=levels, cmap=cmaps,norm=norm,zorder=1,transform=ccrs.PlateCarree())
ax.scatter(ds.station_lon, ds.station_lat,c='k',s=10,marker='o',transform=ccrs.PlateCarree())
cbar = fig.colorbar(cf, ax=ax,orientation="vertical",aspect=25,pad=0.01,shrink=0.7)
cbar.set_label('mm')
cbar.set_ticklabels(levels)
ax.set_title('Annual average summer precipitation')
# 南海
pos1 = ax.get_position()
pos2 = [pos1.x1 - ((pos1.x1 - pos1.x0) / 7), pos1.y0, pos1.width / 7,pos1.height / 5]
proj_n = ccrs.LambertConformal(central_latitude=90, central_longitude=115)
ax_n = fig.add_axes(pos2, projection=proj_n)
ax_n.set_extent([105, 120, 0, 20], ccrs.PlateCarree())
ax_n.add_geometries(Reader(os.path.join(SHP, 'cnhimap.shp')).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=1)
ax_n.add_feature(cfeature.OCEAN.with_scale('50m'))
ax_n.add_feature(cfeature.LAND.with_scale('50m'))
ax_n.add_feature(cfeature.RIVERS.with_scale('50m'))
ax_n.add_feature(cfeature.LAKES.with_scale('50m'))
# 白化
clip = maskout.shp2clip(cf,ax,shpfile=os.path.join(SHP, 'country1.shp'),region='China',proj=proj)
fig.savefig('F:/Rpython/lp11/plot.png', dpi=600, bbox_inches='tight')
plt.show()
原文地址:https://www.cnblogs.com/ljwgis/p/12981101.html
- 冒用数字签名的对抗:亟需加强的签名审核
- Twitter开源云环境时间序列数据断层检测工具BreakoutDetection
- 用Python的长短期记忆神经网络进行时间序列预测
- 【问底】许鹏:使用Spark+Cassandra打造高性能数据分析平台(一)
- 隐秘通讯与跳板?C&C服务器究竟是怎么一回事
- 灵活布置、可二次开发的乌云公开漏洞及知识库搜索
- 干货 | 2014年我国大数据发展分析报告
- 这个恶意软件“奇葩”的反虚拟机技巧
- Android漏洞CVE-2015-3825分析及exploit实战:从Crash到劫持PC
- VaultPasswordView:可用于查看windows Vault密码的工具
- 逆向工厂(二):静态分析技术
- 打开文件夹就运行?COM劫持利用新姿势
- Java集合总览
- 常见面试题之ListView的复用及如何优化
- JavaScript 教程
- JavaScript 编辑工具
- JavaScript 与HTML
- JavaScript 与Java
- JavaScript 数据结构
- JavaScript 基本数据类型
- JavaScript 特殊数据类型
- JavaScript 运算符
- JavaScript typeof 运算符
- JavaScript 表达式
- JavaScript 类型转换
- JavaScript 基本语法
- JavaScript 注释
- Javascript 基本处理流程
- Javascript 选择结构
- Javascript if 语句
- Javascript if 语句的嵌套
- Javascript switch 语句
- Javascript 循环结构
- Javascript 循环结构实例
- Javascript 跳转语句
- Javascript 控制语句总结
- Javascript 函数介绍
- Javascript 函数的定义
- Javascript 函数调用
- Javascript 几种特殊的函数
- JavaScript 内置函数简介
- Javascript eval() 函数
- Javascript isFinite() 函数
- Javascript isNaN() 函数
- parseInt() 与 parseFloat()
- escape() 与 unescape()
- Javascript 字符串介绍
- Javascript length属性
- javascript 字符串函数
- Javascript 日期对象简介
- Javascript 日期对象用途
- Date 对象属性和方法
- Javascript 数组是什么
- Javascript 创建数组
- Javascript 数组赋值与取值
- Javascript 数组属性和方法