# -*- coding: utf-8 -*-
"""
Created on Wed Jun 18 14:12:50 2025
@author: 25636
"""
import netCDF4 as nc
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
import matplotlib.ticker as mticker
import matplotlib.colors as mcolors
# 读取NetCDF文件
file_path = r'C:/Users/25636/Desktop/bi/ziji/zhengduanfx/first/ziliao/200307.nc'
data = nc.Dataset(file_path)
# 提取变量数据
u = data.variables['u'][:] # 纬向风 (m/s)
v = data.variables['v'][:] # 经向风 (m/s)
q = data.variables['q'][:] # 比湿 (kg/kg)
lats = data.variables['latitude'][:]
lons = data.variables['longitude'][:]
levels = data.variables['level'][:]
# 假设选择第一个时次和850hPa层
time_idx = 0
level_idx = np.where(levels == 850)[0][0] # 找到850hPa对应的索引
# 提取特定层和时次的数据
u_850 = u[time_idx, level_idx, :, :]
v_850 = v[time_idx, level_idx, :, :]
q_850 = q[time_idx, level_idx, :, :]
# 计算水汽通量 (kg/(m·s))
qu = u_850 * q_850 * 1000 # 乘以1000将kg/kg转换为g/kg
qv = v_850 * q_850 * 1000
# 计算水汽通量大小用于着色
qv_magnitude = np.sqrt(qu**2 + qv**2)
# 计算水汽通量散度
earth_radius = 6371000 # 地球半径 (m)
dlon = np.deg2rad(np.diff(lons)[0]) # 经度间隔 (弧度)
dlat = np.deg2rad(np.diff(lats)[0]) # 纬度间隔 (弧度)
# 创建经纬度网格
lon_grid, lat_grid = np.meshgrid(lons, lats)
# 计算网格间距 (米)
dx = earth_radius * np.cos(np.deg2rad(lat_grid)) * dlon
dy = earth_radius * dlat
# 计算水汽通量散度
dqu_dx = np.gradient(qu, axis=1) / dx
dqv_dy = np.gradient(qv, axis=0) / dy
div = dqu_dx + dqv_dy # 单位: g/(kg·s·m)
# 创建图形
fig = plt.figure(figsize=(15, 15), dpi=120)
plt.rcParams.update({'font.size': 12, 'font.family': 'sans-serif'})
# =============== 水汽通量图 ===============
ax1 = fig.add_subplot(211, projection=ccrs.PlateCarree())
ax1.set_title('850 hPa Water Vapor Flux Vector', fontsize=16, pad=10, weight='bold')
# 添加地理特征
ax1.add_feature(cfeature.COASTLINE, linewidth=1.2)
ax1.add_feature(cfeature.BORDERS, linestyle='-', linewidth=0.8, alpha=0.7)
ax1.add_feature(cfeature.LAKES, alpha=0.5, edgecolor='black', facecolor='lightblue')
ax1.add_feature(cfeature.RIVERS, linewidth=0.8, edgecolor='blue')
ax1.add_feature(cfeature.OCEAN, facecolor='lightblue', alpha=0.2)
# 设置经纬度刻度和标签
gl = ax1.gridlines(draw_labels=True, linestyle='--', alpha=0.5)
gl.xlocator = mticker.FixedLocator(np.arange(0, 180, 10))
gl.ylocator = mticker.FixedLocator(np.arange(0, 90, 10))
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
gl.xlabel_style = {'size': 10}
gl.ylabel_style = {'size': 10}
gl.xlabels_top = False
gl.ylabels_right = False
# 计算合适的箭头密度
skip_x = max(1, int(lons.size / 25))
skip_y = max(1, int(lats.size / 25))
skip = max(skip_x, skip_y)
# 绘制带颜色的水汽通量箭头
q = ax1.quiver(
lon_grid[::skip, ::skip],
lat_grid[::skip, ::skip],
qu[::skip, ::skip],
qv[::skip, ::skip],
qv_magnitude[::skip, ::skip],
cmap='viridis',
scale=500,
width=0.003,
headwidth=3,
headlength=4,
headaxislength=3.5,
transform=ccrs.PlateCarree()
)
# 添加颜色条
cbar = plt.colorbar(q, ax=ax1, shrink=0.7, pad=0.02)
cbar.set_label('Water Vapor Flux Magnitude (g/kg·m/s)', fontsize=10)
# 添加参考矢量
ax1.quiverkey(q, 0.85, 0.05, 100, '100 g/kg·m/s',
labelpos='S', coordinates='axes', color='red',
fontproperties={'size': 9, 'weight': 'bold'})
# =============== 水汽通量散度图 ===============
ax2 = fig.add_subplot(212, projection=ccrs.PlateCarree())
ax2.set_title('850 hPa Water Vapor Flux Divergence', fontsize=16, pad=10, weight='bold')
# 添加地理特征
ax2.add_feature(cfeature.COASTLINE, linewidth=1.2)
ax2.add_feature(cfeature.BORDERS, linestyle='-', linewidth=0.8, alpha=0.7)
ax2.add_feature(cfeature.LAKES, alpha=0.5, edgecolor='black', facecolor='lightblue')
ax2.add_feature(cfeature.RIVERS, linewidth=0.8, edgecolor='blue')
ax2.add_feature(cfeature.OCEAN, facecolor='lightblue', alpha=0.2)
# 设置经纬度刻度和标签
gl2 = ax2.gridlines(draw_labels=True, linestyle='--', alpha=0.5)
gl2.xlocator = mticker.FixedLocator(np.arange(0, 180, 10))
gl2.ylocator = mticker.FixedLocator(np.arange(0, 90, 10))
gl2.xformatter = LONGITUDE_FORMATTER
gl2.yformatter = LATITUDE_FORMATTER
gl2.xlabel_style = {'size': 10}
gl2.ylabel_style = {'size': 10}
gl2.xlabels_top = False
gl2.ylabels_right = False
# 创建自定义的离散色标
div_cmap = plt.get_cmap('bwr_r')
levels = np.linspace(-5, 5, 21)
norm = mcolors.BoundaryNorm(levels, div_cmap.N)
# 绘制散度填色图
contour = ax2.contourf(
lon_grid,
lat_grid,
div * 1e7, # 缩放显示单位
cmap=div_cmap,
norm=norm,
levels=levels,
extend='both',
transform=ccrs.PlateCarree()
)
# 添加等值线
contour_lines = ax2.contour(
lon_grid,
lat_grid,
div * 1e7,
levels=np.arange(-4, 5, 1),
colors='black',
linewidths=0.5,
transform=ccrs.PlateCarree()
)
ax2.clabel(contour_lines, fmt='%1.0f', fontsize=8)
# 添加颜色条
cbar_div = plt.colorbar(contour, ax=ax2, shrink=0.7, pad=0.02,
extendfrac='auto', ticks=np.arange(-5, 6, 1))
cbar_div.set_label('Divergence (10$^{-7}$ g/kg·s·m)', fontsize=10)
# 添加重要地理标签
def add_region_labels(ax):
"""在图上添加重要区域标签"""
ax.text(115, 35, 'China', transform=ccrs.PlateCarree(),
fontsize=12, weight='bold', ha='center', color='darkblue')
ax.text(135, 35, 'Japan', transform=ccrs.PlateCarree(),
fontsize=10, weight='bold', ha='center', color='darkblue')
ax.text(125, 45, 'Sea of Japan', transform=ccrs.PlateCarree(),
fontsize=9, ha='center', color='darkblue', alpha=0.8)
add_region_labels(ax1)
add_region_labels(ax2)
# 添加日期信息
time_var = data.variables['time']
time_unit = time_var.units
date_str = nc.num2date(time_var[time_idx], time_unit).strftime('%Y-%m-%d %H:%M UTC')
fig.text(0.15, 0.95, f'Analysis Time: {date_str}', fontsize=12, weight='bold')
# 添加数据来源说明
fig.text(0.5, 0.01, 'Data Source: ECMWF Reanalysis | Visualization: Python Cartopy',
ha='center', fontsize=9, alpha=0.7)
# 调整布局并保存
plt.tight_layout(rect=[0, 0.03, 1, 0.95]) # 为页眉页脚留空间
plt.savefig('water_vapor_flux_analysis.png', dpi=300, bbox_inches='tight')
plt.savefig('water_vapor_flux_analysis.pdf', bbox_inches='tight') # 矢量格式
plt.show()
修改代码，使输出的图片更加美观专业
本文将同调球链接的签名不变量理论应用于同调边界链接的覆盖链接，通过同调边界链接的模式和Seifert矩阵，导出了计算覆盖链接签名不变量的公式。利用该公式，我们构造出了一些与 ribbon 链接正变异体但与边界链接不相容的融合边界链接，并证明对于任何有限集合的模式，存在不与集合内模式对应的同调边界链接不相容的情况。
扫一扫
2423
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
