简介：潮汐模型驱动程序（TMD）3.0版配套MATLAB代码包，面向从事海洋潮汐模拟、水文预报及地球物理研究的高校师生与科研人员，用于解决区域潮汐场建模、潮位预测与模型驱动计算等实际问题。压缩包共197个文件，约93.55MB，以113个png图示、38个m脚本、17个md说明文档和13个html页面为主，另含nc数据文件、gif动画及xml、bib等配置与参考文献文件，覆盖模型文档、格式说明与示例数据。已有345人学习下载。资源内含TMD模型文件格式说明、预测文档与EOT教程，配合MATLAB脚本可直接复现潮汐驱动流程，帮助读者理解模型输入输出结构、掌握潮汐预测实现方法，并借助图示与动画直观核对结果，适合作为潮汐建模入门与二次开发的参考素材。
1. 潮汐模型驱动程序 TMD 3.0：一个被低估的潮汐调和分析工具
如果你手头有一年以上的水位观测序列，想分离出八大分潮的振幅和相位，又不想从零推导最小二乘的观测方程，那 TMD（Tidal Model Driver）大概率是你绕不开的一个工具。它最早由俄勒冈州立大学相关团队维护，核心思路是把全球或区域潮汐模型的调和常数预先存成二进制格式，再通过一套驱动程序按经纬度和时间插值出任意时刻的潮位。3.0 版把接口整理得更规整，附带 MATLAB 调用代码，意味着你不需要啃 Fortran 源码，只要准备好模型文件和站位坐标，就能在 MATLAB 里直接算出潮汐预报。
它适合谁？做近海工程水位校核的、跑海洋模式需要开边界潮汐驱动的、做验潮站数据质量控制的，以及写论文需要批量算潮汐调和常数的。不适合谁？只想画一张潮汐曲线图、对精度没要求的人，用现成的在线预报就够了。这篇笔记按「TMD 3.0 到底怎么跑起来 → 模型数据怎么配 → 调和分析怎么落地 → 哪里容易翻车」的顺序讲，代码能直接抄，参数会逐个说明。
2. 把 TMD 3.0 在 MATLAB 里跑通：从解压到第一条潮位曲线
2.1 先搞清楚 TMD 3.0 的目录结构和依赖关系
拿到压缩包后别急着双击运行，先看清楚里面有什么。TMD 3.0 的典型结构分三块：一是核心驱动代码，通常是 Fortran 或编译好的可执行文件，负责读取模型数据并插值；二是模型数据文件，扩展名常见为.nc或自定义二进制格式，存放各分潮的振幅和相位网格；三是 MATLAB 封装层，提供tmd_tide_pred、tmd_constit这类函数，把底层调用包装成 MATLAB 能识别的接口。
这里有个关键认知：TMD 本身不是潮汐模型，它只是「驱动程序」。模型数据才是真正的潮汐信息载体。所以你会看到两种文件要分开管理——驱动代码放一个目录，模型数据放另一个目录，MATLAB 路径只添加驱动代码所在目录，模型数据路径通过参数传入。混在一起放不是不行，但换模型时容易乱。
依赖方面，MATLAB 版本不需要太新，R2016b 之后的版本基本都能跑，因为用到的就是基础的矩阵运算和文件读写。如果你用的是 Windows 且核心是编译好的 exe，注意系统架构要匹配，32 位 MATLAB 调不动 64 位 exe，这是血泪经验。Linux 下则通常需要自己编译 Fortran 源码，gfortran是常见选择。
2.2 最小可运行示例：给定经纬度和时间算潮位
下面这段代码是我一般会先跑的「冒烟测试」，目的是确认驱动能读到模型、能返回数值、数值量级合理。假设模型数据文件叫TPXO9_atlas_v5.nc，站位设在东海某点。
% 添加 TMD 驱动代码路径（按实际解压位置修改） addpath(genpath('D:\TMD3.0\driver')); % 模型数据文件路径 model_file = 'D:\TMD3.0\data\TPXO9_atlas_v5.nc'; % 站位经纬度（东经、北纬为正） lon = 122.5; lat = 30.0; % 时间范围：2024-01-01 00:00 到 2024-01-02 00:00，每小时一个点 t_start = datenum(2024,1,1,0,0,0); t_end = datenum(2024,1,2,0,0,0); t_vec = t_start:1/24:t_end; % 调用潮汐预报函数 [h, constituents] = tmd_tide_pred(model_file, t_vec, lat, lon, 'z'); % 输出前几个值看看量级 disp(h(1:5)); disp(constituents);逻辑说明：tmd_tide_pred是封装好的预报函数，第一个参数是模型文件，第二个是 MATLAB datenum 格式的时间向量，第三、四个是纬度和经度，最后一个'z'表示要水位（有些实现用'h'表示高度，具体看函数签名）。返回值h是潮位序列，单位通常是米；constituents返回参与计算的分潮列表。
参数说明：时间向量必须是 datenum 格式，不是 datetime，这是最容易翻车的地方。经纬度顺序是「先纬后经」，和很多地理函数的习惯相反，写反了会插值到错误网格甚至报错。如果你发现返回值全是 NaN，先检查经纬度是否落在模型覆盖范围内——全球模型一般没问题，区域模型就要看边界了。
2.3 模型数据文件的加载与站位插值逻辑
TMD 的插值不是简单的最近邻，它会在模型网格上做双线性插值，再对每个分潮的振幅和相位分别处理。振幅直接插值，相位则要小心跨 0/360 度的跳变，好的实现会做相位解缠。你不需要自己写这部分，但要知道它存在，因为当站位靠近网格边界或水深变化剧烈处，插值误差会明显增大。
模型文件加载是一次性的，第一次调用会读进内存，后续调用同一模型会快很多。如果你要批量算几百个站位，建议把模型读一次然后循环站位，而不是每个站位都重新加载。我见过有人写循环时把模型文件路径放在循环内反复读，结果跑一晚上没跑完，这就是没注意缓存机制。
另外，模型文件里的分潮数量决定了你能算哪些分潮。常见全球模型包含 M2、S2、N2、K2、K1、O1、P1、Q1 八个主要分潮，有些还包含 M4、MS4 等浅水分潮。如果你的研究需要更多分潮，得换模型或自己补充。MATLAB 封装层一般会返回实际参与计算的分潮名，拿到后先扫一眼，确认没有缺你需要的那个。
3. 潮汐调和分析的落地：从水位序列反算分潮参数
3.1 调和分析的基本方程与 TMD 的角色区分
调和分析的本质是：已知水位序列，反推各分潮的振幅和相位。方程是线性的，每个时刻的水位等于各分潮贡献之和加上残差。TMD 在这里扮演两个角色——一是提供预报值用于对比，二是它的模型数据本身就可以作为先验，帮你判断反算结果是否合理。
常见做法是：用 TMD 预报一段时间的潮位，再把这个预报值当作「真值」，去检验你自己的调和分析代码是否正确。如果你从观测数据反算出的 M2 振幅和 TMD 预报的 M2 振幅差了一个量级，那要么观测数据有问题，要么你的分析代码有 bug。这个交叉验证思路非常实用，比单纯看残差平方和直观得多。
调和分析的观测方程可以写成矩阵形式：A * x = b，其中A是设计矩阵，每列对应一个分潮的 cos 和 sin 项，x是待求系数，b是观测水位。最小二乘解就是x = (A'A)^(-1) A'b。MATLAB 里用反斜杠运算符就能算，不需要自己写求逆。
3.2 用 MATLAB 实现最小二乘调和分析并对比 TMD 结果
下面这段代码演示从合成水位序列反算分潮参数，再和 TMD 预报对比。合成数据用 TMD 生成，所以理论上反算结果应该非常接近。
% 生成合成水位序列（用 TMD 预报作为“观测”） t_vec = datenum(2024,1,1,0,0,0):1/24:datenum(2024,2,1,0,0,0); [h_obs, ~] = tmd_tide_pred(model_file, t_vec, lat, lon, 'z'); % 定义分潮频率（周期，单位小时） constits = {'M2','S2','N2','K2','K1','O1','P1','Q1'}; periods = [12.4206, 12.0000, 12.6583, 11.9672, ... 23.9345, 25.8193, 24.0659, 26.8684]; % 构建设计矩阵 N = length(t_vec); A = ones(N, 1); % 常数项（平均水位） for k = 1:length(constits) omega = 2*pi / periods(k); % 角频率，单位 1/小时 t_hours = (t_vec - t_vec(1)) * 24; % 转成小时 A = [A, cos(omega*t_hours(:)), sin(omega*t_hours(:))]; end % 最小二乘求解 coef = A \ h_obs(:); % 提取振幅和相位 amp = zeros(length(constits),1); pha = zeros(length(constits),1); for k = 1:length(constits) c = coef(1 + 2*k - 1); s = coef(1 + 2*k); amp(k) = sqrt(c^2 + s^2); pha(k) = atan2(s, c) * 180/pi; % 转成度 end % 打印结果 for k = 1:length(constits) fprintf('%s: 振幅=%.4f m, 相位=%.2f 度\n', constits{k}, amp(k), pha(k)); end逻辑说明：设计矩阵第一列是常数项，吸收平均水位；之后每个分潮占两列，分别是 cos 和 sin。时间转成「从起始时刻算起的小时数」是为了避免大数导致的数值精度问题。A \ h_obs直接给出最小二乘解，MATLAB 内部用 QR 分解，比显式求逆稳定。
参数说明：分潮周期必须用正确值，M2 是 12.4206 小时，不是 12 小时整，差这一点会导致相位随时间漂移。相位用atan2而不是atan，因为要区分象限。如果你发现某个分潮振幅异常大，先检查该分潮周期是否写错，或者观测序列长度是否足够——太短的序列无法分离相近频率的分潮，比如 K1 和 P1 周期接近，至少需要一个月以上数据才能分开。
3.3 结果验证：残差分析和与 TMD 预报的交叉对比
算完调和常数后，别急着出图。先做两件事：一是看残差序列是否还有明显周期信号，二是把反算的调和常数代回预报方程，和原始序列对比。如果残差里还能看到明显的半日周期，说明你漏了某个分潮，或者分潮周期设错了。
残差分析用 MATLAB 几行就能做：
h_fit = A * coef; % 拟合值 resid = h_obs(:) - h_fit; % 残差 fprintf('残差标准差: %.4f m\n', std(resid)); fprintf('残差最大值: %.4f m\n', max(abs(resid))); % 画图检查 figure; subplot(2,1,1); plot(t_vec, h_obs, 'b', t_vec, h_fit, 'r--'); legend('原始','拟合'); datetick('x','mm-dd'); subplot(2,1,2); plot(t_vec, resid); datetick('x','mm-dd'); title('残差序列');如果残差标准差在厘米量级，说明拟合不错；如果到了分米量级，要么数据有噪声，要么分潮没选全。交叉对比 TMD 预报值时，注意 TMD 用的是模型网格插值，而你反算用的是观测序列，两者本身就有差异，不要期望完全一致。合理的差异范围取决于站位处模型精度，一般 M2 振幅差在几厘米内算正常。
4. 避坑与排查：TMD 3.0 跑不起来时先看这几条
4.1 现象：调用函数报「未定义函数或变量」
原因：MATLAB 路径没添加对，或者添加了但函数名拼写和实际文件不一致。TMD 3.0 的 MATLAB 封装层函数名在不同版本里可能有细微差别，比如tmd_tide_pred和tide_pred都出现过。
解决：先用which tmd_tide_pred确认 MATLAB 能不能找到这个函数。如果返回空，检查addpath的路径是否指向了包含.m文件的目录，而不是上一级目录。用genpath递归添加子目录更保险。如果函数名确实不对，去解压目录里搜.m文件，看实际名字是什么。
4.2 现象：返回的潮位全是 NaN 或零
原因：经纬度超出模型覆盖范围，或者经纬度顺序写反了。区域模型尤其容易出这个问题，全球模型一般不会。
解决：先确认模型覆盖范围，通常模型文件里会有说明或属性。把经纬度限制在覆盖范围内再试。另外检查调用时的参数顺序，TMD 的 MATLAB 封装一般是「先纬后经」，和interp2的习惯相反。如果拿不准，用两个明显不同的经纬度各跑一次，看返回值有没有变化，没变化就是没插值到正确位置。
4.3 现象：时间向量用 datetime 格式导致报错
原因：TMD 的 MATLAB 封装层大多按 datenum 设计，datetime 是后来引入的格式，不兼容。
解决：用datenum转换。如果原始数据是 datetime，用datenum(dt)转一下。反过来，如果函数返回的时间是 datenum，你想用 datetime 画图，用datetime(datenum_val, 'ConvertFrom', 'datenum')。这个转换在批量处理时容易忘，建议在脚本开头统一转好。
4.4 现象：调和分析结果里某个分潮振幅特别大或特别小
原因：分潮周期设错，或者观测序列太短导致分潮之间无法分离。K1 和 P1、S2 和 K2 这两对是重灾区。
解决：核对分潮周期表，M2=12.4206、S2=12.0000、K1=23.9345、O1=25.8193，这些值不能凭记忆写。序列长度至少覆盖 30 天，最好 60 天以上。如果数据就是短，考虑固定某些分潮的振幅相位关系（比如用模型值约束），而不是全部自由反算。
4.5 现象：Windows 下调用 exe 核心时闪退或无输出
原因：系统架构不匹配，或者缺少运行时库。32 位 MATLAB 调 64 位 exe 会直接失败，反之亦然。
解决：确认 MATLAB 和 exe 的架构一致，用computer('arch')看 MATLAB 架构。如果 exe 是编译好的，检查是否需要安装 Fortran 运行时库。实在搞不定就用纯 MATLAB 实现，TMD 的插值逻辑不复杂，自己写一个双线性插值加调和预报也就几十行，可控性还更高。
5. 进阶技巧：用 TMD 做多站位批量预报与边界条件生成
当你需要给海洋模式提供开边界潮汐驱动时，单个站位预报就不够了，得沿边界生成一系列站位的潮位时间序列。这时候效率是关键。我一般会先把模型读一次，然后循环所有边界点，每个点只做插值和调和求和，不再重复读文件。
具体做法是：把边界点经纬度存成两列矩阵，时间向量只建一次，然后循环调用预报函数。如果 TMD 的封装层不支持批量传入多个站位，就自己包一层循环。注意每次调用返回的潮位是列向量，循环结束后拼成矩阵，每列一个站位，存成 NetCDF 或文本文件供模式读取。
另一个技巧是分潮筛选。有些应用只关心 M2 和 S2，那就在预报时只算这两个分潮，速度会快很多。TMD 的封装层一般允许传入分潮列表参数，查一下函数签名，把不需要的分潮去掉。我见过有人为了算 M2 把八个分潮全算一遍，白白多花七倍时间。
验证批量结果时，抽几个点单独跑一次，对比批量输出是否一致。不一致通常是循环里变量没清干净，或者经纬度索引写错了。这个检查花不了几分钟，但能避免把错误边界条件喂给模式后跑几天才发现问题。
最后说个习惯：每次换模型数据文件，先跑一遍单点冒烟测试，确认量级和分潮列表正常，再上批量。模型文件版本不同，分潮命名和网格分辨率都可能变，直接套用旧脚本容易翻车。希望帮到你。
本文还有配套的精品资源，点击获取
