企业官网建设流程全解析

手里只有一景 MODIS L1B 数据，任务是要交一张 1 km 分辨率的海表温度分布图，装的软件是 ENVI，没有专门的 SST 产品处理链——这大概是做海洋、气象、环境、生态方向的人最常撞上的一种处境。直接去下现成的 SST 产品当然省事，但产品的空间分辨率、时间窗、掩膜策略、算法版本全都由别人定，你要做锋面识别、近岸变化分析或者跟实测数据对齐，往往就被这些"既定条件"卡住。所以真正需要自己动手的时候，绕不开的一件事就是：在 ENVI 里用劈窗算法，把 MODIS 两个热红外通道的亮温差，一步步换算成海表温度。

这篇东西讲的就是这条链路怎么走通：劈窗算法到底在补偿什么物理量，为什么它比单通道更省事，MODIS 数据里哪几个 SDS 才是真正要动的，DN 怎么定标成辐亮度、怎么反演成亮温，系数怎么挑、符号怎么防错，GLT 和蝴蝶结怎么处理，云和陆地掩膜怎么做，最后怎么验证自己算出来的数不是一堆好看的假数据。适合已经会用 ENVI 基本操作、能看懂 HDF 结构，但还没完整跑过一遍热红外反演流程的人；也适合跑过一遍但是结果总觉得不对、想找找问题出在哪的人。文中给出的系数、偏移量都是示例量级，正式使用必须按自己的区域重新标定，这一点我会反复强调。

1. 先想清楚劈窗算法在补什么，再动手点鼠标

很多人上手第一件事是打开 ENVI 找工具，其实真正决定结果对不对的，是脑子里的物理图像清不清楚。劈窗算法不是"两个波段一减就出温度"的魔法，它背后有一条相当朴素的辐射传输逻辑，理解了这条逻辑，后面所有参数选择、异常排查都有依据。

1.1 11 µm 与 12 µm 的亮温差里装着水汽信息

海面向上发射的热辐射，在穿过大气时会被水汽连续吸收削弱，11 µm 附近的吸收弱一些，12 µm 附近强一些。两个通道看到的是同一片海面，但"看得深浅"不同：11 µm 通道透射率更高，看到的海面温度更接近真实值；12 µm 通道被削弱得更多，亮温偏低。两者的差值随着大气水汽含量、观测天顶角、温度廓线的变化而变化，于是这个差值本身就成了一个"大气校正项"的代理变量。

把通道 i 的辐射传输方程写开，形式是这样的：

L_i = ε_i · B_i(Ts) · τ_i + L_i↑ + (1 - ε_i) · τ_i · L_i↓

其中 ε 是海面发射率（11 µm 约 0.99，12 µm 约 0.988），τ 是整层透射率，L↑ 是大气上行辐射，L↓ 是下行辐射被海面反射回来的部分。对 Planck 函数在 300 K 附近做线性化，忽略发射率差异带来的二阶项，两个通道相减之后，大气上行辐射和大部分廓线项会互相抵消，剩下的主导项正比于 (Ts - Ta)·(1-τ31)/τ31 - (Ts - Ta)·(1-τ32)/τ32，也就是正比于水汽造成的差异削弱。这就是为什么最终的公式可以简化成一条直线：

SST = T31 + A · (T31 - T32) + B

A 的量级通常在 1.5 到 3.0 之间，B 在 -2 K 到 +2 K 之间。系数看着简单，但它其实是对整层大气状态的一种"平均化封装"，所以换区域、换季节、换水汽条件，系数就得跟着变。这一点是后面所有坑的源头。

顺带说一个直觉判断：开阔洋面上 (T31 - T32) 一般是正值，干空气时 0.5 K 左右，湿空气时能到 2 到 3 K。如果你算出来的差值大面积是负的，先别怀疑公式，八成是云、陆地或者数据索引错了。

1.2 三条技术路线的取舍，为什么业务化偏爱劈窗

做海温反演可选的路线不止一条，把三条主流方案的代价摆在一起看，选择就很清楚了。

| 路线 | 需要的输入 | 计算代价 | 主要短板 | | 单通道法 | 精确的水汽廓线、气溶胶、发射率 | 中等 | 廓线误差直接进结果，1 g/cm² 的水汽偏差可能带来 1 K 以上误差 | | 劈窗法 | 两个热红外通道亮温、视角 | 很低，可逐像元 | 系数区域适用性、薄卷云敏感 | | 辐射传输正演查表 | 廓线库、逐像元匹配 | 高，构建 LUT 成本大 | 廓线时空代表性不足时反而更差 |

劈窗法的妙处在于它把"绝对定标"的问题转化成了"相对差异"的问题：你不需要知道大气里到底有多少水汽，只要两个通道对水汽的响应差异是稳定的，这个差异就能被系数吸收掉。对于一天几十景、要做时间序列的业务场景，这个代价优势是压倒性的。

ENVI 这边的情况也值得说清楚。ENVI 里针对可见光短波红外的辐射传输校正工具（大气校正模块那一套）覆盖范围到 2.5 µm 左右，热红外不在它的处理区间里；针对高光谱热红外的温度发射率分离工具要求输入发射率光谱和大气廓线，逐景建模成本高。所以对多光谱的 MODIS 来说，自己用 Band Math 搭一条劈窗链路，反而是最直接、最可控的做法。ENVI 的价值主要体现在几何定位工具、可视化试算、掩膜叠加和 IDL 批处理这几块，算法本身还是得自己写。

1.3 这套流程能用到什么程度，什么情况下别硬用

劈窗算法在开阔水体上表现最好，尤其是水温 0 到 32 °C、天顶角小于 45°、大气比较稳定的情形。它不适用的场景同样要心里有数：近岸浑浊水体和河口区域，悬浮物和溶解有机物会改变 11 µm 与 12 µm 的发射率比值，劈窗的核心假设被破坏；薄卷云边缘，云检测漏掉的像元会被系统性拉低 2 K 以上；沙尘、火山气溶胶、强烈逆温条件下，两个通道的差异不再单调对应水汽；太阳耀斑区在白天会引入镜面反射的短波泄漏。

还有一个经常被忽略的点：卫星反演的是"表皮温度"（skin temperature），浮标测的是水下几十厘米到几米的"水体温度"（bulk temperature），两者在夜间平静海面相差约 -0.1 到 -0.3 K，白天有日照时差异更大且方向不定。验证的时候不扣这个偏差，你会误以为算法本身有系统误差。

2. MODIS L1B 里真正需要动的只有几个 SDS

一景 MODIS 一级产品打开后，SDS 列表长得让人眼晕，但这条链路要用的其实就三类：发射波段的辐亮度、地理定位与观测几何、云掩膜。把这三个文件的分工理清楚，后面的操作会顺很多。

2.1 MOD02、MOD03、MOD35 各自负责什么

发射波段数据文件提供的是辐亮度量化值（DN），里面 EV_1KM_Emissive 这一个 SDS 就包含了 1 km 分辨率的全部热红外通道，共 16 个波段，从第 20 波段到第 36 波段，中间跳过了第 26 波段。这个"跳过"是后面索引错位的根源，务必记住。

几何定位文件提供逐像元的经纬度、传感器天顶角和方位角、太阳天顶角，1 km 通道的经纬度数组和发射波段是同一套网格，可以直接对齐。要注意的是传感器天顶角在文件里是 int16 存储、缩放因子 0.01，直接读原始值会得到 3000 这种数，必须先乘 0.01 才是度数。

云掩膜产品提供逐像元的云检测结果和置信度，空间分辨率是 1 km（聚合后的产品），需要重采样或者按最近邻对齐到发射波段网格上。

三个文件的时间戳必须严格对应同一景（同一颗星、同一个 5 分钟段），拿 Terra 的定位文件去配 Aqua 的辐亮度，出来的图会整片错位，而且是那种"看着像那么回事、细看全是歪的"的错。

2.2 定标系数必须从 SDS 属性里读，别手抄

辐亮度定标的缩放因子和偏移量是写在 SDS 属性里的数组，逐波段一个值，而且不同数据版本、不同处理批次之间有细微差别。手抄别人的数值，短期看不出问题，做长时间序列时会引入难以追溯的系统偏差。

用 IDL 把属性直接读出来是最稳的做法：

pro dump_emissive_attrs, fname compile_opt idl2 sd_id = hdf_sd_start(fname, /read) idx = hdf_sd_nametoindex(sd_id, 'EV_1KM_Emissive') sds_id = hdf_sd_select(sd_id, idx) hdf_sd_getinfo, sds_id, name=nm, dims=dm print, nm, dm hdf_sd_attrinfo, sds_id, 'radiance_scales', data=scales hdf_sd_attrinfo, sds_id, 'radiance_offsets', data=offsets hdf_sd_attrinfo, sds_id, 'band_names', data=bnames print, 'band_names: ', bnames print, 'scales 31/32: ', scales[10:11] print, 'offsets 31/32: ', offsets[10:11] hdf_sd_endaccess, sds_id hdf_sd_end, sd_id end

打印出来的 band_names 长这样："20,21,22,23,24,25,27,28,29,30,31,32,33,34,35,36"。数一下位置：第 1 个是 20，第 11 个是 31，第 12 个是 32。所以在 IDL 的 0 基索引里，第 31 波段是索引 10，第 32 波段是索引 11——不是 30 和 31。这个错位我见过太多次，症状是算出来的海温整体偏 3 到 8 K，而且空间形态还挺像样，特别容易蒙混过关。

偏移量的数值往往很大（示例中第 31 波段约为 1577.34），这是正常的，它对应的是空间视场的 DN 水平，不是错误。缩放因子的量级在 1e-5 到 1e-4 之间，具体以属性为准。

2.3 把两个通道单独提出来

在 ENVI 里把 16 个波段拆开有两条路。图形界面走 Spectral Subset 或者 Layer Stacking，选第 11 和第 12 个波段；脚本化的话，Band Math 支持对多波段文件按索引取值，写b1[10]和b1[11]就能直接引用，省掉一次中间文件落盘。批处理场景下后者更省事，也更少出"文件名对应错了"这类低级事故。

提出来之后建议立刻做一次目视检查：两个通道的影像纹理应该高度相似，海岸线、岛屿轮廓、云的分布完全一致，唯一的差别是灰度略有不同。如果两个波段的纹理对不上，说明取错索引或者行列方向反了，这时候停下来改，比后面查半天错要划算得多。

3. 从 DN 到亮温：定标和普朗克反演里的细节账

这一段是纯计算，但恰恰是最容易出隐性错误的地方。定标公式写错一个括号、波长取错一位小数，结果可能只差零点几 K，肉眼完全看不出来，但跟实测一比就原形毕露。

3.1 辐射定标的公式与 Band Math 写法

MODIS 发射波段的定标是线性的：

L = (DN - radiance_offset) × radiance_scale

单位是 W/(m²·sr·µm)。用示例数值写进 ENVI 的 Band Math，第 31 波段就是这样一行：

(b1 - 1577.34) * 7.304e-05

第 32 波段把偏移和缩放换成对应的值即可。这里有三个必须注意的点。

第一，输出类型一定要显式指定为 32 位浮点。ENVI 的 Band Math 默认可能沿用输入的数据类型，输入是 uint16 的话，小数部分直接被截掉，算出来的辐亮度全是整数，后面的亮温会整体偏 1 K 以上，而且偏得很均匀，非常难发现。

第二，填充值要先处理。发射波段的无效像元用 65535 之类的值填充，直接参与减法会得到一个巨大的负数，再进对数函数就报错或者出 NaN。稳妥的写法是先算掩膜，或者在表达式里加判断，让无效像元输出成 0 或者 NaN。

第三，单位别混。有些资料给的定标公式里辐亮度单位是 W/(m²·sr·cm⁻¹)，跟 µm 差一个波数换算，混用会差出好几个数量级，这种错误一般一眼能看出来，反而是上面两种"差一点点"的错误更危险。

3.2 亮温反演的表达式与有效波长取值

辐亮度转亮温用的是 Planck 函数反解：

T = c2 / (λ · ln(c1 / (λ^5 · L) + 1))

其中 c1 = 1.191042e8 W·µm⁴/(m²·sr)，c2 = 14387.69 µm·K，λ 取通道的有效波长，单位 µm。第 31 波段取 11.03，第 32 波段取 12.02。把 λ^5 预先算成常数写进表达式，ENVI 的运算会快一些：

14387.69 / (11.03 * alog(1.191042e8 / (b1 * 1.6327e5) + 1))

这里 1.6327e5 就是 11.03 的五次方（精确值约 163265）。第 32 波段的表达式把 11.03 换成 12.02、1.6327e5 换成 2.5091e5 即可。

波长取值这件事值得多说一句。严格来说应该用通道响应的积分有效波长，而不是中心波长，两者差别在 0.01 µm 量级，对亮温的影响大约 0.05 到 0.1 K。对于劈窗算法 0.5 K 量级的目标精度，这个量级的误差可以接受，但如果你要做多传感器交叉比对，最好统一用同一套有效波长，别一个用中心波长一个用积分波长，凭空造出 0.1 K 的系统偏差。

3.3 抽几个像元手工核对，这一步别省

表达式跑完之后，随机挑几个像元，把 DN、辐亮度、亮温三个数手工算一遍，跟软件结果对一下。举个例子：某像元第 31 波段 DN 为 12900，偏移 1577.34，缩放 7.304e-05，则 L = (12900 - 1577.34) × 7.304e-05 ≈ 0.8267。代入反解式：

λ^5·L = 1.63265e5 × 0.8267 ≈ 1.3497e5 c1 / (λ^5·L) = 1.191042e8 / 1.3497e5 ≈ 882.5 ln(883.5) ≈ 6.784 T = 14387.69 / (11.03 × 6.784) ≈ 192.3

等一下，这个结果明显不对，192 K 不是海水。问题出在辐亮度太小了：0.83 W/(m²·sr·µm) 对应的亮温确实只有 190 K 上下，说明这个 DN 本身就偏低，要么是取到了冷云顶，要么是定标参数用错了波段。正常 300 K 海面的 11 µm 辐亮度在 9 W/(m²·sr·µm) 左右，对应 DN 大约是 9 / 7.304e-05 + 1577 ≈ 124800。如果你读到的 DN 只有一万多，先检查是不是把反射波段的缩放因子用到了发射波段上——这两套参数的来源 SDS 不同（一个是 EV_1KM_RefSB，一个是 EV_1KM_Emissive），复制粘贴的时候极容易串。

这个"手算发现量级不对"的过程，就是我做这套流程时固定保留的自检动作。它比看统计直方图有效得多，因为直方图只能告诉你分布形态，手工核对能直接抓出参数级别错误。

反演完之后看一眼亮温的数值范围：开阔洋面白天 290 到 305 K，夜间 288 到 302 K，极区冷季可以低到 271 K 附近（冰水混合）。如果整景亮温落在 250 K 以下，那不是算法问题，是定标或索引错了。同时检查 T31 - T32 的分布，正常水体上应该以正值为主，均值 1 K 上下。

4. 劈窗公式落地：系数选择比公式本身更考验判断

公式只有一行，真正花时间的是系数。系数选错的代价远大于任何一个操作细节的错误，因为它带来的是整片区域的系统性偏移，而且不会在图上显形。

4.1 线性劈窗与带初猜的 NLSST 差在哪

最基础的形式是纯线性：

SST = a + b·T31 + c·(T31 - T32)

它的隐含假设是"大气状态在水汽维度上可以用一个线性响应描述"，在中等水汽条件下够用。改进版本引入一个温度初猜（first guess）做交叉项：

SST = a + b·T31 + c·(T31 - T32)·Tguess + d·(secθ - 1)

交叉项的作用是让校正量随温度和水汽的耦合关系变化，在水汽跨度大的区域（比如从干冷的副热带高压边缘到湿热的赤道辐合带）比纯线性稳一些。代价是你必须先有一个初猜场——通常拿气候态海温或者再分析资料插值，这一步反而引入了新的不确定性。

我的建议是：做单景、特定季节的快速反演，先用线性式；做跨季节、跨纬度带的长时序，再上交叉项。别一上来就追求复杂模型，先把线性式的误差摸清楚，你会更清楚误差到底来自哪里。

4.2 系数选取的分区思路与一个示例表

系数的本质是对大气廓线的积分近似，所以合理的分区维度是水汽含量和温度层结，而这两个量大致可以用纬度带加季节来代理。下面这张表纯粹是量级示例，用来说明形式和各系数的相对大小，绝对不能直接拿去生产使用：

| 系数集 | 适用情形 | a (K) | b | c (K/K) | d (K) | | A | 低纬暖水、水汽充沛 | -0.8 | 1.002 | 2.6 | 1.2 | | B | 中纬夏季、水汽中等 | -1.2 | 1.004 | 2.1 | 0.9 | | C | 中纬冬季、水汽偏低 | -0.5 | 1.001 | 1.6 | 0.6 |

从这张表能看出两个规律：水汽越充沛，ΔT 的校正权重 c 越大，因为此时 12 µm 被削弱得更厉害，需要更强的补偿；观测天顶角越大，视角项系数 d 越重要，因为斜路径上的大气厚度按 secθ 增长。

这里有一个必须警惕的符号陷阱：不同文献对 ΔT 的定义方向不一致，有的写 (T31 - T32)，有的写 (T32 - T31)。如果你套用了一组按 (T32 - T31) 定义的系数，却按 (T31 - T32) 代入，校正项会整体反号，暖水区偏低、冷水区偏高，误差轻松上到 5 K 以上。拿到任何一组系数，第一件事是确认它的 ΔT 定义和单位（K 还是 °C）。

正式使用前的标定路径有三条：一是用本地探空廓线做辐射传输正演，拟合出适配本区域的系数；二是用匹配的浮标或走航实测数据做最小二乘回归；三是用官方发布的逐日或逐月系数网格，按像元插值。三条路各有代价，第二条最贴近"你这片海域"，但需要足够多的匹配样本；第三条最省事，但会引入外部产品的误差。

4.3 视角项与逐像元的完整表达式

视角项一定要逐像元加，不能用一个景的平均天顶角代替。一景 MODIS 的刈幅宽度超过 2000 km，边缘像元的天顶角能到 60° 以上，secθ - 1 达到 1.0，如果用星下点的 1.0 去乘，边缘区域的校正量直接少一半，误差能到 1 K 以上。

用 MOD03 的传感器天顶角（已乘 0.01 转成度）作为第三个输入波段，完整表达式是这样：

b1 + 2.1 * (b1 - b2) - 1.2 - 0.9 * (1.0 / cos(b3 * 0.0174533) - 1.0)

其中 b1 是亮温 31、b2 是亮温 32、b3 是天顶角度数。0.0174533 是度转弧度的系数，ENVI 的 cos 函数吃的是弧度，这个转换漏掉的话，30° 会被当成 30 弧度，余弦值乱跳，结果会出现莫名其妙的条带。

如果做的是大范围拼接，建议把视角项单独算成一个波段存下来，后面无论是换系数还是换模型，都不用重算。这种"中间量落盘"的习惯在批处理里特别值钱，能省下大量重复计算时间。

5. 几何定位：GLT 和蝴蝶结，绕不开的两件事

辐射量算对了，不代表图能用。MODIS 1 km 数据是原始扫描几何，没有经过投影，直接出图会看到明显的旋转和扭曲，跟矢量边界完全套不上。定位这一步做不干净，后面的掩膜、验证、与实测匹配全都白搭。

5.1 MOD03 的经纬度怎么用

定位文件里的 Latitude 和 Longitude 是 float32，已经是十进制度，不需要缩放，和发射波段同为 1 km 网格，逐像元一一对应。这两个数组就是构建 GLT（地理查找表）的全部输入。

需要提醒的是这两个数组的维度顺序。有的工具读出来是 [行, 列]，有的读出来是 [列, 行]，搞反了以后 GLT 会把整幅图转 90°，症状是海岸线方向明显不对但图像本身还挺"正常"，第一次遇到很容易懵。判断方法很简单：取第一行所有经度，如果是单调递增的，说明这一维是列方向；如果是单调递减或者来回摆动，说明维度搞反了。

传感器天顶角、方位角、太阳天顶角这些量都是 int16 加 0.01 缩放，需要在读取时统一转换。太阳天顶角对海温反演本身不直接参与计算，但可以用来做昼夜判断和耀斑区标记，顺手存一份不吃亏。

5.2 专用工具和手工建 GLT 的取舍

ENVI 在几何校正菜单下有针对 MODIS 一级产品的专用定位工具，可以直接读取辐亮度文件和定位文件，输出带地理坐标的栅格，里面还带了蝴蝶结校正的选项。这是最省事的路，尤其是做单景快速出图的时候，几分钟就能搞定。

但专用工具在批处理上有两个限制：一是它对文件命名和内部结构有约定，遇到重命名过的文件或者跨版本数据可能识别失败；二是它的重采样和输出投影选项相对固定，如果你想严格控制输出网格（比如必须对齐到某个已有的 0.01° 网格），就不太够用。这时候手工建 GLT 更自由：用 Build GLT 工具，X 和 Y 分别指定经度和纬度波段，设定输出投影和像元大小，生成 GLT 文件，再用 Apply GLT 把每个波段重采样过去。多波段的情况下，GLT 只需要建一次，反复套用。

我一般是这样分工：探索阶段用专用工具，图能出来、能看就行；进入正式生产，一律走手工 GLT，因为网格可控、可复现、可版本管理。

5.3 蝴蝶结效应和重采样带来的边缘误差

MODIS 的 1 km 扫描是每 10 行一组，在扫描边缘由于视角变化，相邻扫描带之间会有 2 到 4 行的重叠，形成所谓"蝴蝶结"效应。不做处理的话，拼接处会出现明显的重复或者错位条带，在海岸线附近特别扎眼。

处理方式有两种：一是用工具自带的蝴蝶结校正选项，二是在建 GLT 之后用重采样把重叠区"压平"。前者的代价是会改变原始像元的位置关系，后者的代价是重叠区被重复采样、权重不均。对海温这种空间连续性强的量，我倾向用前者，因为温度场本身平滑，位置准确性比辐射保真度更重要。

重采样方法的选择也有讲究。最近邻保真但会产生块状边界，在温度梯度大的锋面处容易出现阶梯；双线性平滑但会在云边缘把云和海的温度混在一起，制造出并不存在的"暖舌"。我的做法是先做云掩膜再做重采样，并且对掩膜做一次 1 到 2 像元的膨胀，让重采样时云边缘的污染区域被提前切掉。这个方法牺牲了一点有效覆盖，换来的是边界处可信度的明显提升。

6. 云、陆地、异常值：掩膜决定了最终图能不能用

这一步在流程里最不起眼，却最影响成果的可用性。剩下没被掩掉的杂散像元，会在图上变成孤立的冷点或热点，做统计的时候把均值拉偏，做锋面识别的时候制造假信号。

6.1 云掩膜字节的正确读法

云掩膜产品的核心是一个逐像元的字节数组，第一个字节的低两位编码了云的判定结果，大致是：00 表示有云，01 表示不确定，10 表示可能是晴空，11 表示确定晴空。只保留"确定晴空"是保守做法，覆盖会少一些；保留"确定晴空 + 可能晴空"是常规做法，代价是引入少量薄云污染。

薄卷云是这套流程里最阴险的污染源。它在大尺度云图上几乎看不见，对可见光的反射很弱，但热红外上会显著拉低亮温，而且两个通道被拉低的幅度不同，ΔT 反而可能变大，导致劈窗校正在错误的方向上又加了一笔。经验判断：如果某片区域的海温比周边低 2 到 3 K，形状不规则、边界模糊，而云掩膜显示"可能晴空"，多半就是薄卷云。

还有一个常被忽略的细节：云掩膜产品和辐亮度产品的成像几何是同一套，但存储分辨率和聚合方式可能不同，直接叠加会错位。正确做法是把云掩膜重采样到发射波段网格上，用最近邻，不要用双线性——掩膜是分类量，插值出来的中间值没有物理意义。

6.2 陆地掩膜与近岸像元的处理策略

海温反演的成果里出现陆地高温是最常见的低级错误之一。陆地掩膜的数据源有几个选择：定位文件里的海陆标记、可见光波段的反射率阈值、归一化植被指数。用反射率阈值最直接：取第 2 波段（近红外）反射率大于某个阈值的判为陆地，因为水体在近红外几乎全吸收，反射率很低，水陆反差极大。

近岸像元要单独讨论。混合像元里只要掺进一点陆地，发射率就偏离水体假设，劈窗的误差会成倍放大。常规做法是沿海岸线向内陆方向做 1 到 2 个像元的缓冲区，把缓冲区内的像元全部剔掉。这个操作会损失一部分近岸覆盖，但换来的是数据可信——做近岸应用的人往往更在意"这个数能不能用"，而不是"覆盖率高不高"。

浑浊水体是另一个麻烦。高悬沙浓度会让 11 µm 和 12 µm 的发射率差异发生变化，劈窗的 ΔT-水汽关系被扭曲。这种情况没有简单的掩膜办法，只能靠辅助的悬浮物产品做质量标记，或者干脆在成果里附一张"近岸低置信度"标识。

6.3 几条能自动跑的异常自检规则

人工逐景检查不现实，写几条规则让程序自己筛比较实际：

| 自检规则 | 判据示例 | 处理动作 | | ΔT 为负 | T31 - T32 < 0 | 标记为可疑，通常为云或陆地 | | ΔT 过大 | T31 - T32 > 5 K | 标记，多为薄卷云或异常气溶胶 | | 亮温越界 | T31 < 265 K 或 > 315 K | 剔除 | | 孤立点 | 3×3 窗口内与该像元温差 > 3 K | 剔除或替换为中值 | | 视角过大 | 天顶角 > 55° | 标记，精度不保证 |

这些阈值不是金科玉律，得按区域调。低纬海域的亮温下限可以放宽到 270 K，高纬冷季就得放到 265 K 以下。规则跑完之后建议把被剔除的像元单独存一层，做一次目视检查，你会发现很多规则的边界其实是需要微调的，而这个过程本身就是在积累对本地数据的认识。

3×3 中值滤波这一步我做得很克制。它对去掉孤立噪点很有效，但也会磨平真实的温度锋面，尤其是那种宽度只有两三个像元的窄锋。折中方案是只对"被标记为异常"的像元做中值替换，正常像元保持原值不动，这样既清了噪点，又保住了结构。

7. 怎么判断自己算出来的数是真的

图好看和数对是两件事。海温场的空间形态受海岸线、洋流、水深调制，只要定位没错，随便一套系数画出来的图都"看着挺像"，所以必须用独立数据来验证。

7.1 与实测数据匹配的具体做法

拿浮标、走航或者锚系观测做匹配，关键在于匹配规则要严。我用的规则是：时间窗控制在 ±3 小时以内（超过这个窗口，日变化引起的偏差可能超过 0.5 K），空间上取目标像元周围 3×3 窗口内的有效像元中值（避免单个像元的噪声），并且要求这 9 个像元里至少有 5 个是有效晴空像元，同时剔除距离云掩膜边界 1 个像元以内的样本。

匹配完之后不要急着算均方根误差，先画散点图看形态。如果散点沿 1:1 线两侧均匀分布，那是随机误差，可能是系数不够贴合；如果整体偏移一个常数，那是系统偏差，优先怀疑表皮-水体温度差、系数常数项、或者定标参数；如果散点在高温端偏离更大，说明系数的温度依赖性没抓住，得换带交叉项的模型。

那个"整体偏移"里最容易被误判为算法错误的是表皮温度效应。卫星看到的是微米级厚度的表皮，浮标测的是水体次表层，两者夜间差 -0.1 到 -0.3 K，白天在弱风强日照条件下可以差到 -0.5 K 甚至更多。做验证统计时如果不区分昼夜，这个偏差会被摊平成一个说不清来源的常数项，让你在系数上白折腾很久。

7.2 误差预算：知道每一项贡献多少

心里有一张误差预算表，排查的时候会快很多：

| 误差来源 | 量级 (K) | 性质 | | 传感器噪声 | 0.05 | 随机，星下点 | | 星上定标不确定度 | 0.1 到 0.2 | 系统性 | | 劈窗系数区域不匹配 | 0.3 到 0.6 | 系统性，最大项 | | 水汽与温度廓线偏差 | 0.2 到 0.4 | 部分被 ΔT 抵消 | | 薄卷云残留 | 0.5 到 3.0 | 单像元，重尾 | | 大视角临边效应 | 0.