HJ-2B/IRS热红外数据交叉定标：基于双差法与高原湖泊的精度提升实践

1. 项目概述：为什么HJ-2B/IRS的辐射定标需要“二次校准”？

在定量遥感领域，数据质量是生命线。无论是反演地表温度监测城市热岛，还是分析火灾烟雾评估环境影响，我们依赖的底层数据都必须准确可靠。这颗于2020年发射的HJ-2B卫星，其搭载的红外扫描仪（IRS）为我们提供了宝贵的96米分辨率热红外数据。然而，就像任何精密仪器一样，它在太空极端环境中运行一段时间后，其“测量标尺”——即辐射定标系数——可能会出现微小的漂移或偏差。

官方机构（CRESDA）每年会通过星上黑体进行一次两点定标，发布官方定标系数。但这个方法有个固有局限：它只在两个特定温度点（通常是高温和低温黑体）进行校准，对于这两个点之外的温度，其精度依赖外推，而星上黑体自身的性能也可能随时间衰减。这就好比只用冰点和沸点来校准一支温度计，虽然这两点很准，但中间的温度读数就可能存在误差。我们在实际使用HJ-2B/IRS数据进行沙尘监测等研究时，也确实发现其辐射定标精度存在提升空间，这直接制约了数据的深入应用。

因此，我们需要一种不依赖于卫星自身硬件、能够频繁进行且精度更高的“外部校准”方法。这就是交叉定标的核心思想：找一个在轨运行稳定、精度公认极高的“尺子”（参考传感器），让它和我们的“待测尺子”（目标传感器）在相同时间、看着相同的地物，比较两者的读数，从而反推出待测传感器更准确的定标系数。本文要分享的，就是我们团队针对HJ-2B/IRS热红外通道，基于双差法，在高原湖泊这一理想实验场，开展的一次精细化辐射交叉定标实践。整个过程，我们以Terra/MODIS这颗“黄金标准”传感器为参考，最终将定标精度提升到了优于1K的水平。

2. 核心思路与技术选型：为什么是“双差法”与“高原湖泊”？

2.1 双差法：从“直接比对”到“误差抵消”的思维跃迁

最简单的交叉定标思路，是直接比较两个传感器对同一目标观测到的亮度温度或辐射值。但这会引入一个根本性问题：MODIS和IRS的光谱响应函数不同。即使看同一个温度完全均匀的湖面，由于它们“感受”热辐射的波长范围有细微差别，观测值本身就会有差异。这个差异不是传感器不准，而是物理特性使然。

传统方法需要借助复杂的大气辐射传输模型来模拟和修正这个光谱差异，但模型本身有误差，输入的大气廓线数据也有不确定性。这些误差会混入最终结果，让我们分不清到底是模型不准，还是传感器真的需要校准。

双差法的精妙之处，就在于它巧妙地抵消了这部分公共误差。它的核心步骤可以概括为“一模拟，二求差”：

1. 模拟差：利用辐射传输模型（我们用的是MODTRAN）和实时大气廓线数据，模拟在特定地表、大气条件下，MODIS和IRS“理论上”应该观测到的亮度温度差值（ΔT_simu）。这个差值纯粹是由两者光谱响应函数不同造成的。
2. 观测差：从实际卫星影像中，提取匹配像元，计算MODIS和IRS“实际上”观测到的亮度温度差值（ΔT_obs）。
3. 求取真差：将“观测差”减去“模拟差”，得到的结果（ΔT_calib）就剥离了光谱响应和大气效应的影响，理论上只包含两个传感器之间的相对定标偏差。由于MODIS的定标精度极高（不确定度<0.38%），这个偏差主要就反映了IRS需要被校正的量。

提示：双差法的本质是一种“差分消元”思想。它不追求绝对精确地模拟每一个物理过程，而是通过计算差值，让模型误差、大气输入误差等对两个传感器影响相近的部分在相减过程中被大幅削弱，从而更纯粹地暴露传感器自身的性能差异。

2.2 高原湖泊：近乎完美的天然定标场

选择青藏高原上的青海湖、扎陵湖、鄂陵湖、色林错和纳木错作为定标场，是经过深思熟虑的，主要基于其四大优势：

1. 空间均匀性好：湖泊水面在百米尺度上可以近似看作一个温度均一、发射率高的“灰体”，避免了像陆地那样因植被、土壤类型复杂带来的像元内温度剧烈变化。
2. 热容量大，温度稳定：水体热惯性大，温度变化缓慢。即使在卫星过境前后有1小时左右的时间差（如MODIS和IRS），湖面温度也不会发生剧烈变化，保证了时空匹配的有效性。
3. 大气条件洁净：高海拔地区空气稀薄，气溶胶和水汽含量低。大气对热红外辐射的吸收和散射作用弱，大大简化了辐射传输过程的复杂性，降低了大气校正的不确定性。
4. 角度依赖性弱：水面在热红外波段的发射率随观测角度变化较小，这意味着即使MODIS和IRS的观测角度不同，其观测到的辐射值差异也主要来自大气路径长度的不同，而非地表本身的方向性效应，便于进行角度校正。

这些特性使得高原湖泊成为国际公认的热红外传感器辐射定标和验证的优选场地。我们的任务，就是在这个“纯净”的实验室里，设计一套严谨的实验，把IRS的“尺子”校准得更准。

3. 数据准备与预处理：一切精确度的起点

巧妇难为无米之炊，高精度的定标始于高质量、匹配良好的数据。这一步的严谨程度直接决定了后续所有分析的可靠性。

3.1 数据源清单与关键参数

我们构建了一个由目标传感器、参考传感器、验证传感器和辅助数据组成的完整数据链条：

• 目标传感器：HJ-2B/IRS。我们使用其L1级数据，重点关注Band 8 (10.5–11.4 μm)和Band 9 (11.5–12.5 μm)两个长波红外通道。原始数据为数字量化值（DN），需要通过定标系数转换为辐射亮度。
• 参考传感器：Terra/MODIS。选择它是因为其热红外通道（Band 31: 10.78–11.28 μm, Band 32: 11.77–12.27 μm）具有经过近二十年验证的、卓越的定标稳定性和精度（NEΔT < 0.05 K）。我们下载了MOD021KM（表观辐射产品）、MOD03_L1A（地理定位和观测几何产品）和MOD04_3K（气溶胶产品）。
• 验证传感器：Landsat-8/TIRS。这是一把独立的“尺子”，用于交叉验证我们新定标系数的效果。其热红外数据本身也经过严格定标，可作为第三方基准。
• 辅助数据：
  • 大气廓线：欧洲中期天气预报中心的ERA5再分析数据。我们使用“压力层逐小时数据”，空间分辨率约0.25°，能提供过境时刻研究区域上空的温度、气压、湿度、臭氧垂直分布，是驱动辐射传输模型的关键输入。
  • 实测数据：在太湖开展同步观测试验时，使用CE312热红外辐射计测量水面辐射亮度，用于绝对精度验证。
  • 其他：Himawari-9静止卫星数据用于分析时间匹配误差，NCEP再分析数据用于不确定性分析中的大气廓线敏感性测试。

3.2 影像对筛选与匹配：苛刻的时空窗口

不是任意两景同时期的影像都能用来做交叉定标。我们制定了严格的筛选标准，在2023年共筛选出4组有效的IRS-MODIS影像对：

筛选逻辑解析：

1. 时间匹配：要求IRS和MODIS过境时间差小于1小时，以最小化因太阳高度角变化和地表温度日变化带来的影响。
2. 空间匹配：选择完全无云的晴朗天气，确保传感器观测的是纯净的湖面信号。
3. 几何匹配：优先选择两者都接近星下点（观测天顶角≤10°）的影像。但对于MODIS这种大幅宽传感器，星下点影像很少。因此，对于数据集1和4，我们接受了MODIS约20°的观测角，但后续会通过辐射传输模型对其进行严格的天顶角辐射校正。
4. 区域匹配：在影像上精确勾绘湖泊边界，确保在两个传感器的影像上选取的是同一片水体区域。

3.3 定标点筛选：剔除“坏点”，保留“好点”

即便在湖面上，由于近岸浅水区、岛屿、礁石或水下地形的影响，某些像元的温度也可能与主体水域有差异。我们必须筛选出真正均匀的像元作为定标点。

我们的策略是“CV阈值初筛 + RANSAC算法去噪”：

1. 窗口统计：由于IRS和MODIS分辨率不同（96m vs. 1km），我们不进行重采样（避免引入插值误差），而是直接匹配对应尺度的像元窗口。对于MODIS，选取3x3像元窗口（约3km x 3km）；对于IRS，选取10x10像元窗口（约960m x 960m）。计算每个窗口内像元的平均值（μ）和标准差（σ）。
2. 变异系数（CV）阈值筛选：计算每个窗口的变异系数CV = (σ / μ) * 100%。我们设定CV < 1%作为阈值，只保留那些像元值波动非常小的窗口，确保空间均匀性。
3. RANSAC算法剔除离群点：即使通过了CV筛选，由于湖泊形状不规则，某些窗口可能仍包含部分岸线或异常温度区。我们采用随机采样一致性算法。其原理是：随机选取一部分点拟合一个线性模型（辐射值关系），然后计算所有点到该模型的距离，将距离小于设定阈值的点视为“内点”，大于阈值的视为“外点”。经过多次迭代，保留内点最多的模型，并剔除所有外点。这种方法比普通最小二乘更稳健，能有效抵抗离群值的干扰。

经过这套组合拳，我们得到了数百个高质量、空间均匀的定标点对，为后续计算打下了坚实基础。

4. 核心算法实现：双差法定标系数求解全流程

有了干净的数据，接下来就是核心的计算过程。这个过程环环相扣，每一步都需要谨慎处理。

4.1 第一步：MODIS非星下点观测的辐射校正

对于数据集1和4，MODIS是以约20°的天顶角观测湖面的。观测角度越大，辐射穿过的大气路径越长，衰减越严重，测到的辐射值就越低。如果直接用这个斜视观测值与IRS的星下点观测值比较，会引入系统误差。

我们的校正方法是“模拟构建校正曲线”：

1. 输入真实场景：将过境时刻的ERA5大气廓线、湖面温度估计值输入MODTRAN模型。
2. 模拟不同角度辐射：固定地表和大气条件，让MODTRAN模拟从星下点（0°）到实际观测角（20°）之间，不同天顶角下MODIS传感器入瞳处的辐射亮度。
3. 计算辐射比：根据公式R(θ) = [L_obs(θ) - L_nadir] / L_nadir，计算每个角度下的TOA辐射差异比率。如图5所示，辐射值随角度增大而减小，且两个通道趋势一致。
4. 校正：利用公式L_nadir_corr = L_obs(θ) × [1 - R(θ)]，将斜视观测辐射值L_obs(θ)校正到等效的星下点辐射值L_nadir_corr。经过此步骤，我们成功将非星下点数据“归一化”到星下点条件，扩大了可用数据量。若不进行此校正，仅角度差异就可能引入0.1-0.3 K的亮度温度误差。

4.2 第二步：构建模拟亮度温度差函数

这是双差法的精髓所在，目的是量化纯粹由光谱响应差异引起的理论观测差。

1. 设置温度梯度：考虑到高原湖泊的实际温度范围，我们在MODTRAN中设置湖面温度从220 K到320 K，以1 K为步长变化。
2. 逐温度模拟：对每一个设定的湖面温度，结合过境时刻的ERA5大气廓线，运行MODTRAN模拟得到对应的大气顶部单波长辐射亮度L_TOA(λ)。
3. 光谱卷积：将模拟出的L_TOA(λ)分别与MODIS和IRS的光谱响应函数进行卷积积分，得到每个传感器在该温度下的等效辐射亮度L_simu。公式为：L_simu = ∫ [L_TOA(λ) * f(λ) dλ] / ∫ f(λ) dλ。这一步至关重要，它体现了传感器“看到”的辐射是经过其光谱响应加权平均后的结果。
4. 辐射转温度：利用普朗克公式的反函数，将L_simu_MODIS和L_simu_IRS分别转换为模拟亮度温度T_simu_MODIS和T_simu_IRS。
5. 计算并拟合差值函数：计算每一个温度点上的模拟亮度温度差ΔT_simu = T_simu_MODIS - T_simu_IRS。以MODIS观测亮度温度归一化后的黑体温度（作为横轴），以ΔT_simu为纵轴，拟合一条三次多项式曲线。这条曲线ΔT_simu = f(T_BlackBody)就是模拟亮度温度差函数。它告诉我们，在特定的大气和地表条件下，对于任何一个给定的黑体温度，MODIS和IRS“理论上”应该相差多少度。

4.3 第三步：应用双差法求解定标系数

现在，我们将理论模拟与真实观测结合起来。

1. 计算观测差：从筛选出的定标点中，读取经过角度校正后的MODIS辐射L_MODIS_obs，并转换为亮度温度T_MODIS_obs。同时，读取对应IRS窗口的平均DN值DN_IRS。
2. 计算理论IRS亮度温度：对于每一个定标点，根据其T_MODIS_obs找到对应的黑体温度，然后从第二步拟合的函数中，插值得到该黑体温度下的ΔT_simu。根据双差法原理，IRS“理论上”应该观测到的亮度温度为：T_IRS_theoretical = T_MODIS_obs - ΔT_simu。这个T_IRS_theoretical已经扣除了光谱响应差异的影响。
3. 温度转辐射：再次利用普朗克公式，将T_IRS_theoretical转换为理论辐射亮度L_IRS_theoretical。
4. 线性回归求解系数：现在，我们有了多组数据：DN_IRS(x) 和L_IRS_theoretical(y)。对这些数据点进行线性回归拟合，拟合方程L_IRS_theoretical = Gain * DN_IRS + Offset中的斜率Gain和截距Offset，就是我们最终求得的拟合定标系数。

注意：这里我们通常选择高、低温两个场景（例如夏季温暖湖泊和冬季寒冷湖泊）的数据点一起进行拟合，以覆盖更宽的温度动态范围，模拟星上黑体两点定标的原理，但数据点更多、更连续，从而得到更稳健的系数。

5. 结果验证与不确定性分析：新系数到底好在哪？

我们分别计算了2023年上半年和下半年两个时间段的FCCs。为了验证其有效性，我们将其与官方发布的OCCs进行了多角度对比。

5.1 与官方系数的直接对比

单纯比较Gain和Offset的数值大小没有物理意义，因为它们是协同工作的。更科学的做法是，在相同的DN值输入下，比较使用FCCs和OCCs计算出的亮度温度差异。

我们选取了从260 K到320 K这个常用的温度范围。结果显示：

• Band 8：FCCs与OCCs计算出的亮度温度差异很小（< 0.25 K），且基本不随温度变化。
• Band 9：差异相对较大，最大可达2.82 K，并且表现出明显的温度依赖性——在低温端差异更大，高温端差异减小。这说明官方系数在Band 9，尤其是低温区域，可能存在一定的非线性误差，而我们的FCCs通过更连续的温度采样，可能更好地刻画了传感器的响应特性。

5.2 第三方传感器交叉验证：Landsat-8/TIRS

我们选取了与HJ-2B过境时间接近的Landsat-8/TIRS影像（同样是高质量热红外数据），作为独立的“裁判”。分别使用FCCs和OCCs将IRS的DN值转换为亮度温度，然后与TIRS观测的亮度温度进行比对。

验证结果令人振奋：

• 2023年5月24日（色林错）：
  • Band 8: OCCs偏差(Bias)为-0.17 K，FCCs提升至-0.04 K。
  • Band 9: OCCs偏差高达0.96 K，FCCs将其显著降低至0.25 K。
• 2023年10月21日（青海湖）：
  • Band 8: OCCs偏差为0.45 K，FCCs为0.37 K。
  • Band 9: OCCs偏差巨大，达2.71 K，而FCCs将其控制在0.74 K。

在两个验证场景中，FCCs的均方根误差也全面低于OCCs。这强有力地证明，我们通过双差法得到的定标系数，在不同时间和地点都表现出了优于官方系数的精度，特别是在Band 9的提升尤为显著。

5.3 基于地面实测数据的绝对验证

为了进行更严格的绝对精度检验，我们于2023年10月17日IRS过境太湖时，使用CE312热红外辐射计进行了同步水面辐射测量。这是一种“自下而上”的验证方法：

1. 测量：获取水面辐射亮度L_CE312。
2. 大气校正：利用同步的ERA5大气廓线数据，通过MODTRAN模型模拟大气透过率和上行辐射，并结合IRS的光谱响应函数，计算得到传感器入瞳处“应该”接收到的理论辐射值L_true。
3. 对比：分别使用FCCs和OCCs，将IRS影像上对应点的DN值转换为辐射亮度，再换算为亮度温度。将这个计算值与基于地面测量推导出的L_true对应的亮度温度进行比较。

实测验证结果： 在太湖的5个采样点，FCCs计算出的亮度温度与地面实测推导值的偏差均小于OCCs。

• Band 8：OCCs平均偏差为-1.32 K，FCCs改善至-1.04 K。
• Band 9：OCCs平均偏差高达-2.43 K，而FCCs将其大幅降低至-0.92 K。

地面实测验证从绝对物理量的角度证实，FCCs显著提升了IRS热红外数据的辐射精度。

5.4 不确定性来源量化分析

任何定标方法都有其误差来源。我们对本方法涉及的七大不确定性进行了逐一量化分析，以评估结果的可靠度：

将以上各独立不确定度分量按平方和开根号的方式进行合成，得到总不确定度。对于2023年下半年得到的FCCs，Band 8和Band 9的总不确定度分别低于0.937%和1.257%。对于一个300 K的目标，这相当于亮度温度的校准精度分别优于0.64 K和0.93 K。这个精度水平完全满足大多数热红外定量遥感应用的需求。

6. 实操心得与常见问题排查

经过整个项目从设计到实现的全过程，我总结了一些对于后来者非常有价值的实操经验和可能遇到的坑。

6.1 关键步骤中的“魔鬼细节”

1. 光谱响应函数卷积：这是双差法中最容易出错的一步。务必确保你使用的光谱响应函数数据是准确的、波长间隔是均匀的。在编写卷积代码时，要检查积分区间是否完全覆盖了传感器的有效波段范围，避免在波段边缘出现截断误差。建议将卷积结果与文献中的典型值进行交叉检查。
2. 大气廓线时间匹配：ERA5是逐小时数据，而卫星过境是精确到分钟的。我们采用“最近邻”策略，即选择过境时刻之后最近的那个整点数据。虽然这引入了最多59分钟的时间差，但在高原稳定大气条件下，其影响在可接受范围内。对于对流活跃地区，可能需要考虑时空插值。
3. RANSAC参数设置：迭代次数、内点距离阈值和内点数量要求需要根据数据情况调试。迭代次数太少可能找不到最优模型；距离阈值太大会纳入过多离群点，太小则可能剔除有效点。一个实用的技巧是，先用全部数据画散点图，观察主体趋势，根据趋势线的残差分布来设定初始距离阈值。
4. 角度校正的适用性：我们的角度校正模型基于特定的大气和地表条件模拟。它对于清洁大气下的水体目标效果很好。但如果将此法应用于气溶胶含量高或地表类型复杂的区域，需要重新模拟建立校正关系，不可直接套用。

6.2 可能遇到的问题与解决方案

• 问题一：双差法处理后，ΔT_calib 仍然很大且不稳定。

  • 排查思路：
    1. 检查时空匹配：确认IRS和MODIS影像的配准精度是否足够高？时间差是否真的在1小时内？云检测是否彻底，是否有薄云污染？
    2. 检查均匀性筛选：CV阈值是否合适？RANSAC是否真的剔除了离群点？可以可视化筛选前后的定标点空间分布，看是否剔除了靠近岸边的点。
    3. 检查大气廓线：ERA5数据在湖泊上空是否准确？可以对比过境时刻附近的无线电探空仪数据（如果有）进行验证。
    4. 检查光谱响应函数：确认使用的光谱响应函数文件是否对应正确的卫星和传感器版本。

• 问题二：新定标系数（FCCs）在部分影像上效果很好，但在另一些影像上计算出的辐射值明显异常。

  • 排查思路：
    1. 季节与温度范围：检查异常影像的拍摄季节和地表温度。如果定标时只用了夏季高温数据，那么用于冬季低温场景时，外推可能导致误差。务必确保用于拟合FCCs的高低温定标点能覆盖目标应用的温度范围。
    2. 观测几何：检查异常影像的卫星观测天顶角和方位角。我们的方法假设地表为朗伯体，但实际存在一定的方向性。对于观测角度很大的影像，即使进行了大气路径校正，地表的方向性发射率效应也可能引入误差。
    3. 传感器状态：查询卫星工程数据，确认在异常影像获取时段，IRS传感器是否处于特殊工作模式或是否有已知的仪器异常。

• 问题三：与地面实测数据验证时，偏差仍然超过1K。

  • 排查思路：
    1. 实测数据代表性：CE312测量的是“点”尺度辐射，而卫星像元是“面”尺度。需要确保测量点位于大面积均匀水域中心，且测量期间水面平静无风，避免太阳耀斑。测量应在卫星过境前后半小时内完成，并取多次测量的平均值。
    2. 大气校正误差：这是最大的误差来源。确保用于MODTRAN模拟的大气廓线、气溶胶光学厚度等参数与实测时刻、地点尽可能匹配。如果可能，使用实测的大气参数。
    3. 光谱匹配因子：CE312和IRS的光谱响应不同，必须计算精确的星地光谱匹配因子进行校正。忽略这一步或计算不准确会直接导致偏差。

6.3 对业务化运行的思考

本研究验证了双差法用于HJ-2B/IRS辐射定标的有效性。若想将其业务化，用于长期监测传感器性能衰减，还需要建立一个自动化流程：

1. 自动化数据获取与预处理：编写脚本自动查询、下载匹配的IRS、MODIS、ERA5数据。
2. 稳定可靠的定标场库：建立全球多个稳定均匀目标（如深海、沙漠、南极冰盖）的数据库，避免单一区域天气不佳导致定标中断。
3. 趋势分析与预警：定期（如每月）计算FCCs，并与历史序列、官方OCCs进行对比。绘制关键参数（如Gain值）随时间变化的曲线，一旦发现明显的趋势性漂移，即可向用户发布定标更新或数据使用建议。

最后想说的是，辐射定标是一项既需要深厚物理功底，又需要严谨工程实践的工作。双差法为我们提供了一把锐利的“手术刀”，能够更精细地剥离出传感器自身的信号。但再好的方法，也离不开对数据质量一丝不苟的把控和对每一个误差源的深刻理解。这次高原湖泊的定标实践，让我们对HJ-2B/IRS的数据特性有了更深的把握，也为其在生态环境监测、防灾减灾等领域的精准应用增添了一份信心。未来，我们计划探索将深度学习等方法引入定标流程，进一步优化在定标样本有限情况下的模型构建，并尝试将这套方法推广到更多国产卫星传感器的在轨性能监测中。