CN106355600A - 基于非线性结构张量的自适应类圆型结构元素构造形态学算子的方法 - Google Patents

Abstract

这里涉及的是基于非线性结构张量的自适应类圆型结构元素构造形态学算子的方法。针对线性结构张量容易估计出不准确，甚至错误的梯度信息的缺点，首先采用非线性结构张量代替线性结构张量，计算非线性结构张量的特征值和特征向量，由此定义边缘强度测量和角点强度测量两个测量参数，明确区分边缘区域、光滑区域和拐角等图像特征，从而定义自适应椭圆型结构元素的长、短半轴及方向，然后为了使得到的自适应结构元素更加接近理想椭圆形状，对椭圆结构元素和原图像进行了插值处理，最后采用提出的自适应类圆形结构元素重新定义了自适应腐蚀、膨胀、开运算、闭运算和击中击不中变换等形态学算子。验证实验结果表明，本方法与其他方法在结构自适应、角点保护、滤波和目标提取等方面相比，都具有很好的效果。

Description

基于非线性结构张量的自适应类圆型结构元素构造形态学算 子的方法

技术领域

这里涉及的是基于非线性结构张量的自适应类圆型结构元素构造形态学算子的方法研究。首先计算非线性结构张量的特征值和特征向量，并用其定义角点强度测量和边缘强度测量，从而准确估计图像的局部各向异性特征；然后利用角点和边缘强度测量定义椭圆结构元素的长轴、短轴和方向的计算公式；同时对得到的椭圆结构元素及原始图像进行插值处理；最后利用自适应类圆型结构元素构造自适应腐蚀、膨胀、开、闭及击中击不中变换等形态学算子。

背景技术

数学形态学起初被应用于图像处理中时，所用到的结构元素通常是根据一些先验知识人为选择的，因此其在整个图像空间中固定不变，这样非常容易导致不理想的处理结果。所以为了更好地将数学形态学引入日益广泛的图像处理应用中，针对不同的应用领域、不同的图像局部特征来设计自适应的结构元素显得尤为必要。

目前存在多种构造自适应结构元素的方法，可以总结如下三类：完全基于图像中连通性区域的方法，基于灰度和空间信息的组合加权的方法，和基于结构的方法。第一类方法中，如一般自适应邻域，其结构元素没有尺寸和形状限制，唯一的参数m∈Z⁺仅控制结构元素的各向同性，所以该方法完全不关注空间关系。基于灰度和空间信息的组合加权的方法为空间约束提供了方法，典型的如形态学变形虫，它应用加权灰度距离去限制变形虫的形状，但是对于形态学变形虫而言，设置一个合适的权重参数λ是一项具有挑战性的工作，而且这种方法对图像的对比度变化十分敏感。与形态学变形虫方法类似，基于显著性距离变换的显著性自适应结构元素会跨越边缘，但是这种方法不能很好的处理图像中的角点。基于结构的方法，如基于线性结构张量的椭圆型自适应结构元素，依据各向异性比率其形状从线型到盘型不断变化。然而线性结构张量采用高斯平滑技术，其容易模糊图像的细节特征，尤其是接近于图像中两个不同的特征区域的边缘处。

发明内容

由上述不同的方法可以看出，椭圆型自适应结构元素由于充分利用了图像的各向异性信息，很值得被深入研究。因此这里提出了一种新的构造自适应结构元素的方法，采用非线性结构张量代替线性结构张量，并由其特征值定义了角点测量参数和边缘测量参数，以获得更加准确的图像方向及各向异性信息；另外，为了使计算得到的结构元素更接近于理想椭圆形状，这里对椭圆结构元素和原始图像进行了插值处理。最终利用构造的自适应类圆形结构元素重新定义了基本形态学算子。

1、非线性结构张量

定义图像表达式为I(x，y)，则非线性结构张量即采用非线性各向异性扩散过程代替线性结构张量中的高斯平滑操作，如式(1)。

$\frac{\∂u_{i,j}}{\∂t}=div\left(g\right(\underset{k,l=1}{\overset{2}{\Σ}}\▿u_{k,l}\▿u_{k,l}^T)\▿u_{k,l}),(i,j=1,2)---\left(1\right)$

其中，u_i，j为初始结构张量的任一分量。t的意义在离散时域里就是迭代滤波次数，迭代的次数t越大，利用的周围邻域的像素就越多；g(·)是一种随像素张量越大而平滑力递减的平滑滤波函数，由于在式(1)中，是一个矩阵，所以平滑函数g(·)作用后可改写为如式(2)的形式：

$g(\underset{k,l=1}{\overset{2}{\Σ}}\▿u_{kl}\▿u_{kl}^T)=Vdiag\left(g_1\right({\mathrm{\λ}}_1),g_2({\mathrm{\λ}}_2\left)\right)V^T=g_1\left({\mathrm{\λ}}_1\right){\mathrm{\θ}}_1{\mathrm{\θ}}_1^T+g_2\left({\mathrm{\λ}}_2\right){\mathrm{\θ}}_2{\mathrm{\θ}}_2^T---\left(2\right)$

其中，λ₁≥λ₂为ST的特征值，θ₁和θ₂为一对正交基，与其对应的特征向量平行。θ₁指向局部数据最大的跳变方向，θ₂则指向边缘方向。为了沿着图像边缘平滑矩阵图像，使在扩散过程中沿着边缘的扩散强度大于垂直边缘的扩散强度，要求扩散系数g₁(λ₁)、g₂(λ₂)与相应的特征值成反比，且满足g₁(λ₁)＞g₂(λ₂)，如此实现了各通道的各向异性扩散，这里使用如式(3)、(4)所示的形式：

g₁(λ₁)＝1/(1+λ₁+λ₂)^p (3)

g₂(λ₂)＝1/(1+λ₁+λ₂)^q (4)

其中，p、q为控制参数，为了沿着跳变最小的方向平滑张量场数据，要求p＞q。这里设定p＝1.5，q＝0.5。

公式(1)的迭代过程可由(5)表示：

$u_{i,j}^{n+1}=u_{i,j}^n+\mathrm{\Δ}t\·\left(div\right(g\left(\underset{k,l=1}{\overset{2}{\Σ}}\▿u_{k,l}\▿u_{k,l}^T\right)\·\▿u_{i,j}^n\left)\right)---\left(5\right)$

迭代终止条件如式(6)：

$\underset{i,j=1}{\overset{2}{\&Sigma;}}|u_{i,j}^{n+1}-u_{i,j}^n|<\mathrm{\&epsiv;}---\left(6\right)$

其中是迭代过程中非结构张量的任一分量的第n层数据，为第n+1层数据，Δt为步长，ε为预先设置的误差限，g(·)根据及时更新。

2、测量参数

由式(5)、(6)可得到迭代终止后的非线性结构张量，其特征值η₁、η₂和特征向量v₁、v₂的计算公式如式(7)、(8)所示：

${\mathrm{\&eta;}}_{1,2}=\frac{1}{2}(u_{11}+u_{22}\&PlusMinus;\sqrt{{(u_{11}-u_{22})}^2+4u_{12}^2})---\left(7\right)$

$\left\{\begin{array}{c}{v}_1={(2u_{12},u_{22}-u_{11}+\sqrt{{(u_{22}-u_{11})}^2+4u_{12}^2})}^T\\ v_2=v_1^{\&perp;}\end{array}\right.---\left(8\right)$

分析非线性结构张量的特征值η₁、η₂可以得到一些图像的局部结构信息。为了用量化指标测量感兴趣区域，因此定义了两种测量参数：

1)拐角强度测量m_c：

$m_c=\frac{4{\mathrm{\&eta;}}_1{\mathrm{\&eta;}}_2}{{\mathrm{\&eta;}}_1+{\mathrm{\&eta;}}_2}---\left(9\right)$

2)边缘强度测量m_e：

m_e＝η₁-η₂ (10)

对图像中的任一像素点，当m_e≈0、m_c≈0时，表明对图像中的任一像素点，在该点沿任何方向的灰度变化都很小，对应的是图像的光滑区域；当m_e＞＞0、m_c≈0时，表明图像沿特征向量v₁方向的灰度变化率，远大于垂直于此方向的灰度变化率，对应的是图像的边缘区域；当m_e≈0、m_c＞＞0时，表明图像在该像素点沿任何方向的灰度变化都很大，对应的是图像的角点。

3、基于非线性结构张量的自适应椭圆型结构元素

自适应椭圆型结构元素表示为其中，a(x，y)为椭圆型结构元素的长半轴，b(x，y)为其短半轴，表示椭圆型结构元素方向。

在具体设计自适应椭圆型结构元素时应充分考虑其形状、尺寸和方向随着图像的局部各向异性特征而不断变化：

1)在拐角处，长半轴a(x，y)和短半轴b(x，y)应保持很小；

2)在边缘处，a(x，y)应满足很大，而b(x，y)很小或趋近于0，结构元素形状接近线型，方向与梯度方向垂直；

3)在平滑区域，a(x，y)和b(x，y)都很大。结构元素形状趋于圆盘。

因此，这里重新定义椭圆长轴、短轴计算公式如式(11)、(12)所示。

$a(x,y)=r\&CenterDot;\{1-\exp (-\frac{C_m}{{\left(m_c\right(x,y)/{\mathrm{\&beta;}}_1)}^m})\}---\left(11\right)$

$b(x,y)=a(x,y)\&CenterDot;\{1-\exp (-\frac{C_m}{{\left(m_e\right(x,y)/{\mathrm{\&beta;}}_2)}^m})\}---\left(12\right)$

其中，r表示允许的最大半径；β₁和β₂为归一化参数，这里分别取角点测量参数m_c和边缘测量参数m_e的最大值的75％。m、Gm为常量，其关系可由公式(13)计算得出。这里取m＝1.5、G_m＝0.7627。

1-exp(-C_m)(1+mC_m)＝0 (13)

椭圆结构元素的方向定义如式(14)所示。

其中，v_2，x(x，y)、v_2，y(x，y)表示特征向量v₂的分量。

4、提高分辨率

给定输入图像I及参数r，可以得到每个像素点的椭圆参数值a、b和且都是连续数据，而在构造椭圆结构元素时，需要对参数数值进行离散化，直接四舍五入化数据会造成很大误差，且得到的椭圆结构元素边界不够平滑，所以这里提出了提高离散椭圆结构元素分辨率的方法，这也意味着不得不提高输入图像的分辨率。

1)图像分辨率

提高图像分辨率简单而有效的方法是进行内插，但通常的内插方法(如临近插值、双线性插值)容易丢失图像的纹理特征信息。为克服这一缺点，提出了随机分形插值方法。随机分形插值的物理依据是“自相似性”，这是分形的基本特征，它的典型数学模型是分数Brown运动，能充分反应图像的统计纹理特征。所以这里采用随机分形插值的方法来提高图像的分辨率。

随机分形插值的具体过程如下：

①定义原始图像为I(x，y)，内插后的图像为I_H(i，j)，根据分形插值模型，利用式(15)、(16)分别计算像素灰度正态分布的标准差δ和分形参数H：

δ＝E|I_H(t+1)-I_H(t)| (15)

log E|I_H(t+Δt)-I_H(t)|²-2H log||Δt||＝2logδ (16)

②当i，j均为奇数时：

I_H(i，j)＝I(x，y) (17)

③当i，j均为偶数时：

$\begin{array}{c}{I}_H(i,j)=\frac{1}{4}\{I_H(i-1,j-1)+I_H(i+1,j-1)+I_H(i-1,j+1)\\ +I_H(i+1,j+1)\}+\sqrt{1-2^{2H-2}}|\left|\mathrm{\&Delta;}t\right||\&CenterDot;H\&CenterDot;\mathrm{\&delta;}\&CenterDot;Gauss\end{array}---\left(18\right)$

④当i，j为一奇一偶时：

$\begin{array}{c}{I}_H(i,j)=\frac{1}{4}\{I_H(i,j-1)+I_H(i-1,j)+I_H(i+1,j)+I_H(i,j+1)\}\\ +2^{-H/2}\sqrt{1-2^{2H-2}}\left|\right|\mathrm{\&Delta;}t\left|\right|\&CenterDot;H\&CenterDot;\mathrm{\&delta;}\&CenterDot;Gauss\end{array}---\left(19\right)$

其中，E表示数学期望，Gauss为高斯随机变量，服从N(0，1)分布，||Δt||代表两插值点间的距离。

上述步骤可以重复迭代直到所需的空间分辨率为止。由此推导出当迭代k次后图像I_H中的像素点{2^kx-(2^k-1)，2^ky-(2^k-1)}处的灰度值即为原始图像I对应像素点(x，y)的灰度值。

2)椭圆结构元素分辨率

经过公式(11)、(12)、(14)计算得到椭圆参数后，需要构造离散的椭圆矩阵，通常情况下，离散的椭圆是由一系列整数型坐标点(x，y)组成的集合X，且所有属于集合X的点(x，y)满足式(20)所示的椭圆方程的限制，椭圆的原点为(0，0)。

但为了提高离散椭圆的分辨率，这里取(x，y)＝(t/2^k，t/2^k)，若满足公式(20)，则点(x，y)属于集合X，需要强调的是，这里点(x，y)并非整数型数据点。其中-r×2^k≤t≤r×2^k为整数型数据，r为椭圆结构元素允许的最大半径。由此，将得到的集合X转换为矩阵形式，即得到插值后的离散椭圆结构元素。

5、自适应腐蚀、膨胀、开运算和闭运算

定义输入图像I(x，y)，其定义域为D。插值后图像为I_H(i，j)，其定义域为D′。I(x，y)中任一像素点L(x，y)∈D处的自适应结构元素为ESE{2^kx-(2^k-1)，2^ky-(2^k-1)}为SE(x，y)插值后的自适应结构元素，像素点L′{2^kx-(2^k-1)，2^ky-(2^k-1)}∈D′。I_H(i，j)中除L′外其他像素点处没有定义结构元素，所以在进行形态学操作时，只对I_H(i，j)中像素点L′进行处理，由此定义自适应腐蚀、膨胀开、闭运算如公式(21)-(24)所示：

IO(x，y)＝γ(I_H)(i，j)＝(δ(ε))(I_H)(i，j) (23)

IC(x，y)＝Ψ(I_H)(i，j)＝(ε(δ))(I_H)(i，j) (24)

其中，i＝2^kx-(2^k-1)，j＝2^ky-(2^k-1)，(x，y)∈D，(i，j)∈D′，ε、δ分别表示腐蚀、膨胀运算，γ、Ψ分别表示开、闭运算，ESE^c(i，j)为反射邻域。

6、自适应击中击不中变换(Hit-or-Miss Transform)

与最基本的形态学变换(腐蚀、膨胀、开、闭)相比，击中击不中变换涉及一个结构元素对(B_FG，B_BG)：B_FG与研究对象匹配，B_BG与背景匹配，其结果是能够提取图像中所有满足与B_FG相同形状的目标像素。然而，传统击中击不中变换使用固定不变的结构元素，只能提取出相同尺寸、相同形状的目标。为此提出使用自适应的结构元素对，使能够提取出具有不同尺寸的类圆形目标。

1)基于倒角距离变换(Chamfer Distance Transform)的自适应结构元素对为提取出类圆形目标，定义自适应结构元素的形状为圆形，其尺寸根据CDT的结果进行自适应调整。其中，考虑到CDT的偏差和计算效率，CDT的模版大小选择的5×5。

具体计算自适应结构元素的尺寸算法如下：

①利用公式(10)中的非线性结构张量的边缘强度测量m_e来判断输入图像的边缘像素点，如式(25)所示：

$\left\{\begin{array}{cc}f(x,y)=m_e(x,y),& \begin{array}{cc}if& m_e(x,y)>delta\end{array}\\ f(x,y)=0,& other\end{array}\right.---\left(25\right)$

②对①中得到的f(x，y)进行非极大值抑制，令属于局部极大值的像素点赋值为1，非局部极大值点赋值为0，由此可以得到单像素宽的边缘图像f_edge。

③对f_edge进行CDT，得到的图像每一点的灰度值为该点到达最近边缘点的距离。将此距离作为结构元素B_FG的半径，记为R_HMT。

根据得到的半径R_HMT即可定义自适应结构元素B_FG如式(26)所示：

B_FG(x，y)＝{|(x，y)-(p，q)|＜R_HMT(x，y)}，(p，q)∈D (26)

其中，(x，y)、(p，q)分别为任意两个像素点；

因为在定义自适应结构元素对(B_FG，B_BG)时，需要满足B_FG和B_BG具有共同的原点，并且互不连通，所以定义自适应结构元素B_BG如式(27)所示：

B_BG(x，y)＝{R_HMT(x，y)+u＜|(x，y)-(p，q)|＜R_HMT(x，y)+u+b}，(p，q)∈D (27)

其中，u为不确定区域半径，b为B_BG半径，两者在整幅图像中保持固定不变。

附图说明

图1椭圆型自适应结构元素方法与本方法的腐蚀结果对比(a)原图像(b-g)椭圆型自适应结构元素方法结果，r_w＝4，16(top to right)and M＝4，8and 16(left to right)(h-j)本方法结果，其中，r分别为＝4，8and 16。

图2椭圆型自适应结构元素方法、显著性自适应结构元素方法和本方法在角点保护方面的结果对比(a)原图像(b)显著性自适应结构元素方法的膨胀结果(c)椭圆型自适应结构元素方法的膨胀结果(d)本方法腐蚀结果(e)本方法膨胀结果(f)本方法开运算结果

图3椭圆型自适应结构元素方法与本方法在滤波方面的结果对比(a)原图像(b)椭圆型自适应结构元素方法的腐蚀结果(c)本方法的腐蚀结果

图4传统HMT的目标提取结果、基于显著性距离变换的自适应HMT的目标提取结果和这里提出的自适应HMT的目标提取结果对比(a)原图像(b)传统HMT目标提取结果(c)基于显著性距离变换的自适应HMT的目标提取结果(d)这里提出的自适应HMT的目标提取结果

图5这里提出的自适应HMT在不同感兴趣区域半径u下的目标提取结果(a)原图像(b)u＝2(c)u＝3(d)u＝4(e)u＝5(f)u＝6

图6自适应结构元素对(B_FG，B_BG)的感兴趣区域半径u与目标的关系(a)u取值很小时情况(b)u取值很大时情况

具体实施方式

本方法共进行4组不同方面的实验，包括结构自适应、角点保护、噪声滤波和目标提取。

图1是椭圆型自适应结构元素方法与本方法的腐蚀结果对比。(a)图是原始图像，(b)-(g)图是椭圆型自适应结构元素方法，其中(b)-(d)图是r_w为4，M分别为4，8，16，(e)-(g)图是r_w为16，M分别为4，8，16。而(h)-(i)图是本方法腐蚀结果，其中r分别为4，8，16。从图中的对比结果可以看出，椭圆型自适应结构元素的方法在处理具有相同宽度和对比度的环状区域时，处理效果不尽相同，越接近中心，腐蚀的越细，这也说明在进行实验时，r_w与M需要保持相对的平衡；而从(h)-(i)图可以看出，本方法在处理这一区域时，效果相同。

图2是将椭圆型自适应结构元素方法、显著性自适应结构元素方法和本方法在角点保护方面的结果进行对比，其中(a)图是原图像，(b)图是显著性自适应结构元素方法的膨胀结果，(c)图是椭圆型自适应结构元素方法的膨胀结果，(d)-(f)分别是本方法的腐蚀、膨胀和开运算结果结果。从这些结果可以看出，前两种方法与本方法相比，都不能很好的保护角点，而且基于本方法的开运算基本与原始图像保持一致。

图3是将本方法应用到冶金晶界图像的滤波方面，并与椭圆型自适应结构元素的方法进行对比。实验中，本方法采用自适应的长轴允许最大半径r，具体如下：

①将图像按照连通区域的面积大小分为三个区域：白色噪声区域、晶界区域和背景区域，并利用公式(12)计算自适应结构元素的方向

②对于属于晶界或背景区域的任一像素点(x₁，y₁)，r(x₁，y₁)被定义为在方向上，与像素点(x₁，y₁)属于同一区域的像素点个数。

③对于属于白色噪声区域中的任意像素点(x₂，y₂)，r(x₂，y₂)被定义为在方向和与垂直的方向上，与像素点(x₂，y₂)属于同一区域的像素点个数的最大值。

第一列是三幅原始图像，第二列是椭圆型自适应结构元素的腐蚀结果，第三列是本方法的腐蚀结果。从第二列和第三列的结果可以看出，前者并没有将白色噪声点滤除干净，而且模糊了晶界区域，而本文方法在很好的保护了晶界的同时，有效滤除白色噪声。

图4和图5是采用HMT的目标提取结果。在图4中，将传统HMT、基于显著性距离变换的自适应HMT和本方法进行对比。可以看出，传统HMT只能提取出同一尺寸的目标，基于显著性距离变换的自适应HMT虽然提取出了不同尺寸的目标，但对于低亮度的目标并不能被提取出来，而本方法能够将所有不同尺寸、不同亮度的目标全部提取出来。在图5中，考虑了不确定区域半径u的尺寸对检测结果影响，(a)图为原始图像，(b)-(f)分别是u为2到6的HMT结果。很明显可以看出随着u的尺寸增加，提取出的目标越多，提取出的目标数量的统计结果可见表1，这是因为在u很小时，B_BG不能匹配背景，而在u较大时，才能够匹配背景，如图6所示。然而，u的尺寸并不是越大越好，因为u较大时，有可能匹配的并不是背景，而是相邻的目标。

表1图5中目标检测的统计结果

不确定区域u的半径	u＝2	u＝3	u＝4	u＝5	u＝6
						统计结果	29.90％	60.82％	83.51％	86.60％	84.54％

Claims (1)

1.基于非线性结构张量的自适应类圆型结构元素构造形态学算子的方法，包括以下步骤：

A.非线性结构张量

输入图像为I(x，y)的非线性结构张量可表示为：

$\frac{\&part;u_{i,j}}{\&part;t}=div\left(g\right(\underset{k,l=1}{\overset{2}{\mathrm{\&Sigma;}}}\&dtri;u_{k,l}\&dtri;u_{k,l}^T)\&dtri;u_{k,l}),(i,j=1,2)---\left(1\right)$

其中，u_i，j为初始结构张量的任一分量，I_x、I_y分别为图像I对水平和垂直方向的导数，为梯度算子，div表示散度算子，t的意义在离散时域里就是迭代滤波次数，g(·)是一种随像素张量越大而平滑力递减的平滑滤波函数，由于在式(1)中，是一个矩阵，所以平滑滤波函数g(·)作用后，可改写为如式(2)的形式：

$g(\underset{k,l=1}{\overset{2}{\mathrm{\&Sigma;}}}\&dtri;u_{kl}\&dtri;u_{kl}^T)=Vdiag\left(g_1\right({\mathrm{\&lambda;}}_1),g_2({\mathrm{\&lambda;}}_2\left)\right)V^T=g_1\left({\mathrm{\&lambda;}}_1\right){\mathrm{\&theta;}}_1{\mathrm{\&theta;}}_1^T+g_2\left({\mathrm{\&lambda;}}_2\right){\mathrm{\&theta;}}_2{\mathrm{\&theta;}}_2^T---\left(2\right)$

其中，λ₁≥λ₂为ST的特征值，θ₁和θ₂为一对正交基，g₁(λ₁)、g₂(λ₂)分别如式(3)、(4)所示：

g₁(λ₁)＝1/(1+λ₁+λ₂)^1.5 (3)

g₂(λ₂)＝1/(1+λ₁+λ₂)^0.5 (4)

公式(1)的迭代过程可由式(5)表示：

$u_{i,j}^{n+1}=u_{i,j}^n+\mathrm{\&Delta;}t\&CenterDot;\left(div\right(g\left(\underset{k,l=1}{\overset{2}{\mathrm{\&Sigma;}}}\&dtri;u_{k,l}\&dtri;u_{k,l}^T\right)\&CenterDot;\&dtri;u_{i,j}^n\left)\right),(i,j=1,2)---\left(5\right)$

迭代终止条件如式(6)：

$\underset{i,j=1}{\overset{2}{\mathrm{\&Sigma;}}}\left|u_{i,j}^{n+1}-u_{i,j}^n\right|<\mathrm{\&epsiv;}---\left(6\right)$

其中是迭代过程中非线性结构张量的任一分量的第n层数据，为第n+1层数据，Δt为步长，ε为预先设置的误差限，g(·)根据及时更新；

B.测量参数

由式(5)、(6)可得到满足迭代终止条件的的非线性结构张量，其特征值η₁、η₂和特征向量v₁、v₂的计算公式如式(7)、(8)所示：

${\mathrm{\&eta;}}_{1,2}=\frac{1}{2}(u_{11}+u_{22}\&PlusMinus;\sqrt{{(u_{11}-u_{22})}^2+4u_{12}^2})---\left(7\right)$

$\left\{\begin{array}{c}{v}_1={(2u_{12},u_{22}-u_{11}+\sqrt{{(u_{22}-u_{11})}^2+4u_{12}^2})}^T\\ v_2=v_1^{\&perp;}\end{array}\right.---\left(8\right)$

分析非线性结构张量的特征值η₁、η₂可以得到图像的一些局部结构信息，为了用量化指标测量感兴趣区域，定义了两个测量参数，如式(9)(10)：

①拐角强度测量m_c：

$m_c=\frac{4{\mathrm{\&eta;}}_1{\mathrm{\&eta;}}_2}{{\mathrm{\&eta;}}_1+{\mathrm{\&eta;}}_2}---\left(9\right)$

②边缘强度测量m_e：

m_e＝η₁-η₂ (10)

分析这两个测量参数，可以得出：

①当m_e≈0，m_c≈0时，表明对图像中的任一像素点，在该点沿任何方向的灰度变化都很小，对应的是图像的光滑区域；

②当m_e＞＞0，m_c≈0时，表明图像沿特征向量v₁方向的灰度变化率，远大于垂直于此方向的灰度变化率，对应的是图像的边缘区域；

③当m_e≈0，m_c＞＞0时，表明图像在该像素点沿任何方向的灰度变化都很大，对应的是图像的角点；

C.基于非线性结构张量的自适应椭圆结构元素的构建

自适应椭圆结构元素表示为其中，a(x，y)为椭圆型结构元素的长半轴，b(x，y)为其短半轴，表示椭圆型结构元素方向，三者定义如式(11)、(12)、(13)：

$a(x,y)=r\&CenterDot;\{1-\exp (-\frac{C_m}{{\left(m_c\right(x,y)/{\mathrm{\&beta;}}_1)}^m})\}---\left(11\right)$

$b(x,y)=a(x,y)\&CenterDot;\{1-\exp (-\frac{C_m}{{\left(m_e\right(x,y)/{\mathrm{\&beta;}}_2)}^m})\}---\left(12\right)$

其中，r表示允许的最大半径，β₁和β₂为归一化参数，分别取角点测量参数m_c和边缘测量参数m_e的最大值的75％，m、C_m为常量，这里取m＝1.5，C_m＝0.7627，v_2，x(x，y)、v_2，y(x，y)表示特征向量v₂的x，y分量；

D.采用随机分形插值法提高图像分辨率

随机分形插值的具体过程如下：

①定义原始图像为I(x，y)，内插后的图像为I_H(i，j)，根据分形插值模型，利用式(14)、(15)分别计算像素灰度正态分布的标准差δ和分形参数H：

δ＝E|I_H(t+1)-I_H(t)| (14)

log E|I_H(t+Δt)-I_H(t)|²-2H log||Δt||＝2logδ (15)

②当i，j均为奇数时：

I_H(i，j)＝I(x，y) (16)

③当i，j均为偶数时：

$\begin{array}{c}{I}_H(i,j)=\frac{1}{4}\{I_H(i-1,j-1)+I_H(i+1,j-1)+I_H(i-1,j+1)\}\\ +I_H(i+1,j+1)\}+\sqrt{1-2^{2H-2}}|\left|\mathrm{\&Delta;}t\right||\&CenterDot;H\&CenterDot;\mathrm{\&delta;}\&CenterDot;Gauss\end{array}---\left(17\right)$

④当i，j为一奇一偶时：

$\begin{array}{c}{I}_H(i,j)=\frac{1}{4}\{I_H(i,j-1)+I_H(i-1,j)+I_H(i+1,j)+I_H(i,j+1)\}\\ +2^{-H/2}\sqrt{1-2^{2H-2}}\left|\right|\mathrm{\&Delta;}t\left|\right|\&CenterDot;H\&CenterDot;\mathrm{\&delta;}\&CenterDot;Gauss\end{array}---\left(18\right)$

其中，E表示数学期望，Gauss为高斯随机变量，服从N(0，1)分布，||Δt||代表两插值点间的距离；

上述步骤可以重复迭代直到达到所需的空间分辨率为止，由此推导出当迭代k次后图像I_H中的像素点{2^kx-(2^k-1)，2^ky-(2^k-1)}处的灰度值即为原始图像I对应像素点(x，y)的灰度值；

E.提高自适应椭圆结构元素的分辨率

经过公式(11)、(12)、(13)计算得到椭圆参数后，需要构造离散的椭圆矩阵，通常情况下，离散的椭圆是由一系列整数型坐标(x，y)组成的集合X，且这些点(x，y)满足式(19)所示的椭圆方程的限制，椭圆的原点为(0，0)；

然而，为了提高离散椭圆的分辨率，这里取(x，y)＝(t/2^k，t/2^k)，若满足公式(19)，则像素点(x，y)属于集合X，其中-r×2^k≤t≤r×2^k为整数型数据，r为椭圆结构元素允许的最大半径；由此，将得到的集合X转换为矩阵形式，即得到插值后的离散椭圆结构元素；

F.自适应腐蚀、膨胀、开运算和闭运算的构造

输入图像I(x，y)，其定义域为D，插值后图像为I_H(i，j)，其定义域为D′，I(x，y)中任一像素点L(x，y)∈D处的自适应结构元素为ESE{2^kx-(2^k-1)，2^ky-(2^k-1)}为SE(x，y)插值后的自适应结构元素，像素点L′{2^kx-(2^k-1)，2^ky-(2^k-1)}∈D′，由此定义自适应腐蚀、膨胀、开、闭运算如公式(20)-(23)所示：

IO(x，y)＝γ(I_H)(i，j)＝(δ(ε))(I_H)(i，j) (22)

IC(x，y)＝ψ(I_H)(i，j)＝(ε(δ))(I_H)(i，j) (23)

其中，i＝2^kx-(2^k-1)，j＝2^ky-(2^k-1)，(x，y)∈D，(i，j)∈D′，ε、δ分别表示腐蚀、膨胀运算，γ、ψ分别表示开、闭运算，ESE^c(i，j)为ESE(i，j)反射邻域；

G.自适应击中击不中(Hit-or-Miss)变换

击中击不中变换涉及一个结构元素对(B_FG，B_BG)，为了提取类圆形目标，这里定义自适应结构元素的形状为圆形，尺寸计算过程如下：

①利用公式(10)中的非线性结构张量的边缘强度测量m_e来判断输入图像的边缘像素点，如式(24)所示，其中delta为经验阈值；

$\left\{\begin{array}{cc}f(x,y)=m_e(x,y),& \begin{array}{cc}if& m_e(x,y)>delta\end{array}\\ f(x,y)=0,& other\end{array}\right.---\left(24\right)$

②对①中得到的f(x，y)进行非极大值抑制，然后令属于局部极大值的像素点赋值为1，非局部极大值点赋值为0，由此可以得到单像素宽的边缘图像f_edge；

③对f_edge进行倒角距离变换(Chamfer Distance Transform)，得到的图像每一点的灰度值为该点到达最近边缘点的距离，将此距离作为结构元素B_FG的半径，记为R_HMT；

根据得到的半径R_HMT即可定义自适应结构元素B_FG如式(25)所示：

B_FG(x，y)＝{|(x，y)-(p，q)|＜R_HMT(x，y)}，(p，q)∈D (25)

其中，(x，y)、(p，q)均为任意两个点；

因为在定义自适应结构元素对(B_FG，B_BG)时，需要满足B_FG和B_BG具有共同的原点，并且互不连通，所以定义自适应结构元素B_BG如式(26)所示：

B_BG(x，y)＝{R_HMT(x，y)+u＜|(x，y)-(p，q)|＜R_HMT(x，y)+u+b}，(p，q)∈D (26)

其中，u为不确定区域半径，b为B_BG半径，两者在整幅图像中保持固定不变。