動機
有個同事吧，寫論文，讓我?guī)兔Τ鰝€圖，就寫了個代碼，然后我的博客好久沒更新了，就順便貼上來了！
很多人感興趣風速的箭頭怎樣畫，可能這種圖使用 NCL 非常容易，很多沒用過代碼的小朋友，就有點犯怵，怕 python 畫起來很困難。但是不然，看完我的代碼，就會發(fā)現(xiàn)很簡單，并且也可以批量，同時還能自定義國界等shp文件，這對于發(fā)sci等國際論文很重要，因為有時候內(nèi)置的國界是有問題的。
數(shù)據(jù)
本次博客使用的數(shù)據(jù)為 ERA5 hourly data on pressure levels from 1940 to present數(shù)據(jù)，數(shù)據(jù)的下載方式及注冊賬號，我在前面的博客中都寫過，詳細可參考以下兩篇博客：
以下為我們數(shù)據(jù)介紹界面和需要下載的變量：
數(shù)據(jù)介紹地址：https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-pressure-levels?tab=overview
數(shù)據(jù)選擇界面
代碼
廢話不多說，直接上代碼。文章來源：http://www.zghlxwxcb.cn/news/detail-857146.html
導入包
import xarray as xr
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
import geopandas as gpd
# 設置全局字體為新羅馬
plt.rcParams['font.family'] = 'serif'
plt.rcParams['font.serif'] = ['Times New Roman']
# plt.rcParams['font.serif'] = ['SimSun']
# 設置全局字體權(quán)重為normal
plt.rcParams['font.weight'] = 'normal'
# 設置全局字體大小
matplotlib.rcParams['font.size'] = 19 # 設置全局字體大小為12
畫水汽通量散度圖
# 加載shapefile
gdf = gpd.read_file(r'./shp/Pronvience.shp')
# 使用geopandas讀取地理數(shù)據(jù)，這里我們手動創(chuàng)建一個GeoDataFrame
gdf_point = gpd.GeoDataFrame({
'City': ['Mingfeng Station', 'Kalasai Station'],
'Latitude': [37.5,37],
'Longitude': [80,81]
}, geometry=gpd.points_from_xy([80,81], [37.5,37]))
# 載入數(shù)據(jù)
data_path = r'./20170731_case.nc' # 替換為您的文件路徑
ds = xr.open_dataset(data_path)
time = '2017-07-30T22:00:00'
# level_hPa = 700
# for level_hPa in [200,500,700,850]:
for level_hPa in [600]:
# 選擇特定時間和氣壓層
ds_selected = ds.sel(time= time, level=level_hPa) # 示例：2022年1月1日0時，850hPa
# 獲取數(shù)據(jù)變量
u = ds_selected['u'] # 東西向風速
v = ds_selected['v'] # 南北向風速
q = ds_selected['q'] # 比濕
# 獲取經(jīng)度和緯度，假設這些是坐標維度
longitude = u.longitude
latitude = u.latitude
# 計算水汽通量
qu = q * u # 東西向水汽通量
qv = q * v # 南北向水汽通量
# 計算水汽通量散度 單位為
div_q = (qu.differentiate('longitude') + qv.differentiate('latitude'))* 10
# 打印結(jié)果
# print(div_q)
# 創(chuàng)建圖形和軸對象
fig, ax = plt.subplots(figsize=(6, 6),dpi=500) # 圖形尺寸為10x6英寸
# 可視化散度結(jié)果
contour = div_q.plot(add_colorbar=False, cmap="RdBu_r", vmin=-1, vmax=1) # 使用黑色線條繪制20個等級的等高線
#
# 在ax上繪制等高線圖
div_q.plot.contour(levels=25, colors='black',linewidths=0.6)
# 添加顏色條
fig.colorbar(contour, ax=ax, label='Water Vapor Flux Divergence (g/cm2/s)')
# 使用quiver函數(shù)需要確保數(shù)據(jù)的間隔，這里我們每隔5個點取樣
Q = ax.quiver(longitude[::5], latitude[::5], u[::5, ::5], v[::5, ::5], scale=300,color="red")
# 繪制shapefile
gdf.plot(ax=ax, color='none', edgecolor='green',linewidths=0.7) # 無填充，黑色邊界
# gdf_point.plot(ax=ax, color='red') # 標記紐約的位置
# 繪制點
ax.scatter(gdf_point['Longitude'], gdf_point['Latitude'], color='red', s=100)
# 標注城市名稱
for x, y, city in zip(gdf_point['Longitude'], gdf_point['Latitude'], gdf_point['City']):
ax.text(x, y, ' ' + city, verticalalignment='center', fontsize=15)
# 設置經(jīng)緯度范圍
ax.set_xlim(75, 90)
ax.set_ylim(30, 45)
ax.set_xlabel('Longitude')
ax.set_ylabel('Latitude')
ax.set_title('') # 清除標題
# 添加標題在圖片正下方
# fig.suptitle('{}hPa {}'.format( level_hPa,time.replace("T"," ") ), y=-0.01,va='bottom')
# 調(diào)整布局以避免重疊和裁剪
fig.tight_layout()
plt.savefig("./{}hPa {}.jpg".format( level_hPa,time.replace(":","") ), dpi=500)
plt.show()
水汽通量圖
# 加載shapefile
gdf = gpd.read_file(r'./shp/Pronvience.shp')
# 載入數(shù)據(jù)
data_path = r'./20170731_case.nc' # 替換為您的文件路徑
ds = xr.open_dataset(data_path)
time = '2017-07-30T22:00:00'
for level_hPa in [200,500,600,700,850]:
# 選擇特定時間和氣壓層
ds_selected = ds.sel(time= time, level=level_hPa) # 示例：2022年1月1日0時，850hPa
# 獲取數(shù)據(jù)變量
u = ds_selected['u'] # 東西向風速
v = ds_selected['v'] # 南北向風速
q = ds_selected['q'] # 比濕
# 獲取經(jīng)度和緯度，假設這些是坐標維度
longitude = u.longitude
latitude = u.latitude
# 計算水汽通量
qu = q * u * 100 # 東西向水汽通量
qv = q * v * 100 # 南北向水汽通量
wvf = np.sqrt(qu**2 + qv**2)
# 計算水汽通量散度 單位為
# div_q = (qu.differentiate('longitude') + qv.differentiate('latitude'))* 10
# 打印結(jié)果
# print(div_q)
# 創(chuàng)建圖形和軸對象
fig, ax = plt.subplots(figsize=(6, 6),dpi=400) # 圖形尺寸為10x6英寸
# 可視化散度結(jié)果
contour = wvf.plot(add_colorbar=False, cmap="RdBu_r", vmin=0, vmax=10) # 使用黑色線條繪制20個等級的等高線
#
# 在ax上繪制等高線圖
wvf.plot.contour(levels=25, colors='black',linewidths=0.6)
# 添加顏色條
fig.colorbar(contour, ax=ax, label='Water Vapor Flux(g/cm/s)')
# 使用quiver函數(shù)需要確保數(shù)據(jù)的間隔，這里我們每隔5個點取樣
Q = ax.quiver(longitude[::5], latitude[::5], u[::5, ::5], v[::5, ::5], scale=300,color="red")
# 繪制shapefile
gdf.plot(ax=ax, color='none', edgecolor='green',linewidths=0.7) # 無填充，黑色邊界
# 設置經(jīng)緯度范圍
ax.set_xlim(75, 90)
ax.set_ylim(30, 45)
ax.set_xlabel('Longitude')
ax.set_ylabel('Latitude')
ax.set_title('') # 清除標題
# 添加標題在圖片正下方
# fig.suptitle('{}hPa {}'.format( level_hPa,time.replace("T"," ") ), y=-0.01,va='bottom')
# 調(diào)整布局以避免重疊和裁剪
fig.tight_layout()
plt.savefig("./WVF_{}hPa {}.jpg".format( level_hPa,time.replace(":","") ), dpi=500)
plt.show()
結(jié)果圖
文章來源地址http://www.zghlxwxcb.cn/news/detail-857146.html
到了這里，關于python ERA5 畫水汽通量散度圖地圖：風速風向矢量圖、疊加等高線、色彩分級、添加shp文件、添加位置點及備注的文章就介紹完了。如果您還想了解更多內(nèi)容，請在右上角搜索TOY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關文章，希望大家以后多多支持TOY模板網(wǎng)！
