import numpy as np
import xarray as xr
import cartopy.crs as ccrs
from matplotlib import pyplot as plt
data = xr.open_dataset(r"D:\考试\Geo-500hPa-98-20-NH-1deg.nc",
drop_variables=["time_bnds"])
#drop_variables 除去多余的变量
data.hgt
f = data
print(f.dims)#维数
print('-----------------')
print(f.coords)#坐标
print('-----------------')
print(f.attrs)#资料详细来源
help(f)
import datetime as dt
#处理 datetime64 格式数据使用,可以解析时间格式
f.hgt.time
#查看 f.sst.time
f.hgt.time.dt.month
#将 f.hgt.time 解析为月份2
# 索引 1991-2020 年之间 1 月的数据 500hPa
z = f.hgt.loc[f.time.dt.month.isin([1])].loc['1991-01-01':'2020-02-01',
500]
z.loc['2008-01-01']
ave_hgt = np.array(z).reshape(30, 73, 144).mean(0)
#求平均位势高度场
(即气候态)
ave_hgt
ave_hgt.shape
#维数
from cartopy.util import add_cyclic_point
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER,
LATITUDE_FORMATTER
import matplotlib.ticker as mticker
import cartopy.feature as cfeature
import matplotlib
import cartopy.mpl.ticker as cticker
#刻度
#防止中文出错
matplotlib.rcParams['font.sans-serif'] = ['SimHei']
#matplotlib.rcParams['font.family']='sans-serif'
plt.rcParams['axes.unicode_minus'] = False
#除去 0 和 360 处空白
chgt, cycle_lon = add_cyclic_point(ave_hgt, coord=lon)
LON, LAT = np.meshgrid(cycle_lon, lat)
def contour_map(fig, img_extent, spec):
fig.set_extent(img_extent, crs=ccrs.PlateCarree())
fig.add_feature(cfeature.COASTLINE.with_scale('50m'))
#添加海岸线
fig.add_feature(cfeature.LAKES, alpha=0.5)
#添加湖泊
fig.set_xticks(np.arange(leftlon, rightlon + spec, spec),
crs=ccrs.PlateCarree())
fig.set_yticks(np.arange(lowerlat, upperlat + spec, spec),
crs=ccrs.PlateCarree())
lon_formatter = cticker.LongitudeFormatter()
lat_formatter = cticker.LatitudeFormatter()
fig.xaxis.set_major_formatter(lon_formatter)
fig.yaxis.set_major_formatter(lat_formatter)
leftlon, rightlon, lowerlat, upperlat = (-180, 180, -90, 90)
img_extent = [leftlon, rightlon, lowerlat, upperlat]
#经纬度范围
# 生成画布
fig = plt.figure(figsize=(12, 5))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
# 绘制等值线图
contour = ax.contourf(LON, LAT,3
chgt/10.0,levels=np.arange(500,600,6),linewidths=0.90,cmap='jet',tran
sform=ccrs.PlateCarree())
contour1 = ax.contour(LON, LAT,
chgt/10.,colors='black',levels=np.arange(500,600,6),linewidths=0.9)
contour_map(ax, img_extent, 30)
#ax.clabel(contour)
#添加等值线的值
ax.set_title('1991-2020 年平均 1 月 500hPa 位势高度')
#设置标签
plt.xticks(fontsize=13)#设置标签大小
plt.yticks(fontsize=13)
plt.colorbar(contour)#色标
plt.savefig('D:\data\DQXHYC\sx01\图 1.jpg')
plt.show()
z2008 = f.hgt.loc[f.time.dt.month.isin([1])].loc['2008-01-01', 500]
z2008
dep=z2008-ave_hgt#求距平
dep
from cartopy.util import add_cyclic_point
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER,
LATITUDE_FORMATTER
import matplotlib.ticker as mticker
import cartopy.feature as cfeature
import cartopy.mpl.ticker as cticker
#刻度
#除去 0 和 360 处空白
cdep, cycle_lon = add_cyclic_point(dep, coord=lon)
LON, LAT = np.meshgrid(cycle_lon, lat)
def contour_map(fig, img_extent, spec):
fig.set_extent(img_extent, crs=ccrs.PlateCarree())
fig.add_feature(cfeature.COASTLINE.with_scale('50m'))
#添加海岸
线
fig.add_feature(cfeature.LAKES, alpha=0.5)
#添加湖泊
fig.set_xticks(np.arange(leftlon, rightlon + spec, spec),
crs=ccrs.PlateCarree())
fig.set_yticks(np.arange(lowerlat, upperlat + spec, spec),
crs=ccrs.PlateCarree())
lon_formatter = cticker.LongitudeFormatter()
lat_formatter = cticker.LatitudeFormatter()
fig.xaxis.set_major_formatter(lon_formatter)
fig.yaxis.set_major_formatter(lat_formatter)4
leftlon, rightlon, lowerlat, upperlat = (-180, 180, -90, 90)
img_extent = [leftlon, rightlon, lowerlat, upperlat]
#经纬度范围
# 生成画布
fig = plt.figure(figsize=(12, 5))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
# 绘制等值线图
contour = ax.contourf(LON, LAT,
cdep/10.0,levels=np.arange(-20,20,2),linewidths=0.90,cmap='seismic',
transform=ccrs.PlateCarree())
contour1 = ax.contour(LON, LAT,
cdep/10.,colors='black',levels=np.arange(-20,20,2),linewidths=0.9)
contour_map(ax, img_extent, 30)
#ax.clabel(contour)
#添加等值线的值
ax.set_title('2008 年 1 月 500hPa 位势高度距平(即相对于气候态的偏差)')
#设置标签
plt.xticks(fontsize=13)#设置标签大小
plt.yticks(fontsize=13)
plt.colorbar(contour)#色标
plt.savefig('D:\data\DQXHYC\sx01\图 2.jpg')
plt.show()
ave_lon = np.array(z2008).reshape(73,144).mean(1)
#纬圈平均值
type(ave_lon)
#%%
(z2008).shape
#%%
ave_lons=np.array([ave_lon]*144).T
ave_lons#转置矩阵,使其维度与 z2008 相同
#%%
wp=z2008-ave_lons#求纬偏 wp
#%% md
## 绘制环流纬偏图
from cartopy.util import add_cyclic_point
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER,
LATITUDE_FORMATTER
import matplotlib.ticker as mticker
import cartopy.feature as cfeature
import cartopy.mpl.ticker as cticker
#刻度5
#除去 0 和 360 处空白
cwp, cycle_lon = add_cyclic_point(wp, coord=lon)
LON, LAT = np.meshgrid(cycle_lon, lat)
def contour_map(fig, img_extent, spec):
fig.set_extent(img_extent, crs=ccrs.PlateCarree())
fig.add_feature(cfeature.COASTLINE.with_scale('50m'))
#添加海岸 线
fig.add_feature(cfeature.LAKES, alpha=0.5)
#添加湖泊
fig.set_xticks(np.arange(leftlon, rightlon + spec, spec),
crs=ccrs.PlateCarree())
fig.set_yticks(np.arange(lowerlat, upperlat + spec, spec),
crs=ccrs.PlateCarree())
lon_formatter = cticker.LongitudeFormatter()
lat_formatter = cticker.LatitudeFormatter()
fig.xaxis.set_major_formatter(lon_formatter)
fig.yaxis.set_major_formatter(lat_formatter)
leftlon, rightlon, lowerlat, upperlat = (-180, 180, -90, 90)
img_extent = [leftlon, rightlon, lowerlat, upperlat]
#经纬度范围
# 生成画布
fig = plt.figure(figsize=(12, 5))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
# 绘制等值线图
contour = ax.contourf(LON, LAT,
cwp/10.0,levels=np.arange(-30,30,3),linewidths=1.0,cmap='RdBu_r',
transform=ccrs.PlateCarree())
contour1 = ax.contour(LON, LAT,
cwp/10.,colors='black',levels=np.arange(-30,30,3),linewidths=1.0)
contour_map(ax, img_extent, 30)
#ax.clabel(contour)
#添加等值线的值
ax.set_title('2008 年 1 月 500hPa 位势高度纬偏值')
#设置标签
plt.xticks(fontsize=13)#设置标签大小
plt.yticks(fontsize=13)
plt.colorbar(contour)#色标
#position = fig.add_axes([0.91, 0.28, 0.02, 0.45])
# 坐标+长宽
#cbar = fig.colorbar(contour, cax=position, orientation='vertical')
# 水平 horizontal
plt.savefig('D:\data\DQXHYC\sx01\图 3.jpg')
plt.show()