gcmfaces工具箱：Matlab/Octave处理立方球网格海洋模式数据指南

简介：gcmfaces 是一款面向 Matlab 与 Octave 的开源工具箱，专为全球气候模型（GCM）海洋环流数据处理而设计。它帮助科研人员高效读取、管理、可视化和计算大规模分块网格数据，支持物理量诊断与并行加速，尤其适合海洋科学与气候模拟方向的研究者。资源包含完整源码与配套资料，共 317 个文件；以 296 个 .m 函数/脚本为主体，覆盖数据读取、网格操作与输出处理等核心功能，另有 7 个 .rst、4 个 PDF 文档作为使用指南和 API 参考，以及示例数据、构建脚本等辅助文件，整体压缩包约 3MB，轻量且可直接使用。目前已有 118 人学习下载。借助该资源，用户可系统掌握 face 分块数据的读写与可视化、常用物理量计算、并行处理等基本流程，内置 diags_set_A.m、process2nctiles.m 等示例也能辅助快速上手，方便开展 GCM 海洋环流的诊断分析与后续二次开发。

1. gcmfaces：让 Matlab 和 Octave 理解六个立方球面的海洋模式数据

拿到一个命名里带“zip”的 gcmfaces 资料包时，你多半正准备处理全球海洋环流模型（比如 MITgcm）在立方球网格上的输出。这类输出不是我们习惯的等经纬度矩阵，而是按六个立方球面（cubed-sphere face）分别保存的二维场。若直接用普通 Matlab 数组去算全球平均、插值或画图，六个面之间的接缝和极区畸变会把结果污染成一条条锯齿。gcmfaces 就是专为这个场景出现的 Matlab/Octave 工具箱：它把六个面的数据封装成对象，在加减乘除、mean、sum、插值等常用操作背后替你维护面索引和网格几何，最终把“处理全球模式场”的代码压缩到和单块数组一样短。适合做模式诊断、气候统计和卫星产品验证的工程师与科研人员。

2. 安装 gcmfaces：在 Matlab 与 Octave 下把工具箱变成可调用模块

很多教程把安装简单说成“解压后加路径”，放到自己机器上却报错，大多是因为没有处理根目录下的子模块、mex 目录和 Octave 依赖。这个工具箱本质是纯 m 脚本与少量辅助函数的组合，不依赖优化工具箱或机器人工具箱，安装步骤比想象中少，但路径管理必须干净。

2.1 解压后先看 Contents.m 和 Examples 目录

先建一个固定目录，再解压。Linux/macOS 下：

mkdir -p ~/software unzip gcmfaces_matlab_octave_toolbox.zip -d ~/software/gcmfaces

Windows 用户直接解压到D:\software\gcmfaces，不要放在“下载”目录里，因为后续savepath会要求拥有写入权限。进入解压目录后，先确认两样东西：根目录有没有Contents.m，以及有没有Examples文件夹。Contents.m是 Matlab/Octave 识别“工具箱”的标记，有了它，执行help gcmfaces才能显示模块列表；Examples里的脚本则是比文档更贴近真实调用方式的入口。

下表总结了安装阶段最需要留意的几类文件：

提示：如果根目录没有Contents.m，不要急着下载另一个版本。先看顶层有哪些目录，再通过which逐个验证目标文件是否进入路径，很多时候只是解压层级多了一层。

2.2 用 startup.m 固化 addpath，避免每次重新配

我习惯把所有自装工具箱统一挂到用户级startup.m里，而不是每次都敲一段addpath。在用户主目录下建一个：

% startup.m —— 用户级启动配置 toolboxRoot = fullfile(getenv('HOME'), 'software', 'gcmfaces'); if exist(toolboxRoot, 'dir') addpath(genpath(toolboxRoot)); end

参数说明：genpath会递归加入该目录下所有子目录，对 gcmfaces 这种带examples、utils、grid多级目录的结构最省事；getenv('HOME')在 Windows 上可能为空，你可以换成fullfile('D:', 'software', 'gcmfaces')。写入后执行savepath，把当前路径保存到pathdef.m。Matlab 在 R2023b 之后的安装教程里，很多都要求用setenv或userpath管理用户目录，但startup.m始终是通用做法；Octave 5.x 以后同样支持这个文件。

2.3 Octave 下核对 netcdf 依赖与 mex 目录

如果你是 Octave 用户，安装后第一件事不是跑数据，而是确认 netcdf 包可用：

pkg list pkg load netcdf

pkg list会列出已安装包，若 netcdf 未安装，用pkg install -forge netcdf从 Octave Forge 拉取。这一步常见报错是编译器版本与包不匹配，Octave 会提示用mkoctfile -p检查编译器。老版本 gcmfaces 还可能在mex子目录里提供编译好的二进制，Matlab 版和 Octave 版不能混用。解压后若看到mex/octave和mex/matlab两套文件，请按当前解释器只把对应目录加进路径。同时加入两个目录会在运行时出现“Invalid MEX-file”的报错，因为平台特定符号表对不上。

2.4 一个最小自检脚本：确认 gcmfaces_global 可运行

路径配好不等于能跑，我用这个脚本做快速验证：

function ok = check_gcmfaces() ok = false; assert(~isempty(which('gcmfaces_global')), 'gcmfaces_global 不在搜索路径'); fprintf('gcmfaces 主入口：%s\n', which('gcmfaces_global')); gcmfaces_global; ok = true; fprintf('全局初始化完成\n'); end

逻辑说明：which('gcmfaces_global')能找到文件，说明gcmfaces_global.m所在的目录确实在搜索路径中；随后调用gcmfaces_global触发工具箱的全局状态初始化。gcmfaces 大量脚本依赖一个全局网格对象作为“当前网格”，这个函数就是负责把它准备好的。若这里失败，后面的read_bin、interp都会跟着报错，所以不要继续配置别的模块。

3. gcmfaces 的数据组织：六个面的对象如何被 Matlab/Octave 所理解

安装只是第一步，真正让人卡住的是理解“对象”这个抽象层。gcmfaces 没有把六个面硬拼成一个大矩阵，而是保留每个面自己的二维索引，同时把六个面放进一个统一类型里。

3.1 立方球网格为什么拆成六个面

立方球（cubed sphere）把球面投影到包裹它的立方体上，每个立方体面变成一个水平块。全球模式用这种网格后，北极和南极不再是容易畸变的“单点”，每条纬线也都不再是等经纬度意义上的直线。六个面之间的接缝处，相邻点的坐标和距离由网格文件单独给出，普通数组无法表达这种拓扑关系。gcmfaces 的解决办法是让每个变量都带着六个面块，配合网格文件里的距离和面积信息操作。

3.2 先看一个 gcmfaces 对象的 size 与 fields

在没读文档前，最安全的做法是先用size观察输出：

gcmfaces_global; global mygrid; disp(size(mygrid)); class(mygrid) fieldnames(mygrid)

如果你这时候看到的是类似[48 48 1 6]的尺寸，说明该发行版把六个面作为最后一维；如果看到[288 192]，说明它已经被某种方式展平。两种表示都见过，写法完全取决于版本。因此不要背“某对象一定是几维”，而是每次拿到数据后先打印size和class。gcmfaces 类内部通常有一个faces字段，存六块面数据，字段名可能随版本变化，直接fieldnames(obj)查更可靠。

3.3 用 help 读取当前版本的真结构

我常在项目里写一段启动时就执行的帮助命令：

help gcmfaces help gcmfaces/mean open gcmfaces.m

help gcmfaces显示类说明；help gcmfaces/mean显示重载方法说明；open gcmfaces.m打开类定义，让你确认它的保存属性。参数含义：gcmfaces/mean这样的“类/方法”语法在 Matlab 和 Octave 中都能定位方法，比搜索文档更快。加上了open这一步，你还能看到该版本里faces到底是cell还是普通数组，这决定后面访问单面时用obj.faces{i}还是obj.faces(:,:,i)。

3.4 纯数组模拟：六个面在内存里的排布

为了不绕进版本差异，可以用普通 cell 模拟六个面的存储逻辑：

nFaces = 6; nx = 48; ny = 48; facecell = cell(nFaces, 1); for f = 1:nFaces facecell{f} = rand(ny, nx); % 每个面一块独立数组 end

这段代码展示的是 gcmfaces 最核心的存储思想：每个面用独立的(ny,nx)矩阵，六个面之间不共享索引。实际工具箱会把这种 cell 封装成类，并在矩阵运算时判断哪些索引恰好落在面接缝上。接缝处的点是两个面共用的，直接循环求平均时会重复计数，gcmfaces 则会在生成对象时记录接缝掩膜，运算后统一做一次“合并相邻面”的处理。下面的表格列出了最常用的重载操作：

代码与表格配合起来，是初学者最容易上手的路径：先把六个面当成若干独立数组，理解cell层面的循环，再去看工具箱的重载方法。

4. 用 gcmfaces 做全局平均、插值和分区统计的参数细节

gcmfaces 被引入项目，多数是为了解决“全局平均”和“换网格”这两类高频操作。这个库的价值在于能把六块面当成一个整体来算，但前提是参数含义正确。下面几个常用场景来自我实际处理 MITgcm 输出时的经验。

4.1 mean 与 sum 的维度语义，以及逐面退化的坑

直接调用mean(gcmfacesObj)，很多版本返回的不是一个标量，而是六个面上各自的结果。这是最容易踩的坑。先执行：

help gcmfaces/mean mean(gcmfacesObj)

help会告诉你当前版本是否支持'all'或维度参数。如果你的版本返回 6 个值，就需要显式传维度。垂直维和时间维一般不是面维度，调用mean(obj, 3)时工具箱会把三维数组的前两维视为水平面，第三维作为汇总维；但不同发行版对维度的约定并不完全一致，建议先在只有两个时间层的测试文件上跑一次，确认结果维度后再批量处理。

4.2 做全球面积加权平均

立方球网格各面网格面积不均，简单平均没有意义。常见做法是构造面积权重对象：

global mygrid % 假设 mygrid 属于当前版本的网格结构 % 若字段叫 rA 直接用；若叫 areaCell，则替换 area = mygrid.rA; globalMean = sum(data .* area) / sum(area);

参数说明：rA在 MITgcm 中代表网格单元面积，量纲可能是平方米；data为 gcmfaces 对象时，.*会触发逐面元素相乘，sum汇总全部六个面。结果是否接近物理平均，可以用全球海表温度的一个已知参考值粗验：若偏差远大于 0.5℃，检查面序和掩膜是否对齐。另一种更稳妥的做法是先用fieldnames(mygrid)确认面积字段名，再参与计算，避免硬编码字段名造成的项目间迁移问题。

4.3 插值到规则经纬网格的常用接口

要从立方球面转换到 1°×1° 规则网格，常见调用是：

lon = 0.5:1:359.5; lat = -89.5:1:89.5; [LON, LAT] = meshgrid(lon, lat); out = interp(gcmfacesObj, LON, LAT, 'method', 'linear');

各参数含义：LON、LAT为插值目标坐标矩阵；method选linear或nearest，气候场上linear平滑但会抹掉锋面，nearest适合对比离散点位。如果interp在当前版本中未定义，则执行which interp，看是否被其他工具箱遮蔽。常见的遮蔽来源是映射工具箱或用户自定义的interp.m，用which -all interp列出全部同名函数，再决定是否把 gcmfaces 的目录前置到路径开头。

4.4 用掩膜对象做分区统计

分区均值直接利用重载的关系运算：

landMask = data > -1e10; % 根据数据有效范围生成陆海掩膜 basinMean = sum(data .* oceanMask) / sum(oceanMask);

这里oceanMask是通过逻辑比较生成的 gcmfaces 对象，data .* oceanMask会把陆地点清零，再求和得到海洋区域总和。需要注意：不要把掩膜对象和普通double数组混用，混合运算可能先把对象转成 cell，从而报出“未定义 operator .*”的错误。遇到这种报错，先检查class(mask)，再把双方都改成 gcmfaces 对象。

5. gcmfaces 可视化与 NetCDF 输出：把六个面拼成能放进论文的图

处理海洋模式数据，最终要落到图和文件。可视化这一步也是检查数据是否读对的快捷方式。

5.1 先画单个面，确认方向与投影

用imagesc看单面：

for f = 1:6 subplot(2, 3, f); imagesc(data.faces{f}); % 字段名以 fieldnames 为准 axis xy; colorbar; end

关键点是axis xy。很多模式面数据 Y 轴向上增长，不加这一句，图会上下颠倒，让你误判海流向。先逐个面看，确认陆海边界与模式文档描述一致，再做拼接图。

5.2 用拼接方式快速观察全球场

六面单独画已经能定位大部分问题，但论文需要全球拼接图。常见做法是把六个面按立方球展开顺序拼成三行两列：

globalRows = [1 2 3; 4 5 6]; canvas = zeros(size(data.faces{1}, 1) * 3, size(data.faces{1}, 2) * 2); for f = 1:6 [r, c] = ind2sub([2 3], f); r0 = (r - 1) * ny + 1; c0 = (c - 1) * nx + 1; canvas(r0:r0+ny-1, c0:c0+nx-1) = data.faces{f}; end imagesc(canvas);

这是一段示意代码，真正的拼接顺序要按网格文件里的face索引来排，不同模式对外观顺序不同。若六个面的拼接图和文献的立方球展示图不一致，不要怀疑图，先看网格编号。另一个经常出错的地方是imagesc默认会把 NaN 显示成最小值，导致接缝处出现黑色线条；可以先用set(gca,'Color',[0.8 0.8 0.8])把底色设为灰，再接缝会显示成中性色，不容易误导颜色映射。

5.3 NetCDF 输出前的维度整理与字段约定

给下游气象/海洋团队提交 NetCDF 前，gcmfaces 对象本身不能直接写入，一般先插值到规则经纬网格。用 Matlab 内置函数可以快速落盘：

nccreate('out.nc', 'temp', 'Dimensions', {'lon', 360, 'lat', 180}); ncwrite('out.nc', 'temp', temp2D);

nccreate第一参数文件名，第二参数变量名，Dimensions定义两个坐标轴的长度；ncwrite写数据时顺序要与nccreate一致。常见约定如下表：

Octave 下更直接的做法是ncwrite前用pkg load netcdf，并检查netcdf.open是否正常。转换后的规则网格文件能被 NCO、CDO 正常读取，gcmfaces 就不再参与后续工作。写多时间层时，把time放在最后一位可以避免与nccreate的写入顺序冲突。

6. 把 gcmfaces 放进批量作业与无界面环境：三个值得养成的习惯

到这里，安装、数据组织、常用运算和输出都已经能覆盖日常需求。最后说几个在真实项目里最值得养成的习惯，尤其是当你要跑 20 年逐日输出的时候。

第一个习惯是不要save整个 gcmfaces 对象。六个面的大头数据用save确实能存，但类对象的序列化格式与版本强相关，同一份.mat在 Matlab 和 Octave 之间来回读，经常会遇到“类定义不一致”。我一般只保存最原始的二进制文件或来源数组，另写一个重建脚本，每次运行时重新生成 gcmfaces 对象。

第二个习惯是在无界面环境里用非交互模式跑批处理。Matlab R2019a 以后可以用：

matlab -batch "run('process_daily.m')"

Octave 则用：

octave --no-gui --eval "run('process_daily.m')"

参数说明：-batch运行完会自动退出，不会弹桌面；--no-gui屏蔽 Octave 的图形界面。两者的共同问题是启动都要读一遍startup.m，因此路径配置必须在启动时就固定好。

第三个习惯是保留一个“对照诊断量”。我第一次把 gcmfaces 接入项目时，先算了一个已经知道答案的场——全球平均海表温度。模式官方文档给出的参考值是 18.4℃ 左右，自己的脚本算出的偏差要小于 0.01℃。这个检查不是可选的，它能同时暴露面序错误、面积权重错误和掩膜不匹配三类问题。极区附近偶尔出现NaN是正常的，但如果接缝处出现规则条纹，就要回到第 3 章的size检查，看六个面是否被当成一个整体在运算。

最后一个小技巧：把“读数据 → 重建网格 → 输出规则网格”的步骤封装成一个函数，保证每次调用都从同一个入口进去，这样模式输出版本更新时，只需要改这一个函数。gcmfaces 不是一个让你全文照抄的对象，它值得你为它写一层属于自己的薄封装。

本文还有配套的精品资源，点击获取