气动热辐射效应校正方法

1.本发明属于图像处理技术领域，具体涉及一种气动热辐射效应校正方法。


背景技术：

2.当前，高速飞行器是航空航天技术领域的一个重要方向，而气动光学是空气动力学和光学相结合的一个交叉学科，也是飞行器发展的重要支撑方向。在气动光学校正这个领域里，气动光学热辐射校正恰好是应用中的一个重要方向。
3.退化图像的热辐射效应主要来自高温光学窗口自身辐射，是一种加性空间非均匀性噪声模型。高速飞行器在大气层内飞行时，光学图像由于气动加热而处于严重的气动热环境中，从而产生强烈的热辐射效应，进而影响目标的探测和识别。此外，带有光学成像探测制导系统的高速飞行器在大气层内飞行时，其光学头罩与来流之间发生剧烈的相互作用，引起周围大气显著变化。此外由于空气粘性的作用，与光学窗口表面相接触的气流将受到阻滞，使得气流速度降低，在窗口表面附近形成边界层。边界层内具有很大速度梯度的各层会产生强烈的摩擦，气流的动能会不可逆地转变成热能，造成窗口壁面温度的升高。光学窗口由于被气动加热而处于严重的气动热环境中，因此光学图像产生了热辐射背景噪声，使成像饱和以及信噪比劣化，降低了光学成像探测系统的信噪比和图像质量。
4.近年来，气动热辐射校正得到了一定的发展，但是至今未从气动热辐射强度中心和范围为突破点进行热辐射校正。此外，现有方法并未考虑强度中心和范围，使得较多高频干扰目标像素点参与到热辐射效应的估计中，导致低频热辐射效应层估计不够精确；与此同时，现有的方法是将同图像尺寸大小的零矩阵作为校正模型的初始值，加入到校正模型中，不断迭代和优化，整个过程需要求解多个参数，并且需要进行大量的迭代才能得到校正后的图像。


技术实现要素：

5.针对现有技术中存在的问题，本发明提供一种气动热辐射效应校正方法，解决气动热辐射校正精确差以及校正计算复杂的问题。
6.为实现上述目的，本发明提供了一种气动热辐射效应校正方法，包括以下步骤：
7.s1、获取气动热辐射退化图像s，通过图像分解算法分解为平滑层p和纹理细节层w；
8.s2、对平滑层p进行像素阈值标记处理，确定热辐射强度中心；
9.s3、获得平滑层p的梯度图，通过设定两个不同大小的邻域窗口遍历梯度图，确定强度中心范围；
10.s4、对平滑层p进行滤波处理，再根据强度中心范围内的像素点对平滑层p进行高斯曲面拟合，获得初始热辐射效应层的高斯分布曲面参数；
11.s5、对纹理细节层w添加帧波正则化约束，再结合高斯分布曲面b的正则化约束，建立基于范数最小化的图像校正模型；
12.s6、对图像校正模型进行split bregman迭代求解，得到校正后的清晰图像f。
13.进一步地，步骤s2中，对平滑层p进行像素阈值标记处理，确定热辐射强度中心具体为：
14.通过设置一定的灰度值阈值，在平滑层p图像中进行标记以确定标记范围，再通过图像中标记范围的像素，根据重心公式确定热辐射强度中心(x0,y0)：
[0015][0016]
其中，(xi,yi)表示标记像素的坐标，n是标记像素的总个数。
[0017]
进一步地，步骤s3中，通过设定两个不同大小的邻域窗口遍历梯度图，确定强度中心范围具体为：
[0018]
设置两个大小不同的窗口，从强度中心向图像边界遍历，当同时满足以下三个条件：
[0019]
(a)小窗口内像素间梯度的最大值大于设置的阈值；
[0020]
(b)大窗口内像素间梯度的最大值大于设置的阈值；
[0021]
(c)大窗口内梯度的最大值大于或者等于小窗口内梯度的最大值；
[0022]
则停止遍历，输出对应像素坐标位置(xi，yi)，以此来确定气动热辐射层的目标区域。
[0023]
进一步地，两个窗口的大小分别为3*3和5*5。
[0024]
进一步地，通过sobel算子获得平滑层p的梯度图。
[0025]
进一步地，步骤s4中，根据强度中心范围内的像素点对平滑层p进行高斯曲面拟合，获得初始热辐射效应层的高斯分布曲面参数具体为：
[0026]
对强度中心范围内的像素点进行等间隔或等比例选取，将选取的像素点和强度中心位置的像素点带入到高斯曲面函数中进行求解和拟合，从而获得初始高斯分布曲面的参数。
[0027]
进一步地，将选取的像素点和强度中心位置的像素点带入到高斯曲面函数中进行求解和拟合包括：
[0028]
高斯曲面函数为：
[0029][0030]
其中，z为获取的采样像素点坐标(x,y)的灰度值，(μ1，μ2)为强度中心坐标，σ为高斯函数的均方差，为高斯分布的幅值；
[0031]
将高斯曲面函数两边同时取对数，得：
[0032][0033]
再将上式进行转换，得：
[0034]
lnz＝k0+k1x+k2y+k3x2+k4y2[0035]
其中，
[0036]
高斯曲面的拟合，即求q＝min∑(k0+k1x+k2y+k3x2+k4y
2-lnz)2，由最小值条件可知分别对k0、k1、k2、k3、k4求偏导，得：
[0037][0038]
将上式化简，如下所示：
[0039][0040]
将上式转换为矩阵的形式，有：
[0041]
ak＝b
[0042]
其中，k＝(k0、k1、k2、k3、k4)是高斯曲面的拟合参数，此为线性方程的解问题，得：
[0043]
k＝a-1b[0044]
结合高斯曲面函数，此时拟合的曲面方程为：
[0045]
b＝exp(k0+k1x+k2y+k3x2+k4y2)。
[0046]
进一步地，对平滑层p进行3*3大小的中值滤波。
[0047]
进一步地，步骤s5中，对纹理细节层w添加帧波正则化约束，再结合高斯分布曲面b的正则化约束，建立基于范数最小化的图像校正模型具体为：
[0048]
对纹理细节层w添加帧波正则化约束‖wf‖1，再结合初始高斯分布曲面的正则化约束建立基于l1范数和l2范数的图像校正模型，图像校正模型表达式为：
[0049][0050]
其中，α、β、γ为惩罚函数，||.||1为l1范数，为l2范数，b表示高斯分布曲面，
表示梯度算子。
[0051]
进一步地，步骤s6中，对图像校正模型进行split bregman迭代求解，得到校正后的清晰图像f具体为：
[0052]
f的迭代求解过程如下：
[0053][0054]
令d2→
wf，并引入split bregman辅助变量b2，此时上式表示为：
[0055][0056]
使用split bregman迭代求解f和d2，首先对x求解过程如下：
[0057][0058]
对于清晰图像f的迭代求解公式进行更新表示：
[0059][0060]
其中k为迭代次数，且k为大于0的整数，而d2和b2的更新如下:
[0061][0062][0063]
而关于偏置场b的迭代求解过程如下：
[0064][0065]
令d3→
b，并引入split bregman辅助变量b3，此时上式表示为：
[0066][0067]
将初始高斯拟合曲面b作为初始值代入上式，即上式表示为：
[0068][0069]
使用split bregman迭代求解b和d3如下：
[0070][0071]
对于上式进行更新可表示：
[0072]
[0073]
故d3和b3的求解公式为：
[0074][0075]
对其求解得d3和b3如下：
[0076][0077][0078]
本发明与现有技术相比，具有以下优点及有益效果：
[0079]
本发明具有拟合曲面所需参数少、可以在校正中更好地保边缘、以及运算迭代次数少等特点，能够对气动热辐射领域中热辐射图像有着比较好的校正效果，满足在热辐射等干扰下对图像进行快速校正的要求。
附图说明
[0080]
图1为气动热辐射效应校正方法流程图；
[0081]
图2为实施例的气动热辐射效应校正示意图；
[0082]
图3为曲面拟合流程图；
[0083]
图4(a)是气动热辐射退化图像；
[0084]
图4(b)是图像分层后的纹理细节层；
[0085]
图4(c)是图像分层后的平滑层；
[0086]
图4(d)是平滑层中标记的热辐射强度中心；
[0087]
图4(e)是气动热辐射偏置场；
[0088]
图4(f)是校正后的清晰图像。
具体实施方式
[0089]
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。
[0090]
本发明针对现有气动热辐射相关领域的图像校正的改进需求，通过选取平滑层的强度中心和强度范围求取拟合曲面参数，并将拟合曲面设置为校正模型热辐射层的迭代初始值。此外，还在校正模型中加入正则化约束，这样可以对气动热辐射效应达到一个快速校正的效果。
[0091]
本发明的基于高斯分布和帧波正则化约束的气动热辐射效应校正方法，包括以下步骤：获取退化图像，对其进行分层处理，得到平滑层图像和纹理细节层图像；对平滑层图像进行像素阈值标记处理，获得热辐射强度中心；在确定热辐射强度中心后，对平滑层图像进行梯度计算操作，通过设定两个不同大小的邻域窗口遍历平滑层图像，其中大窗口内的梯度差值大于小窗口内的梯度差值，以此确定强度中心范围；对平滑层进行滤波处理，根据确定的强度中心范围对范围内的像素点等间距选取，再进行高斯曲面拟合，获得初始高斯分布曲面；对纹理细节层添加帧波正则化约束，再结合高斯分布曲面的正则化约束，建立基
于范数最小化的图像校正模型；对图像校正模型进行split bregman迭代求解，从而获得校正后的清晰图像。
[0092]
本发明提供的基于高斯分布曲面和帧波正则化约束的气动热辐射效应校正方法，拟合曲面所需参数少、可以在校正中更好地保边缘、以及运算迭代次数少，如图1所示，具体包括以下步骤：
[0093]
s1、从图像库中获取气动热辐射退化图像s，通过图像分解算法分解为平滑层p和纹理细节层w；
[0094]
s2、对平滑层p进行像素阈值标记处理，确定热辐射强度中心；
[0095]
s3、求出平滑层p的梯度图，确定两个不同大小的邻域窗口，其中大窗口内的梯度差值大于小窗口内的梯度差值，根据热辐射效应区域的梯度特性，以此确定强度中心范围；
[0096]
s4、首先对平滑层p进行滤波处理，再对s3确定的热辐射效应目标区域根据区域内的像素点等比例抽取，进行高斯曲面拟合，获得初始热辐射效应层b0的高斯分布曲面参数；
[0097]
s5、令潜在清晰图像f的初值f0＝w，即等于纹理细节层。为了在校正后图像中较好地保持图像边缘和一些重要细节特征，故对纹理细节层w添加帧波正则化约束‖wf‖1，再结合高斯分布曲面b的正则化约束，建立基于范数最小化的图像校正模型；
[0098]
s6、对校正模型进行split bregman迭代求解，得到校正后的清晰图像f。
[0099]
进一步地，本发明的步骤s1中进行图像分层的方法具体为：
[0100]
采用图像分层算法，将气动热辐射退化图像分解成平滑层和纹理细节层，其中，气动热辐射退化图像记为s，图像分层模型公式为：
[0101]
s＝p+w
ꢀꢀꢀ
(1)
[0102]
其中p表示平滑层图像，w表示纹理细节层图像。
[0103]
同时，气动热辐射退化模型对应的公式为：
[0104]
s＝f+b+n
ꢀꢀꢀ
(2)
[0105]
其中f表示清晰图像，b表示气动热辐射偏置场，n表示噪声。
[0106]
进一步地，本发明的步骤s2的具体方法为：
[0107]
对步骤s1中获得的平滑层p进行热辐射强度中心标记操作，通过设置一定的灰度值阈值并进行标记，从而在平滑层p图像中确定大致范围，再通过图像中标记范围的像素，根据重心公式：
[0108][0109]
其中(xi,yi)表示标记像素的坐标，n是标记像素的个数，再通过求取其重心从而确定热辐射强度中心位置。
[0110]
进一步地，本发明的步骤s3的具体方法为：
[0111]
通过sobel算子获得平滑层p的梯度图。以强度中心为中心点，设置两个大小不同的窗口，分别为3*3和5*5，从强度中心向图像边界遍历，当满足(1)小窗口内像素与像素间梯度的最大值大于设置的阈值；(2)大窗口内像素的梯度的最大值也大于阈值；(3)大窗口
内梯度的最大值大于或者等于小窗口内的梯度最大值，停止遍历，输出对应像素坐标位置，以此来确定气动热辐射层的目标区域。
[0112]
进一步地，本发明的步骤s4的具体方法为：
[0113]
对平滑层p进行3*3大小的中值滤波，过滤掉图像中干扰的噪声像素点，根据平滑层p的梯度图可知，平滑层p近似为一个高斯曲面，故可利用二维高斯曲面来对其进行拟合。所以，根据步骤s3中确定的热辐射效应目标区域，对目标区域内的像素进行间隔为4，像素坐标等比例取整的选取像素点，再根据区域内选取的像素点和强度中心位置的像素点带入到高斯曲面函数中进行求解和拟合，从而获得初始高斯分布曲面的参数。具体计算过程如下：
[0114]
二维高斯曲面函数为：
[0115][0116]
其中，z为获取的采样像素点坐标(x,y)的灰度值，(μ1，μ2)为步骤s2中确定的强度中心坐标，σ为高斯函数的均方差，为高斯分布的幅值。将公式(4)两边同时取对数，得：
[0117][0118]
再将公式(5)进行转换，可得：
[0119]
lnz＝k0+k1x+k2y+k3x2+k4y2ꢀꢀꢀ
(6)
[0120]
其中，
[0121]
高斯曲面的拟合，即求q＝min∑(k0+k1x+k2y+k3x2+k4y
2-lnz)2，由最小值条件可知分别对k0、k1、k2、k3、k4求偏导，得：
[0122][0123]
将公式(7)进一步化简，如下所示：
[0124][0125]
将公式(8)转换为矩阵的形式，有：
[0126]
ak＝b
ꢀꢀꢀ
(9)
[0127]
其中，k＝(k0、k1、k2、k3、k4)是高斯曲面的拟合参数，此为线性方程的解问题，可得：
[0128]
k＝a-1bꢀꢀꢀ
(10)
[0129]
结合公式(4)，此时拟合的曲面方程为：
[0130]
b＝exp(k0+k1x+k2y+k3x2+k4y2)。
ꢀꢀꢀ
(11)
[0131]
进一步地，本发明的步骤s5的具体方法为：
[0132]
由于帧波具有对称性，光滑性，以及更小的支撑域等特性，这些性质可以更好的表达图像。帧波将大小为m
×
n的图像用一个n维向量(n＝mn)来表示，并用矩阵记号w和w
t
表示帧波中的分解和重建。由于，基于帧波正则化方法能比传统的正则化技术获得更好的解，所以为了获得退化图像y中的热辐射层b和清晰图像f，对纹理细节层w添加帧波正则化约束‖wf‖1。再结合初始高斯分布曲面的正则化约束设计了一个基于l1范数和l2范数的图像校正模型，其校正模型表达式为：
[0133][0134]
其中，α、β、γ为惩罚函数，||.||1为l1范数，为l2范数，b表示高斯分布曲面，表示梯度算子。
[0135]
进一步地，本发明的步骤s6的具体方法为：
[0136]
首先，对图像校正模型进行split bregman迭代求解，然后，通过对初始清晰图像和初始热辐射层分别固定求解。
[0137]
关于清晰图像f的迭代求解过程具体如下：
[0138][0139]
令d2→
wf，并引入split bregman辅助变量b2，此时公式(13)可表示为：
[0140][0141]
使用split bregman迭代求解f和d2，首先对x求解过程如下：
[0142]
[0143]
对于公式(13)进行更新可表示：
[0144][0145]
其中k为迭代次数，且k为大于0的整数，而d2和b2的更新如下：
[0146][0147][0148]
而关于偏置场b的迭代求解过程具体如下：
[0149][0150]
令d3→
b，并引入split bregman辅助变量b3，此时公式(19)可表示为：
[0151][0152]
由于在公式(11)中已经得到初始高斯拟合曲面b，故在交替迭代中b只需将其作为初始值代入，即公式(20)可表示为：
[0153][0154]
使用split bregman迭代求解b和d3如下：
[0155][0156]
对于公式(22)进行更新可表示：
[0157][0158]
故d3和b3的求解公式为：
[0159][0160]
对其求解得d3和b3如下：
[0161][0162][0163]
具体地，如图2所示，本发明实施例的基于高斯分布和帧波正则化约束的气动热辐射效应校正方法，包括如下步骤：
[0164]
s1、采用图像分层算法，将气动热辐射退化图像，如图4(a)所示，分解成图4(b)所示的平滑层和图4(c)所示的纹理细节层，其中，气动热辐射退化图像记为s，图像分层模型
公式为：
[0165]
s＝p+w
ꢀꢀꢀ
(1)
[0166]
其中p表示平滑层图像，w表示纹理细节层图像。
[0167]
同时，气动热辐射退化模型对应的公式为：
[0168]
s＝f+b+n
ꢀꢀꢀ
(2)
[0169]
其中f表示清晰图像，b表示气动热辐射偏置场，n表示噪声。
[0170]
s2、对步骤s1中获得的平滑层p进行热辐射强度中心标记操作，通过设置一定的灰度值阈值并进行标记，从而在平滑层p图像中确定大致范围，再通过图像中标记范围的像素，根据重心公式：
[0171][0172]
其中(xi,yi)表示标记像素的坐标，n是标记像素的个数，再通过求取其重心从而确定热辐射强度中心位置，如图4(d)所示。
[0173]
s3、通过sobel算子获得平滑层p的梯度图。以强度中心为中心点，设置两个大小不同的窗口，分别为3*3和5*5，从强度中心向图像边界遍历，当满足(1)小窗口内像素与像素间梯度的最大值大于设置的阈值；(2)大窗口内像素的梯度的最大值也大于阈值；(3)大窗口内梯度的最大值大于或者等于小窗口内的梯度最大值，停止遍历，输出对应像素坐标位置，以此来确定气动热辐射层的目标区域。
[0174]
s4、对平滑层p进行3*3大小的中值滤波，过滤掉图像中干扰的噪声像素点，根据平滑层p的梯度图可知，平滑层p近似为一个高斯曲面，故可利用二维高斯曲面来对其进行拟合。所以，根据步骤s3中确定的热辐射效应目标区域，对目标区域内的像素进行间隔为4，像素坐标等比例取整的选取像素点，再根据区域内选取的像素点和强度中心位置的像素点带入到高斯曲面函数中进行求解和拟合，从而获得初始高斯分布曲面的参数，如图3所示。具体计算过程如下：
[0175]
二维高斯曲面函数为：
[0176][0177]
其中，z为获取的采样像素点坐标(x,y)的灰度值，(μ1，μ2)为步骤s2中确定的强度中心坐标，σ为高斯函数的均方差，为高斯分布的幅值。将公式(4)两边同时取对数，得：
[0178][0179]
再将公式(5)进行转换，可得：
[0180]
lnz＝k0+k1x+k2y+k3x2+k4y2ꢀꢀꢀ
(6)
[0181]
其中，
[0182]
高斯曲面的拟合，即求q＝min∑(k0+k1x+k2y+k3x2+k4y
2-lnz)2，由最小值条件可知分别对k0、k1、k2、k3、k4求偏导，得：
[0183][0184]
将公式(7)进一步化简，如下所示：
[0185][0186]
将公式(8)转换为矩阵的形式，有：
[0187]
ak＝b
ꢀꢀꢀ
(9)
[0188]
其中，k＝(k0、k1、k2、k3、k4)是高斯曲面的拟合参数，此为线性方程的解问题，可得：
[0189]
k＝a-1bꢀꢀꢀ
(10)
[0190]
结合公式(4)，此时拟合的曲面方程为：
[0191]
b＝exp(k0+k1x+k2y+k3x2+k4y2)。
ꢀꢀꢀ
(11)
[0192]
s5、由于帧波具有对称性，光滑性，以及更小的支撑域等特性，这些性质可以更好的表达图像。帧波将大小为m
×
n的图像用一个n维向量(n＝mn)来表示，并用矩阵记号w和w
t
表示帧波中的分解和重建。由于，基于帧波正则化方法能比传统的正则化技术获得更好的解，所以为了获得退化图像s中的热辐射层b和清晰图像f，对纹理细节层w添加帧波正则化约束‖wf‖1。再结合初始高斯分布曲面的正则化约束设计了一个基于l1范数和l2范数的图像校正模型，其校正模型表达式为：
[0193][0194]
其中，α、β、γ为惩罚函数，||.||1为l1范数，为l2范数，b表示高斯分布曲面，表示梯度算子。
[0195]
s6、首先，对图像校正模型进行split bregman迭代求解，然后，通过对初始清晰图像和初始热辐射层分别固定求解。
[0196]
关于清晰图像f的迭代求解过程具体如下：
[0197][0198]
令d2→
wf，并引入split bregman辅助变量b2，此时公式(13)可表示为：
[0199][0200]
使用split bregman迭代求解f和d2，首先对x求解过程如下：
[0201][0202]
对于公式(13)进行更新可表示：
[0203][0204]
其中k为迭代次数，且k为大于0的整数，而d2和b2的更新如下:
[0205][0206][0207]
而关于偏置场b的迭代求解过程具体如下：
[0208][0209]
由于在公式(11)中已经得到初始高斯拟合曲面b，故令lnd3→
lnb，并引入split bregman辅助变量b3，在交替迭代中b只需将其作为初始值代入，此时公式(19)可表示为：
[0210][0211]
将公式(20)进一步化简后，可表示为：
[0212][0213]
使用split bregman迭代求解b和d3如下：
[0214][0215]
对于公式(22)进行更新可表示：
[0216][0217]
故d3和b3的求解公式为：
[0218][0219]
对其求解得d3和b3如下：
[0220][0221][0222]
热辐射退化图像校正输出偏置场b和清晰图像f如图4(e)和图4(f)所示。
[0223]
需要指出，根据实施的需要，可将本技术中描述的各个步骤/部件拆分为更多步骤/部件，也可将两个或多个步骤/部件或者步骤/部件的部分操作组合成新的步骤/部件，以实现本发明的目的。
[0224]
本领域的技术人员容易理解，以上仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。