基于高光谱成像的纤维图像高光检测与消除方法与流程

1.本发明涉及计算机辅助图像处理领域，具体是基于高光谱成像的纤维图像高光检测与消除方法。


背景技术：

2.高光谱图像中的每个像素具有连续的波长光谱，信息量巨大，是光谱分析的一大优势，极大地提高了测量分析的效率。目前，显微高光谱成像技术以光谱成像的形式采集微米级物体的高光谱图像，是一种非接触且非破坏的微米级物体颜色测量方法。然而，在高光谱成像系统在采集纺织品高光谱图像过程中，由于样品特性、相机与样本之间的相对布置、样品几何形状、照明光源与观测角度等的影响，少部分光线直接从纤维表面反射，产生高光现象。纤维表面因受高光影响而导致颜色信息丢失，增加了颜色测量的误差。
3.为了消除高光谱图像存在的高光，众多学者对此进行了广泛研究。bochko等人提出了一种用于检测图像中光谱和漫反射部分的概率主成分分析方法，不需要光源信息。fu等人提出了高光谱图像的正交子空间投影模型，这种方法需要光源的先验信息。然而，高光谱成像系统实际采集得到的高光谱图像，光源的先验信息并非保持稳定一致，同时高光消除的效率有待提高。
4.为有效恢复高光谱图像高光区域的漫反射光谱信息，通常期望能够从采集到的高光谱数据中消除镜面反射光谱信息，提高高光谱成像系统采集图像的准确性。因此，快速且有效地消除高光谱数据中的镜面反射分量是一个亟待解决的问题。


技术实现要素：

5.本发明要解决的技术问题是提供基于高光谱成像的纤维图像高光检测与消除方法，用以对高光谱图像采用结合局部空间明度差异性和高光约束的近似稀疏表示方法及用上述方法得到无高光的高光谱图像，消除高光影响，还原高光谱图像中受镜面反射影响的光谱信息。
6.为了解决上述技术问题，本发明提供基于高光谱成像的纤维图像高光检测与消除方法，具体过程包括：
7.s01、利用显微高光谱成像仪采集原始高光谱图像，然后在计算机中进行图像预处理；
8.s02、基于局部空间明度差异性，检测预处理后的高光谱图像中的高光区域，输出高光约束掩膜b；
9.s03、基于高光约束的近似稀疏表示方法，对预处理后的高光谱图像中的高光区域的混合光谱信号重复迭代分离为镜面反射分量、漫反射分量和高斯噪声分量，最终保留漫反射分量；
10.s04、将预处理后的高光谱图像的高光区域的混合光谱信号代替为漫反射分量，得到无高光的高光谱图像。
11.作为本发明的基于高光谱成像的纤维图像高光检测与消除方法的改进：
12.步骤s01中所述图像预处理按下式(1)处理获得预处理后的高光谱图像：
[0013][0014]
其中：w，h表示像元空间坐标，λ表示光谱波段；i
c(w,h；λ)
表示像素相对反射率单位归一化后的值，i
s(w,h；λ)
表示原始高光谱图像中的信号响应值，i
b(w；λ)
表示显微高光谱成像仪采集关闭系统镜头后的暗信号响应值，i
w(w；λ)
表示显微高光谱成像仪采集的标准白板数据的信号响应值，k为标准白板数据的校正系数。
[0015]
作为本发明的基于高光谱成像的纤维图像高光检测与消除方法的进一步改进：
[0016]
步骤s02的具体过程为：
[0017]
步骤s0201、对所述预处理后的高光谱图像按式(2)提取明度特征：
[0018][0019]
其中r(λ)为光谱反射率，s(λ)的标准光源的相对光谱功率，为cie1964标准三刺激值函数，光谱分辨率为10nm，y为高光谱图像明度值，y取值范围为0～100之间，k为比例系数；
[0020]
步骤s0202、提取图像像元的明度差异性特征
[0021]
局部空间的明度差异性方法描述如下：
[0022]
l
out
(q)＝σ
p∈ω
g
s
(p,q)g
y
(p,q,σ
y
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
[0023]
其中，l
out
(q)表示局部领域内中心像元与领域内其他像元的空间明度差异性，ω为领域内的像元集合，q为局部空间内的中心像元，p为领域内的像元，σ
y
为领域内的自适应方差，多层空间近邻度函数g
s
和明度差异性函数定义g
y
如下：
[0024][0025]
其中，m为子窗口内的像元数，y
p
为子窗口内的像元明度值，y
q
为子窗口内的中心像元明度值；
[0026][0027]
其中x，y为像元p，q的空间维坐标，σ
s
为空间领域的方差，z为子窗口的中心像元与领域内像元之间的近邻度之和；
[0028]
s0203、基于全图像统计分析：
[0029]
遍历所述预处理后的高光谱图像中所有空间像元的明度差异性l
out
，其中当l
out
(w,h)≥l
threshold
(w,h)时，x(w,h)为高光区域，l
out
(w,h)<l
threshold
(w,h)时，x(w,h)为非高光区域，从而获得具有高光区域标记的高光约束掩膜x(w,h)表示高光约束掩膜b的元素：
[0030][0031]
其中l
threshold
为判断是否为高光像元的阈值，mean(l
out
)为l
out
的平均值，std(l
out
)为l
out
的标准偏差，为缩放因子。
[0032]
作为本发明的基于高光谱成像的纤维图像高光检测与消除方法的进一步改进：
[0033]
步骤s03的具体过程为：
[0034]
s0301、将所述预处理后的高光谱图像转换为高光谱图像矩阵其中，w和h为预处理后的高光谱图像空间维的长度和宽度，λ为预处理后的高光谱图像的光谱维，e＝w*h；
[0035]
将步骤s0203得到的高光约束掩膜b在光谱维以相同值进行扩充，扩充的光谱维波段数为λ，并转换为高光约束掩膜矩阵高光约束掩膜矩阵w与高光谱图像矩阵a的矩阵大小相同；
[0036]
s0302、将高光谱图像矩阵a分解为低秩矩阵表示所述漫反射分量、稀疏矩阵表示所述镜面反射分量、高斯噪声矩阵e表示所述高斯噪声分量，将步骤s0301所得的高光约束掩膜矩阵w约束稀疏矩阵s，得到高光约束的稀疏矩阵，模型如下所示：
[0037][0038]
其中：rank(l)为低秩矩阵的秩函数；||
·
||0表示l0范数；
⊙
为哈达玛积，表示两个矩阵相应元素的相乘；α为稀疏矩阵的正则化参数，β为高斯噪声矩阵的正则化参数；
[0039]
s0303、高光谱图像的高光消除问题转化为求解核范数与l1范数组合的最小化：
[0040][0041]
其中，||
·
||
*
为矩阵的核范数，σ
i
(l)为矩阵l的第i个奇异值；||
·
||1为l1范数，e
i
表示高斯噪声矩阵e中的元素；s
i
表示稀疏矩阵s中的元素；
[0042]
将式(8)改写为增广拉格朗日形式：
[0043][0044]
其中<
·
,
·
>表示内积运算；y和y
s
为拉格朗日乘子；μ为惩罚系数；s

为松弛变量，s
‑
s

＝0，||
·
||
f
定义为frobenius范数，指向量各元素的平方和再计算平方根；
[0045]
s0304、通过交替方向乘子法重复迭代计算式(9)，计算误差ε，当误差ε小于预设阈值或迭代次数大于预设最大迭代次数，停止迭代，从而获得低秩矩阵l。
[0046]
作为本发明的基于高光谱成像的纤维图像高光检测与消除方法的进一步改进：
[0047]
步骤s04的过程具体为：
[0048]
将步骤s0304获得的所述低秩矩阵l转换为低秩高光谱图像l’；
[0049]
然后将步骤s01所述预处理后的高光谱图像中的高光区域数据替换为相同空间区域的低秩高光谱图像l’，最终在计算机上输出无高光的高光谱图像。
[0050]
作为本发明的基于高光谱成像的纤维图像高光检测与消除方法的进一步改进：
[0051]
步骤s0304中所述交替方向乘子法的具体过程为：
[0052]
(1)、输入所述高光谱图像矩阵a、高光约束掩膜矩阵w；
[0053]
(2)、初始化参数，包括：惩罚系数μ、误差阈值ε、迭代次数初始值k、迭代参数ρ、拉格朗日乘子y、y
s
和最大迭代次数k
max
；
[0054]
(3)、固定稀疏矩阵s、高斯噪声矩阵e和拉格朗日乘子y，估计低秩矩阵l：
[0055][0056]
令t
k
＝x y
k
/μ
k
‑
s
k
‑
e
k
，将式(10)改写为：
[0057][0058]
t
k
的奇异值分解为t
k
＝udv
t
，其中u和v均为单位正交阵，u称为左奇异矩阵，v称为右奇异矩阵，d为对角矩阵，d的特征为除了主对角线上的元素以外全为0，主对角线上的每个元素都称为奇异值；
[0059]
l
k 1
近似解为：l
k 1
＝us
η
(d)v
t
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)
[0060]
其中：
[0061]
其中，σ
i
为对角矩阵d的第i个奇异值，η＝1/μ
k
；
[0062]
(4)、固定低秩矩阵l、稀疏矩阵s和拉格朗日乘子y，估计高斯噪声矩阵e：
[0063][0064]
高斯噪声矩阵e近似求解：
[0065][0066]
其中g
δ
(
·
)定义为：
[0067]
g
δ
(x
λ
)＝sgn(x
λ
)max(|x
λ
|
‑
δ
λ
,0)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(16)
[0068]
(5)、用基于高光约束的近似稀疏表示方法表示稀疏矩阵s＝a
‑
l，固定低秩矩阵l、松弛变量s

、拉格朗日乘子y和y
s
，估计稀疏矩阵s为：
[0069]
[0070]
依据软阈值函数计算l1范数最小化问题：
[0071][0072]
(6)、更新松弛变量s

：
[0073][0074]
(7)、更新拉格朗日乘子y：
[0075]
y
k 1
＝y
k
μ
k
(a
‑
l
k 1
‑
s
k 1
‑
e
k 1
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(20)
[0076]
(8)、更新拉格朗日乘子y
s
：
[0077]
y
s
＝y
s
μ
k
(s
k 1
‑
s

)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(21)
[0078]
(9)、更新惩罚系数μ：
[0079]
μ
k 1
＝ρ
μ
ꢀꢀꢀ
(22)
[0080]
(10)、计算误差ε：
[0081][0082]
(11)、更新迭代次数k：
[0083]
k＝k 1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(24)
[0084]
(12)、重复步骤(3)～(11)，当误差ε小于预设阈值或迭代次数k大于预设最大迭代次数kmax，停止迭代，从而获得所述低秩矩阵l。
[0085]
本发明的有益效果主要体现在：
[0086]
1)、本发明融合图像局部领域的多层空间邻近度特征和明度差异性特征的局部空间明度差异性方法进行高光检测，对纤维表面高光区域的检测更客观，能快速检测高光谱图像中的高光区域；
[0087]
2)、本发明提出了高光约束的近似稀疏表示方法用于高光消除，对高光消除区域进行了区域约束，抑制算法对非高光区域的影响，有效分离高光区域混合光谱反射率的的镜面反射分量而不影响非高光区域，保留高光谱图像的空间特征和漫反射光谱信息；对表示镜面反射的稀疏矩阵进行了非负约束，使得算法不过度消除高光而导致高光区域变暗；
[0088]
3)、本发明采用近似稀疏表示方法，不仅提高了高光消除效果，同时减少了计算复杂性和计算时间，更适用于实际应用。
附图说明
[0089]
下面结合附图对本发明的具体实施方式作进一步详细说明。
[0090]
图1为本发明基于高光谱成像的纤维图像高光检测与消除方法的流程图；
[0091]
图2为图1中的高光检测的流程示意图；
[0092]
图3为图1中的高光消除的流程示意图；
[0093]
图4为本发明实验1羊毛纤维高光谱图像高光消除前后的效果图；
[0094]
图5为本发明实验2羊毛纤维高光谱图像高光消除前后的光谱反射率的曲线图；
[0095]
图6为本发明实验3高光约束的近似稀疏表示方法与高光约束的鲁棒主成分分析方法的高光消除迭代速度对比曲线图。
具体实施方式
[0096]
下面结合具体实施例对本发明进行进一步描述，但本发明的保护范围并不仅限于此。本发明所述方法的应用范围为高光谱成像系统获取图像时，由于受物体表面高光反射影响而在图像表面产生高光现象的高光谱图像，主要应用为基于高光谱成像的颜色测量方法所采集的高光谱图像。
[0097]
实施例1、基于高光谱成像的纤维图像高光检测与消除方法，如图1
‑
3所示，包括对高光谱图像进行高光检测和高光消除。本实施例所述图像是指基于显微高光谱成像仪采集的羊毛纤维高光谱图像。实施例中的高光检测与消除方法基于如下假设：高光谱图像中存在稀疏分布的高光区域；高光区域的混合光谱反射率由漫反射分量、镜面反射分量和高斯噪声分量组成；基于物理理论，漫反射分量、镜面反射分量具有非负性。具体地流程如下：
[0098]
步骤1、基于显微高光谱成像仪采集的羊毛纤维高光谱图像的高光检测流程，如图2所示，包括：显微高光谱成像仪采集的羊毛纤维的原始高光谱图像，将采集到的原始高光谱图像进行预处理，得到预处理后的高光谱图像，然后根据本发明提出的基于局部空间明度差异性的高光检测方法，在计算机中对预处理后的高光谱图像进行高光检测，包括提取高光谱数据明度特征、提取局部空间明度差异性、图像差异性统计，输出高光约束掩膜b和有高光区域标记的灰度图像，具体过程为：
[0099]
步骤1.1原始高光谱图像进行预处理：
[0100]
利用显微高光谱成像仪采集羊毛纤维的原始高光谱图像，在图像采集期间，噪声在图像传感器上累积，使得传感器响应值超出样本的真实强度；其次照明的不均匀性以及镜头表面上的灰尘可能会导致各种图像伪影，预处理的目的是消除上述影响，使得处理后的数据独立于成像系统和测量条件。对原始高光谱数据根据标准白板数据和暗光谱数据在不同坐标以及不同光谱通道的响应值进行校正，从而获得预处理后的高光谱图像，具体算法如式(1)所示：
[0101][0102]
式中：w，h表示像元空间坐标，λ表示光谱波段；i
c(w,h；λ)
表示像素相对反射率单位归一化后的值，i
s(w,h；λ)
表示样品原始高光谱图像中的信号响应值，i
b(w；λ)
表示显微高光谱成像仪采集关闭系统镜头后的暗信号响应值，i
w(w；λ)
表示显微高光谱成像仪采集的标准白板数据的信号响应值，k为标准白板数据的校正系数。
[0103]
步骤1.2、提取明度特征
[0104]
在计算机中，对预处理后的高光谱图像根据式(2)提取明度特征：
[0105][0106]
其中r(λ)为光谱反射率，s(λ)的标准光源的相对光谱功率，为cie1964标准三刺激值函数，光谱分辨率为10nm，y为高光谱图像明度值，y取值范围为0～100之间，y值越大表示明度越高，k为比例系数，保证计算所得的y的最大值为100。
[0107]
步骤1.3、融合图像局部领域的多层空间邻近度特征和明度差异性特征的局部空间明度差异性方法进行高光检测
[0108]
大多数高光检测方法中只单独考虑了每个像元的明度特征，而没有考虑图像中相邻像元之间存在的空间相关性。在高光谱图像空间中，相邻像元归属于同一种类别的可能性很大，即它们的光谱和颜色信息在一定程度上具有高度的相关性，当局部区域受到高光影响，导致该区域的光谱和颜色信息在一定程度上存在较大的差异性。因此，本发明提出融合图像局部领域的多层空间邻近度特征和明度差异性特征的局部空间明度差异性方法进行高光检测。
[0109]
局部空间明度差异性方法描述如下：
[0110]
l
out
(q)＝∑
p∈ω
g
s
(p,q)g
y
(p,q,σ
y
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
[0111]
其中l
out
(q)表示局部领域内中心像元与领域内其他像元的空间明度差异性，值越大，其差异性越大，ω为领域内的像元集合，q为局部空间内的中心像元，p为领域内的像元，σ
y
为领域内的自适应方差，g
s
为多层空间近邻度函数，g
y
为明度差异性函数；
[0112]
(1)高光区域的特征可以定义为高光中心区域的明度高于领域内像元的明度且局部空间差异性远高于图像的平均局部空间差异性。因此，本发明提出明度差异性定义如下：
[0113][0114]
其中m为子窗口内的像元数，y
p
为子窗口内的像元明度值，y
q
为子窗口内的中心像元明度值。当像元q与p的明度差异越大且q的明度大于p的明度时，g
y
越大，中心像元为高光像元的可能性上升；反之，当y
q
<y
p
时，g
y
为越小，即表示领域像元p比中心像元q的明度高，中心像元为高光像元的可能性下降。
[0115]
(2)当纤维表面存在高光时，高光区域的光谱和颜色信息受到较大的影响，导致子窗口内像元之间的相关性下降。根据中心像元q与领域内的像元p之间空间关系建立多层空间近邻度g
s
，多层空间近邻度描述如下：
[0116][0117]
其中x，y为像元p，q的空间维坐标，σ
s
为空间领域的方差，本实例中σ
s
＝3，z为子窗
口的中心像元与领域内像元之间的近邻度之和。
[0118]
步骤1.4、基于全图像统计分析：
[0119]
对步骤1.3获得的全图像的局部空间明度差异性l
out
统计分析图像明度差异性，根据统计结果将较大差异值标记为高光区域，得到二维的高光约束掩膜b，b的空间大小为w*h，w和h分别为空间维坐标，x(w,h)表示b的元素，当高光约束掩膜b中的元素x(w,h)等于0时表示该像元为漫反射像元，当x(w,h)等于1表示该像元为高光像元：
[0120][0121]
其中l
out
为式(3)计算的局部空间明度差异性，l
threshold
为判断是否为高光像元的阈值。mean(l
out
)为l
out
的平均值，std(l
out
)为l
out
的标准偏差，为缩放因子，通过控制缩放因子的大小，动态调整高光区域的范围。
[0122]
遍历预处理后的高光谱图像中所有空间像元的明度差异性，其中当l
out
(w,h)≥l
threshold
(w,h)时，x(w,h)为高光区域，l
out
(w,h)<l
threshold
(w,h)时，x(w,h)为非高光区域，输出高光约束掩膜b和具有高光区域标记的灰度图像。
[0123]
基于局部空间明度差异性的方法，同时融合了多层空间近邻度和明度差异性的特征，改进了高光检测的方式，降低了非高光现象引起的明度变化影响，对纤维表面高光区域的检测更客观。
[0124]
步骤2、基于显微高光谱成像仪采集的羊毛纤维高光谱图像的高光消除：提出基于高光约束的近似稀疏表示方法，对预处理后的高光谱图像中高光区域内含有镜面反射的混合光谱反射率信号重复迭代分离为镜面反射分量、漫反射分量和高斯噪声分量，最后保留漫反射分量；如图3所示，高光消除流程包括：将预处理后的高光谱图像转换为二维的高光谱图像矩阵；基于高光约束的近似稀疏表示方法重复迭代：将高光谱图像矩阵分解为低秩矩阵、稀疏矩阵和高斯噪声矩阵、计算误差并判断收敛条件；然后将重复迭代结果中的二维的低秩矩阵转换为三维的低秩高光谱图像；以及用低秩高光谱图像代替预处理后的高光谱图像中的高光区域得到无高光的高光谱图像，具体过程为：
[0125]
步骤2.1、将步骤1.1预处理后的高光谱图像转换为二维的高光谱图像矩阵其中，w和h为预处理后的高光谱图像空间维的长度和宽度，λ为预处理后的高光谱图像的光谱维，e为空间维w和h的乘积，即e＝w*h。
[0126]
将步骤1.4得到的高光约束掩膜在光谱维λ进行扩充，并转换二维的高光约束掩膜矩阵高光约束掩膜矩阵的矩阵大小与高光谱图像矩阵相同。
[0127]
步骤2.2、将高光谱图像矩阵a分解为三部分：低秩矩阵l，稀疏矩阵s和高斯噪声矩阵e，低秩矩阵l表示漫反射分量，稀疏矩阵s表示镜面反射分量，高斯噪声矩阵e表示高斯噪声分量，低秩矩阵l、稀疏矩阵s和高斯噪声矩阵e的矩阵大小与高光谱图像矩阵a相同，模型如下所示：
[0128]
[0129]
其中：式中为二维的高光谱图像矩阵，为低秩矩阵，l的取值范围约束在[0,1]，rank(l)为低秩矩阵的秩函数，定义为矩阵中非0奇异值的个数，约束矩阵的低秩性；为稀疏矩阵，为二维的高光约束掩膜矩阵，用于约束稀疏矩阵s得到高光约束的稀疏矩阵，||
·
||0表示l0范数，定义为矩阵中非零元素的个数，l0范数描述矩阵的稀疏性；
⊙
为哈达玛积，表示两矩阵相应元素的相乘；α为稀疏矩阵的正则化参数，β为高斯噪声矩阵的正则化参数。
[0130]
步骤2.3、由于式(7)是非凸问题，采用核范数和l1范数分别对rank(l)和l0范数进行凸松弛，因此，高光谱图像的高光消除问题转化为求解核范数与l1范数组合的最小化：
[0131][0132]
式中||
·
||
*
为矩阵的核范数，定义为矩阵奇异值之和，即σ
i
(l)为矩阵l的第i个奇异值；||
·
||1为l1范数，定义为矩阵元素之和，即e
i
表示高斯噪声矩阵e中的元素；s
i
表示稀疏矩阵s中的元素；为高光约束掩膜矩阵，作为稀疏矩阵s的高光约束条件；α为稀疏矩阵的正则化参数，β为高斯噪声矩阵的正则化参数。
[0133]
本发明将式(8)改写为增广拉格朗日形式，解决等式约束下的优化问题：
[0134][0135]
式中<
·
,
·
>表示内积运算；y和y
s
为拉格朗日乘子；μ为惩罚系数；α为稀疏矩阵的正则化参数，β为高斯噪声矩阵的正则化参数。s

为松弛变量，目的是将稀疏矩阵的非负不平等约束转换为平等约束，保证s
‑
s

＝0。
[0136]
步骤2.4、稀疏矩阵s和低秩矩阵l计算是独立子问题，本发明采用交替方向乘子法分步重复迭代计算式(9)中的每个变量，并计算误差ε，以下为交替方向乘子法的具体求解过程：
[0137]
(1)、输入高光谱图像矩阵a，高光约束掩膜矩阵w；
[0138]
(2)、参数初始化：惩罚系数μ＝1.25，误差阈值ε＝10
‑
5，迭代次数k的初始值为1，迭代参数ρ＝1.05，y＝max(e,λ)
‑
0.5
，y
s
＝max(e,λ)
‑
0.5
，最大迭代次数k
max
＝10
^4
；
[0139]
(3)、固定稀疏矩阵s、高斯噪声矩阵e和拉格朗日乘子y，估计低秩矩阵l：
[0140][0141]
令t
k
＝x y
k
/μ
k
‑
s
k
‑
e
k
，将式(10)改写为：
[0142][0143]
t
k
的奇异值分解为t
k
＝udv
t
，其中u和v均为单位正交阵，即有uu
t
＝i和vv
t
＝i，u称为左奇异矩阵，v称为右奇异矩阵，d为对角矩阵，d的特征为除了主对角线上的元素以外全
为0，主对角线上的每个元素都称为奇异值。||
·
||
f
定义为frobenius范数，指向量各元素的平方和再计算平方根。
[0144]
采用软阈值函数近似计算核范数最小化问题，因此l
k 1
近似解为：
[0145]
l
k 1
＝us
η
(d)v
t
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)
[0146]
其中：
[0147][0148]
其中σ
i
为对角矩阵d的第i个奇异值，η＝1/μ
k
。
[0149]
(4)、固定低秩矩阵l、稀疏矩阵s和拉格朗日乘子y，估计高斯噪声矩阵e：
[0150][0151]
高斯噪声矩阵e近似求解：
[0152][0153]
其中g
δ
(
·
)定义为：
[0154]
g
δ
(x
λ
)＝sgn(x
λ
)max(|x
λ
|
‑
δ
λ
,0)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(16)
[0155]
(5)、为了提高收敛速度，因此本发明提出了基于高光约束的近似稀疏表示方法将a
‑
l表示稀疏矩阵s，近似求解稀疏矩阵s，即s＝a
‑
l，固定低秩矩阵l、松弛变量s

、拉格朗日乘子y和y
s
，估计近似的稀疏矩阵s为：
[0156][0157]
依据软阈值函数计算l1范数最小化问题：
[0158][0159]
(6)、更新松弛变量s

：
[0160][0161]
(7)、更新拉格朗日乘子y：
[0162]
y
k 1
＝y
k
μ
k
(a
‑
l
k 1
‑
s
k 1
‑
e
k 1
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(20)
[0163]
(8)、更新拉格朗日乘子y
s
：
[0164]
y
s
＝y
s
μ
k
(s
k 1
‑
s

)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(21)
[0165]
(9)、更新惩罚系数μ：
[0166]
μ
k 1
＝ρμ
k
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(22)
[0167]
(10)、计算误差ε：
[0168][0169]
计算误差ε是为了达到误差阈值后停止迭代，即表明高光消除完成。
[0170]
(11)、更新迭代次数k：
[0171]
k＝k 1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(24)
[0172]
(12)、重复步骤(3)～(11)，当误差ε小于预设阈值或迭代次数k大于预设最大迭代次数k
max
，停止迭代，从而获得低秩矩阵l。
[0173]
步骤2.5、将二维的低秩矩阵l转换为三维的低秩高光谱图像l’；
[0174]
步骤2.6、将预处理后的高光谱图像中的高光区域数据替换为相同空间区域的三维的低秩高光谱图像l’，最终在计算机上输出无高光的高光谱图像。
[0175]
实验1：
[0176]
利用显微高光谱成像仪采集的羊毛纤维的原始高光谱图像，采集三张具有不同高光度的原始图像如附图4(a)所示，羊毛纤维的原始高光谱图像的显示表明高光随机分布在羊毛纤维表面上。
[0177]
然后按实施例1中步骤1所述，将采集到的原始高光谱图像进行高光检测，分别对应地获得三张具有高光区域标记的高光检测结果图像，如附图4(b)所示，其中白色部分，为标记的高光检测结果；图4(b)的检测结果显示本发明的算法能够对羊毛纤维表面存在的高光进行有效检测，排除非高光区域，同时保留纤维的空间细节，为高光消除方法提供有效的区域约束。
[0178]
再按实施例1中步骤2所述，对具有高光区域标记的高光检测结果图像进行高光消除，输出的无高光的羊毛纤维图像，如附图4(c)为高光谱图像高光消除处理后的无高光图像；通过对图像中的高光区域进行区域约束和非负约束，防止了算法对非高光区域的影响以及高光的过度消除，保留漫反射分量，使得修复的效果更合理。
[0179]
实验2：
[0180]
利用显微高光谱成像仪采集羊毛纤维的原始高光谱图像，按实施例1的步骤1.1原始高光谱图像进行预处理，获得的预处理后的高光谱图像，如附图5(a)所示，其中以a、b和c区域为例，a和c区域为高光区域，区域中的光谱发射率存在镜面反射分量，b区域为漫反射区域，光谱发射率不包含镜面反射分量。将预处理后的高光谱图像按实施例1所述经高光检测和高光消除处理后输出无高光的羊毛纤维图像，然后分别对预处理后的高光谱图像和无高光的羊毛纤维图像获取a、b和c区域的光谱反射率，得到a、b和c区域在高光消除前后的反射率曲线，如附图5(b)所示，其中实线为高光消除前的反射率曲线，虚线为高光消除后的反射率曲线。
[0181]
高光消除前后的反射率曲线显示：高光消除方法仅对羊毛纤维表面的高光区域进行高光消除，图中b区域作为漫反射区域，其在高光消除前后的光谱反射率差异微小，而在a和c的高光区域，高光消除前后光谱反射率变化较大。高光消除算法对高光区域有较大的改善，对于漫反射区域影响微小。因此，算法通过高光约束仅对含有镜面反射的高光区域消除镜面反射分量，保留漫反射分量。采用对数光谱相似性式(rmsle)、光谱角相似性(sam)以及cie de2000(2:1:1)色差公式对高光消除前后的反射率曲线差异进行分析，结果如表1，a和
c区域在高光消除前后光谱和色差有明显变化，而区域b的差异可以忽略不计。
[0182]
表1羊毛纤维高光消除前与高光消除后的光谱相似性和色差对比
[0183][0184]
实验3：
[0185]
将本发明提出的高光约束的近似稀疏表示方法(scsa)与高光约束的鲁棒主成分分析方法(sc
‑
rpca)在收敛速度上进行比较：
[0186]
利用显微高光谱成像仪采集羊毛纤维的原始高光谱图像，按实施例1的步骤1.1对原始高光谱图像进行预处理，获得的预处理后的高光谱图像按实施例1所述经高光检测、高光消除处理后输出高光消除后的高光谱图像，同时将获得的预处理后的高光谱图像经过鲁棒主成分分析方法(rpca)获得高光消除后的高光谱图像，然后采用峰值信噪比(psnr)评估了高光消除方法的性能，即采用的psnr计算公式将预处理后的高光谱图像与高光消除后的高光谱图像作为输入量，psnr值作为输出量，统计结果如图6所示，psnr值越低，表示高光谱图像中的镜面反射分量消除越多。随着迭代次数增加时，psnr下降幅度逐渐减小，高光消除的迭代误差越逼近误差阈值，表明镜面反射分量越少。
[0187]
相比鲁棒主成分分析方法(sc
‑
rpca)的高光消除方法，本发明提出的高光约束的稀疏近似表示方法(scsa)能够快速分离高光分量，峰值信噪比(psnr)快速下降，高光区域图像更接近真实的无高光图像，说明稀疏近似表示方法提高了高光消除性能，同时降低了计算时间。本发明提高了收敛速度，更适用于实际应用。
[0188]
与现有技术相比，本发明提出了基于局部空间明度差异性方法，检测高光谱图像的高光区域，输出结果作为高光消除的约束条件；提出了高光约束的近似稀疏表示方法将高光区域的混合光谱分离为镜面反射分量和漫反射分量和高斯噪声分量；本发明不需要预先确定有光源的光谱分布信息，能够根据图像特征进行高光消除；本发明提出的方法在不影响高光消除性能的基础上能够提高收敛速度，适用于实际应用。
[0189]
最后，还需要注意的是，以上列举的仅是本发明的若干个具体实施例。显然，本发明不限于以上实施例，还可以有许多变形。本领域的普通技术人员能从本发明公开的内容直接导出或联想到的所有变形，均应认为是本发明的保护范围。