Non-uniformity correction algorithm of space-borne long-wave infrared focal plane arrays based on improved Gamma curve

理论研究表明，红外探测器的电信号响应输出和入射红外辐射功率是连续单调递增函数关系，可表示为：

式中：${{R = }}{\tau _{eff}}{t_{{\rm{int}}}}\left[ {\pi {{\cos }^4}\theta /4{F^2}{\rm{ + }}1} \right]A$，$ {\tau _{eff}} $为光学系统透过率， ${{{t}}_{\rm{int} }}$为积分时间，θ为像元相对出射光轴的偏轴角，F为成像光学系统的F数，A为像元有效面积；$ {\lambda _{\rm{1}}} $、$ {\lambda _{\rm{2}}} $分别为光学系统截止波长的上下限；η为像元量子效率；D为像元暗电荷，$L{\rm{(}}\lambda {{,T)}}$为物空间单元在温度T下的光谱光子辐射功率，根据普朗克定律，其有如下函数关系：

式中：h为普朗克常量；c为光速。

由公式(1)、(2)可见，黑体温度、红外辐射强度和探测器响应的关系可用图1(a)所示曲线表示，红外探测器响应与光谱光子辐射是非线性对应关系，仅在较窄的特定动态范围内可近似为准线性关系^[7]。而在动态范围较大时，IRFPA探测元的输出响应信号与入射辐照度之间的非线性特性表现突出，表现为强非线性关系^[8]。进一步研究表明，实际探测元对红外辐射信号响应都有一个由截止到饱和的过程，响应曲线一般经历由截止到缓慢上升，准线性增长直到饱和，即响应曲线服从“S”型非线性分布^[9-10]。归纳而言，可将探测单元的曲线分为低辐射截止饱和区、低辐射非线性区、线性工作区、高辐射非线性区和高辐射截止饱和区五个区域^[11]。

Figure 1.  (a) Schematic diagram of the original response characteristic curves of the IRFPA detector; (b) Schematic diagram of response characteristic curves after correcting by traditional linear correction algorithm

从图1(b)可以看出，由于传统的星载非均匀性校正算法是建立在IRFPA探测元的响应是线性定常的假设之上的，经过非均匀性校正后，线性区域内的像元通过斜率和偏置的校正可达到一致，校正后的图像可以很好地保证线性区内的成像质量，然而在非线性区内，像元间的非均匀性难以消除，从而产生校正误差。

典型的“S”型曲线是伽马曲线^[12]，也称为通用Logistic曲线，最早由比利时科学家P.E.Vehulst导出其数学模型，它是一个广泛使用而且灵活的函数，该函数及其改进形式已被应用于政治学、种群动态、林木生长等多领域的模型建立^[13]。可通过函数参数配置拟合典型的“S”型曲线，曲线的拐点及拐点处的增长速率均可控制，其函数形式如下式所示：

式中：y_i和x_i为函数的变量和自变量；A_i为低端渐近线数值；B_i为高端渐近线数值；C_i为控制最大增长位置（拐点）；D_i为控制拐点速度。然而，将Gamma曲线直接用于红外探测器响应曲线拟合存在一定局限性，由于Gamma函数关于拐点对称，会导致不能够结合实际的探测器响应特性对强非线性区域进行增强描述。笔者课题组对大量实验室和在轨数据的统计特性进行分析，结果表明在高辐射非线性区，红外探测器响应的非线性程度明显强于低辐射非线性区。可见，传统的Gamma模型的拐点对称特性无法准确模拟实际红外探测器响应在高辐射区和低辐射区非线性程度的不同。因此，需要找到一种改进的“S”型曲线模拟红外探测器的响应模型。归纳而言，该模型应满足如下条件：

（1）函数不能是拐点对称的，从而可以分别描述低端温度和高端温度区域的非线性，且非线性程控可通过参数控制。

（2）通过调整模型参数可控制函数的斜率、截距、非线性区域范围、值域和定义域，从而定量化描述红外探测器输出增益、偏置、非线性工作区间、输出图像灰度值、输入动态范围。

（3）模型的函数结构形式简单，易于通过非线性压缩方法实现线性域和非线性域之间的映射和逆映射。

经过综合分析比较，在经典Gamma曲线4参数模型的基础上，增加参数t表征瞬时增长率（密度制约参数），公式（3）函数模型演化成如下非线性响应曲线：

式中：y_i和x_i分别表示探测像元i的辐照输入和响应输出；新增的参数t的取值范围为[0, ∞]，用于控制曲线在截止区和饱和区的非线性程度，使其具有表征像元非线性响应特性更好的适应性，同时可控制最大增长出现在低端还是高端；系数A_i的物理意义是探测器响应的最低值，表征探测器响应偏置，系数(A_i，B_i)共同决定探测器单元响应的动态范围；系数C_i控制探测器响应的非线性区起点位置；D_i决定曲线的增长速度，物理意义是探测器的广义增益。

传统的星载红外非均匀性校正算法建立在探测器的线性或分段线性响应模型基础之上，其基本原理是通过定期辐射定标获取密集温度点的定标图像数据，利用线性校正算法抑制像元间的响应不一致性。文中拟采用改进伽马曲线模型更真实地模拟探测器的实际响应特性，因此，若要借助星载线性化校正算法实现校正参数的估计和更新，则必须首先完成对非线性响应模型的线性化。图2所示的流程图给出了文中算法实现过程中的信息流示意，其中，如图2(a)和(b)所示，算法首先运用非线性映射（取对数运算）方法处理红外探测器的输出图像，将图像响应特性线性化；而后，如图2(c)所示，对线性化后的图像实施基于线性校正算法的非均匀性校正，获取均匀的响应图像；最后，如图2(d)所示，通过对线性化图像的校正图像进行逆-非线性映射（取指数运算）还原出红外探测器输出图像。

Figure 2.  (a) Algorithm flow step 1: get raw output of detector unit; (b) Algorithm step 2: non-linear mapping; (c) Algorithm step 3: linear non-uniformity correction; (d) Algorithm step 4: inverse non-linear mapping

由于IRFPA探测元响应具有时间漂移特性，因此需要通过周期性的黑体定标更新校正参数矩阵，用于消除模型参数时移引起的非均匀性校正残差。然而，黑体定标会导致星载红外成像停机，从而影响卫星在轨红外成像的重访能力和应急成像能力，因此，笔者希望尽可能降低每次黑体定标的时长。文中算法的创新性还在于可大幅降低黑体定标时采样温度点的数量，从而达到降低黑体定标时长的目的。

传统黑体定标方法无差别的设置密集定标温度点获取探测元的响应曲线特征，然而在线性区设置多个定标温度点对校正精度并无有效提升。文中算法采用提取IRFPA探测元响应曲线特征点的方式确定黑体定标温度点，特征点定义为曲线变化率最大的位置，即曲线的两个拐点，其物理意义分别是低辐射非线性区至线性区、线性区至高辐射非线性区的工作温度转换点。通过该方法，配合所建立的改进伽马曲线模型，采用更少的定标温度点即可实现对探测器响应曲线的准确拟合，同时，考虑到校正算法对黑体定标温度点附近的工作温度范围具有更好的校正效果，该方式还能够进一步消除建模误差，从而提升校正精度。

通过数学推导的方式来详细描述文中算法实现过程。文中算法可归纳为四个步骤，首先，对图像进行线性化变换；其次，在变换域采用线性校正算法完成非均匀性校正；同时，提取伽马曲线特征点作为下一节拍黑体定标温度点；最后，通过逆-非线性变换还原出实际的校正输出图像。图3所示的算法原理图直观地反映了文中算法完整的实施流程。

Figure 3.  Schematic diagram of the proposed algorithm

首先，通过非线性压缩方法对公式(4)两端进行移项及取对数运算，完成对IRFPA探测元输出信号的线性化：

由此可以看出，对探测元响应近似满足S形曲线的输出图像采用公式(5)进行线性化处理后，非线性响应${{{y}}_{{i}}}$转化为线性响应的${{{y'}}_{{i}}}$。

其次，在变换域对探测元进行线性非均匀性校正，校正输出${{{\hat Y}}_{{i}}}$可表示为：

式中：${{{X}}_{{i}}} = {\left[ {{{{{y'}}}_{{i}}}{{,}}1} \right]^{\rm T}}$为观测矩阵；${{{W}}_{{i}}} = {\left[ {{{{{\hat g}}}_{{i}}}{{,}}{{{{\hat o}}}_{{i}}}} \right]^{\rm T}}$为校正参数矩阵，其中${{{\hat g}}_{{i}}}$和${{{\hat o}}_{{i}}}$为校正系数的增益和偏移参数，${{{\hat g}}_{{i}}}$和${{{\hat o}}_{{i}}}$由高低温黑体辐射定标图像经公式(5)非线性变换后由以下公式获取：

式中：$\overline {{{DN}}} {{(}}{T_1}{{)}}$和$\overline {{{DN}}} {{(}}{T_{{2}}}{{)}}$分别表示在${{{T}}_1}$和${{{T}}_2}$温度点获取的黑体定标图像经公式(5)进行非线性化压缩后图像的均值，相应的${{D}}{{{N}}_i}{{(}}{T_1}{{)}}$和${{D}}{{{N}}_i}{\rm{(}}{T_{\rm{2}}}{\rm{)}}$为对应的第i个像元在${{{T}}_1}$和${{{T}}_2}$温度点获取的数值。将${{{\hat g}}_{{i}}}$和${{{\hat o}}_{{i}}}$代入公式(6)即可获得线性化图像的校正输出${{{\hat Y}}_{{i}}}$。

同时，获取伽马曲线的特征作为黑体定标温度点。令伽马曲线二阶导数为0，可得下一节拍高低温黑体定标温度点${{{T}}_{\rm{H}}}$和${{{T}}_{\rm{L}}}$如下式所示：

最后，利用公式(6)中获得的线性化图像的校正输出${{{\hat Y}}_{{i}}}{{(k)}}$可看做${{{y'}}_{{i}}}{{(k)}}$的估计值，因此，将${{{\hat Y}}_{{i}}}{{(k)}}$代入公式(5)，并对等式各项进行逆-非线性变换，采用公式(11)即还原出实际的图像输出：

为了检验星载长波红外成像通道在轨大动态范围工作条件下文中算法对空间非均匀性的抑制效果，将文中算法与传统的基于探测器线性响应模型的多点定标非均匀性校正算法进行对比，多点定标算法的黑体定标温度点选取10个，采用分段线性法校正原始图像。同时，还给出了典型的基于曲线拟合的非均匀性校正算法^[14]作为对比。测试数据分别采用了实验室定标黑体图像序列数据和基于在轨实际非均匀性特性提取的半实景图像序列数据，通过性能指标定量化分析和视觉效果两方面验证文中算法针对国产化大阵列红外探测器强非线性非均匀性校正的效果。

校正效果定量化性能评价指标采用国家标准规定的非均匀性指标^[15]${{{U}}_r}$和图像粗糙度指标^[2]$\; \rho $。其中，指标${{{U}}_r}$为红外探测器各像元响应均方根偏差与平均响应的百分比，表达式如下：

式中：$ \bar R $和${R_{{{ij}}}}$分别表示红外探测器输出图像的全局平均值和i像元j时刻响应输出值。

指标$\; \rho $用于衡量给定图像的粗糙程度，其表达式为：

式中：${{{h}}_{\rm{1}}}$为水平差分模板[1,−1]；${{{h}}_{\rm{2}}} = {{{h}}_{\rm{1}}}^{\rm T}$表示垂直差分模板；$ {\left\| \cdot \right\|_1} $表示${{{L}}_{\rm{1}}}$范数。指标${{{U}}_r}$和$\; \rho $的取值越小，表示图像的响应一致性越好，表明校正后图像的残留非均匀性越小，校正效果越好。

采用实验室黑体辐射定标系统采集不同黑体温度条件下的定标图像，分别采用文中算法和传统星载多点校正算法完成辐射校正，并对比两种校正算法在大动态范围信号输入条件下的校正性能。为充分模拟空间环境，辐射定标在真空低温环境下进行，真空罐真空度小于10⁻³ Pa，热沉温度小于100 K，定标黑体模拟240~340 K范围相机入瞳处的辐射亮度。校正参数采用270 K和300 K黑体图像进行计算，待校正图像选取动态范围内240、275、305、340 K共4个温度点，校正前后图像灰度值见图4，其中图4(a)为校正前4个温度点图像灰度分布，图4(b)为采用文中算法的校正结果，图4(c)为采用传统算法的校正结果，量化对比结果如表1所示。

Figure 4.  (a) Image gray distribution before correction; (b) Result for NUC of the proposed algorithm; (c) Result for NUC of typical curve fitting based algorithm;(d) Result for NUC of multi-point correction algorithm

采用非均匀性指标${{{U}}_r}$和图像粗糙度指标$\; \rho $量化比较两种算法对黑体图像的校正结果如表1所示。

从黑体定标实验验证结果可以看出，传统的星载多点校正算法由于模型失配的原因，校正效果不理想，残留分均匀性较大。而经过文中和参考文献[14]所述的曲线拟合算法校正后，在中等温度区域，两种算法均能得到较好的校正效果；然而，参考文献[14]应用场景为地面热成像，工作动态范围相对较低，在天基红外遥感应用领域具有一定局限性，表现为探测器响应非线性程度更强的高度温度区域校正效果不佳。文中算法通过改进的5参数S曲线模型，以其更好的自由度实现了对于探测器高温区域强非线性的准确拟合，从而具有更好的校正效果，表现了文中算法对探测器响应强非线性的良好适应性。

将卫星红外成像通道在轨实际获取的非均匀性特性叠加在一幅典型的标准长波红外遥感图像上，分别采用文中算法、参考文献[14]算法和传统星载多点校正算法完成辐射校正，并对比分析校正效果，校正效果如图5所示。从视觉效果的对比而言，相对于传统星载多点校正算法，文中和参考文献[14]所述的曲线拟合算法对表现为条纹图案的非均匀性噪声抑制效果更好，校正后的图像场景细节更加平滑清晰。文中算法的优势还在于拟合探测器响应曲线所需的黑体定标温度点数量相对更少，从而可以大幅降低星载遥感红外相机停机定标时间。

Figure 5.  (a) Semi-real infrared remote sensing image before correction; (b) Result for NUC of the proposed algorithm; (c) Result for NUC of curve fitting based algorithm; (d) Result for NUC of multi-point correction algorithm

此外，再次借助参数${{{U}}_r}$和$ \;\rho $量化对比上述两种校正算法对实景红外图像的校正效果。表2所示的实验结果表明经文中算法校正后的红外图像具有相对较小的粗糙度和残余非均匀性，校正效果更好。

星载红外探测器的非线性响应特性是当前天基红外遥感非均匀性校正中亟待解决的技术难题。传统的星载非均匀性校正算法能在一定程度上消除红外探测器非均匀性，然而，传统算法是基于探测元的线性响应模型导出的，天基应用要求红外探测器工作在大动态范围，探测元实际响应的非线性势必会引起较大的校正误差。为此，提出一种基于动态伽马函数的星载红外非均匀性校正算法，该算法采用更接近于探测器实际响应特性的改进伽马曲线模型以拓宽信号响应的动态范围，然后利用非线性压缩运算将校正过程映射为基于线性化模型的校正参数估算，同时根据模型参数动态确定黑体定标温度点，最终获得准确的校正输出。针对黑体图像和实景图像的实验结果表明，文中算法具有动态范围大、校正精度高和黑体定标效率高的优点。