1. 项目缘起从一张“光秃秃”的流场图说起如果你用Python的Cartopy、Matplotlib画过南海的海流矢量图大概率会经历过这样一个阶段费了九牛二虎之力终于把经向流、纬向流的数据读进来用quiver或者quiverkey画出了箭头地图边界也设置好了标题、坐标轴一应俱全。满心欢喜地保存图片发给同事或者放进报告里结果对方第一句话就问“这箭头跑多快单位是什么” 你才猛然发现图上缺了一个最关键的元素——比例尺Scale或者说图例Legend用来指明“一个单位长度的箭头代表多大的流速”。没有比例尺的流场图就像没有刻度的尺子再精美的箭头也失去了定量的意义。读者无法判断南海黑潮的流速是0.5米/秒还是1.5米/秒也无法比较不同区域流动的强弱。这个看似微小的细节恰恰是科学可视化从“能看”到“专业”的关键一步。本次要解决的就是如何为南海年平均海流矢量图添加一个清晰、准确且美观的比例尺。我们手头的数据通常是NetCDF格式的网格化数据包含经度、纬度、纬向流速U东西方向和经向流速V南北方向四个变量。我们的目标不是简单地调用plt.quiver而是要在南海区域地图上绘制能反映真实地理空间关系的流场箭头并附上一个指明物理量大小和方向的参考箭头。这个需求在海洋学、气象学领域非常普遍但Matplotlib默认的quiverkey在复杂地图投影上定位和缩放时常常会遇到箭头变形、位置难以控制等问题。接下来我将分享一套经过实战检验的、稳定可靠的解决方案。2. 核心工具链搭建Cartopy与Matplotlib的深度协同在开始画图之前搭建一个正确且高效的环境是重中之重。很多人卡在第一步不是因为代码多难而是库的版本不匹配或者依赖没装全。这里我强烈推荐使用conda来管理环境它能很好地处理地理绘图库复杂的二进制依赖。2.1 环境配置与避坑指南首先创建一个专属的绘图环境conda create -n ocean_plot python3.9 conda activate ocean_plot接下来安装核心库。这里有一个关键坑点不要直接用pip install cartopy。Cartopy依赖GEOS和PROJ等C库用pip安装极易编译失败。务必使用conda的conda-forge频道安装它能自动处理好所有系统依赖。conda install -c conda-forge cartopy matplotlib netCDF4 xarray numpy这条命令一次性安装了我们的核心工具链Cartopy地理绘图的神器负责地图投影、海岸线、边界等。没有它在Python里画专业地图会非常痛苦。Matplotlib绘图的基石Cartopy构建在其之上。NetCDF4 / Xarray用于读取和处理海洋模式常用的NetCDF格式数据。Xarray提供了更友好、类似pandas的接口强烈推荐。NumPy数值计算基础。安装完成后强烈建议在Python交互环境里快速验证一下Cartopy的关键功能是否能正常导入和运行特别是涉及投影的部分可以避免在写了大段代码后才发现环境有问题。2.2 数据准备与读取策略假设我们有一个名为south_china_sea_current_annual_mean.nc的NetCDF文件。使用Xarray来打开它比直接用NetCDF4库要直观得多import xarray as xr # 打开数据集 ds xr.open_dataset(south_china_sea_current_annual_mean.nc) # 查看数据结构 print(ds)典型的输出会显示数据变量和坐标例如Dimensions: (longitude: 144, latitude: 73, depth: 1, time: 1) Coordinates: * longitude (longitude) float32 100.0 100.25 100.5 ... 134.75 135.0 * latitude (latitude) float32 0.0 0.25 0.5 ... 18.0 18.25 18.5 * depth (depth) float32 5.0 * time (time) datetime64[ns] 2020-07-02 Data variables: uo (time, depth, latitude, longitude) float32 ... vo (time, depth, latitude, longitude) float32 ...这里uo和vo分别代表纬向和经向流速。我们的目标是绘制表层比如第一个深度层的年平均流场。由于时间维可能也只有一层年平均我们需要将其数据提取出来并通常转换为NumPy数组供Matplotlib使用import numpy as np # 提取表层索引0年平均数据并压缩掉大小为1的维度 u ds[uo].isel(depth0, time0).values # 形状变为 (latitude, longitude) v ds[vo].isel(depth0, time0).values lon ds[longitude].values lat ds[latitude].values # 创建经纬度的网格矩阵对于quiver是必需的 lon2d, lat2d np.meshgrid(lon, lat)一个重要经验海洋模式数据常常包含缺失值FillValue在NetCDF中通常用一个大数如1e20或9.96921e36表示。直接绘图会导致箭头异常。务必在绘图前处理# 假设缺失值标记为 1e20 missing_value 1e20 u np.ma.masked_values(u, missing_value) v np.ma.masked_values(v, missing_value)使用NumPy的掩码数组masked array可以优雅地处理这个问题绘图时这些位置会被自动忽略。3. 绘制基础南海流场地图有了干净的数据我们就可以开始绘制地图了。这一步的重点是选择合适的投影和设置美观的底图。3.1 地图投影与区域裁剪对于南海区域PlateCarree等经纬度投影虽然简单但会导致区域形状在高纬度略有拉伸。Mercator墨卡托投影在低纬度地区变形小更适合南海。这里我们使用PlateCarree作为数据坐标系并用Mercator来显示地图这是Cartopy的常规操作。import cartopy.crs as ccrs import cartopy.feature as cfeature import matplotlib.pyplot as plt # 创建图形和坐标系 # 使用PlateCarree表示数据本身的坐标系经度/纬度 # 使用Mercator作为地图的显示投影 fig plt.figure(figsize(12, 10)) ax plt.axes(projectionccrs.Mercator(central_longitude115.0)) # 设置地图范围覆盖整个南海主要区域 ax.set_extent([105, 125, 0, 25], crsccrs.PlateCarree()) # 添加地理特征 ax.add_feature(cfeature.LAND, colorlightgray) ax.add_feature(cfeature.OCEAN, colorazure) ax.add_feature(cfeature.COASTLINE, linewidth0.8) ax.add_feature(cfeature.BORDERS, linestyle:, linewidth0.5) # 添加经纬度网格线 gl ax.gridlines(draw_labelsTrue, linewidth0.5, colorgray, alpha0.5, linestyle--) gl.top_labels False # 关闭顶部标签 gl.right_labels False # 关闭右侧标签 gl.xlabel_style {size: 10} gl.ylabel_style {size: 10}3.2 绘制海流矢量箭头现在在准备好的地图上绘制流场。这里直接使用Matplotlib的quiver函数但关键在于要指定transform参数告诉Cartopy这些箭头所在的位置是基于什么坐标系的。我们的lon2d和lat2d是经纬度所以对应ccrs.PlateCarree()。# 绘制流场箭头 # 为了图面清晰通常需要对数据进行稀疏化处理比如每隔3个点取一个 stride 3 Q ax.quiver(lon2d[::stride, ::stride], lat2d[::stride, ::stride], u[::stride, ::stride], v[::stride, ::stride], transformccrs.PlateCarree(), # 声明箭头位置的坐标系 anglesuv, # 让箭头方向由(U,V)分量直接决定 scale_unitsinches, # 比例尺单位与后续scale参数配合 scaleNone, # 先不设置自动缩放便于后续手动控制 colordarkblue, width0.002, # 箭头柄的宽度 headwidth3, # 箭头头部的宽度倍数 headlength4, # 箭头头部的长度倍数 headaxislength3.5) # 箭头头部到转折点的长度 ax.set_title(南海年平均海流, fontsize16, pad20)此时一张带有南海地形和流场箭头的图就初步完成了。但是图上还没有告诉我们箭头到底代表多快的流速。这就是接下来要解决的核心问题。4. 为流场图添加定量的比例尺添加比例尺是本文的重中之重。Matplotlib提供了quiverkey函数但在Cartopy变换过的地图上直接使用很容易出现比例尺箭头显示异常比如变得巨大或看不见或者位置错乱。其根本原因在于quiverkey默认在“数据坐标系”下工作而我们的地图经过了投影变换屏幕上的一个“单位长度”对应的真实地理长度在不同位置是不同的。4.1 传统quiverkey方法的局限与调试我们先看看直接使用quiverkey会发生什么# 尝试在图形坐标 (0.05, 0.05) 处添加比例尺代表 0.5 m/s qk ax.quiverkey(Q, 0.05, 0.05, 0.5, r$0.5 \, m/s$, labelposE, coordinatesaxes) plt.show()你会发现这个比例尺箭头可能大小完全不对或者根本看不到。因为coordinatesaxes使用的是轴坐标0到1但箭头的缩放比例scale与主quiver图关联而主图的箭头长度经过了地图投影的扭曲计算导致这个参考箭头无法正确映射。4.2 稳健的比例尺绘制方案自定义参考箭头经过多次踩坑我总结出一个稳定可靠的方法不在原quiver对象上添加quiverkey而是单独在同一个地图投影下绘制一个“虚拟”的参考箭头。这个箭头的U、V分量就是我们想要代表的流速值将其绘制在地图上的一个固定位置比如右下角空白海域。这个方法的精髓在于这个参考箭头和主流场箭头处于完全相同的绘图环境下相同的transform相同的scale参数因此它们的长度比例关系是绝对正确的。步骤一确定参考流速和绘制位置首先我们需要决定比例尺代表多大的流速。可以计算主流场数据的典型值比如95%分位数避免用最大值导致比例尺过小。假设我们决定用0.5 m/s。然后在南海地图的右下角选取一个点作为参考箭头的起点例如经度122°E纬度8°N。步骤二绘制参考箭头# 定义参考箭头的起点经纬度 ref_lon, ref_lat 122.0, 8.0 # 定义参考箭头代表的流速分量 (U, V)。这里我们画一个指向正东的箭头代表0.5 m/s ref_u, ref_v 0.5, 0.0 # 在相同的地图投影和坐标系下单独绘制这个参考箭头 # 注意这里的scale参数必须与主quiver图完全一致如果主图用了scaleNone这里也要用scaleNone。 # 但更常见的做法是在主quiver中指定一个scale值使箭头长度适中。 # 我们先为主quiver设置一个合适的scale # ... 回到主quiver绘制部分修改scale参数 ... # Q ax.quiver(..., scale50, ...) # 例如scale50这个值需要根据你的数据尝试调整 # 假设主quiver的scale已设为50那么参考箭头也使用相同的scale ref_quiver ax.quiver(ref_lon, ref_lat, ref_u, ref_v, transformccrs.PlateCarree(), anglesuv, scale_unitsinches, scale50, # 必须与主图一致 colorred, # 用不同颜色突出显示 width0.003, headwidth4, zorder5) # 确保比例尺在最上层 # 在参考箭头旁边添加文本标签 # 使用text函数并指定文本的坐标系也为地图坐标系 ax.text(ref_lon 0.5, ref_lat, r$0.5 \, m/s$, transformccrs.PlateCarree(), verticalalignmentcenter, fontsize12, colorred, bboxdict(boxstyleround,pad0.2, facecolorwhite, alpha0.8))步骤三确定主图的scale值上面的代码有个前提主quiver的scale参数需要设置一个合适的值。scale的含义是scale值越大箭头长度越短。它是一个经验值需要通过反复尝试来获得最佳的视觉密度。一个实用的方法是# 在主quiver绘制前估算一个初始scale # 计算流速的模的平均值或中位数作为参考 speed np.sqrt(u**2 v**2) median_speed np.ma.median(speed) # 使用掩码数组的中位数 print(f流速中位数: {median_speed} m/s) # 初始scale可以设为 50 / median_speed 之类的经验公式然后手动调整 initial_scale 30 # 从30开始尝试 Q ax.quiver(..., scaleinitial_scale, ...)运行代码查看箭头是否过密或过疏然后调整scale值直到图面清晰美观。记住这个最终的scale值在绘制参考箭头时必须使用完全相同的值。4.3 进阶创建更专业的比例尺图例上述方法画出了一个正确的比例尺但美观度可以进一步提升。我们可以模仿专业绘图软件绘制一个带框和说明的图例。思路是在地图坐标系外的轴坐标系axescoordinates里用FancyArrowPatch或普通的arrow画一个箭头并配上文字和边框。但要注意这个箭头的长度需要与地图上的箭头有“视觉上”的对应关系。一个更工程化的方法是获取主quiver中某个代表性箭头的屏幕显示长度单位英寸或点。在轴坐标系中绘制一个具有相同屏幕长度的箭头。这涉及到坐标变换实现起来稍复杂。一个更简单的妥协方案是我们放弃“绝对精确的屏幕长度匹配”而是追求“视觉上的和谐与信息明确”。我们可以直接在轴坐标系画一个看起来大小合适的箭头并用文字明确标注其代表的流速值。虽然不够完美但在大多数科研报告中已完全够用。import matplotlib.patches as mpatches from matplotlib.offsetbox import AnchoredOffsetbox, TextArea, HPacker # 在轴坐标系 (0.85, 0.08) 处创建一个 anchored offsetbox 作为图例 # 这个位置是 (x, y)基于整个图的左下角为(0,0)右上角为(1,1) legend_x, legend_y 0.85, 0.08 # 创建一个水平排列的容器 ob AnchoredOffsetbox( loclower center, # 锚点位置这里用lower center但会被bbox_to_anchor覆盖 childNone, pad0.4, frameonTrue, bbox_to_anchor(legend_x, legend_y), bbox_transformax.transAxes, # 关键使用轴坐标系变换 borderpad0.5, ) # 创建图例内容一个箭头 文字 # 使用 FancyArrowPatch 在容器内部坐标系画箭头 # 容器内部坐标系是相对于容器左下角的 (0,0) 到 (1,1) arrow mpatches.FancyArrowPatch((0.1, 0.5), (0.6, 0.5), arrowstyle-, mutation_scale15, linewidth2, colorred) text TextArea(r$0.5 \, \mathrm{m \, s^{-1}}$, textpropsdict(size12)) # 将箭头和文字水平打包 packer HPacker(children[arrow, text], sep10, # 间距 pad0, packTrue) ob.child packer ax.add_artist(ob)这种方法将比例尺图例固定在图框的右下角不受地图缩放和平移的影响并且样式统一看起来非常专业。5. 流程优化与完整代码整合将以上所有步骤整合并加入一些优化如颜色映射表示流速大小、更好的图形布局等就得到了一套完整的绘图流水线。5.1 完整脚本示例import xarray as xr import numpy as np import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeature import matplotlib.patches as mpatches from matplotlib.offsetbox import AnchoredOffsetbox, TextArea, HPacker # 1. 数据读取与预处理 ds xr.open_dataset(south_china_sea_current_annual_mean.nc) u ds[uo].isel(depth0, time0).values v ds[vo].isel(depth0, time0).values lon ds[longitude].values lat ds[latitude].values # 处理缺失值 missing_value 1e20 u np.ma.masked_values(u, missing_value) v np.ma.masked_values(v, missing_value) # 创建网格 lon2d, lat2d np.meshgrid(lon, lat) # 计算流速大小用于颜色映射 speed np.ma.sqrt(u**2 v**2) # 2. 创建地图 fig plt.figure(figsize(14, 10)) ax plt.axes(projectionccrs.Mercator(central_longitude115.0)) ax.set_extent([105, 125, 0, 25], crsccrs.PlateCarree()) # 添加底图特征 ax.add_feature(cfeature.LAND, colorwheat, edgecolork, linewidth0.5) ax.add_feature(cfeature.OCEAN, colorlightcyan) ax.add_feature(cfeature.COASTLINE, linewidth0.8) ax.add_feature(cfeature.BORDERS, linestyle:, linewidth0.5, alpha0.7) # 添加网格线 gl ax.gridlines(draw_labelsTrue, linewidth0.5, colorgray, alpha0.5, linestyle:) gl.top_labels False gl.right_labels False # 3. 绘制彩色流场图 stride 3 # 使用颜色表示流速大小箭头表示方向 scatter ax.scatter(lon2d[::stride, ::stride], lat2d[::stride, ::stride], cspeed[::stride, ::stride], cmapRdYlBu_r, # 红-黄-蓝反转后蓝色表示高速 s5, # 点的大小 transformccrs.PlateCarree(), alpha0.7) # 绘制箭头 Q ax.quiver(lon2d[::stride, ::stride], lat2d[::stride, ::stride], u[::stride, ::stride], v[::stride, ::stride], transformccrs.PlateCarree(), anglesuv, scale_unitsinches, scale40, # 经过调试的scale值使箭头长度适中 colork, # 箭头用黑色与彩色点形成对比 width0.0025, headwidth3.5, headlength4.5, headaxislength4.0, zorder3) # 添加颜色条 cbar plt.colorbar(scatter, axax, orientationvertical, pad0.05, shrink0.7) cbar.set_label(流速 (m/s), fontsize12) # 4. 添加自定义比例尺图例轴坐标系方法 legend_x, legend_y 0.82, 0.12 ob AnchoredOffsetbox( loclower center, pad0.4, frameonTrue, bbox_to_anchor(legend_x, legend_y), bbox_transformax.transAxes, borderpad0.8, framealpha0.9, # 图例框透明度 edgecolorgray ) # 在图例框内绘制箭头和文字 arrow mpatches.FancyArrowPatch((0.15, 0.5), (0.65, 0.5), arrowstyle-, mutation_scale20, # 控制箭头大小 linewidth2.5, colordarkred) text TextArea(r$\mathbf{0.5 \, m \, s^{-1}}$, textpropsdict(size11, weightbold)) packer HPacker(children[arrow, text], sep12, pad0, packTrue) ob.child packer ax.add_artist(ob) # 5. 添加标题和修饰 ax.set_title(南海年平均海流场颜色表示流速大小, fontsize16, fontweightbold, pad20) # 6. 保存与显示 plt.tight_layout() plt.savefig(south_china_sea_current_with_scale.png, dpi300, bbox_inchestight) plt.show()5.2 关键调试经验与注意事项scale参数的调试这是影响图面美观度的最关键参数。没有固定公式需要根据你的数据范围和图形大小手动调整。建议写一个循环或使用交互模式plt.ion()快速尝试几个值如20, 30, 40, 50, 60观察箭头密度。稀疏化采样 (stride)高分辨率数据如果不做稀疏化quiver会绘制数十万个箭头导致图形文件巨大、渲染极慢且视觉上是一片混乱的“毛毯”。stride值需要根据你的数据分辨率经纬网格数和图幅大小决定通常尝试2, 3, 4, 5直到箭头清晰可辨又不至于太稀疏。颜色映射 (cmap)表示流速大小时选择合适的色带很重要。viridis,plasma是通用的连续色带。海洋学中常用RdYlBu_r红黄蓝反转表示冷暖或高低。避免使用如jet这类虽然鲜艳但感知不均匀的色带。图形保存保存高分辨率图片时使用bbox_inchestight可以自动裁剪掉图形周围的白边让图片更紧凑。dpi300足以满足大部分出版要求。性能问题如果数据量极大绘制和保存可能会很慢。除了稀疏化还可以考虑将底图海岸线、陆地等保存为背景图片或者使用更底层的绘图方法进行优化。对于超大规模数据专业的可视化软件如Paraview可能是更好的选择。通过以上步骤你不仅能为南海海流图加上一个标准的比例尺更能掌握一套在复杂地图投影下进行科学矢量可视化的完整方法。这套方法同样适用于风场、应力场等任何矢量数据的制图。记住好的可视化是科研沟通的桥梁而比例尺就是这座桥上不可或缺的刻度。
)
