1. 项目概述：从数据文件到气象图景的桥梁

如果你正在处理气象、海洋或气候数据，那么netCDF文件对你来说一定不陌生。这种自描述、跨平台的数据格式，几乎是地球科学领域的“通用语言”。然而，从拿到一个.nc文件，到最终生成一张能清晰表达科学含义的图表，中间的路往往并不平坦。数据维度复杂、坐标系统抽象、绘图库选择困难，每一步都可能让新手感到无从下手。

这个系列，我们就来系统地解决这个问题。作为开篇，我们将聚焦最基础也最关键的一步：如何用Python高效地读取netCDF文件，并完成一次标准的可视化绘图。我不会只给你一段干巴巴的代码，而是会带你理解每一步背后的“为什么”，分享我在处理上百个气象数据集时踩过的坑和总结的技巧。无论你是大气科学、海洋学的研究生，还是从事环境数据分析的工程师，掌握这套从数据到图形的标准化流程，都能让你的工作效率大幅提升。

我们将使用目前公认的、处理这类多维网格数据最得心应手的库——xarray。它就像是专门为netCDF数据定制的“Pandas”，让你能用类似处理表格数据的直观方式，去操作具有经度、纬度、时间、高度等多维度的复杂数据。配合matplotlib和cartopy进行绘图，你将能轻松复现论文中的各类气象要素空间分布图、时间序列图。

2. 核心工具链解析：为什么是xarray+Cartopy？

在开始写代码之前，我们得先搞清楚手里的“兵器”。气象数据可视化不是简单的plt.plot()，其特殊性在于数据本身的多维性和地理投影的需求。盲目组合工具只会事倍功半。

2.1 xarray：为多维网格数据而生

你可能会问，用numpy直接读取netCDF4库不行吗？当然可以，但你会立刻陷入维度管理的泥潭。一个典型的气象再分析数据（如ERA5）可能包含(time: 365, latitude: 721, longitude: 1440)三个维度。当你只想提取北京地区夏季的平均温度时，你需要手动计算经纬度索引、处理时间切片、管理缺失值……代码会变得冗长且难以阅读。

xarray的核心优势在于引入了带标签的数组（DataArray）和数据集（Dataset）。每个维度都有明确的名称（如lat,lon,time）和坐标值。这意味着你可以用ds.sel(lat=39.9, lon=116.4, method=‘nearest’)这样的语义化方式提取数据，而不是data[250, 600]这样令人困惑的魔法数字。它自动处理了netCDF文件中的变量、属性和坐标，让你能专注于科学分析本身，而不是数据索引的算术。

注意：xarray底层依赖于netCDF4或h5netcdf等库来执行实际的I/O操作。通常直接pip install xarray netcdf4即可，xarray会自动选择可用的后端。

2.2 Cartopy：让地图投影不再是噩梦

气象图十有八九需要画在地图上。matplotlib自带的Basemap工具包已经停止维护，而Cartopy是其现代替代品，也是目前业界的标准选择。它的强大之处在于：

1. 丰富的投影系统：从全球常用的PlateCarree（等经纬度）、Robinson，到区域性的LambertConformal、PolarStereographic，Cartopy都提供了原生支持。这对于需要保持面积、方向或距离特性的分析图至关重要。
2. 便捷的地理特征添加：海岸线、国界、河流、湖泊等地理特征，往往只需一行代码就能添加，并且能自动适配你设置的投影。
3. 与matplotlib无缝集成：Cartopy通过定义新的GeoAxes子类来扩展matplotlib，这意味着你熟悉的几乎所有matplotlib绘图函数和样式设置方法都能继续使用。

将xarray和Cartopy结合，前者负责高效、语义化的数据操作，后者负责专业、准确的地理可视化，构成了处理气象绘图任务的“黄金搭档”。

2.3 环境搭建与包管理建议

为了避免版本冲突，强烈建议使用conda来管理你的科学计算环境，特别是Cartopy的安装，因为它有一些地理数据库的依赖。

# 创建一个新的conda环境 conda create -n meteorology-viz python=3.9 conda activate meteorology-viz # 通过conda-forge频道安装所有核心包（最稳定） conda install -c conda-forge xarray netcdf4 matplotlib cartopy jupyter # 可选但推荐的包：用于更美观的配色和进度条 conda install -c conda-forge cmocean tqdm

使用conda-forge频道能确保所有包的二进制依赖兼容。如果你习惯用pip，在Linux或macOS上安装Cartopy可能会遇到GEOS、PROJ等库的编译问题，conda帮你省去了这些麻烦。

3. 实战第一步：深度解析netCDF文件结构

在画图之前，我们必须先“读懂”数据。很多绘图错误都源于对数据结构的误解。让我们用一个实际的例子来演练。

假设我们下载了一个ERA5再分析数据的月平均海平面气压（msl）文件era5_monthly_msl_2022.nc。第一步不是直接绘图，而是彻底检查它。

import xarray as xr # 使用xarray打开netCDF文件 file_path = 'era5_monthly_msl_2022.nc' ds = xr.open_dataset(file_path) print(ds)

执行这段代码后，你会看到一个结构化的输出，这是理解数据的钥匙。输出通常包含以下几个部分：

Variables (数据变量)：这是文件的核心，比如msl。你需要关注它的dimensions，例如(time: 12, latitude: 721, longitude: 1440)，这告诉你它是一个包含12个时间点、全球高分辨率格点的三维数组。还要看它的units属性，这里是Pa，这关系到绘图时颜色栏的标注。

Coordinates (坐标)：定义了每个维度的具体数值。

• longitude: 从0到359.75，间隔0.25度。这里有一个关键点：很多数据集的经度范围是0-360°，而Cartopy绘图时通常期望-180°到180°。我们后续需要处理。
• latitude: 从90到-90，间隔0.25度。
• time: 12个时间点，格式为datetime64[ns]。

Attributes (全局属性)：描述了数据集的来源、版本、历史等元数据，对于论文绘图时撰写图注非常重要。

实操心得：open_datasetvsload_dataset

• xr.open_dataset()：这是懒加载。它只读取元数据和数据结构，而不将庞大的数据数组立即读入内存。这对于动辄几个GB的气象数据来说是默认且推荐的方式，只有在实际用到数据（如切片、计算、绘图）时，相应的部分才会被加载。
• xr.load_dataset()：这是立即加载。它会将所有数据一次性读入内存。除非你确认文件很小，否则不要轻易使用，否则可能导致内存溢出。

更进一步的检查可以这样做：

# 查看某个具体变量的详细信息 print(ds['msl']) # 查看经纬度坐标的具体值 print(ds.longitude.values[:10]) # 查看前10个经度值 print(ds.latitude.values[:10]) # 查看前10个纬度值 # 检查是否有缺失值及其填充方式 print(ds['msl'].encoding) # 查看编码信息，如scale_factor, add_offset, _FillValue

4. 数据预处理：绘图前的关键清洗与转换

直接从数据集中取数据来画，很可能得到一张奇怪或错误的图。以下是几个必须检查的预处理步骤。

4.1 经度坐标转换：从0-360°到-180-180°

大多数全球绘图库（包括Cartopy）默认使用-180°到180°的经度坐标系（本初子午线在中间）。而很多模式输出或再分析数据（如ERA5, CMIP6）为了计算方便，使用0°到360°的经度坐标系（本初子午线在边缘）。

如果你用0-360°的数据直接绘图，你会发现地图被“切”了一刀，格林威治0度经线两侧的数据不连续。解决方法是用assign_coords和roll函数：

# 方法：将经度从0-360调整到-180-180 # 首先，给经度坐标重新赋值 ds = ds.assign_coords(longitude=(((ds.longitude + 180) % 360) - 180)) # 然后，按新的经度顺序对数据重新排序 ds = ds.sortby('longitude')

为什么这么做？(ds.longitude + 180) % 360)这个操作将经度范围平移到180-540，再对360取模，得到0-360的范围，最后减去180，就得到了-180到180的范围。sortby确保数据数组按照新的经度值升序排列，这是Cartopy正确渲染所必需的。

4.2 单位换算与变量选择

数据的单位可能不适合直接展示。例如，海平面气压原始单位是帕斯卡（Pa），但在气象学中常用百帕（hPa）或毫巴（mb）。1 hPa = 100 Pa。

# 将海平面气压从Pa转换为hPa if ds['msl'].units == 'Pa': ds['msl'] = ds['msl'] / 100.0 ds['msl'].attrs['units'] = 'hPa' # 记得更新属性

同时，你可能只需要某个时间点或某个月平均的数据：

# 选择2022年7月的数据（假设time坐标是datetime类型） data_july = ds['msl'].sel(time='2022-07-01', method='nearest') # 或者计算2022年的年平均 data_annual_mean = ds['msl'].mean(dim='time')

4.3 处理缺失值与无效值

NetCDF文件通常会用特定的_FillValue或missing_value属性来标记无效数据。xarray在读取时通常会将其转换为NaN。但有时需要手动处理：

# 检查是否存在填充值属性 fill_value = ds['msl'].encoding.get('_FillValue', None) if fill_value is not None: # xarray通常已自动处理，但可再次确认 ds['msl'] = ds['msl'].where(ds['msl'] != fill_value)

5. 核心绘图实战：绘制全球海平面气压分布图

现在，数据已经准备就绪，让我们绘制一张专业的全球海平面气压分布填色图。我们将一步步拆解，并解释每个参数的意义。

5.1 创建地图投影与子图

首先，导入必要的库并创建画布和地理坐标轴。

import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeature import numpy as np # 1. 创建图形和坐标轴，并指定投影 # 这里使用PlateCarree投影（等经纬度），这是最常用的全球投影 fig = plt.figure(figsize=(14, 8)) # 设置一个较大的画布 ax = plt.axes(projection=ccrs.PlateCarree()) # 关键：创建GeoAxes # 2. 添加地理特征，让地图更易读 ax.add_feature(cfeature.COASTLINE, linewidth=0.8) # 海岸线 ax.add_feature(cfeature.BORDERS, linewidth=0.5, linestyle=':') # 国界线（虚线） ax.add_feature(cfeature.OCEAN, color='lightblue', alpha=0.3) # 海洋填充（浅蓝色，半透明） ax.add_feature(cfeature.LAND, color='lightgray', alpha=0.3) # 陆地填充（浅灰色，半透明） # 3. 设置网格线 gl = ax.gridlines(draw_labels=True, linewidth=0.5, color='gray', alpha=0.5, linestyle='--') gl.top_labels = False # 不显示顶部标签 gl.right_labels = False # 不显示右侧标签 gl.xlabel_style = {'size': 10} gl.ylabel_style = {'size': 10}

关键参数解读：

• projection=ccrs.PlateCarree()：这告诉Cartopy我们使用等经纬度投影。这里有一个极易混淆的点：PlateCarree既是数据本身的坐标系（我们数据是等经纬度格点），也是地图的投影方式。对于其他投影（如兰伯特投影），数据可能需要被转换。
• ax.gridlines(draw_labels=True)：自动绘制并标注经纬度网格线。通过top_labels和right_labels控制标签位置，避免重叠。

5.2 绘制数据填色图

接下来，将我们处理好的数据画到地图上。

# 假设我们绘制2022年7月的平均海平面气压 plot_data = data_july # 从4.2节获得的数据 # 4. 绘制填色图 # `transform=ccrs.PlateCarree()` 至关重要！它声明了数据本身的坐标系。 im = ax.contourf(plot_data.longitude, plot_data.latitude, plot_data, levels=60, # 颜色分级数，越多越平滑 cmap='RdBu_r', # 使用红蓝渐变色，_r表示反转 transform=ccrs.PlateCarree()) # 5. 添加等值线（可选，使高低压中心更明显） contour = ax.contour(plot_data.longitude, plot_data.latitude, plot_data, levels=15, # 等值线数量 colors='black', linewidths=0.5, transform=ccrs.PlateCarree()) ax.clabel(contour, inline=True, fontsize=8, fmt='%1.0f') # 在等值线上标注数值 # 6. 添加颜色栏 cbar = plt.colorbar(im, ax=ax, orientation='horizontal', pad=0.05, shrink=0.8) cbar.set_label('Sea Level Pressure (hPa)', fontsize=12)

核心技巧：transform参数这是Cartopy绘图中最容易出错的地方。transform参数用于告诉Cartopy你的数据是在什么地理坐标系下。我们的数据是规则的经纬度网格，所以transform=ccrs.PlateCarree()。无论地图的投影（projection）设置成什么（比如Robinson），只要数据是经纬度的，这个transform参数就不变。Cartopy会在内部负责将数据从PlateCarree坐标系转换到目标投影上。如果忘记设置或设置错误，数据可能会被画在错误的地理位置上。

5.3 美化与输出

最后，添加标题，调整布局并保存图片。

# 7. 添加标题 ax.set_title('Global Sea Level Pressure - July 2022 (ERA5)', fontsize=16, pad=20) # 8. 调整图形布局，确保所有元素都能显示 plt.tight_layout() # 9. 保存图片（高分辨率，适合出版物） output_path = 'global_msl_july_2022.png' plt.savefig(output_path, dpi=300, bbox_inches='tight') print(f'图表已保存至：{output_path}') # 10. 显示图表（在Jupyter Notebook或脚本中） plt.show()

至此，一张标准的全球气象要素分布图就完成了。它包含了清晰的地理背景、规整的网格、直观的色标和等值线，以及可供引用的标题和文件来源。

6. 进阶技巧与常见问题排查

掌握了基础绘图后，你可能会遇到一些更复杂的需求或棘手的问题。这里分享一些实战中积累的经验。

6.1 绘制区域子图

研究往往聚焦于特定区域，如东亚、北大西洋。Cartopy可以方便地设置图形范围。

# 在创建坐标轴时，通过`map_extent`参数设置区域 [西经, 东经, 南纬, 北纬] extent = [70, 140, 10, 60] # 东亚区域 fig, ax = plt.subplots(figsize=(10, 8), subplot_kw={'projection': ccrs.PlateCarree()}) ax.set_extent(extent, crs=ccrs.PlateCarree()) # 关键：设置图形显示范围 # 然后在此ax上绘图，步骤同上... # 添加地理特征时，可以只添加分辨率更高的特征，如河流 ax.add_feature(cfeature.COASTLINE) ax.add_feature(cfeature.RIVERS, linewidth=0.5, edgecolor='blue')

6.2 处理高分辨率数据与性能优化

全球高分辨率数据（如0.25°）在绘图时，contourf可能会很慢。可以考虑以下策略：

1. 数据裁剪：在绘图前，先用sel或isel将数据裁剪到目标区域，大幅减少数据量。

   regional_data = plot_data.sel(latitude=slice(10, 60), longitude=slice(70, 140))

2. 使用pcolormesh替代contourf：对于单纯的填色图，pcolormesh速度更快，但它不进行插值，是直接渲染网格。

   im = ax.pcolormesh(regional_data.longitude, regional_data.latitude, regional_data, cmap='RdBu_r', shading='auto', # ‘auto’自适应网格渲染 transform=ccrs.PlateCarree())

3. 降低绘图分辨率：对于快速预览，可以对数据先进行粗化处理。

   coarse_data = plot_data.coarsen(latitude=2, longitude=2, boundary='trim').mean()

6.3 常见错误与解决方案速查表

6.4 配色方案的选择

配色不仅关乎美观，更影响信息的准确传达。避免使用彩虹色系（jet），因为它在感知上不均匀，可能误导对数据梯度的判断。气象海洋学领域有专门的配色库cmocean。

import cmocean # 使用cmocean中针对热力学变量的配色 im = ax.contourf(..., cmap=cmocean.cm.thermal) # 温度场 # 或针对差异的配色 im = ax.contourf(..., cmap=cmocean.cm.balance) # 正负差异场（如异常图）

如果无法安装cmocean，matplotlib内置的viridis、plasma（序列数据）、RdBu_r、coolwarm（发散数据）也是很好的选择。

从打开一个陌生的netCDF文件，到生成一张可用于分析或展示的专业气象图表，这个过程需要理解数据、工具和地理可视化的基本原理。xarray让你从繁琐的维度索引中解放出来，Cartopy则为你提供了绘制专业地图的武器。核心在于记住数据坐标系（transform）与地图投影（projection）的区别，并养成良好的数据检查习惯。在接下来的系列文章中，我们会探讨更多主题，比如绘制垂直剖面图、风矢量和流线图、以及制作动画等。