CN103473759B - 一种WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法 - Google Patents

Abstract

本发明公开了一种WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法。首先提出一种WKPCA算法，对高维特征空间的各特征向量进行特征向量角匹配(FAM)加权，削弱或排除CRF区域内病态或异常的特征数据干扰，更精确的提取CRF主成分；在此基础上，定义一种同质度概念和计算方法，通过nCRF特征向量在中心主成分的投影，计算环境-中心的同质性；最后基于该同质度对nCRF中各抑制量校正，使得同质区域相互抑制量大、异质区域抑制量小或不相互抑制、同时尽可能削弱轮廓元素自抑制作用，从而提高抑制作用的准确率。因此本发明提出的模型能够更全面的检测环境-中心的差异，降低噪声干扰、更精确的抑制纹理细节、提高轮廓响应强度和完整性。

Description

一种WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法

技术领域

本发明属于一种基于视觉建模的复杂场景下微光图像显著轮廓提取方法，特别是一种WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法。

背景技术

轮廓提取在夜视图像理解分析方面发挥着重要作用。目前夜视目标探测识别方面的应用大部分是针对自然场景的，因此微光图像中包含了大量的自然纹理（例如树和草），传统的边缘检测算子的作用结果保留大量非轮廓的边缘成分（canny算子）。而且微光图像本身噪声干扰强，如何针对微光图像特征，去除这些由纹理和噪声所产生的局部非兴趣边缘，并且保持轮廓的完整性是夜视图像轮廓检测主要面临的问题。

针对复杂场景的轮廓提取问题提出了诸多解决方法，其中基于生物视觉机理的非经典感受野模型的轮廓提取在高质量可见光图像中获得了显著效果。视皮层（V1）神经元感受野（CRF）的大外周（非经典感受野nCRF）对CRF起调制作用，这种调制主要是抑制性的，能够实现同质区域相互抑制，而使得孤立的边缘要比群体边缘更为显著。基于侧抑制区的仿生模型，较好地去除了背景纹理产生的边缘。Grigorescu等人（ContourdetectionbasedonNonclassicalReceptiveFieldinhibition;ImprovedContourDetectionbyNon-classicalReceptiveFieldInhibition;Contourandboundarydetectionimprovedbysurroundsuppressionoftextureedges）利用nCRF的抑制特性进行轮廓检测，利用环境对中心的方向抑制，减少了环境纹理的影响，并提出各向异性抑制和各向同性抑制模型。Giuseppe等人（ABiologicallyMotivatedMultiresolutionApproachtoContourDetection）基于贝叶斯降噪和环境抑制技术提出一种生物启发的多分辨率轮廓检测技术。Ursino等人（Amodelofcontourextractionincludingmultiplescales,flexibleinhibitionandattention）引入层级注意机制，考虑不同尺度下背景抑制的轮廓提取。这些抑制模型能够一定程度的抑制背景纹理，提取出显著轮廓；但没有对中心环境进行特征差异分析，不能很好的解决异质成份互抑制和轮廓元素自抑制问题，可能导致轮廓响应微弱和轮廓间断。

为减小共线轮廓抑制的作用，桑农等人（Contourdetectionbasedoninhibitionofprimaryvisualcortex）根据nCRF刺激方位与CRF刺激方位差异，对抑制作用加权，建立了基于侧抑制区的蝶形模型。Zeng等人（Center-surroundinteractionwithadaptiveinhibition:Acomputationalmodelforcontourdetection;Contourdetectionbasedonanon-classicalreceptivefieldmodelwithbutterfly-shapedinhibitionsubregions）提出一种双尺度轮廓提取方法以及一种改进的方向选择抑制模型，该模型采用蝶形区域计算环境抑制。Giuseppe等人（Animprovedmodelforsurroundsuppressionbysteerablefiltersandmultilevelinhibitionwithapplicationtocontourdetection）利用可变滤波器和多层级抑制提出一种环境抑制模型。方向差异加权的抑制模型一定程度上减少了轮廓自抑制作用，但存在异向异质区域相互抑制的情况。而且相对于高质量的可见光图像，微光图像噪声干扰严重，轮廓局部方向特征不突出，环境抑制不准确。一方面高噪声、低对比度导致微光图像背景纹理无法完全抑制；另一方面抑制作用削弱轮廓强度，轮廓受到周边纹理的抑制，容易出现断裂，影响后续目标识别。

发明内容

本发明的目的在于提供WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法，该方法能够解决单一方位差加权抑制的不足，结合生物视觉机理，针对微光图像多维特征差异，从复杂场景下微光图像中有效抑制噪声和纹理，提取显著轮廓。

实现本发明目的的技术解决方案为：首先在高维特征空间研究一种特征向量角匹配(FAMFeaturevectorAngleMatching)加权KPCA(WKPCAWeightingKPCA)方法计算环境-中心同质度(HDHomogeneityDegree)：该方法抑制了CRF像元中的病态或异常数据，再对其进行主成分分析，提高主成分准确率；将nCRF各像元在中心主成分投影计算其与CRF中心同质度。其次基于该环境-中心多特征差异对nCRF抑制校正，构建一种同质度校正抑制模型：该模型在同质元素互抑制的同时，避免异质像元间相互作用及轮廓成份自抑制，实现微光图像中的噪声、纹理抑制，以及复杂自然场景下夜视图像的显著轮廓背景分离。

本发明与现有技术相比，其显著优点：

（1）本发明针对微光图像多维特征分析，首次提出一种WKPCA算法，在高维特征空间基于特征向量角匹配准则对CRF协方差矩阵加权，有效削弱异常特征数据和噪声干扰，提高CRF主成分分析的准确性。

（2）本发明创新性的定义一种同质度概念及其计算方法，采用WKPCA算法进行CRF与nCRF区域同质性分析，在高维特征空间nCRF各像元特征向量投影到CRF中心主成分，计算两者同质度。该方法能够更全面准确的表针环境-中心的差异。

（3）本发明提出一种WKPCA同质度校正nCRF抑制模型，根据nCRF各像元与CRF中心同质度对其抑制量校正。该模型更有效排除噪声干扰、削弱背景细节的同时，更精确的提高同质成份的抑制程度，降低异质像元、轮廓元素间的相互抑制，提高轮廓整体响应及减少轮廓间断的可能性。

附图说明

图1是nCRF距离衰减函数。

图2是WKPCA同质度校正nCRF抑制微光图像显著轮廓提取方法流程图。

图3是多维特征空间非经典感受野结构。

图4是以真实目标轮廓像元为中心的nCRF结构放大图，及nCRF中各像元特征向量与CRF中心均值向量FAM权值三维图。

图5是PCA、KPCA、WKPCA算法对高斯数据的分类效果比较。

图6是PCA、KPCA、WKPCA算法同质性分析结果比较。

图7是背景、目标、轮廓中心的WKPCA同质度校正抑制分析。行向从上到下分别是原图、WKPCA同质度值输出、Gabor能量极值输出以及有效抑制输出（同质度值与Gabor能量极值乘积）。

图8是各模型对复杂自然场景下微光图像轮廓提取效果比较。

具体实施方式

下面结合附图对本申请做进一步说明。

1、nCRF环境抑制：

二维Gabor函数能有效地描述视皮层简单细胞的感受野剖面，通过奇偶对简单感受野滤波器的反应模（Gabor能量），能很好地模拟典型复杂细胞的基本特性。这些复杂细胞可以看成局部方位能量算子，用复杂细胞活动的最大值可以对图形边与线进行准确定位，因此本发明通过Gabor能量来模拟复杂细胞的响应。二维Gabor滤波器表示如下。

其中θ是CRF的偏好方位；是相差；由和分别表示奇偶Gabor滤波器；长宽比γ决定高斯椭圆的离心率；λ是波长；高斯标准差σ是决定了CRF作用面积；σ/λ表示空间频率带宽。

简单细胞的响应为Gabor函数与输入图像I的卷积。由正交简单细胞响应定义CRF刺激的神经元响应E_σ(x,y;θ_i)。

$E_{\mathrm{\σ}}^2(x,y;{\mathrm{\θ}}_i)=r_{\mathrm{\λ\σ}0}^2(x,y;{\mathrm{\θ}}_2)+r_{\mathrm{\λ\σ\π}/2}^2(x,y;{\mathrm{\θ}}_i)---\left(3\right)$

其中Gabor能量被分为N_θ个响应方位：θ_i=iπ/N_θ,i=0,1,…,N_θ-1。根据[9]设置参数：N_θ=12，γ=0.5，σ/λ=0.56，选择σ作为自变量。CRF的响应输出为各方位Gabor能量极值。

E_σ(x,y)=max{E_σ(x,y;θ_i)|i=0,1,…,N_θ-1}(4)

生物研究表明V1区细胞的nCRF主要是抑制性的，这种环境抑制作用旨在分离同质区域，突显轮廓。在nCRF模型中采用距离衰减函数ω_σ(x,y)描述环境抑制的空间结构。ω_σ(x,y)定义为半波校正和L1规范的同轴高斯差函数DOG_σ,ρσ(x,y)。

${\mathrm{DOG}}_{\mathrm{\σ},\mathrm{\ρ\σ}}(x,y)=\frac{1}{\sqrt{2\mathrm{\π}{\left(\mathrm{\ρ\σ}\right)}^2}}\exp (-\frac{x^2+y^2}{2{\left(\mathrm{\ρ\σ}\right)}^2})-\frac{1}{\sqrt{2\mathrm{\π}{\mathrm{\σ}}^2}}\exp (-\frac{x^2+y^2}{2{\mathrm{\σ}}^2})---\left(5\right)$

ω_σ(x,y)=|DOG_σ,ρσ(x,y)|⁺/‖|DOG_σ,ρσ(x,y)|⁺‖₁

${\left|z\right|}^{+}=\left\{\begin{array}{cc}0& z<0\\ z& z\&GreaterEqual;0\end{array}\right.---\left(6\right)$

抑制形式t_σ(x,y)定义为环形背景中各像元的距离衰变局部方向能量。

t_σ(x,y)=ω_σ(x,y)*E_σ(x,y)(7)

如附图1所示，|DOG|⁺确保ω_σ(x,y)只作用于环境区域。σ和ρσ分别决定了CRF面积S_c与nCRF面积S_nc。生物实验表明nCRF直径能够达到CRF直径的2-5倍，因此取ρ=4。通过DOG_σ,kσ(x,y)=0，可以得到CRF半径r_c和nCRF半径ρr_c。

$r_c=\sqrt{2\ln \mathrm{\&rho;}/{\mathrm{\&rho;}}^2-1}\&CenterDot;\mathrm{\&rho;\&sigma;}---\left(8\right)$

nCRF环境抑制输出即为：

R_σ(x,y)=|aE_σ(x,y)-bt_σ(x,y)|⁺

(9)

2、同质度校正抑制

如附图2所示，为提高环境-中心抑制的准确性，使得与中心同质的区域提供较大抑制量、异质的区域不对中心抑制，以及轮廓区域自抑制减少，本发明在nCRF模型基础上，对各环境像元与中心进行同质性分析，提出一种同质度计算方法，对环境抑制成份校正，实现强噪声干扰下微光图像背景纹理抑制和显著轮廓提取。

为计算环境各像元与中心同质度δ_w(x,y)，在多维特征空间，拟采用KPCA方法对CRF中心像元进行主成分分析，并将nCRF中各像元在CRF主成分投影计算两者差异性。但由于CRF中心存在病态或异常数据(噪声数据、背景中夹杂的目标数据、目标中夹杂的背景数据等)，导致主成分提取不准确，因此本发明提出一种特征向量角匹配(FAM)加权KPCA方法(WKPCA)，抑制病态或异常数据，提高主成分分析准确率。

2.1多维特征分析

上述环形nCRF结构在各方位上同等抑制，会导致异质像元及轮廓元素之间相互抑制，容易削弱轮廓响应及破坏轮廓完整性。蝶形抑制模型就是采用一种环境-中心方位差加权方法，在抑制纹理细节的同时减少共线轮廓元素间的相互作用，其对高质量可见光图像轮廓提取性能优于各向同性、各向异性模型，但是仍存在异向异质像元互抑制作用。而且相比于高质量可见光图像，微光图像中边缘受噪声干扰Gabor方向能量不准确甚至无明显方向性，因此只具有方向差加权的蝶形抑制无法很好提取微光图像轮廓。

生物实验表明当CRF和nCRF内的图形特征一致时抑制作用最强，而两者存在差异时抑制作用减弱或消失，这种相关特征包括方位、空间频率和对比度等。因此为弥补单一方向差加权抑制的不足，本发明采用多维特征分析，在特征空间更全面准确的计算nCRF与CRF差异。为有效表针微光图像特征，同时减少计算冗余，选取各像素5×5邻域内的4个旋转不变的统计特征以及方向作为像素的特征集，其中4个统计特征包括空间频率、对比度、灰度共生统计熵和灰度共生统计相关性。方向表针梯度信息；空间频率和对比度反映局部灰度统计分布；当像元间距设为1或2时，灰度共生统计特性能够描述局部高频信息，而且不同于Gabor滤波器采用高斯函数对局部区域加权，灰度共生统计特性假设邻域窗口内为均匀分布，因而对局部噪声不敏感。因此融合多种特征形成更具表针性的特征集。

定义fx＝[O_c,SF_c,CON_c,cooE_c,cooCOR_c]^T、fy＝[O_nc,SF_nc,CON_nc,cooE_nc,cooCOR_nc]^T分别表示CRF和nCRF内各像元特征向量。这里O_c、O_nc、SF_c、SF_nc、CON_c、CON_nc、cooE_c、cooE_nc、cooCOR_c和cooCOR_nc分别表示CRF和nCRF内各像元归一化的偏好方向、空间频率、对比度及灰度共生统计熵和相关性。

2.2同质性分析

为计算特征空间中nCRF各环境像元与CRF中心的差异，可以采用主成分分析(PCA)算法，即计算中心像元特征集的主成分，并将环境像元特征向量在主成分投影计算两者同质度。线性PCA算法要求数据服从高斯分布，该条件对真实场景很难满足；KPCA算法利用各维特征间的相关性，将线性空间特征信号映射到高维特征空间中进行差异检测，可以挖掘各特征间的非线性信息，更好的区分轮廓和背景。

KPCA

KPCA算法以检测点为中心，计算其CRF和nCRF数据高维特征空间差异。设子空间维度为P，定义P×N_c矩阵表示CRF多维特征数据，样本fx_i为P维特征向量fx＝[O_c,SF_c,CON_c,cooE_c,cooCOR_c]^T，N_c为CRF内像元数目；P×N_nc矩阵表示nCRF多维特征数据，样本fy_j为P维特征向量fy＝[O_nc,SF_nc,CON_nc,cooE_nc,cooCOR_nc]^T，N_nc为nCRF内像元数目。利用非线性映射函数φ将FX、FY映射到高维特征空间后，C_φ定义为CRF高维特征空间协方差：

$C_{\mathrm{\&phi;}}=\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}(\mathrm{\&phi;}\left({\mathrm{fx}}_i\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}}){(\mathrm{\&phi;}\left({\mathrm{fx}}_i\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})}^T---\left(10\right)$

其中V_φ为C_φ非零特征值Λ_φ对应的特征向量。

nCRF各像元与CRF中心同质度δ(r)定义为

${\mathrm{\&delta;}}^{-1}\left(r\right)={(\mathrm{\&phi;}\left(r\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})}^T{V}_{\mathrm{\&phi;}}{V}_{\mathrm{\&phi;}}^T(\mathrm{\&phi;}\left(r\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})---\left(11\right)$

其中r∈FY为nCRF各像元特征向量，即在高维特征空间计算环境像元特征在中心主成分的投影特性。

WKPCA

KPCA能够更有效的提取和利用数据的非线性特征，但KPCA法也存在一定的缺陷，即若实际中心数据为病态分布(噪声干扰)、背景成份中目标点数目或目标成份中背景点数目较多时，高维特征空间协方差矩阵不能完全描述目标或背景数据。由式10可看出：CRF相关矩阵C_φ中各像元的权重相等，但若中心数据为病态分布、中心组成为背景数据目标点参杂或目标数据背景点参杂，C_φ就不能完全描述CRF数据分布，主成分分析不准确。所以高维特征空间检测点邻域的CRF协方差矩阵C_φ中各像素均需引入权重因子，抑制或降低目标数据中的背景点或异常点、背景数据中的目标点或异常点对C_φ的影响，C_φ就能更准确的描述目标、背景数据分布特性。

因此本发明提出了一种WKPCA算法，该方法在高维特征空间以C_φ中各像元特征向量φ(fx_i)与数据中心向量均值μ_φ的特征角匹配(FAM)为准则，对C_φ中各像元特征向量引入对应的权重因子。即若像元特征向量与数据中心向量均值的夹角较小，表示两者特征信息相似，则该像元获得较大权值，反之亦然，如此去除或抑制C_φ中的异常数据。高维特征空间中加权CRF协方差矩阵C_φw表示为：

$C_{\mathrm{\&phi;w}}=\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{w}_{{\mathrm{fx}}_i}(\mathrm{\&phi;}\left({\mathrm{fx}}_i\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}}){(\mathrm{\&phi;}\left({\mathrm{fx}}_i\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})}^T/\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{w}_{{\mathrm{fx}}_i}---\left(12\right)$

式中为各像元的FAM权重因子：

$w_{{\mathrm{fx}}_i}=\cos (\mathrm{\&phi;}\left({\mathrm{fx}}_i\right),{\mathrm{\&mu;}}_{\mathrm{\&phi;}})=\frac{{\mathrm{\&phi;}}^T\left({\mathrm{fx}}_i\right){\mathrm{\&mu;}}_{\mathrm{\&phi;}}}{\left|\right|{\mathrm{\&phi;}}^T\left({\mathrm{fx}}_i\right)\left|\right|\left|\right|{\mathrm{\&mu;}}_{\mathrm{\&phi;}}\left|\right|}=\frac{\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_i,{\mathrm{fx}}_m)}{\sqrt{k({\mathrm{fx}}_i,{\mathrm{fx}}_i)}\sqrt{\frac{1}{N_c^2}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}\underset{n=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,{\mathrm{fx}}_n)}}---\left(13\right)$

i＝1,…,N_c，其中利用了核函数性质：k(fx,fy)＝<φ(fx),φ(fy)>。

WKPCA算法即计算高维特征空间中C_φw的特征向量V_φw及nCRF内各像元特征向量φ(r)在V_φw上的投影<V_φw,φ(r)>。

定义高维特征空间加权数据 ${\{{\mathrm{\&phi;}}_w\left({\mathrm{fx}}_i\right)=\sqrt{w_{{\mathrm{fx}}_i}/w}(\mathrm{\&phi;}\left({\mathrm{fx}}_i\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})\}}_{i=1}^{N_c},$ 由式12可推出：

$C_{\mathrm{\&phi;w}}=\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{\mathrm{\&phi;}}_w\left({\mathrm{fx}}_i\right){\mathrm{\&phi;}}_w{\left({\mathrm{fx}}_i\right)}^T---\left(14\right)$

为C_φw非零特征值构成的对角矩阵，为各特征值对应的特征向量。由式14可知，每个特征向量都在φ_w(FX)的度量空间内，因此可以表示成的线性组合：

$V_{\mathrm{\&phi;w}}^l=\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i^l{\mathrm{\&phi;}}_w\left({\mathrm{fx}}_i\right)={\mathrm{\&phi;}}_w\left(\mathrm{FX}\right){\mathrm{\&alpha;}}^l---\left(15\right)$

其中 ${\mathrm{\&alpha;}}^l={[{\mathrm{\&alpha;}}_1^l,{\mathrm{\&alpha;}}_2^l,...,{\mathrm{\&alpha;}}_{N_c}^l]}^T.$

将14式和15式代入可得：

$N_c{\mathrm{\&lambda;}}_w^l{\mathrm{\&phi;}}_w\left(\mathrm{FX}\right){\mathrm{\&alpha;}}^l={\mathrm{\&phi;}}_w\left(\mathrm{FX}\right){\mathrm{\&phi;}}_w^T\left(\mathrm{FX}\right){\mathrm{\&phi;}}_w\left(\mathrm{FX}\right){\mathrm{\&alpha;}}^l---\left(16\right)$

将16式两边左乘且可得则核矩阵K_w的特征值为协方差矩阵C_φw的N_c倍，对应的特征向量利用核函数性质将K_w表示为：

$\begin{array}{cc}{K}_w=D_w{N}_c{\mathrm{\&Lambda;}}_{\mathrm{\&phi;w}}{D}_w^T={\left[{\mathrm{\&phi;}}_w^T\left({\mathrm{fx}}_i\right){\mathrm{\&phi;}}_w\left({\mathrm{fx}}_j\right)\right]}_{N_c\&times;N_c}& i,j=1,...,N_c\end{array}$

$=\begin{array}{cc}{\left[\frac{\sqrt{w_{{\mathrm{fx}}_i}{w}_{{\mathrm{fx}}_j}}}{w}{(\mathrm{\&phi;}\left({\mathrm{fx}}_i\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})}^T(\mathrm{\&phi;}\left({\mathrm{fx}}_j\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})\right]}_{N_c\&times;N_c}& i,j=1,...,N_c---\left(17\right)\end{array}$

$=\begin{array}{cc}{\left[\frac{\sqrt{w_{{\mathrm{fx}}_i}{w}_{{\mathrm{fx}}_j}}}{w}(k({\mathrm{fx}}_i,{\mathrm{fx}}_j)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_i,{\mathrm{fx}}_m)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,{\mathrm{fx}}_j)+\frac{1}{N_c^2}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}\underset{n=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,{\mathrm{fx}}_n))\right]}_{N_c\&times;N_c}& i,j=1,...,N_c\end{array}$

由15式可得，nCRF像元特征向量φ(r)在高维特征空间加权特征向量V_φw上的投影为：

$<V_{\mathrm{\&phi;w}},\mathrm{\&phi;}\left(r\right)>=D_w^T{\mathrm{\&phi;}}_w^T\left(\mathrm{FX}\right)\mathrm{\&phi;}\left(r\right)=D_w^T{k}_w(Z,r)---\left(18\right)$

式中 $\begin{array}{cc}{k}_w(Z,r)={\mathrm{\&phi;}}_w^T\left(\mathrm{FX}\right)\mathrm{\&phi;}\left(r\right)={\left[\sqrt{w_{{\mathrm{fx}}_i}/w}(k({\mathrm{fx}}_i,r)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,r))\right]}_{N_c\&times;1}& i=1,...,N_c\end{array}$

高维特征空间nCRF各像元与中心CRF同质度δ_w(r)表示为：

${\mathrm{\&delta;}}_w^{-1}\left(r\right)={(\mathrm{\&phi;}\left(r\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})}^T{V}_{\mathrm{\&phi;w}}{V}_{\mathrm{\&phi;w}}^T(\mathrm{\&phi;}\left(r\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})$

$={(\mathrm{\&phi;}\left(r\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})}^T{\mathrm{\&phi;}}_w\left(\mathrm{FX}\right)D_w{D}_w^T{\mathrm{\&phi;}}_w^T\left(\mathrm{FX}\right)(\mathrm{\&phi;}\left(r\right)-{\mathrm{\&mu;}}_{\mathrm{\&phi;}})---\left(19\right)$

$={(k_w(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{k}_w(Z,{\mathrm{fx}}_i))}^T{D}_w{D}_w^T(k_w(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{k}_w(Z,{\mathrm{fx}}_i))$

因此选择合适的核函数k构造正定的核矩阵K_w，计算其特征向量D_w，即可得出各环境像元与中心的同质度δ_w(x,y)。通过多次实验选用径向基核函数它能够表征特征向量间能量差异；同时采用FAM函数能够表征特征向量曲线的形状差异，如此WKPCA算法既考虑了像元间特征能量差异，也考虑了像元间特征曲线形状差异，更全面的区分多维特征数据。

2.3同质度校正抑制

仍然考虑抑制量随距离衰减特性，结合同质度δ_w(x,y)和距离衰减函数ω_σ(x,y)对nCRF中各像元抑制量校正，将nCRF模型中环境抑制T_σ(x,y)修订如下：

T_σ(x,y)=(δ_w(x,y)×ω_σ(x,y))*E_σ(x,y)

(20)

由于WKPCA更准确计算CRF中心主成分，目标与背景能够有效区分。对于检测点为背景内部像元情况，主要受与其同质区域的抑制，nCRF中有效抑制像元的数目和同质度值大，由于自然纹理Gabor方向能量极值大且分布杂乱，即使受距离衰减函数调制，CRF中心检测点仍获得较大抑制量而刺激响应被严重削弱；同样目标内部检测点主要受与其同质区域的抑制，有效抑制像元的数目取决于nCRF中同质区域大小(即目标大小)，目标区域同质度值大，其Gabor方向能量极值主要分布于目标内部和边缘，若目标较大则nCRF中同质区大，CRF中心目标像元获得较大抑制量而响应弱，若目标较小则nCRF中同质区小，CRF中心目标像元不受抑制，保留整体小目标；而对于显著目标轮廓检测点，由于WKPCA去除了异常干扰数据，CRF主成分与目标或背景特征吻合程度高，nCRF中目标内部的虚假边缘或背景纹理同质度值高，由于轮廓像元的多维特征区别于目标和背景、成像空间分辨率约束导致目标边缘产生混淆像元，其光谱特征与中心主成分产生偏差，轮廓像元本身同质度值较低，而其Gabor方向能量极值主要集中于边缘，因此中心轮廓像元主要受源自同一轮廓的像元抑制，有效抑制像元数目小且同质度值较小，加之距离衰减调制，轮廓像元几乎不受抑制而响应强。终上所述，WKPCA同质度校正抑制能够有效抑制环境纹理，获得较大轮廓响应，且最大程度减少轮廓自抑制，保留完整轮廓。

下面结合实施例对本发明做进一步说明。

参照流程图2，本发明的具体实现步骤如下：

步骤1，输入微光图像。

步骤2，设置非经典感受野模型参数：nCRF和CRF区域半径比ρ，高斯标准差σ，核函数尺度σ_k，调制系数a∈[0,1]、b∈[0,1]。根据式1及参数计算CRF和nCRF区域圆形半径r_c和ρr_c。

$r_c=\sqrt{2\ln \mathrm{\&rho;}/{\mathrm{\&rho;}}^2-1}\&CenterDot;\mathrm{\&rho;\&sigma;}---\left(21\right)$

步骤3，逐行扫描各个像素点，计算各像素点Gabor方向能量极值E_σ(x,y)；由式2计算nCRF区域内的距离加权函数ω_σ(x,y)。

ω_σ(x,y)=|DOG_σ,ρσ(x,y)|⁺/‖|DOG_σ,ρσ(x,y)|⁺‖₁(22)

${\mathrm{DOG}}_{\mathrm{\&sigma;},\mathrm{\&rho;\&sigma;}}(x,y)=\frac{1}{\sqrt{2\mathrm{\&pi;}{\left(\mathrm{\&rho;\&sigma;}\right)}^2}}\exp (-\frac{x^2+y^2}{2{\left(\mathrm{\&rho;\&sigma;}\right)}^2})-\frac{1}{\sqrt{2\mathrm{\&pi;}{\mathrm{\&sigma;}}^2}}\exp (-\frac{x^2+y^2}{2{\mathrm{\&sigma;}}^2})---\left(23\right)$

步骤4，结合微光图像特性，计算图像各像素点特征向量，包括方向、空间频率、对比度、灰度共生统计熵和灰度共生统计相关性。

构建各像素点CRF区域内多维特征数据集其中fx＝[O_c,SF_c,CON_c,cooE_c,cooCOR_c]^T，N_c为CRF内像元数目，O_c、SF_c、CON_c、cooE_c和cooCOR_c为CRF区域内归一化的方向、空间频率、对比度、灰度共生统计熵和灰度共生统计相关性；nCRF区域内多维特征数据集其中fy＝[O_nc,SF_nc,CON_nc,cooE_nc,cooCOR_nc]^T，N_nc为nCRF内像元数目，O_nc、SF_nc、CON_nc、cooE_nc和cooCOR_nc为nCRF区域内归一化的方向、空间频率、对比度、灰度共生统计熵和灰度共生统计相关性。

步骤5，基于WKPCA算法，采用径向基核函数，由式4计算CRF各像元的特征向量角匹配FAM值i＝1,…,N_c；根据FAM值对高维特征空间CRF各像元特征向量加权，计算加权后CRF多维特征数据集主成分，以及nCRF中各像元特征向量在CRF主成分上的投影，即由式6计算图像各像素点的nCRF区域内各像元与CRF中心同质度δ_w(r)。其中r∈FY，D_w为核矩阵K_w的特征向量。

$w_{{\mathrm{fx}}_i}=\frac{\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_i,{\mathrm{fx}}_m)}{\sqrt{k({\mathrm{fx}}_i,{\mathrm{fx}}_i)}\sqrt{\frac{1}{N_c^2}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}\underset{n=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,{\mathrm{fx}}_n)}}---\left(24\right)$

$k(\mathrm{fx},\mathrm{fy})=\exp [(-{\left|\right|\mathrm{fx}-\mathrm{fy}\left|\right|}^2)/{2\mathrm{\&sigma;}}_k^2]---\left(25\right)$

${\mathrm{\&delta;}}_w^{-1}\left(r\right)={(k_w(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{k}_w(Z,{\mathrm{fx}}_i))}^T{D}_w{D}_w^T(k_w(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{k}_w(Z,{\mathrm{fx}}_i))---\left(26\right)$

$K_w=\begin{array}{cc}{\left[\frac{\sqrt{w_{{\mathrm{fx}}_i}{w}_{{\mathrm{fx}}_j}}}{w}(k({\mathrm{fx}}_i,{\mathrm{fx}}_j)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_i,{\mathrm{fx}}_m)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,{\mathrm{fx}}_j)+\frac{1}{N_c^2}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}\underset{n=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,{\mathrm{fx}}_n))\right]}_{N_c\&times;N_c}& i,j=1,...,N_c---\left(27\right)\end{array}$

$\begin{array}{cc}{k}_w(Z,r)={\mathrm{\&phi;}}_w^T\left(\mathrm{FX}\right)\mathrm{\&phi;}\left(r\right)={\left[\sqrt{w_{{\mathrm{fx}}_i}/w}(k({\mathrm{fx}}_i,r)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,r))\right]}_{N_c\&times;1}& i=1,...,N_c---\left(28\right)\end{array}$

步骤6，结合同质度δ_w(x,y)和距离衰减函数ω_σ(x,y)对nCRF中各像元抑制量校正，由式9和式10计算图像各像素点的环境抑制量T_σ(x,y)及抑制输出R_σ(x,y)。其中δ_w(x,y)即为δ_w(r)。

T_σ(x,y)=(δ_w(x,y)×ω_σ(x,y))*E_σ(x,y)(29)

R_σ(x,y)=|aE_σ(x,y)-bT_σ(x,y)|⁺(30)

步骤7，对输出R_σ(x,y)进行非极大值抑制，设置灰度阈值和轮廓长度阈值，对非极大值抑制结果进行二值化，得到最终的轮廓输出。

本发明的效果可以通过以下仿真结果进一步说明：

1、由图4可见，结构图中CRF窗口中主要为目标像元，其均值特征向量近似目标特性。因此FAM三维图中CRF窗口和nCRF环形窗口内目标像元的FAM值较大(白色)，背景像元的FAM值较小(黑色)。因此FAM权重因子能够有效削弱CRF协方差矩阵中的异常数据。实验中设置σ＝2.0，CRF窗口半径5，nCRF窗口半径20。

2、由图5可见，对模拟数据仿真，测试比较了PCA、KPCA和WKPCA算法性能。相比于线性PCA，非线性KPCA一定程度的区分目标和背景，但是无法进一步对不同分布目标分类；而WKPCA不仅能够分离目标和背景，而且能够有效区分不同分布的目标。实验中设置σ_k＝0.5。

3、如图6所示，WKPCA算法性能优势进一步表现在真实微光图像数据分析中。CRF中心以目标为主，KPCA比PCA主成分分析准确(a)。但是受异常数据(噪声、目标中的背景数据等)干扰，KPCA主成分产生偏差，导致目标、背景内部像元的同质度均较低；而目标轮廓像元多维统计特征包含目标和背景信息，其特征向量易与干扰主成分近似，导致目标轮廓的同质度较高(b)。WKPCA由于去除了异常数据的干扰，中心主成分计算准确，因而目标内部同质度较高；目标轮廓的灰度统计特征包含目标与背景信息，特征向量与中心主成分不一致而同质度稍低；背景同质度较低(c)。实验中设置σ＝2.0,σ_k＝0.5。

采用式31计算PCA同质度，式中μ为CRF中心像元特征均值向量，v为采用PCA算法选定的主成分；采用式32计算KPCA同质度，w为KPCA算法选定的主成分；采用式19计算WKPCA同质度。

δ^-1(r)＝(r-μ)^Tvv^T(r-μ)(31)

${\mathrm{\&delta;}}^{-1}\left(r\right)={(k(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k(Z,{\mathrm{fx}}_i))}^T{\mathrm{ww}}^T(k(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k(Z,{\mathrm{fx}}_i))$

(32)

其中 $k(Z,r)={\mathrm{\&phi;}}^T\left(\mathrm{FX}\right)\mathrm{\&phi;}\left(r\right)=\begin{array}{cc}{[k({\mathrm{fx}}_i,r)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\mathrm{\&Sigma;}}}k({\mathrm{fx}}_m,r)]}_{N_c\&times;1}& i=1,...,N_c.\end{array}$

4、如图7所示，针对自然场景中微光图像验证WKPCA同质度校正抑制的合理有效性。对于背景纹理(列1-2)nCRF与CRF同质度高(行2)，Gabor方向能量极值杂乱(行3)，导致nCRF大部分像元对中心抑制，有效抑制像元数目大、各像元抑制值大(行4)，因此背景中心响应弱；对于目标内部(列3-4)nCRF中同一目标部分与CRF同质度高(行2)，Gabor方向能量极值集中于边缘与目标内部(行3)，导致nCRF中目标像元均以较大值对中心抑制，有效抑制像元数目取决于目标大小、抑制值大(行4)，因此目标中心响应弱；轮廓部分(列5-6)受WKPCA作用，nCRF中与CRF主成分一致的目标部分同质度高(行2)，轮廓部分同质度稍低，Gabor方向能量极值主要集中于边缘(行3)，目标内部的虚假边缘(列5行3)、背景纹理(列6行3)的Gabor能量较低，因此轮廓检测像元主要受源自同一轮廓的像元抑制，有效抑制像元数目小，且同质度值较小(行4)，加之距离校正，各像元抑制值小，因此轮廓中心响应强。实验中设置σ＝2.0,σ_k＝0.5。

5、如图8所示，利用本申请方法提取微光图像中显著轮廓。预处理中采用低通滤波器降噪；后处理中进行非极大值抑制和二值化排除虚假边缘，以选择较高调制响应及较大长度的边缘作为轮廓输出。即使经过低通滤波仍存在噪声干扰，蝶形抑制残留了大量纹理和虚假边缘；线性PCA同质度校正抑制由于考虑了多特征，降低噪声干扰；KPCA同质度校正抑制采用非线性映射，提高了中心主成分准确度，更有效的区分目标和背景，进一步削弱了部分纹理和虚假边缘响应；WKPCA同质度校正抑制去除中心样本中的异常数据，最大程度的抑制同质区域响应，保留轮廓信息。实验中设置σ＝2.0,a＝1.0,b＝0.8,σ_k＝0.5。

Claims (3)

1.一种WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法，其特征在于步骤如下：

步骤1，输入微光图像；

步骤2，设置非经典感受野模型参数：nCRF和CRF区域半径比ρ，高斯标准差σ，核函数尺度σ_k，调制系数a∈[0,1]、b∈[0,1]，根据式1及参数计算CRF和nCRF区域圆形半径r_c和ρr_c；

$r_c=\sqrt{2ln\mathrm{\&rho;}/{\mathrm{\&rho;}}^2-1}\&CenterDot;\mathrm{\&rho;}\mathrm{\&sigma;}---\left(1\right)$

步骤3，逐行扫描各个像素点，计算各像素点Gabor方向能量极值E_σ(x,y)；由式2计算nCRF区域内的距离加权函数(x,y)为像素在图像中的位置；

${\mathrm{\&omega;}}_{{\mathrm{\&sigma;}}_1}(x,y)=|{\mathrm{DOG}}_{\mathrm{\&sigma;},\mathrm{\&rho;}\mathrm{\&sigma;}}(x,y){|}^{+}/|\left|\right|{\mathrm{DOG}}_{\mathrm{\&sigma;},\mathrm{\&rho;}\mathrm{\&sigma;}}(x,y){|}^{+}|{|}_1---\left(2\right)$

${\mathrm{DOG}}_{\mathrm{\&sigma;},\mathrm{\&rho;}\mathrm{\&sigma;}}\left(x,y\right)=\frac{1}{\sqrt{2\mathrm{\&pi;}{\left(\mathrm{\&rho;}\mathrm{\&sigma;}\right)}^2}}\exp \left(-\frac{x^2+y^2}{2{\left(\mathrm{\&rho;}\mathrm{\&sigma;}\right)}^2}\right)-\frac{1}{\sqrt{2{\mathrm{\&pi;\&sigma;}}^2}}\exp \left(-\frac{x^2+y^2}{2{\mathrm{\&sigma;}}^2}\right)---\left(3\right)$

步骤4，计算微光图像各像素点特征向量，构建各像素点CRF区域内多维特征数据集nCRF区域内多维特征数据集CRF区域内特征向量fx_i元素为Gabor方向、空间频率、对比度、灰度共生统计熵和灰度共生统计相关性，fx_i＝[O_c,SF_c,CON_c,cooE_c,cooCOR_c]^T，i＝1,...,N_c，N_c为CRF内像元数目，O_c、SF_c、CON_c、cooE_c和cooCOR_c为CRF区域内归一化的Gabor方向、空间频率、对比度、灰度共生统计熵和灰度共生统计相关性；nCRF区域内特征向量fy_j，j＝1,...,N_nc，fy_j＝[O_nc,SF_nc,CON_nc,cooE_nc,cooCOR_nc]^T，N_nc为nCRF内像元数目，O_nc、SF_nc、CON_nc、cooE_nc和cooCOR_nc为nCRF区域内归一化的Gabor方向、空间频率、对比度、灰度共生统计熵和灰度共生统计相关性；

步骤5，采用WKPCA算法，计算图像各像素点的nCRF区域内各像元与CRF中心同质度δ_w(r)；

步骤6，结合同质度δ_w(x,y)和距离衰减函数对nCRF中各像元抑制量校正，由式9和式10计算图像各像素点的环境抑制量T_σ(x,y)及抑制输出R_σ(x,y)；其中δ_w(x,y)即为δ_w(r)；

$T_{\mathrm{\&sigma;}}(x,y)=\left({\mathrm{\&delta;}}_w\right(x,y)\&times;{\mathrm{\&omega;}}_{{\mathrm{\&sigma;}}_2}(x,y\left)\right)*E_{\mathrm{\&sigma;}}(x,y)---\left(9\right)$

R_σ(x,y)＝|aE_σ(x,y)-bT_σ(x,y)|⁺(10)

步骤7，对输出R_σ(x,y)进行非极大值抑制，设置灰度阈值和轮廓长度阈值，对非极大值抑制结果进行二值化，得到最终的轮廓输出。

2.根据权利要求1所述的WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法，其特征在于：所述步骤5中的WKPCA算法基于KPCA算法，在高维特征空间CRF协方差矩阵中引入特征向量角匹配FAM值对CRF各像元特征向量加权，FAM值为高维特征空间中CRF各像元特征向量与CRF特征数据集向量均值的向量夹角。

3.根据权利要求1或2所述的WKPCA同质度校正nCRF抑制的微光图像显著轮廓提取方法，其特征在于：所述步骤5中的WKPCA算法计算CRF特征数据集主成分的步骤如下：

首先由式4计算CRF各像元的特征向量角匹配FAM值i＝1,...,N_c，并根据FAM值对高维特征空间CRF各像元特征向量加权；再计算加权后CRF多维特征数据集主成分，以及nCRF中各像元特征向量在CRF主成分上的投影即δ_w(r)，式6中r∈FY，D_w为核矩阵K_w的特征向量，k(·)为径向基核函数；

$w_{{\mathrm{fx}}_i}=\frac{\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\&Sigma;}}k\left({\mathrm{fx}}_i,{\mathrm{fx}}_m\right)}{\sqrt{k\left({\mathrm{fx}}_i,{\mathrm{fx}}_i\right)}\sqrt{\frac{1}{N_c^2}\underset{m=1}{\overset{N_c}{\&Sigma;}}\underset{n=1}{\overset{N_c}{\&Sigma;}}k\left({\mathrm{fx}}_m,{\mathrm{fx}}_n\right)}}---\left(4\right)$

$k(fx,fy)=\exp \&lsqb;(-||fx-fy|{|}^2)/2{\mathrm{\&sigma;}}_k^2\&rsqb;---\left(5\right)$

${\mathrm{\&delta;}}_w^{-1}\left(r\right)={\left(k_w\right(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\&Sigma;}}{k}_w(Z,{\mathrm{fx}}_i\left)\right)}^T{D}_w{D}_w^T\left(k_w\right(Z,r)-\frac{1}{N_c}\underset{i=1}{\overset{N_c}{\&Sigma;}}{k}_w(Z,{\mathrm{fx}}_i\left)\right)---\left(6\right)$

$K_w={\&lsqb;\frac{\sqrt{w_{{\mathrm{fx}}_i}{w}_{{\mathrm{fx}}_j}}}{w}\left(k\left({\mathrm{fx}}_i,{\mathrm{fx}}_j\right)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\&Sigma;}}k\left({\mathrm{fx}}_i,{\mathrm{fx}}_m\right)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\&Sigma;}}k\left({\mathrm{fx}}_m,{\mathrm{fx}}_j\right)+\frac{1}{N_c^2}\underset{m=1}{\overset{N_c}{\&Sigma;}}\underset{n=1}{\overset{N_c}{\&Sigma;}}k\left({\mathrm{fx}}_m,{\mathrm{fx}}_n\right)\right)\&rsqb;}_{\begin{array}{cc}{N}_c\&times;N_c& i,j=1,...,N_c\end{array}}---\left(7\right)$

$k_w\left(Z,r\right)={\mathrm{\&phi;}}_w^T\left(FX\right)\mathrm{\&phi;}\left(r\right)={\&lsqb;\sqrt{w_{{\mathrm{fx}}_i}/w}\left(k\left({\mathrm{fx}}_i,r\right)-\frac{1}{N_c}\underset{m=1}{\overset{N_c}{\&Sigma;}}k\left({\mathrm{fx}}_m,r\right)\right)\&rsqb;}_{\begin{array}{cc}{N}_c\&times;1& i=1,...,N_c\end{array}}---\left(8\right).$