# 绘制水汽通量散度填色图
plt1 = ax.contourf(
lon,
lat,
div * 1e7,
zorder=1,
cmap='bwr',
levels=np.arange(-2, 2.1, 0.2),
transform=ccrs.PlateCarree(),
extend='both' # 显示超出范围的数值
)
# 绘制水汽通量矢量
cv = ax.quiver(
lon[::4],
lat[::4],
qu[::4, ::4] * 1e3,
qv[::4, ::4] * 1e3,
scale=300,
color='k',
width=0.0015,
headwidth=5,
zorder=2
)
# 添加矢量图例
font_props = FontProperties()
font_props.set_size(15)
ax.quiverkey(
cv,
0.6,
1.02,
1,
r'10g/(cm*s*hPa)',
labelpos='E',
coordinates='axes',
fontproperties=font_props
)
替换import os
import datetime
import tkinter as tk
from tkinter import ttk, filedialog, messagebox
import matplotlib
matplotlib.use('TkAgg')
from matplotlib.backends.backend_tkagg import FigureCanvasTkAgg
import xarray as xr
import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import cartopy.mpl.ticker as cticker
from cartopy.io.shapereader import Reader
import threading
import queue
import warnings
import matplotlib.colors as mcolors
from metpy.calc import divergence # 添加水汽通量散度计算
from metpy.units import units # 添加单位支持
# 忽略警告信息
warnings.filterwarnings('ignore')
# 解决中文乱码问题
plt.rcParams['font.family'] = 'SimHei'
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
class WeatherVisualizerApp:
def __init__(i, root):
i.root = root
i.root.title("气象数据可视化分析系统 - 水汽通量散度图")
i.root.geometry("1100x850")
# 设置应用图标
try:
i.root.iconbitmap("weather_icon.ico")
except:
pass
# 初始化变量
i.current_figure = None
i.canvas = None
i.output_dir = os.getcwd()
i.plot_queue = queue.Queue()
i.file_vars_info = "请选择NC文件查看变量信息"
i.station_position = [100.7625, 21.973611] # 默认站点位置(昆明)
# 创建主框架
i.main_frame = ttk.Frame(root)
i.main_frame.pack(fill=tk.BOTH, expand=True, padx=10, pady=10)
# 创建输入控件
i.create_input_controls()
# 创建绘图区域
i.create_plot_area()
# 加载地理边界数据
i.load_geographic_data()
# 设置默认值
i.time_var.set(datetime.datetime.now().strftime("%Y%m%d00"))
i.lon_range.set("85-120")
i.lat_range.set("15-35")
i.level_var.set("850") # 水汽分析通常在850hPa或700hPa
# 定期检查队列
i.root.after(100, i.process_queue)
# 状态栏
i.status_var = tk.StringVar(value="就绪")
status_bar = ttk.Label(root, textvariable=i.status_var, relief=tk.SUNKEN, anchor=tk.W)
status_bar.pack(side=tk.BOTTOM, fill=tk.X)
def create_input_controls(i):
"""创建用户输入控件"""
# 控制面板框架
control_frame = ttk.LabelFrame(i.main_frame, text="控制面板")
control_frame.pack(fill=tk.X, pady=(0, 10))
# 文件选择
file_frame = ttk.Frame(control_frame)
file_frame.pack(fill=tk.X, pady=5)
ttk.Label(file_frame, text="NC数据文件:").pack(side=tk.LEFT, padx=(5, 0))
i.file_path = tk.StringVar()
file_entry = ttk.Entry(file_frame, textvariable=i.file_path, width=50)
file_entry.pack(side=tk.LEFT, fill=tk.X, expand=True, padx=5)
ttk.Button(file_frame, text="浏览", command=i.browse_file).pack(side=tk.LEFT, padx=(0, 5))
# 参数设置
param_frame = ttk.Frame(control_frame)
param_frame.pack(fill=tk.X, pady=5)
# 时间设置
time_frame = ttk.Frame(param_frame)
time_frame.pack(side=tk.LEFT, fill=tk.X, expand=True, padx=5)
ttk.Label(time_frame, text="时间(北京时):").pack(anchor=tk.W)
i.time_var = tk.StringVar()
time_entry = ttk.Entry(time_frame, textvariable=i.time_var, width=15)
time_entry.pack(fill=tk.X, pady=2)
ttk.Label(time_frame, text="格式: YYYYMMDDHH", foreground="gray", font=("", 8)).pack(anchor=tk.W)
# 经度设置
lon_frame = ttk.Frame(param_frame)
lon_frame.pack(side=tk.LEFT, fill=tk.X, expand=True, padx=5)
ttk.Label(lon_frame, text="经度范围:").pack(anchor=tk.W)
i.lon_range = tk.StringVar()
lon_entry = ttk.Entry(lon_frame, textvariable=i.lon_range, width=15)
lon_entry.pack(fill=tk.X, pady=2)
ttk.Label(lon_frame, text="格式: min-max", foreground="gray", font=("", 8)).pack(anchor=tk.W)
# 纬度设置
lat_frame = ttk.Frame(param_frame)
lat_frame.pack(side=tk.LEFT, fill=tk.X, expand=True, padx=5)
ttk.Label(lat_frame, text="纬度范围:").pack(anchor=tk.W)
i.lat_range = tk.StringVar()
lat_entry = ttk.Entry(lat_frame, textvariable=i.lat_range, width=15)
lat_entry.pack(fill=tk.X, pady=2)
ttk.Label(lat_frame, text="格式: min-max", foreground="gray", font=("", 8)).pack(anchor=tk.W)
# 高度层设置
level_frame = ttk.Frame(param_frame)
level_frame.pack(side=tk.LEFT, fill=tk.X, expand=True, padx=5)
ttk.Label(level_frame, text="高度层(hPa):").pack(anchor=tk.W)
i.level_var = tk.StringVar()
level_entry = ttk.Entry(level_frame, textvariable=i.level_var, width=15)
level_entry.pack(fill=tk.X, pady=2)
ttk.Label(level_frame, text="水汽分析通常使用700或850hPa", foreground="gray", font=("", 8)).pack(anchor=tk.W)
# 站点位置设置
station_frame = ttk.Frame(control_frame)
station_frame.pack(fill=tk.X, pady=5)
ttk.Label(station_frame, text="站点经度:").pack(side=tk.LEFT, padx=(5, 0))
i.station_lon = tk.StringVar(value="100.76")
ttk.Entry(station_frame, textvariable=i.station_lon, width=8).pack(side=tk.LEFT)
ttk.Label(station_frame, text="纬度:").pack(side=tk.LEFT, padx=(5, 0))
i.station_lat = tk.StringVar(value="21.97")
ttk.Entry(station_frame, textvariable=i.station_lat, width=8).pack(side=tk.LEFT)
ttk.Button(station_frame, text="设置站点位置", command=i.set_station_position).pack(side=tk.LEFT, padx=10)
# 变量映射配置
var_frame = ttk.LabelFrame(control_frame, text="变量映射配置")
var_frame.pack(fill=tk.X, pady=5)
ttk.Label(var_frame, text="比湿:").pack(side=tk.LEFT, padx=(10, 2))
i.var_q = tk.StringVar(value="q")
ttk.Entry(var_frame, textvariable=i.var_q, width=8).pack(side=tk.LEFT)
ttk.Label(var_frame, text="相对湿度:").pack(side=tk.LEFT, padx=(10, 2))
i.var_rh = tk.StringVar(value="rh")
ttk.Entry(var_frame, textvariable=i.var_rh, width=8).pack(side=tk.LEFT)
ttk.Label(var_frame, text="U风:").pack(side=tk.LEFT, padx=(10, 2))
i.var_u = tk.StringVar(value="u")
ttk.Entry(var_frame, textvariable=i.var_u, width=8).pack(side=tk.LEFT)
ttk.Label(var_frame, text="V风:").pack(side=tk.LEFT, padx=(10, 2))
i.var_v = tk.StringVar(value="v")
ttk.Entry(var_frame, textvariable=i.var_v, width=8).pack(side=tk.LEFT)
# 操作按钮
btn_frame = ttk.Frame(control_frame)
btn_frame.pack(fill=tk.X, pady=10)
ttk.Button(btn_frame, text="生成水汽通量散度图", command=i.generate_plot, width=18).pack(side=tk.LEFT, padx=10)
ttk.Button(btn_frame, text="保存图表", command=i.save_image, width=12).pack(side=tk.LEFT, padx=10)
ttk.Button(btn_frame, text="批量导出", command=i.batch_export, width=12).pack(side=tk.LEFT, padx=10)
ttk.Button(btn_frame, text="重置设置", command=i.reset_settings, width=12).pack(side=tk.LEFT, padx=10)
def create_plot_area(i):
"""创建绘图区域"""
# 信息显示区域
info_frame = ttk.LabelFrame(i.main_frame, text="数据信息")
info_frame.pack(fill=tk.X, pady=(0, 10))
i.var_info_text = tk.Text(info_frame, height=6, width=80, wrap=tk.WORD)
i.var_info_scroll = ttk.Scrollbar(info_frame, orient="vertical", command=i.var_info_text.yview)
i.var_info_text.configure(yscrollcommand=i.var_info_scroll.set)
i.var_info_text.pack(side=tk.LEFT, fill=tk.BOTH, expand=True)
i.var_info_scroll.pack(side=tk.RIGHT, fill=tk.Y)
i.var_info_text.insert(tk.END, i.file_vars_info)
i.var_info_text.configure(state="disabled")
copy_btn = ttk.Button(info_frame, text="复制信息", command=i.copy_var_info, width=10)
copy_btn.pack(side=tk.RIGHT, padx=5, pady=2)
# 绘图区域
plot_frame = ttk.LabelFrame(i.main_frame)
plot_frame.pack(fill=tk.BOTH, expand=True)
# 使用PanedWindow允许调整绘图区域大小
i.plot_pane = ttk.PanedWindow(plot_frame, orient=tk.VERTICAL)
i.plot_pane.pack(fill=tk.BOTH, expand=True)
# 绘图区域框架
i.plot_frame = ttk.Frame(i.plot_pane)
i.plot_pane.add(i.plot_frame, weight=1)
def load_geographic_data(i):
"""加载地理边界数据"""
try:
# 尝试加载中国边界和云南边界
china_path = "./china_shp/country.shp"
yunnan_path = "./yunnan/yunnan.shp"
if os.path.exists(china_path):
i.china_outline = cfeature.ShapelyFeature(
Reader(china_path).geometries(),
ccrs.PlateCarree(),
edgecolor='k',
facecolor='none',
)
else:
i.china_outline = cfeature.COASTLINE
if os.path.exists(yunnan_path):
i.yunnan = cfeature.ShapelyFeature(
Reader(yunnan_path).geometries(),
ccrs.PlateCarree(),
edgecolor='k',
facecolor='none',
)
else:
i.yunnan = cfeature.COASTLINE
except Exception as e:
print(f"加载地理数据时出错: {str(e)}")
i.china_outline = cfeature.COASTLINE
i.yunnan = cfeature.COASTLINE
def process_queue(i):
"""处理队列中的任务"""
try:
while not i.plot_queue.empty():
result = i.plot_queue.get_nowait()
if result[0] == "error":
messagebox.showerror("错误", result[1])
i.status_var.set("错误: " + result[1])
elif result[0] == "figure":
fig = result[1]
i.update_plot(fig)
i.status_var.set("绘图完成")
i.current_figure = fig
i.output_dir = os.path.dirname(i.file_path.get())
except queue.Empty:
pass
i.root.after(100, i.process_queue)
def browse_file(i):
"""浏览并选择NC文件"""
file_path = filedialog.askopenfilename(
filetypes=[("NetCDF files", "*.nc"), ("All files", "*.*")]
)
if file_path:
i.file_path.set(file_path)
i.auto_detect_variables(file_path)
i.display_file_vars(file_path)
i.status_var.set("已加载文件: " + os.path.basename(file_path))
def display_file_vars(i, file_path):
"""显示NC文件的变量信息"""
try:
with xr.open_dataset(file_path) as ds:
info = f"文件: {os.path.basename(file_path)}\n"
# 查找时间变量
time_var = None
for candidate in ['time', 'valid_time', 'forecast_time', 't', 'times']:
if candidate in ds:
time_var = candidate
break
if time_var:
time_values = ds[time_var].values
if len(time_values) > 0:
info += f"时间范围: {time_values.min()} 至 {time_values.max()}\n"
else:
info += "时间范围: 无可用数据\n"
else:
info += "时间范围: 未找到时间变量\n"
# 查找经度变量
lon_var = None
for candidate in ['longitude', 'lon', 'x', 'longitudes']:
if candidate in ds:
lon_var = candidate
break
if lon_var:
lon_values = ds[lon_var].values
if len(lon_values) > 0:
info += f"经度范围: {lon_values.min()} 至 {lon_values.max()}\n"
else:
info += "经度范围: 无可用数据\n"
else:
info += "经度范围: 未找到经度变量\n"
# 查找纬度变量
lat_var = None
for candidate in ['latitude', 'lat', 'y', 'latitudes']:
if candidate in ds:
lat_var = candidate
break
if lat_var:
lat_values = ds[lat_var].values
if len(lat_values) > 0:
info += f"纬度范围: {lat_values.min()} 至 {lat_values.max()}\n"
else:
info += "纬度范围: 无可用数据\n"
else:
info += "纬度范围: 未找到纬度变量\n"
# 查找高度层变量
level_var = None
for candidate in ['pressure_level', 'level', 'lev', 'height', 'z']:
if candidate in ds:
level_var = candidate
break
if level_var:
level_values = ds[level_var].values
if len(level_values) > 0:
info += f"高度层: {level_values}\n"
else:
info += "高度层: 无可用数据\n"
else:
info += "高度层: 未找到高度层变量\n"
info += "变量:\n"
for var_name in ds.variables:
var = ds[var_name]
dims = ", ".join(var.dims)
long_name = var.attrs.get('long_name', '')
units = var.attrs.get('units', '')
info += f" - {var_name}: {dims}, {long_name}, 单位: {units}\n"
i.file_vars_info = info
i.var_info_text.configure(state="normal")
i.var_info_text.delete(1.0, tk.END)
i.var_info_text.insert(tk.END, info)
i.var_info_text.configure(state="disabled")
except Exception as e:
error_msg = f"无法读取文件信息: {str(e)}"
i.file_vars_info = error_msg
i.var_info_text.configure(state="normal")
i.var_info_text.delete(1.0, tk.END)
i.var_info_text.insert(tk.END, error_msg)
i.var_info_text.configure(state="disabled")
def auto_detect_variables(i, file_path):
"""尝试自动检测文件中的变量名"""
try:
with xr.open_dataset(file_path) as ds:
mappings = {
'q': ['q', 'shum', 'specific_humidity', 'humidity'],
'rh': ['rh', 'r', 'relhum', 'relative_humidity'],
'u': ['u', 'u_component', 'u_wind'],
'v': ['v', 'v_component', 'v_wind']
}
for var, candidates in mappings.items():
found = next((c for c in candidates if c in ds), None)
if found:
getattr(i, f'var_{var}').set(found)
except Exception as e:
print(f"自动检测变量失败: {str(e)}")
def set_station_position(i):
"""设置站点位置"""
try:
lon = float(i.station_lon.get())
lat = float(i.station_lat.get())
i.station_position = [lon, lat]
i.status_var.set(f"站点位置已设置为: 经度={lon}, 纬度={lat}")
except ValueError:
messagebox.showerror("错误", "请输入有效的经度和纬度数值")
def copy_var_info(i):
"""复制变量信息到剪贴板"""
if i.file_vars_info and i.file_vars_info != "请选择NC文件查看变量信息":
i.root.clipboard_clear()
i.root.clipboard_append(i.file_vars_info)
i.status_var.set("变量信息已复制到剪贴板")
def validate_inputs(i):
"""验证用户输入"""
if not i.file_path.get():
messagebox.showerror("错误", "请选择NC数据文件")
return False
if len(i.time_var.get()) != 10 or not i.time_var.get().isdigit():
messagebox.showerror("错误", "时间格式错误,应为YYYYMMDDHH (10位数字)")
return False
try:
# 验证范围参数
lon_min, lon_max = map(float, i.lon_range.get().split('-'))
lat_min, lat_max = map(float, i.lat_range.get().split('-'))
levels = [int(level) for level in i.level_var.get().split()]
if not levels:
messagebox.showerror("错误", "请输入至少一个高度层")
return False
except Exception as e:
messagebox.showerror("错误", f"参数格式错误: {str(e)}")
return False
return True
def generate_plot(i):
"""生成绘图线程"""
if not i.validate_inputs():
return
i.status_var.set("正在处理数据...")
threading.Thread(target=i.plot_data_thread, daemon=True).start()
def plot_data_thread(i):
"""处理数据并绘图"""
try:
# 获取用户输入
inputs = i.get_user_inputs()
if not inputs:
return
# 打开NC文件
with xr.open_dataset(inputs['file']) as ds:
# 检查必需的变量
required_vars = [inputs['vars']['q'], inputs['vars']['u'], inputs['vars']['v']]
if missing := i.check_missing_vars(ds, required_vars):
i.plot_queue.put(("error", f"缺少变量: {', '.join(missing)}"))
return
# 获取数据
data = i.get_data(ds, inputs)
# 创建图形
fig = i.create_moisture_flux_divergence_map(data, inputs)
# 放入队列
i.plot_queue.put(("figure", fig))
except Exception as e:
i.plot_queue.put(("error", f"绘图失败: {str(e)}"))
def get_user_inputs(i):
"""整理用户输入"""
try:
return {
'file': i.file_path.get(),
'time': i.time_var.get(),
'lon_range': list(map(float, i.lon_range.get().split('-'))),
'lat_range': list(map(float, i.lat_range.get().split('-'))),
'levels': [int(l) for l in i.level_var.get().split()],
'station_pos': i.station_position,
'vars': {
'q': i.var_q.get(),
'rh': i.var_rh.get(),
'u': i.var_u.get(),
'v': i.var_v.get()
}
}
except Exception as e:
i.plot_queue.put(("error", f"输入错误: {str(e)}"))
return None
def check_missing_vars(i, dataset, vars):
"""检查缺失的变量"""
return [var for var in vars if var and var not in dataset]
def get_data(i, dataset, inputs):
"""获取并处理数据"""
# 查找时间变量
time_var_candidate = None
for candidate in ['time', 'valid_time', 'forecast_time', 't', 'times']:
if candidate in dataset:
time_var_candidate = candidate
break
if not time_var_candidate:
raise ValueError("未找到时间变量")
# 时间转换
beijing_time = datetime.datetime.strptime(inputs['time'], "%Y%m%d%H")
world_time = beijing_time - datetime.timedelta(hours=8)
world_time_str = world_time.strftime("%Y-%m-%dT%H")
# 查找纬度变量
lat_var_candidate = None
for candidate in ['latitude', 'lat', 'y', 'latitudes']:
if candidate in dataset.dims:
lat_var_candidate = candidate
break
if not lat_var_candidate:
raise ValueError("未找到纬度变量")
# 查找经度变量
lon_var_candidate = None
for candidate in ['longitude', 'lon', 'x', 'longitudes']:
if candidate in dataset.dims:
lon_var_candidate = candidate
break
if not lon_var_candidate:
raise ValueError("未找到经度变量")
# 查找高度层变量
level_var_candidate = None
for candidate in ['pressure_level', 'level', 'lev', 'height', 'z']:
if candidate in dataset:
level_var_candidate = candidate
break
# 解包范围
lon_min, lon_max = inputs['lon_range']
lat_min, lat_max = inputs['lat_range']
level = inputs['levels'][0] # 暂时只处理第一个高度层
# 获取数据
data = {}
for var_type, var_name in inputs['vars'].items():
if var_name and var_name in dataset:
# 构建选择条件
sel_args = {
time_var_candidate: world_time_str,
lat_var_candidate: slice(lat_max, lat_min), # 注意纬度顺序
lon_var_candidate: slice(lon_min, lon_max)
}
if level_var_candidate:
sel_args[level_var_candidate] = level
# 获取数据子集
data[var_type] = dataset[var_name].sel(**sel_args).sortby(lat_var_candidate)
# 添加额外信息
data['time'] = beijing_time
data['level'] = level
data['lat'] = data[list(data.keys())[0]][lat_var_candidate].values
data['lon'] = data[list(data.keys())[0]][lon_var_candidate].values
data['lat_var'] = lat_var_candidate
data['lon_var'] = lon_var_candidate
return data
def calculate_moisture_flux_divergence(i, data):
"""计算水汽通量散度"""
# 获取比湿和风分量
q = data['q'] * 1000 # 转换为g/kg
u = data['u']
v = data['v']
# 计算水汽通量 (g/kg * m/s)
qu = q * u
qv = q * v
# 添加单位以便metpy计算
qu = qu.values * units('g/kg') * units('m/s')
qv = qv.values * units('g/kg') * units('m/s')
# 获取经纬度网格
lats = data['lat']
lons = data['lon']
try:
# 计算水汽通量散度 (使用metpy库)
div = divergence(qu, qv, lons, lats) * 1e7 # 缩放为10^-7 g/(kg·s·m)
return div
except Exception as e:
# 备选方案:使用numpy梯度计算
print(f"使用metpy计算散度失败: {str(e)}, 使用numpy梯度计算")
dx = np.gradient(lons) * 111000 # 转换为米 (1度≈111km)
dy = np.gradient(lats) * 111000
dqu_dx = np.gradient(qu, dx, axis=1)
dqv_dy = np.gradient(qv, dy, axis=0)
div = (dqu_dx + dqv_dy) * 1e7 # 缩放为10^-7 g/(kg·s·m)
return div
def create_moisture_flux_divergence_map(i, data, inputs):
"""创建水汽通量散度图"""
# 设置图形
fig = plt.figure(figsize=(10, 8), dpi=100)
ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
# 设置地图范围
lon_min, lon_max = inputs['lon_range']
lat_min, lat_max = inputs['lat_range']
ax.set_extent([lon_min, lon_max, lat_min, lat_max], crs=ccrs.PlateCarree())
# 添加地理特征
ax.add_feature(i.china_outline, linewidth=1.5, zorder=2)
ax.add_feature(i.yunnan, linewidth=2.0, zorder=2)
ax.add_feature(cfeature.COASTLINE.with_scale('50m'), linewidth=1.0, zorder=1)
ax.add_feature(cfeature.OCEAN, facecolor='lightblue', zorder=1)
ax.add_feature(cfeature.LAND, facecolor='lightyellow', zorder=1)
# 设置网格和刻度
i.set_map_ticks(ax, lon_min, lon_max, lat_min, lat_max)
# 添加标题
title = f"{data['level']}hPa 水汽通量散度图 {data['time'].strftime('%Y-%m-%d %H:00')} (北京时)"
ax.set_title(title, fontsize=14, loc='left')
# 计算水汽通量散度
div = i.calculate_moisture_flux_divergence(data)
# 创建自定义颜色映射
cmap = plt.cm.RdBu_r # 使用红蓝反转色带
levels = np.linspace(-10, 10, 21) # -10到10之间划分21个级别
# 绘制水汽通量散度
cf = ax.contourf(
data['lon'], data['lat'], div,
levels=levels,
cmap=cmap,
extend='both',
alpha=0.7,
transform=ccrs.PlateCarree(),
zorder=3
)
# 添加零值线(辐合辐散分界线)
contour = ax.contour(
data['lon'], data['lat'], div,
levels=[0],
colors='black',
linewidths=1.5,
linestyles='-',
transform=ccrs.PlateCarree(),
zorder=4
)
plt.clabel(contour, inline=True, fontsize=8, fmt='%1.0f')
# 添加颜色条
cbar = fig.colorbar(cf, ax=ax, orientation='vertical', pad=0.05, aspect=30)
cbar.set_label('水汽通量散度 (10$^{-7}$ g/(kg$\cdot$s$\cdot$m))', fontsize=10)
cbar.ax.tick_params(labelsize=8)
# 绘制风场(箭头)
if 'u' in data and 'v' in data:
u_step = max(1, len(data['lon']) // 15) # 自动计算步长
v_step = max(1, len(data['lat']) // 15)
i.draw_wind_vectors(ax, data['lon'], data['lat'], data['u'], data['v'], u_step, v_step)
# 标注站点位置
i.mark_station(ax, inputs['station_pos'])
# 添加版权信息
plt.figtext(0.95, 0.02, "气象数据可视化分析系统 - 水汽通量散度图 © 2023",
ha='right', va='bottom', fontsize=8, color='gray')
plt.tight_layout()
return fig
def set_map_ticks(i, ax, lon_min, lon_max, lat_min, lat_max):
"""设置地图刻度"""
lon_span = max(5, int((lon_max - lon_min) / 5))
lat_span = max(5, int((lat_max - lat_min) / 5))
ax.set_xticks(np.arange(lon_min, lon_max + lon_span * 0.1, lon_span), crs=ccrs.PlateCarree())
ax.set_yticks(np.arange(lat_min, lat_max + lat_span * 0.1, lat_span), crs=ccrs.PlateCarree())
lon_formatter = cticker.LongitudeFormatter()
lat_formatter = cticker.LatitudeFormatter()
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)
ax.tick_params(axis='both', labelsize=10)
def draw_wind_vectors(i, ax, lon, lat, u, v, u_step, v_step):
"""绘制风矢量箭头"""
try:
# 创建风矢量的网格点
x = lon[::u_step]
y = lat[::v_step]
# 提取对应网格点的风分量
u_sub = u[::v_step, ::u_step]
v_sub = v[::v_step, ::u_step]
# 计算风速大小用于颜色映射
speed = np.sqrt(u_sub ** 2 + v_sub ** 2)
# 绘制风矢量箭头
quiver = ax.quiver(
x, y, u_sub, v_sub, speed,
cmap='cool', # 使用冷色调表示风速
scale=200, # 箭头缩放比例
scale_units='inches', # 缩放单位
width=0.003, # 箭头宽度
transform=ccrs.PlateCarree(),
zorder=5
)
# 添加风速颜色条
cbar = plt.colorbar(quiver, ax=ax, orientation='vertical', pad=0.01, aspect=20)
cbar.set_label('风速 (m/s)', fontsize=9)
cbar.ax.tick_params(labelsize=8)
# 添加参考箭头
ax.quiverkey(quiver, 0.92, 0.95, 10, '10 m/s',
labelpos='E', coordinates='axes', fontproperties={'size': 8})
return quiver
except Exception as e:
print(f"绘制风矢量失败: {str(e)}")
return None
def mark_station(i, ax, station_pos):
"""标注站点位置"""
# 获取站点位置
station_lon, station_lat = station_pos
# 绘制星形标记
ax.plot(
station_lon, station_lat,
marker='*', # 使用星形标记
markersize=12, # 增大标记尺寸
color='red', # 红色标记
transform=ccrs.PlateCarree(), # 添加投影转换
zorder=10, # 更高的zorder确保在最上层
markeredgecolor='gold', # 金色边缘增加可见性
markeredgewidth=0.5 # 边缘宽度
)
# 添加站点标签(带背景框)
ax.text(
station_lon + 0.7, station_lat + 0.7, '站点', # 调整位置避免重叠
fontsize=10, # 增大字体
color='red', # 红色文字
weight='bold', # 加粗
transform=ccrs.PlateCarree(), # 添加投影转换
bbox=dict(
facecolor='white', # 白色背景
alpha=0.85, # 稍高的不透明度
boxstyle='round,pad=0.3', # 圆角矩形,增加内边距
edgecolor='red', # 红色边框
linewidth=0.8 # 边框线宽
),
zorder=10
)
def update_plot(i, fig):
"""更新绘图区域的显示"""
# 清除旧的画布
if i.canvas:
i.canvas.get_tk_widget().destroy()
# 创建新的画布
i.canvas = FigureCanvasTkAgg(fig, master=i.plot_frame)
i.canvas.draw()
i.canvas.get_tk_widget().pack(fill=tk.BOTH, expand=True)
def save_image(i):
"""保存图片"""
if not i.current_figure:
messagebox.showinfo("提示", "请先生成图片")
return
# 获取当前时间作为默认文件名
current_time = datetime.datetime.now().strftime("%Y%m%d%H%M%S")
default_filename = f"水汽通量散度图_{current_time}.png"
# 打开文件保存对话框
file_path = filedialog.asksaveasfilename(
initialdir=i.output_dir,
initialfile=default_filename,
defaultextension=".png",
filetypes=[("PNG files", "*.png"), ("JPEG files", "*.jpg"), ("PDF files", "*.pdf"),
("All files", "*.*")]
)
if file_path:
try:
# 根据文件扩展名选择保存格式
file_ext = os.path.splitext(file_path)[1].lower()
if file_ext == '.pdf':
i.current_figure.savefig(file_path, format='pdf', bbox_inches='tight')
elif file_ext == '.jpg':
i.current_figure.savefig(file_path, format='jpeg', dpi=300, bbox_inches='tight')
else:
i.current_figure.savefig(file_path, dpi=300, bbox_inches='tight')
i.status_var.set(f"图片已保存至: {file_path}")
messagebox.showinfo("成功", f"图片已保存至: {file_path}")
except Exception as e:
messagebox.showerror("错误", f"保存失败: {str(e)}")
i.status_var.set("保存失败")
def batch_export(i):
"""批量导出图片"""
if not i.file_path.get():
messagebox.showerror("错误", "请先选择NC数据文件")
return
# 选择输出目录
output_dir = filedialog.askdirectory(title="选择输出目录")
if not output_dir:
return
# 获取时间范围
try:
with xr.open_dataset(i.file_path.get()) as ds:
# 查找时间变量
time_var_candidate = None
for candidate in ['time', 'valid_time', 'forecast_time', 't', 'times']:
if candidate in ds:
time_var_candidate = candidate
break
if not time_var_candidate:
messagebox.showerror("错误", "文件中没有有效的时间变量")
return
times = ds[time_var_candidate].values
if len(times) == 0:
messagebox.showerror("错误", "文件中没有有效的时间数据")
return
# 显示进度对话框
progress = tk.Toplevel(i.root)
progress.title("导出进度")
progress.geometry("300x150")
progress_label = ttk.Label(progress, text="准备开始导出...")
progress_label.pack(pady=10)
progress_bar = ttk.Progressbar(progress, orient="horizontal", length=250, mode="determinate")
progress_bar.pack(pady=10)
progress_cancel = ttk.Button(progress, text="取消",
command=lambda: setattr(progress, "cancelled", True))
progress_cancel.pack(pady=5)
progress.cancelled = False
progress.update()
# 导出每个时间点的图片
total = len(times)
success_count = 0
for idx, t in enumerate(times):
if progress.cancelled:
break
# 更新进度
progress_label.config(text=f"正在处理 {idx + 1}/{total} ({t})")
progress_bar['value'] = (idx + 1) / total * 100
progress.update()
# 生成图片
try:
# 转换时间格式
time_dt = t.astype('datetime64[s]').astype(datetime.datetime)
beijing_time = time_dt + datetime.timedelta(hours=8)
time_str = beijing_time.strftime("%Y%m%d%H")
i.time_var.set(time_str)
inputs = i.get_user_inputs()
if not inputs:
continue
# 获取数据
data = i.get_data(ds, inputs)
# 创建图形
fig = i.create_moisture_flux_divergence_map(data, inputs)
# 保存图片
filename = f"水汽通量散度图_{time_str}_{inputs['levels'][0]}hPa.png"
fig.savefig(os.path.join(output_dir, filename), dpi=150, bbox_inches='tight')
plt.close(fig)
success_count += 1
except Exception as e:
print(f"处理时间点 {t} 时出错: {str(e)}")
progress.destroy()
messagebox.showinfo("完成", f"已导出 {success_count}/{total} 张水汽通量散度图到目录: {output_dir}")
except Exception as e:
messagebox.showerror("错误", f"批量导出失败: {str(e)}")
i.status_var.set("批量导出失败")
def reset_settings(i):
"""重置所有设置为默认值"""
i.file_path.set("")
i.time_var.set(datetime.datetime.now().strftime("%Y%m%d00"))
i.lon_range.set("85-120")
i.lat_range.set("15-35")
i.level_var.set("850")
i.station_lon.set("100.76")
i.station_lat.set("21.97")
i.var_q.set("q")
i.var_rh.set("rh")
i.var_u.set("u")
i.var_v.set("v")
i.station_position = [100.7625, 21.973611]
i.file_vars_info = "请选择NC文件查看变量信息"
i.var_info_text.configure(state="normal")
i.var_info_text.delete(1.0, tk.END)
i.var_info_text.insert(tk.END, i.file_vars_info)
i.var_info_text.configure(state="disabled")
i.status_var.set("设置已重置")
# 清除当前图形
if i.canvas:
i.canvas.get_tk_widget().destroy()
i.canvas = None
i.current_figure = None
# 主程序入口
if __name__ == "__main__":
root = tk.Tk()
app = WeatherVisualizerApp(root)
root.mainloop()
最新发布