一种基于三次样条的红外热图像盲元补偿算法

摘要

本发明提出了一种基于三次样条的红外热图像盲元补偿算法。首先，通过盲元检测，获得红外热图像的盲元分布图。其次，对于每一个盲元点，分别以红外热图像中与该盲元同行同列的像元的坐标为插值节点，构造行方向和列方向二个方向的三次样条，分别求出这二个三次样条在盲元点的插值，取二个插值的均值为盲元的初次补偿值。经过对红外热图像所有盲元进行初次补偿后，再对盲元进行二次补偿：对于每一个盲元点，分别以红外热图像中与该盲元同行同列的像元和经过初次补偿的其他盲元的坐标为插值节点，构造行方向和列方向二个方向的三次样条，分别求出这二个三次样条在该盲元点的插值，取二个插值的均值为该盲元的二次补偿值。遍历所有红外热图像中所有盲元，得到实现二次盲元补偿的红外热图像。

说明书

技术领域

本发明涉及红外热图像处理领域，具体涉及一种基于红外焦平面热成像系统中红外热图像的盲元补偿新机制的算法研究。

背景技术

目前，红外焦平面阵列(IRFPA)作为新一代探测器已广泛应用于军事，医疗，消防，交通，电力，建筑，电子制造和石化冶金等诸多领域。由红外焦平面阵列(IRFPA)构成的红外热视频动态成像系统的成像质量已经越来越多地影响和制约着相关行业的发展。由于半导体制造工艺水平的限制，使得(IRFPA)器件上难免会存在着某些缺陷。因此，使得所成红外热图像产生盲元(即坏元和响应异常像元)。这严重地影响了红外热视频成像系统的成像质量。因此，对盲元进行有效地补偿显得具有极其重要的意义。

目前，国内外已经提出了许多种盲元补偿方法。但总的来说可以分为两大类：第一类为基于时间序列的补偿法，它是利用序列图像的帧间相关性，从相邻帧获取盲元补偿信息。其优点在于能够很好地保持目标的边缘信息，但缺点是对前后帧的依赖性强；第二类为基于序列图像空间的补偿法，它是借助盲元周围像素的灰度信息对其进行补偿。一般使用相邻元替代，线性插值和中值滤波等方法。该类算法的主要优点在于流程简单，可操作性强，但目前使用的方法往往会使图像的细节部分变得模糊。因此，本专利提出了一种新的基于序列图像空间的补偿方法，这种方法利用三次样条函数插值来实现盲元补偿，不仅能够有效地剔除盲元，而且还能够很好地保持目标的边缘信息，所恢复出的图像具有良好的视觉效果。

盲元补偿是采用盲元周围的有效图像信息或前后帧的图像信息对盲元位置的信息进行预测和替代的过程。在参考国内外的各种检测算法的基础上，根据红外热图像自身的特点和计算数学上的插值理论，本专利提出一种基于三次样条函数插值的红外热图像盲元补偿新机制。

(一)三次样条函数的概念：

三次样条插值法是在插值节点之间，通过构造三次多项式(一个或一组不同阶多项式)来形成一条把所有主干节点连接起来的平滑的三次样条曲线，对插值节点处的函数值进行拟合的过程。

三次样条函数插值定义如下：给定区间[a，b]的一个划分：

Δ：a＝x₁＜x₂＜…＜x_n＝b

和一个函数f在这些分点上的值：

f₁＝f(x₁)，f₂＝f(x₂)，…，f_n＝f(x_n)

若函数S(x)满足下列条件：

(1)S(x_i)＝f_i，i＝1，…，n；

(2)S(x)在每个小区间[x_i-1，x_i]上是一个三次多项式：

S(x)＝a_i+b_ix+c_ix²+d_ix³，i＝2，…，n；

(3)S(x)，S′(x)，S″(x)在[a，b]上连续；

则称这样的三次样条函数S(x)为在区间[a，b]上的三次样条插值函数。

从上面的定义可知，S(x)有4(n-1)个参数需要确定。三次样条的插值要求给出了n个条件，三次样条的连续性要求给出了3(n-2)个条件，总共给出了(4n-6)个条件，尚缺的2个条件由边界条件给出。常用的边界条件有3种：

(1)Ⅰ类边界条件：S′(x₁)＝f′₁，S′(x_n)＝f′_n。

(2)Ⅱ类边界条件：S″(x₁)＝f″₁，S″(x_n)＝f″_n。

(3)Ⅲ类边界条件：S^(k)(x₁+0)＝S^(k)(x_n-0)，k＝1，2。

(二)三次样条插值函数的构造：

记m_i＝S″(x_i)，i＝1，…，n。因为S(x)在区间[x_i-1，x_i]上是一个三次多项式，因此，S″(x)是一次多项式，利用S(x)二阶导数的连续性，根据两点间距离公式，可得：

$S^{\prime \prime }\left(x\right)=m_{i-1}\frac{x-x_i}{x_{i-1}-x_i}+m_i\frac{x-x_{i-1}}{x_i-x_{i-1}},x\in [x_{i-1},x_i]---\left[1\right]$

对S″(x)在区间[x_i-1，x_i]上做二次积分，并利用S(x)的插值条件确定两次积分中出现的常数项，由此可得：

$S\left(x\right)=m_{i-1}\frac{{(x_i-x)}^3}{6h_i}+m_i\frac{{(x-x_{i-1})}^3}{6h_i}+(f_{i-1}-\frac{m_{i-1}}{6}{h}_i^2)\frac{x_i-x}{h_i}+(f_i-\frac{m_i}{6}{h}_i^2)\frac{x-x_{i-1}}{h_i}$

h_i＝x_i-x_i-1，x∈[x_i-1，x_i]    【2】

利用S(x)一阶导函数的连续性，进一步可得：

μ_im_i-1+2m_i+λ_im_i+1＝g_i    【3】

这里

${\mu }_i=\frac{h_i}{h_i+h_{i+1}},$${\lambda }_i=\frac{h_{i+1}}{h_i+h_{i+1}},$$g_i=\frac{6}{h_i+h_{i+1}}(\frac{f_{i+1}-f_i}{h_{i+1}}-\frac{f_i-f_{i-1}}{h_i}),$i＝2，…，n-1

方程组【3】中含有n个未知数m₁，m₂，…，m_n，而方程只有(n-2)个，另外2个方程由边界条件给出：

(1)在Ⅰ类边界条件下，由S′(x₁)＝f′₁和S′(x_n)＝f′_n可以导出以下2个方程：

$2m_1+m_2=\frac{6}{h_2}(\frac{f_2-f_1}{h_2}-f_1^{\prime })=g_1---\left[4\right]$

$m_{n-1}+2m_n=\frac{6}{h_{n-1}}(f_n^{\prime }-\frac{f_n-f_{n-1}}{h_n})=g_n---\left[5\right]$

将方程【4】【5】与方程【3】合在一起，即可得到确定m₁，m₂，…，m_n的线性方程组

(2)在Ⅱ类边界条件下，由于在区间端点处二阶导数m₁＝f″₁和m_n＝f″_n已知，所以方程组【3】中实际上只包含有(n-2)个未知数m₂，m₃，…，m_n-1，因此有：

(3)在Ⅲ类边界条件下，由S″(x₁+0)＝S″(x_n-0)和S′(x₁+0)＝S′(x_n-0)可以得到2个方程

m₁＝m_n                   【8】

λ_nm₂+μ_nm_n-1+2m_n＝g_n    【9】

其中

${\lambda }_n=\frac{h_2}{h_2+h_n},$${\mu }_n=\frac{h_n}{h_2+h_n},$$g_n=\frac{6}{h_2+h_n}(\frac{f_2-f_1}{h_2}-\frac{f_n-f_{n-1}}{h_n})$

将方程【8】【9】和方程【3】合在一起，即可得到确定m₁，m₂，…，m_n的线性方程组

发明内容

本专利提出了一种通过利用构造三次样条插值函数的方法，对红外热图像中的盲元的灰度进行插值，从而实现对整幅红外热图像的盲元补偿。

首先，通过盲元检测生成一幅表征红外热图像盲元原始分布信息的二进制标志图像(见图1)；

其次，以盲元分布图中任意一个盲元点为中心，分别找出该盲元点所在行和所在列中的像元，以其灰度值为元素构造行向量组和列向量组，再分别以其行向量组和列向量组中的元素为插值节点对该中心盲元点的灰度进行三次样条插值，取其行插值灰度和列插值灰度的均值灰度作为该中心盲元点的补偿灰度值。遍历图像中的所有盲元点，用此方法对其进行盲元补偿，从而实现对整幅红外热图像的初次盲元补偿(见图2)；

最后，以初次盲元补偿图像为基础，以图像中任意一个盲元点(已实现过初次盲元补偿的点)为中心，分别找出该盲元点所在行和所在列中的所有像素点(包含原始像元点和实现初次补偿后的盲元点)以其灰度值为元素构造行向量组和列向量组，再分别以其行向量组和列向量组中的元素为插值节点对该中心盲元点的灰度进行三次样条插值补偿，取其行和列插值灰度的均值作为该中心盲元点的补偿灰度值。遍历初次盲元补偿图像中的所有盲元点，用此方法对其进行二次插值盲元补偿，从而实现对整幅红外热图像的盲元补偿(见图3)。具体补偿过程如下：(见图5)

(一)盲元检测

设表示一幅带有盲元的红外热图像，这里i＝1，…，n，j＝1，…，n。通过盲元检测产生一幅表征红外热图像原始盲元分布信息的二进制标志图像矩阵(见图2)。具体盲元检测过程：(见图4)

(二)初次盲元补偿

设表示红外热图像中的一个盲元，这里1≤p≤n，1≤q≤n。又设该盲元所在行和所在列的有效像元分别为和这里K_p和L_q分别表示盲元所在行和所在列的有效像元数(见图1)。我们利用这些有效像元对盲元进行三次样条插值补偿，具体步骤如下：

(1)分别以有效像元和的坐标值为插值节点，按照背景技术介绍的方法，构造三次样条插值函数和

$S_{x_p}^1\left(y\right)={\hat{m}}_{j-1}\frac{{({\hat{y}}_j-y)}^3}{6{\hat{h}}_j}+{\hat{m}}_j\frac{{(y-{\hat{y}}_{j-1})}^3}{6{\hat{h}}_j}+(I_{(x_i,{\hat{y}}_{j-1})}-\frac{{\hat{m}}_{j-1}}{6}{\hat{h}}_j^2)\frac{{\hat{y}}_j-y}{{\hat{h}}_j}$

$+(I_{(x_i,{\hat{y}}_j)}-\frac{{\hat{m}}_j}{6}{\hat{h}}_j^2)\frac{y-{\hat{y}}_{j-1}}{{\hat{h}}_j}---\left[12\right]$

这里j＝2，…，K_p；

$S_{y_q}^1\left(x\right)={\hat{M}}_{i-1}\frac{{({\hat{x}}_i-x)}^3}{6{\hat{H}}_i}+{\hat{M}}_i\frac{{(x-{\hat{x}}_{i-1})}^3}{6{\hat{H}}_i}+(I_{({\hat{x}}_{i-1},y_j)}-\frac{{\hat{M}}_{i-1}}{6}{\hat{H}}_i^2)\frac{{\hat{x}}_i-x}{{\hat{H}}_i}$

$+(I_{({\hat{x}}_i,y_i)}-\frac{{\hat{M}}_i}{6}{\hat{H}}_i^2)\frac{x-{\hat{x}}_{i-1}}{{\hat{H}}_i}---\left[13\right]$

这里i＝2，…，L_q；

(2)将盲元点的坐标值(x_p，y_q)分别代入上述函数，求出点在(x_p，y_q)上的插值和于是，盲元的初次补偿值为：

${I^1}_{(x_p,y_q)}=\frac{S_{x_p}^1\left(y_q\right)+S_{y_q}^1\left(x_p\right)}{2}$

(3)按照上述方法，对红外热图像中的n×n个盲元点进行插值补偿，得到经过初次盲元补偿的红外热图像

(三)二次盲元补偿

设为原始红外热图像的一个盲元，这里1≤p≤n，1≤q≤n，利用初次盲元补偿得到的红外热图像对盲元进行二次补偿，具体步骤如下：

(1)利用初次盲元补偿的结果：根据背景技术介绍的方法，构造三次样条插值函数和其中(x_p，y_q)所在区间的三次样条插值函数分别为：

$S_{x_p}^2\left(y\right)={\tilde{m}}_{q-1}\frac{{(y_{q+1}-y)}^3}{6{\tilde{h}}_{q+1}}+{\tilde{m}}_{q+1}\frac{{(y-y_{q-1})}^3}{6{\tilde{h}}_{q+1}}+({\hat{I}}_{(x_p,y_{q-1})}-\frac{{\tilde{m}}_{q-1}}{6}{\tilde{h}}_{q+1}^2)\frac{y_{q+1}-y}{{\tilde{h}}_{q+1}}$

$+({\hat{I}}_{(x_p,y_{q+1})}-\frac{{\tilde{m}}_{q+1}}{6}{\tilde{h}}_{q+1}^2)\frac{y-y_{q-1}}{{\tilde{h}}_{q+1}}---\left[14\right]$

这里$y\in [y_{q-1},y_{q+1}],{\tilde{h}}_{q+1}=y_{q+1}-y_{q-1},$${\tilde{m}}_{q-1}={S_{x_p}^2}^{\prime \prime }\left(y_{q-1}\right),$${\tilde{m}}_{q+1}={S_{x_p}^2}^{\prime \prime }\left(y_{q+1}\right);$

$S_{y_q}^2\left(x\right)={\tilde{M}}_{P-1}\frac{{(x_{p+1}-x)}^3}{6{\tilde{H}}_{P+1}}+{\tilde{M}}_{p+1}\frac{{(x-x_{p-1})}^3}{6{\tilde{H}}_{p+1}}+({\hat{I}}_{(x_{p-1},y_q)}-\frac{{\tilde{M}}_{p-1}}{6}{\tilde{H}}_{p+1}^2)\frac{x_{p+1}-x}{{\tilde{H}}_{p+1}}$

$+({\hat{I}}_{(x_{p+1},y_q)}-\frac{{\tilde{M}}_{p+1}}{6}{\tilde{H}}_{P+1}^2)\frac{x-x_{p-1}}{{\tilde{H}}_{p+1}}---\left[15\right]$

这里$x\in [x_{p-1},x_{p+1}],{\tilde{H}}_{p+1}=x_{p+1}-x_{p-1},$${\tilde{M}}_{p+1}={S_{y_q}^2}^{\prime \prime }\left(x_{p+1}\right),$${\tilde{M}}_{p-1}={S_{y_q}^2}^{\prime \prime }\left(x_{p-1}\right).$

(2)将盲元的坐标值(x_p，y_q)代入上述函数，求出在(x_p，y_q)上的插值和于是，盲元的二次补偿值为：

$I_{(x_p,y_q)}^2=\frac{S_{x_p}^2\left(y_q\right)+S_{y_q}^2\left(x_p\right)}{2}$

(3)按照上述方法，对红外热图像中的n×n个盲元点进行补偿，得到经过二次补偿的红外热图像

本发明特点

(1)本专利提出了一种利用三次样条插值函数，基于序列红外热图像的空间信息对盲元进行有效补偿的新算法。对红外热图像的盲元补偿分三步实现：①通过盲元检测生成一幅表征红外热图像盲元原始分布信息的二进制标志图像。②以任意一个盲元点为中心，以其所在行和所在列中的像元坐标为插值节点，对盲元点的灰度进行行和列的插值补偿。以行和列插值灰度的均值作为实现初次盲元补偿后盲元点的灰度值，遍历红外热图像中的所有盲元点，对其进行插值补偿。从而实现对整幅红外热图像的初次盲元补偿。③以实现过初次盲元补偿后的盲元分布图为基础。再次，以其中任意一个盲元点为中心，以其行向量和列向量中的原始像元及所有实现过初次盲元补偿后的盲元的坐标为插值节点，对该盲元点进行行和列的插值补偿。以行和列插值灰度的均值作为该盲元点的二次补偿灰度值。遍历初次盲元补偿红外热图像中的所有盲元点，对其灰度进行插值补偿。从而实现对整幅红外热图像的盲元补偿。

(2)三次样条插值是“非局部化”的分段插值，每个插值节点上的数据都要影响到全局。但是，样条节点上的数据误差对离该点较远处的影响逐渐变小，因此样条插值具有良好的数值稳定性。用三次样条曲线去拟合时，其结果要比线性插值估计更接近于图像的真实情况。因此，本专利提出了一种采用更为平滑的三次样条函数进行插值的方法，对盲元点的灰度进行插值补偿。它能够很好地解决线性插值存在的缺陷，有效地剔除红外热图像中的盲元，得到良好的盲元补偿效果。

附图说明

图1.盲元分布图

图2.初次盲元补偿图，其中◎表示像元，○表示盲元

图3二次盲元补偿图，其中◎表示像元；表示实现初次补偿后的盲元

图4.盲元检测流程图

图5.盲元补偿流程图

具体实施方式：

(一)盲元建模

(1)盲元的定义及分类

盲元包括两类：死像元和响应异常像元。

死像元指像元响应率小于平均响应率的1/10的像元，它的输出信号不随输入辐射的改变而改变，它在原始图像上表现为亮点或者暗点。

响应异常像元指其响应值高于平均响应值4倍或者低于平均响应值1/4的像素元。响应异常像元可以正常工作，但是如果不被替代的话会引起较大的残留空间噪声。

(2)盲元的数学模型

①像元响应率

假设探测器(IRFPA)为M×N的红外焦平面阵列。在像元响应率为R(i，j)的红外焦平面阵列中，在一定帧周期和一定动态范围条件下，像元对每单位辐射功率产生的输出信号电压为：

$R(i,j)=\frac{V_s(i,j)}{p}$

式中，i＝1～M，j＝1～N；V_s(i，j)像元对应于辐射功率p的响应电压；p为第(i，j)像元所接收的辐射功率。

②平均响应率

红外焦平面阵列各有效像元响应率的平均值：

$\bar{R}=\frac{1}{m\times n-(d+h)}\underset{i=1}{\overset{m}{\Sigma }}\underset{j=1}{\overset{n}{\Sigma }}R(i,j)$

式中，m和n分别是红外焦平面阵列中像元的行数和列数；d和h分别是死像元数和过热像元数。实际测量中，d和h是经过多次迭代计算得到的。

③盲元率

红外焦平面阵列(IRFPA)的盲元数占总像元数的百分比，由下式表示：

$N_b=\frac{d+h}{m\times n}\times 100\%$

其中，d和h分别是死像元和过热像元数。

(二)盲元检测

红外焦平面阵列(IRFPA)中的盲元与正常探测单元在响应特性上有很大差异。通常采用的盲元检测方法包括：双参考辐射源法和基于场景的检测方法；正常探测单元的温度响应特性曲线在一定的动态范围内是线性的，而盲元的动态范围远离正常的探测单元的动态范围。针对盲元响应特性的这一特点，本文采用了基于两点参考辐射源的盲元自动检测技术。(见图4)具体操作过程如下：

(1)用一低温的均匀辐射源Φ₁照射各探测单元，记录每一探测单元的响应并存储起来；

(2)再用一高温的均匀辐射源Φ₂照射(IRFPA)各探测单元，记录每一探测单元的响应并存储起来；

(3)求出这两组数据的差值：$\Delta V_{(x_i,y_j)}=V_{(x_i,y_j)}\left({\Phi }_2\right)-V_{(x_i.y_j)}\left({\Phi }_1\right);$

(4)将这一差值与预设的比较门限值δ作比较；

(5)如果则认为该像(x_i，y_j)元为正常像元；如果则认为该像元(x_i，y_j)为盲元；

δ的选定对整个盲元的检测非常重要。从理论上说，由红外焦平面阵列(IRFPA)成像的数学模型可得到温度差值为ΔT条件下的红外响应灰度差值δ，即可将此值定为比较门限值。但是，由于该数学模型的建立较复杂，且应用价值不大。在实际应用中一般采用经验值来设定δ。即根据各单元对双参考源的实际响应数据及盲元定义的临界值进行有限次的迭代修正以确定δ。

(三)盲元补偿

(1)初次盲元补偿

①对红外热图像进行盲元检测，这里i＝1，…，n，j＝1，…，n。通过盲元检测产生一幅表征红外热图像原始盲元分布信息的二进制标志图像矩阵i＝1，…，n，j＝1，…，n。(见图2)。具体盲元检测过程(见图4)。

其中：

②取红外热图像中的任意一个盲元1≤p≤n，1≤q≤n。分别以该盲元所在行和所在列的有效像元和的坐标值为插值节点，按照背景技术介绍的方法，构造三次样条插值函数和

③将盲元点的坐标值(x_p，y_q)分别代入上述函数，求出点在(x_p，y_q)上的插值和于是，盲元的初次补偿灰度值为：

${I^1}_{(x_p,y_q)}=\frac{S_{x_p}^1\left(y_q\right)+S_{y_q}^1\left(x_p\right)}{2}$

④按照上述方法，对红外热图像中的所有n×n个盲元点进行插值补偿，得到经过初次盲元补偿的红外热图像

(2)二次盲元补偿

①取原始红外热图像中的一个盲元利用初次盲元补偿的结果：根据背景技术介绍的方法，构造三次样条插值函数和

②将盲元的坐标值(x_p，y_q)代入上述函数，求出在(x_p，y_q)上的插值和于是，盲元的二次补偿值为：

$I_{(x_p,y_q)}^2=\frac{S_{x_p}^2\left(y_q\right)+S_{y_q}^2\left(x_p\right)}{2}$

③按照上述方法，对红外热图像中的所有的n×n个盲元点进行插值补偿，得到经过二次盲元补偿的红外热图像从而最终实现对整幅红外热图像的盲元补偿。