局部特征提取

局部特征提取（通用7篇）

局部特征提取 篇1

摘要：局部保持投影算法仅能保持近邻样本的局部结构,无法保证提取的特征有利于后续分类识别。为此,提出一种半监督保持投影特征提取算法。SPP算法能够利用标记样本所携带的类别信息来约束未标记样本,从而提高样本的可分性;同时,还在目标函数中加入一正则项,避免了因矩阵奇异导致算法无法求解的问题。利用实际高光谱数据进行对比实验,结果表明,用SPP算法进行特征提取后的分类精度较LPP算法有显著提升,验证了它的有效性。

关键词：局部保持投影,特征提取,半监督,高光谱

0 引 言

高光谱遥感技术利用成像光谱仪纳米级的光谱分辨率,以几十或几百个波段同时对地物成像,能够获得地物接近连续的光谱信息,成为高分辨对地观测系统的重要组成部分之一。然而,高光谱传感器为地物分类和识别提供细致光谱特征的同时也带来了大量的冗余数据,给后续分类处理带来了难度。为了去除这些冗余数据,通常需对高光谱数据进行特征提取。

局部保持投影[1](Locality Preserving Projections, LPP)是最近提出的一种线性特征提取算法,它是Laplacian Eigenmap[2] 流形学习算法的线性近似,该算法既克服了主成分分析和Fisher判别分析[3]等线性特征提取算法难以保持原始数据非线性流形的缺点,同时又解决了非线性流形学习方法[4]难以直接映射新样本的问题,因此在模式分类领域得到广泛应用[5,6]。然而,LPP算法仅能保持近邻样本对的局部结构,对原本相互远离的样本,LPP算法并未加以约束。对此,Yang等提出无监督鉴别保持投影算法[7],其在建模的同时考虑局部和非局部散度,从一定程度上克服了LPP存在的问题,但该算法在小样本条件下会出现矩阵奇异,导致算法无法求解。王建国等则利用非局部散度和局部散度之差作为鉴别准则,避免了矩阵奇异问题,然而分类精度有所降低[8]。

上述算法虽对LPP性能有一定改善作用,但其仅能够保证经过特征提取后样本的相对位置能够得到保持,无法保证提取的特征更容易分类。 针对这一缺点,本文提出一种改进算法,该算法在LPP算法的目标函数中加入标记样本信息,通过标记样本所携带的类别信息来约束未标记样本的相对位置;同时,还通过在目标函数中加入正则项以避免矩阵奇异问题。由于这种算法在训练过程中同时用到了标记样本和未标记样本——这是一种介于有监督学习和无监督学习之间的半监督学习思想[9],故将这种改进算法称为半监督保持投影(Semi-supervised Preserving Projections,SPP)算法。实验结果表明,经过改进的SPP算法进行特征提取后,不同类别样本的区分性增强,分类精度比LPP算法有较大提高。

1 局部保持投影算法

局部保持是一种典型的线性流形学习算法,其原理可用谱图理论进行解释。设X=[x1,x2,…,xN]是由m维向量构成的数据集合,x1,x2,…,xN∈Rm。LPP目标是寻找一个投影矩阵W,该矩阵可以将位于高维空间的数据映射至一个低维子空间Rd(d≪m)中,并同时尽可能使数据的局部相对位置保持不变,即:如果样本对在原始高维空间中是相互靠近的,那么经投影矩阵变换至低维子空间后,该数据对仍然保持相互靠近。令数据集在低维子空间Rd中表示为Z=[z1,z2,…,zN],则Z=WTX。

为了使数据的局部相对位置保持不变,LPP算法的目标函数定义如下:

$\min \left(\frac{1}{2}{\displaystyle \sum _{i=1}^{{\rm N}}{\displaystyle \sum _{j=1}^{{\rm N}}\parallel }}\mathit{y}{}_i-\mathit{y}{}_j\parallel {}^2\mathit{S}{}_{ij}\right)(1)$

式中:yi=WTxi,且W=[w1,w2,…,wd];Sij是权重矩阵,采用k近邻法来定义:

$\mathit{S}{}_{ij}=\{\begin{array}{l}\exp \left(-\frac{\parallel \mathit{x}{}_i-\mathit{x}{}_j\parallel {}^2}{t}\right),\mathit{x}{}_i\in {\rm N}{}_k(x{}_j)\\ 0,\text{o}\text{t}\text{h}\text{e}\text{r}\text{w}\text{i}\text{s}\text{e}\end{array}(2)$

式中:t为大于0常数;Nk(xi)表示由xi的k个最近邻样本点组成的集合。

对式(1)进行简单的变量代换,可得到下式:

$\begin{array}{l}\frac{1}{2}{\displaystyle \sum _{i=1}^{{\rm N}}{\displaystyle \sum _{j=1}^{{\rm N}}\parallel }}\mathit{y}{}_i-\mathit{y}{}_j\parallel {}^2\mathit{S}{}_{ij}\\ =\frac{1}{2}{\displaystyle \sum _{i=1}^{{\rm N}}{\displaystyle \sum _{j=1}^{{\rm N}}\parallel }}\mathit{W}{}^{\text{Τ}}\mathit{x}{}_i-\mathit{W}{}^{\text{Τ}}\mathit{x}{}_j\parallel {}^2\mathit{S}{}_{ij}\\ ={\displaystyle \sum _i\mathit{W}}{}^{\text{Τ}}\mathit{x}{}_{i}d{}_{ii}\mathit{x}{}_i^{\text{Τ}}\mathit{W}-{\displaystyle \sum _{i,j}\mathit{W}}{}^{\text{Τ}}\mathit{x}{}_i\mathit{S}{}_{ij}\mathit{x}{}_j^{\text{Τ}}\mathit{W}\\ =\mathit{W}{}^{\text{Τ}}\mathit{X}(\mathit{D}-\mathit{S})\mathit{X}{}^{\text{Τ}}\mathit{W}\\ =\mathit{W}{}^{\text{Τ}}\mathit{X}\mathit{L}\mathit{X}{}^{\text{Τ}}\mathit{W}(3)\end{array}$

式中:D=diag(d11,d22,…,dNN)为对角矩阵,其对角线上的元素的矩阵S中对应的列(或行)元素之和,即$d{}_{ii}={\displaystyle \sum _{j=1}^{{\rm N}}\mathit{S}}{}_{ij}‚\mathit{L}=\mathit{D}-\mathit{S}$,矩阵L为Laplacian矩阵。

2 半监督保持投影算法

局部保持投影算法能够使近邻样本对在低维子空间中仍然保持相互邻近,然而由于其没有利用标记样本所携带的类别信息,无法保证原始数据投影至低维空间后更有利于进行分类判别。为了克服LPP算法的这一缺点,本文提出一种同时利用大量未标记样本和少量标记样本的半监督保持投影(SPP)算法。在SPP算法中,为了利用标记样本的类别信息,在目标函数中引入了类似于Fisher判别思想[10]的类内散布矩阵和类间散布矩阵;并且定义了非局部结构信息,使非近邻样本对的相对结构在低维空间中也能得到保持;同时还在目标函数中增加一正则项,克服了矩阵奇异问题。

2.1 算法描述

假定标记样本集为L={xi,yi}${}_{i=1}^{{\rm N}}$,类别标签yi∈{1,2,…,c},其中c为样本类别总数; Nk为第k类的样本数目;x${}_i^k$为属于第k个类别的第i类样本。令S(w)和S(b)分别代表样本的类内散布矩阵和类间离散度矩阵,其定义分别如式(4)和式(5)所示:

$\begin{array}{l}\mathit{S}{}^{(w)}={\displaystyle \sum _{k=1}^c{\displaystyle \sum _{i=1}^{{\rm N}{}_k}(}}\mathit{x}{}_i^k-\mathit{\mu }{}^k)(\mathit{x}{}_i^k-\mathit{\mu }{}^k){}^{\text{Τ}}(4)\\ \mathit{S}{}^{(b)}={\displaystyle \sum _{k=1}^c{\rm N}}{}_k(\mathit{\mu }{}^k-\mathit{\mu })(\mathit{\mu }{}^k-\mathit{\mu }){}^{\text{Τ}}(5)\end{array}$

式中:μk为第k类样本的均值向量;μ为所有样本的均值向量,其定义如式(6)和式(7)所示:

$\begin{array}{l}\mathit{\mu }{}^k=\frac{1}{{\rm N}{}_k}{\displaystyle \sum _{i=1}^{{\rm N}{}_k}\mathit{x}}{}_i^k(6)\\ \mathit{\mu }=\frac{1}{{\rm N}}{\displaystyle \sum _{i=1}^{{\rm N}}\mathit{x}}{}_i(7)\end{array}$

通过求解式(8)所示的Fisher准则,即可实现同类样本在低维子空间中相互聚拢,非同类样本相互远离。

$J(\mathit{W})=\text{a}\text{r}\text{g}\text{m}\text{a}\text{x}\frac{\mathit{W}{}^{\text{Τ}}\mathit{S}{}^{(b)}\mathit{W}}{\mathit{W}{}^{\text{Τ}}\mathit{S}{}^{(w)}\mathit{W}}(8)$

式(8)的求解方法与式(3)类似,采用Lagrange乘子法求解,因此可将类内散布矩阵和类间散布矩阵加入目标函数式(3)中。最佳投影矩阵W的求解方法等同于求解如式(9)的特征值问题:

$\mathit{S}{}^{(b)}\mathit{W}=\lambda (\alpha \mathit{S}{}^{(w)}+(1-\alpha )\mathit{X}\mathit{L}\mathit{X}{}^{\text{Τ}})\mathit{W}(9)$

利用式(9)求解出投影矩阵W,样本在低维子空间不仅能够保持局部结构;而且还利用标记样本的类别信息,使得同类样本聚集在一起,不同类样本相互远离。参数α控制着标记样本对投影矩阵的影响,由其取值范围为[0,1]:当α=1时,投影矩阵完全由标记样本确定;当α=0时,式(13)退化为LPP算法。

然而,仅保持样本在低维子空间中的局部结构不变并不能达到最佳的效果;对于原本相互远离的样本,应该在其投影到低维空间后也保持相互远离。为此,本文参照式(3)局部结构的定义,又定义了非局部结构的目标函数:

$\begin{array}{l}\max (\frac{1}{2}{\displaystyle \sum _{i=1}^{{\rm N}}{\displaystyle \sum _{j=1}^{{\rm N}}\parallel }}\mathit{y}{}_i-\mathit{y}{}_j\parallel {}^2\tilde{\mathit{S}}{}_{ij})\\ =\max ({\displaystyle \sum _i\mathit{W}}{}^{\text{Τ}}\mathit{x}{}_i\tilde{d}{}_{ii}\mathit{x}{}_i^{\text{Τ}}\mathit{W}-{\displaystyle \sum _{i,j}\mathit{W}}{}^{\text{Τ}}\mathit{x}{}_i\tilde{\mathit{S}}{}_{ij}\mathit{x}{}_j^{\text{Τ}}\mathit{W})\\ =\max (\mathit{W}{}^{\text{Τ}}\mathit{X}\tilde{\mathit{L}}\mathit{X}{}^{\text{Τ}}W)(10)\end{array}$

式中:$\tilde{\mathit{L}}=\tilde{\mathit{D}}-\tilde{\mathit{S}}‚\tilde{\mathit{D}}=\text{d}\text{i}\text{a}\text{g}$$(\tilde{d}{}_{11},\tilde{d}{}_{22},\cdots ,\tilde{d}{}_{{\rm N}{\rm N}})‚\tilde{d}{}_{ii}={\displaystyle \sum _{j=1}^{{\rm N}}\tilde{\mathit{S}}}{}_{ij}$。$\tilde{\mathit{S}}$ij的定义如下:

$\tilde{\mathit{S}}{}_{ij}=\{\begin{array}{l}1-\exp \left(-\frac{\parallel \mathit{x}{}_i-\mathit{x}{}_j\parallel {}^2}{t}\right),\mathit{x}{}_i\notin {\rm N}{}_k(\mathit{x}{}_j)\\ 0,\text{o}\text{t}\text{h}\text{e}\text{r}\text{w}\text{i}\text{s}\text{e}\end{array}(11)$

将式(10)也加入式(9)中,得到式(12):

$(\alpha \mathit{S}{}^{(b)}+(1-\alpha )\mathit{X}\tilde{\mathit{L}}\mathit{X}{}^{\text{Τ}})\mathit{W}=\lambda (\alpha \mathit{S}{}^{(w)}+(1-\alpha )\mathit{X}\mathit{L}\mathit{X}{}^{\text{Τ}})\mathit{W}(12)$

此时,利用式(12)得到投影矩阵W后,样本在低维子空间中不仅能够保持近邻样本的局部结构,并且可以保持非近邻样本对的非局部结构;同时,在式(12)中,标记样本的类别信息也得到利用,这使样本投影到低维空间后,属于同一类别的样本会更加靠近,属于不同类别的样本彼此远离。因此从分类角度考虑,利用该算法进行特征提取比原始的LPP算法更加适合进行后续分类处理。但是,在小样本条件下,式(12)等式右端的S(w)+XLXT可能为奇异矩阵,从而导致该特征值问题无法正常求解。为此,本文依照文献[11]的处理方法,通过添加一正则项从而使其变为非奇异矩阵:

$\begin{array}{l}(\alpha \mathit{S}{}^{(b)}+(1-\alpha )\mathit{X}\tilde{\mathit{L}}\mathit{X}{}^{\text{Τ}})\mathit{W}=\\ \lambda (\alpha \mathit{S}{}^{(w)}+(1-\alpha )\mathit{X}\mathit{L}\mathit{X}{}^{\text{Τ}}+\gamma {\rm I}{}^m)\mathit{W}(13)\end{array}$

最终,本文提出的SPP算法的投影矩阵由式(13)求解,该式是一个特征值求解问题,投影矩阵W为前r个最大特征值对应的特征向量组成的矩阵,r为子空间维数。

2.2 算法步骤

综上所述,SPP算法步骤可以归纳如下:

输入:N个m维原始空间数据组成的N×m矩阵X;

输出: N个d维特征空间数据组成的N×d矩阵Y,m≥d;

(1) 在原始空间中计算样本点之间的距离,找出每个未标记样本点的k个最近样本点;

(2) 计算未标记样本对间的权重矩阵Sij和$\tilde{\mathit{S}}$ij,由此可以分别计算未标记样本的局部结构XLXT和非局部结构X$\tilde{\mathit{L}}$XT;

(3) 对于标记样本,分别计算其类内散度S(w)和类间散布矩阵S(b);

(4) 根据式(13)求解投影矩阵W;

最终,原始高维数据在特征空间中的投影为Y=WTX。

3 实验结果与分析

实验采用机载可见/红外成像光谱仪AVIRIS获取的佛罗里达州肯尼迪中心(KSC)高光谱数据作为本文算法的验证数据(如图1所示)。该数据共有224个波段,实验中去除了受大气水分影响以及低信噪比的波段,保留了155个信噪比较高的波段构成数据集合。该数据集中有13类地物,本文选取数量较多的10类地物作为实验的样本集(见表1)。实验平台为Matlab 7.8,Pentium 4处理器 1.8 GHz,1 GB RAM,操作系统为Windows XP。

由第2节可知,本文提出的SPP算法有k,t,α,γ四个参数待定。由于这四个参数取值并无理论指导,故先通过实验确定这四个参数的最优取值。实验发现,参数k的取值对算法性能影响较小,故本文中直接给出该参数的取值k=5。实验选取表1中每类样本的20%作为训练样本集,其余样本作为测试集,采用K-近邻分类器对特征提取后的数据进行分类,实验结果如图2和图3所示。由实验结果可以看出:参数t,α,γ的取值对算法性能有极大影响,在该实验中,当α=0.8,t=0.1,γ=0.5时,算法SPP可以达到最优的特征提取效果。

为了验证SPP算法的有效性,实验还给出了在不同训练样本数量条件下SPP算法的分类精度,并与PCA,LPP以及文献[10]提出的半监督特征提取算法LGSSDR进行性能比较。实验采用K-近邻分类器(K-NN),使用分类结果精度衡量4种特征提取算法性能。实验选取表1所示样本集中的β%作为训练样本,剩余样本作为测试样本,β的取值范围为5～50。降维后的子空间维数取 [12]10。图4给出了在不同训练样本数量情况下4种特征提取算法在K-NN(k=1)分类器下的分类精度。

由图4的分类结果可以看出,由于PCA和LPP均为无监督特征提取算法,无法利用标记样本带来的监督信息,因此分类精度很低,难以满足实际应用要求。而本文提出的半监督保持投影算法SPP在原有LPP基础上加入了类内散度和类间散度,能够充分利用标记样本的类别信息来约束未标记样本点;使得样本投影到特征空间后,同类样本相互靠近,非同类相互远离,降低了分类难度,因此分类精度LPP算法有显著提高。同时,SPP算法的求解类似于Fisher判别分析,在标记样本数量较少时也会出现奇异矩阵,导致算法无法求解。为此,在目标函数的分母上增加了正则项,保证了算法在所有条件下均能够正常求解。

与LGSSDR算法相比:在标记样本数目较少时,本文提出的SPP算法的分类精度与LGSSDR算法相当;但是,随着标记样本数目的增加,SPP算法逐渐显现出其优势,分类精度高于LGSSDR算法。表2给出了4种特征提取算法在10种地物类型上的分类精度。由该表可以看出,在大多数地物类型中SPP算法的分类精度均高于其他三种算法,特别是在第2,3,5类中,SPP算法的分类精度得到显著提高。

4 结 论

LPP算法是基于样本的局部性进行建模,无法保证其提取出的特征有利于后续分类处理。对这一问题,本文提出了一种半监督保持映射(SPP)的特征提取算法。该算法充分利用标记样本的类别信息来约束未标记样本,使样本经SPP算法提取特征后,分类精度显著提高;同时正则项的加入避免了算法出现矩阵奇异问题。最后,通过实际KSC高光谱数据进行对比实验表明了SPP算法的有效性。

参考文献

[1]HE Xiao-fei,NIYOGI P.Locality preserving projections[J].Advance in Neural Information Processing Systems,2004,16:153-160.

[2]BELKIN M,NIYOGI P.Laplacian eigenmaps for dimen-sionality reduction and data representation[J].NeuralComputation,2003,15(6):1373-1396.

[3]DUDA R,HART P.Pattern classification and scene analy-sis[M].New York:Wiley,1973.

[4]宋欣,叶世伟.基于局部线性逼近的流形学习算法[J].计算机仿真,2008,25(7):86-89.

[5]郭金玉,苑玮琦.基于局部保持投影的掌纹识别[J].光学学报,2008,28(10):1920-1924.

[6]TANG Y,RICHARD R.A study of using locality preser-ving projections for feature extraction in speech recognition[C]//Proceedings of IEEE International Conference on A-coustics,Speech and Signal Processing.[S.l.]:IEEE,2008:159-1572.

[7]王建国,杨万扣,杨静宇.基于保持投影的最大散度差的特征抽取方法[J].模式识别与人工智能,2009,22(4):610-613.

[8]YANG Jian,ZHANG D,YANG Jing-yu,et al.Globallymaximizing locally minimizing unsupervised discriminantprojection with applications to face and palm biometrics[J].IEEE Trans.on Pattern Analysis and Machine Intelli-gence,2007,29(4):650-664.

[9]CHAPELLE O,SCHOLKOPF B.Semi-supervised learning[M].Cambridge,MIT Press:2006.

[10]韦佳,彭宏.基于局部与全局保持的半监督维数约减方法[J].软件学报,2008,19(11):2833-2842.

[11]KUO B C,LI C H,YANG J M.Kernel nonparametricweighted feature extraction for hyperspectral image classi-fication[J].IEEE Trans.on Geoscience and RemoteSensing,2009,47(4):1139-1155.

[12]HONG Z Q,YANG J Y.Optimal discriminant plane for asmall number of samples and design method of classifier onthe plane[J].Pattern Recognition,1991,24(4):317-324.

局部特征提取 篇2

关键词：分形维数,医学超声图像,轮廓提取

1 引言

医学超声图像轮廓提取是医学超声图像处理的关键技术之一, 目的是在存有大量噪声背景的医学超声图像中确定目标腔体的轮廓[1]。图像轮廓的清晰与否对医学诊断具有十分重要的作用。为了做出更准确的诊断, 医生需要定量地知道如心输出量等一些重要的参数信息, 这就希望在医学超声图像中能够较为精确地确定出目标腔体 (如心房和心室等) 的舒张和收缩轮廓, 从而为后续的诊断和治疗提供有效的信息。

然而, 由于受超声图像成像原理等固有因素影响, 图像的轮廓通常比较模糊[2]。经典的图像轮廓提取主要有Sobel算子[3]、Laplace算子[4]、Preivitt算子[5]、Kirsch算子[6]、Robert算子[7]等方法。这些基于微分算子的提取方法对普通图像的轮廓有较好的检测效果, 然而超声图像中存在大量的Speckle噪声, 这使得图像轮廓难以检测。直接将这些方法应用于医学超声图像, 难以获得准确的图像轮廓提取结果。

分形几何中通常采用分形维数的方法来估计一幅图像的复杂程度, 分形维数能够反映一幅图像中更深层次的结构信息。本文将分形维数概念扩展, 介绍一种局部模糊分形维数LFFD (Local Fuzzy Fractal Dimension) , 用来提取图像的轮廓信息。鉴于局部模糊分形维数对噪声不敏感, 本文采用基于局部分形维数特征的方法, 通过提取图像的局部模糊分形维数特征, 对医学超声图像进行轮廓提取。

2 局部模糊分形维数的原理及基于分形维数特征的轮廓提取方法

2.1 分形维数

在欧氏空间中, 人们习惯把空间看成三维的, 面看成二维, 而把线看成一维。也可以稍加推广, 认为点是零维的, 还可以引入高维空间, 但通常人们习惯于整数的维数。为了定量地描述客观事物的复杂程度, 1919年, 数学家F.Hausdorff从测量的角度引入了维数概念, 将维数从整数扩大到分数, 突破了一般拓扑集维数为整数的界限。

分形维数FD (Fractal Dimension) 有许多不同的定义, 计盒维数易于计算, 得到广泛使用。它可以用点覆盖方法来估计二维单色图像分形维数。用点覆盖的方法计算二值图像计盒维数时, 图像F被分割成M个宽度为δ的方块Ciδ。二值图像白点用1表示, 黑点用0表示。

计盒维数定义如下:

其中, Nδ (F) 是用半径不大于δ的盒子覆盖图像F时至少包含1个白点的盒子数目。可由下式求得:

其中, 若象素点p是白点, f (p) 取值为1;若象素点p是黑点, f (p) 取值为0。

实验中, 通过获取一组样本数据 (-logδ, log Nδ (F) ) , 采用最小二乘线形回归方法估计这群样本数据点, 计算直线斜率就可以得到该图形的分形维度。

2.2 模糊分形维数

分形维数只能估计二值图像的复杂程度。然而实际应用中, 元素的特征值常常不只取两个值。例如, 对于灰度图像, 在用点覆盖方法估计分形维数时就需要先二值化。二值化将导致大量信息的丢失。这对于医学图像分析诊断是不合适的。

因此将分形维数概念扩展, 用模糊分形维数FFD (Fuzzy Fractal Dimension) 作为灰度图像的特征。将图像F分割成M个宽度为δ的方块Ciδ。模糊分形维数定义如下:

其中, Nδ (F) 是用半径不大于δ的盒子覆盖图像F时所用到的盒子数目。可由下式求得:

其中, 象素点p的特征值f (p) 反映象素点p与白点的相似程度。对于256级灰度图像, 点的特征值采用点灰度值G (p) 除以256。特征值越大, 点与白点越相似。

2.3 局部模糊分形维数

FFD可以表明灰度图像整体的复杂程度。为了能够提取医学图像的局部特征, 图像点的局部模糊分形维数LFFD (Local Fuzzy Fractal Dimension) 被定义为该象素点周围邻域的FFD。由于FFD表示的是灰度图像整体的复杂程度, LFFD也就反映了灰度点周围邻域图像的复杂程度。局部模糊分形维数LFFD的计算流程如图1所示。

2.4 基于局部分形维数特征的医学超声图像轮廓提取

由于不同的灰度、不同的纹理形成的图像轮廓通常有不同的复杂度, 大量实验表明局部模糊分形维数可以反映图像的轮廓信息。又因为局部模糊分形维数对噪声不敏感, 于是, 实验中采用基于局部模糊分形维数的方法, 通过提取图像的局部模糊分形维数特征, 对医学超声图像进行轮廓提取。

3 实验结果与讨论

实验选取了大量医学超声图像, 图2 (a) 为其中一幅心室声像图。下面将基于局部分形维数特征的轮廓提取方法与传统微分算子方法相比较, 对该心室声像图进行轮廓提取。

图2 (b) 到图2 (e) 依次为采用Laplace算子[3]、Preivitt算子[4]、Kirsch算子[5]、Robert算子[6]对图2 (a) 进行轮廓提取的结果, 可以看出, 对于医学超声图像, 这些方法的检测结果受噪声影响都相当严重, 几乎不能达到轮廓提取的效果。图2 (f) 是基于LFFD方法的检测结果, 它较上述轮廓提取方法对噪声有较好的抗干扰性, 同时图像轮廓得到较好的提取, 检测效果明显优于上述轮廓提取方法。

4 结论

本文采用了基于局部分形维数特征的轮廓提取方法提取医学超声图像的轮廓。实验结果表明, 与经典微分算子轮廓提取方法相比, 基于分形维数特征的轮廓提取对医学超声图像中存在的噪声有较好的抗干扰性, 能够较好地提取图像轮廓信息, 为医学超声图像轮廓提取开辟了一条新的途径。

参考文献

[1]陈文山, 汪天富, 林江莉.基于相似度测量的医学超声图像对比度增强[J].中国医学影像技术.2006 (9) :1432-1434.

[2]陈文山.基于曲面形变模型的心脏超声图像三维分割[J].福建电脑.2011 (9) :96-97.

[3]N.C.Woods.A Sobel Edge Detection Algorithm Based System for Analyzing and Classifying Image Based Spam[J].Journal of Emerging Trends in Computing and Information Sciences.2012 (4) :506-111.

[4]Binghui Wei.The Effective Particles Edge Detection Method Based on Laplace[J].Procedia Engineering.2012 (29) :3096-3099.

[5]Jair Garcia-Lamont.Analogue CMOS prototype vision chip with prewitt edge processing[J].Analog Integrated Circuits and Signal Processing.2012 (3) :507-514.

[6]Shao Ping.Fast Kirsch edge detection based on templates decomposition and integral image[J].Acta Automatica Sinica.2007 (8) :795-800.

局部特征提取 篇3

关键词：低质量指纹,方向场,局部字典

指纹方向场能够清晰地展示指纹的基本形状、结构和脊线方向等全局信息, 同时也能够保留绝大部分的细节信息。因而在实际的应用中, 正确的指纹方向场是指纹后续处理的关键。

对于低质量指纹, 由于存在噪声, 导致识别率显著降低。而在现实指纹采集和我国指纹档案中大部分都是低质量指纹, 低质量指纹对指纹识别带来了严重的干扰, 大大提高了指纹识别的难度。为了从低质量指纹中正确地提取其特征, 冯[1]提出了一种基于指纹结构先验知识 (即全局字典) 的新的算法。该方法由线下学习和线上估计2部分组成, 由于在线下学习中构建的是一个针对全局方向场块的字典, 即全局字典, 其中方向场块数目太多, 导致构建和查找速度慢。对此杨[2]使用局部字典代替全局字典, 同时也提出了基于霍夫变换的姿势估计, 其大大提高了总体速度和精确度, 本文对杨[2]文献进行研究时, 发现其中构建字典大小和选取候选方向场块列表设定的阈值, 这2个参数对提取出来的方向场有较大影响, 因此本文对这2个参数进行了介绍和优化, 同时在挑选原型方向场块时, 使用的是kmean聚类方法[3], 以及最后基于上下文方向场纠正时, 为了更快地重建出方向场块, 本文只考虑了局部最优解。详细步骤可查询文献2。

1 线下学习

线下学习的任务是完成原型指纹方向场块的空间分布估计及局部字典的构建, 为线上的姿势估计和字典查询作准备。主要是通过对真实指纹样本和人工标注的方向场进行学习。

1.1 原型方向场块的生成

首先将指纹样本方向场分割成NP*NP (设置为4) 互不重叠且位置已知的方向场块, 这些样本方向场块称为训练方向场块。然后对训练方向场块进行聚类, 挑选每类中具有代表性的方向场块作为原型方向场块。我们使用收敛速度快的K-m e a n s聚类方法[3], 因为该聚类方法比杨采用足本文的要求。因为每类中心处的方向场块之间的差别较大, 且在每类中最具有代表性, 在此选取每类中心处的方向场块作为原型方向场块。

1.2 原型方向场块的分布估计

原型方向场块的分布情况, 表示的是每一块原型方向场块在指纹图像中分布情况, 这里采用的是无参数方式进行估计的, 原型方向场块在坐标系统中位置处的概率分布[2]用P (ϑi|u) 表示, 通过第i类所有方向场块在u位置出现的数目来表示。

Nu, i为第i类样本方向场块出现在u位置处的数目, Nf为指纹样本总数目。其中P (ϑi|u) 原型方向场块在一个确定位置的概率分布可以看成先验概率, 而原型方向场块出现在任意位置的概率分布[2]作为后验概率, 即P (u|ϑi) , 如 (9) 式所示:

u'表示在同一坐标系统中任意位置, 一般条件下P (u) =P (u') 。

1.3 局部字典的构建

将所有样本根据自身指纹中心点位置建立坐标系, 对于每一个位置u, 构建出一个对应的局部字典Du。为了使构建出来的字典具有多样性, 周围邻域内的方向场块都考虑到该局部字典当中来。先将对应位置的邻域内的所有训练方向场块先都加入到临时数组变量中Tu, 然后构建局部字典Du, 构建方法如文献2一样。

邻域的大小对挑选候选方向场块列表有很大的影响, 在此本文分别选取窗口3*3和5*5进行了实验, 其中3*3大小表示以该位置为中心点, 其周围8邻域的点也加入到该位置局部字典中来, 5*5大小则是除了8邻域还有8邻域相邻的点都加入进来。若窗口过小, 局部字典包含的方向场块数目少, 存在和低质量指纹方向场块相似的候选方向场块概率小, 但是当窗口过大时, 从局部字典查找候选方向场块耗时。经过实验发现, 当窗口取5*5时, 既能够保证局部字典中含有和初始方向场块相似的块, 也能够保证查找速度快。

通过对所有的训练方向场进行上述构建步骤, 得到了一系列不同位置的局部字典。局部字典包含的方向场块越多, 则说明该位置会出现更多类型的方向场块, 反之则越少。

2 线上估计

线上估计的主要任务是从局部字典中挑选合适的方向场块代替低质量指纹初始方向场块, 以达到重建低质量指纹方向场的目的。主要步骤为先估计出低质量指纹图像的初始方向场, 再利用线下学习阶段中的原型方向场块的空间分布和局部字典, 分别对估计出低质量指纹的姿势和挑选候选方向场块, 最后根据上下文确定最终方向场块。

2.1 低质量指纹的初始方向场估计

低质量指纹的初始方向场估计, 将对挑选候选方向场块有较大的影响。若估计出来的初始方向场较好, 则经过后面的纠正, 能够很好地重建出丢失的信息, 反之, 则会为后面的工作带来很大的难度。因此初始方向场估计非常的关键, 最常用的有梯度法[5]、短时傅里叶变换 (STFT) [6]和基于二维傅里叶扩张 (FOMFE) [7]的等方法。在此本文选用了STFT, 其一它不像梯度法那样容易受噪声影响, 其二没有像FOMFE那样复杂。

2.2 指纹姿势估计及字典查询

指纹的姿势表示指纹中心点的位置和指纹的方向。所以对指纹姿势估计即求出上述2个值, 通常在采集指纹时, 手指都是垂直向上的, 所以本文直接将垂直向上作为指纹的方向, 因此本文中所有表述的指纹姿势估计都是指指纹的中心点估计, 姿势估计采用是文献2所述的方法。

将低质量指纹估计出的中心点作为原点, 垂直向上作为y轴方向, 建立一个坐标系, 估计该坐标系下位置u的低质量方向场块时, 只需从局部字典Du中位置u处查询出块候选方向场块。为了让nc个含有更多信息的方向场块, 文献1设定候选方向场块之间的相似度需小于0.8, 当候选方向场块之间相似度越接近0.8时, 则候选方向场块较相似, 很难达到候选方向场块中含有更多信息, 本文将其设置为0.4, 这样既保证了候选方向场块中包含更多信息, 也包含了和低质量指纹方向场块相似度很高的候选方向场块。具体的字典查询方法请参考文献1。

通过上述字典查询后, 我们得到了c (1≤c≤nc) 个和初始方向场相似度高, 且包含更多类型的方向场块。图1 (b) 为 (a) 中方框内的候选方向场块, 图1 (c) 展示了由人工标志的方向场。由图可知, 当初始方向场含有大量噪声时 (图 (a) 左上角) , 候选方向场块列表中仍然包含和图1 (c) 对应位置相似的方向场块。

2.3 挑选方向场块, 基于上下文的方向场纠正

我们将从每个位置u中的候选方向场块列表中, 挑选出一个最优方向场块, 以达到纠正初始方向场块的目的。其中每个位置含有c (1≤c≤nc) 个候选方向场块, 并构成一个方向场块列表Ru, Ru={ru, 1, ru, 2, ..., ru, c}, ru, i是u位置处字典中c个候选方向场块中的第i个方向场块, 为了从候选方向场块列表中找到合适的方向场块, 在此考虑了2个方面:其一, 候选方向场块与初始方向场块的相似度;其二, 候选方向场块与其邻域 (本文考虑左领域) 的相似性 (即兼容性) 。综合上述2方面, 根据马尔科夫模型[8]构建的一个局部能量函数, 当局部能量函数达到最大值时, 上述2方面的和最大, 可以判定该候选方向场块最接近原始方向场块, 局部能量函数定义如下:

Es (ru, i) 表示u位置处初始方向场块Θu与字典中第i个方向场块ru, i的相似度, 得:

ns表示原型方向场块和初始方向场块相对应元素之间差值小于阈值 (阈值为π/18) 的个数。Ec (ru, i) 表示u位置处字典中第i个方向场块与其左领域C (ru, i, ru左) 中重叠部分的兼容性, 计算方式[2]如 (6) 所示:

方向场块之间的兼容性[2]的定义为式 (7) :

其中, 表示2个不同方向场块中N0个重叠区域的方向元素。

为初始方向场块每个位置寻找局部能量最大值时所对应的候选方向场块, 以此代替初始方向场块, 最后将整个方向场重建出来, 并用高斯低通滤波, 是方向场更加连续, 寻找算法如下: (1) 将候选列表中的方向场块和初始方向场块最相似的候选方向场块作为第一个方向场块。 (2) 计算该位置处的能量值 (E) , 将其能量最大值所对应的候选方向场块代替其初始方向场块。 (3) 移到下一个位置, 重复上述步骤, 直至重建出整个指纹方向场。 (4) 滤波。经过上述寻找算法, 我们最终从各自的C个候选方向场块中为各自初始方向场块挑选出局部能量最大所对应的方向场块, 最后达到纠正低质量方向场的目的。

3 实验结果及分析

本文实验所用的指纹数据库都是来源于NIST SD27, NIST SD4, FVC2004以及清华大学提供的一个补充包, 其补充包可通过 (http://ivg.au.tsinghua.edu.cn/) 找到。为了平等的和STFT, FOMFE这2种算法进行对比, 本文直接使用了补充包里面提供的2个经过STFT和FOMFE算法估计的方向场, 并使用了误差评价函数和人工标注的方向场进行对比, 3种方法都是对SD27指纹库中的258张低质量指纹进行的估计, 其中包含高质量 (88) , 低质量 (85) 和丑陋 (85) 3种不同质量的指纹数据, 计算出均方误差, 如表1所示。

在本文中, 在构建局部字典上, 利用5*5邻域的方向场块, 快速地构建出包含更多方向场块的局部字典;最后在寻求最优解上, 本文提出了一种快速寻求最优解的方法, 而且也达到了较好的效果。表1展示了在NIST SD27指纹库中258张低质量指纹进行方向场重建, 并和2种有名算法 (STFT和FOMFE) 的对比实验结果数据, 结果表明本文的估计均方误差小于这2种算法, 所以本文对于低质量指纹方向场能够很好地重建出丢失的信息。

4 结论

本文在构建局部字典处, 使用5*5邻域的信息, 使得构建出来的字典更加具有多样性, 使其能够更好地包含低质量指纹方向场块;并在挑选候选方向场块时将原来的阈值0.8设置为0.4, 这使得候选方向场块列表中包含种类更多。我们统计的结果表明优化后的算法, 能够快速地重构成低质量方向场。

参考文献

[1]Feng J, Zhou J, Jain A K, Orientation fields estimation for latent fingerprint enhancement[J], IEEE Trans on Pattern Analysis Machine Intelligence, 2013 (4) :925-940.

[2]Yang X, Feng J J, Zhou J, Localized Dictionaries Based Orientation Field Estimation for Latent Fingerprint[J].IEEE transactions on pattern analysis and machine intelligence, 2014 (5) :955-969.

[3]Kaufman L and Rousseeuw P J, Finding Groups in Data:An Introduction to Cluster Analysis.New York NY, USA:Wiley-Inc, 1990.

[4]Kaufman L, Rousseeuw P J, Finding Groups in Data:An Introduction to Cluster Analysis[M].1th, New York NY, USA:Wiley-Inc, 2005.

[5]Bazen A M, Gerez S H, Systematic methods for the computation of the directional fields and singular points of fingerprints[J].IEEE Trans On Pattern Analysis and Machine, 2007 (7) :905-919.

[6]Chikkerur S, Cartwright A N, Govindaraju V, Fingerprint Enhancement Using STFT Analysis[J].Pattern Recognition, 2007 (1) :198-211.

[7]Wang Y, Hu J, Phillips D, A Fingerprint Orientation Model Based on 2D Fourier Expansion (FOMFE) and Its Application to Singular-Point Detection and Fingerprint Indexing[J].IEEE Trans On Pattern Analysis and Machine Intelligence, 2007, 29 (4) :573-585.

局部特征提取 篇4

1前言

众所周知, 地理数据库是地理信息系统的核心, Google Maps正是通过其庞大的地理数据库为WebGIS提供强有力的地理信息可视化及查询功能支持的。然而, 在一些特殊的路径规划、行车诱导应用中, 用户可能希望通过自己定义路径查找的方式 (如通过限定所经过的道路类型、旅行时间、两点间距离等特殊要求) 进行路径规划。由于Google Maps API仅提供内置的路径规划算法, 封装了对底层地理数据的访问, 导致了这类特殊应用无法实现。

为支持这类特殊应用, 本文在对Google Maps的路网深入研究的基础上, 通过探讨其开放的路径规划API对象及方法, 探索提取局部区域内交通地理数据的策略及算法, 为这类特殊应用提供地理数据支持。

2 Google Maps API常用对象

为支持JavaScript集成Google Maps地图, API为开发者提供了GMap2、GLargeMapControl、GMapTypeControl等对象及相应的方法实现对地图定位、漫游、缩放三种基本功能的支持, 而对于行车导航、本地搜索等服务, API提供了GDirections、GPolyLine、GMarker、GInfoWindow、GTitleLayerOverlay等对象, 以实现对两点间路径查找、路径显示、位置标示、黄页内容显示和用户新增图层的支持。图1展示了Google Maps API中常用对象, 以及它们所支持的常用方法和事件。

3局部交通道路提取策略分析

3.1道路的类型

要提取交通道路网中各路段数据, 首先必须对其道路类型进行分类。Google Maps中展现的道路类型 (特指公路) 包括了高速路、主干道、支道三类。其中高速路、主干道具有中间分车带或双实线标志, 每个车道只允许单向通行;支道分两类, 一类允许在同一车道下双向行驶, 另一类只允许单向行驶。根据这一特点, 我们按照车道的通行能力将Google Maps中的最小路段抽象划分为三类, 依次以doubleWay、twoWay和oneWay区分。

由于这三种类型的路段寻径方式有着很大的区别, 因此有必要首先进行区分。

3.2 Google Maps公路网寻径特点

通过对Google Maps公路网的深入分析, 我们发现:①Google Maps的路段之间是连通的 (即路网上任意两点 (v1’, v2’) 间总可以找到一条分别以这两点为起点和终点的路径) ;②从宏观上看 (如图2所示) , Google Maps地图被公路网划分成不同的封闭区域, 各个区域以环绕该区域的一条封闭的环形路径为边界;③从微观上看, 由于doubleWay类型道路的存在, 构成一个区域的封闭环形路径边界往往并不是一条单线, 而表现为由三种路段连接而成的边界 (如图3所示, 区域PL1的边界由pw1、pw2、pw3、pw4和pw5所构成, 其中pw2、pw3、pw4是doubleWay类型路段, 分别表现为一条环路) 。

事实上, oneWay和twoWay两类路段非常相似, 在路网中均表现为由多条边构成的一条单线 (如图中pw1、pw5、pw10) , 但寻径方式不同:而doubleWay类型路段分车道双向行驶, 因此在路网中表现为多条边构成的一条环路 (如图中pw2、pw3、pw4、pw6、pw7、pw8、pw9) , 以便与实际交通法规中沿行驶方向靠右车道行驶的交通规则相一致。

3.3 Directions对象的道路探索能力

仔细分析Google Maps所提供的用于路径搜索的函数可以发现, 路径搜索主要是通过Directions对象完成的。通过调用GDirections对象的LoadFromWaypoints () 方法可以实现路径查找功能, 并具有以下特点:①地图上任意选择的两点 (v1, v2) 总能定位到最近道路上的两点 (v1’, v2’) ;②在其”load”事件中总能返回由v01到v02所经过的路径polyLine。进一步, 通过得到的polyLine可以获取其包含的顶点数, 并可获得其中任意节点所在的经、纬度坐标。

从3.2节我们知道, 最小环形区域是构成整个公路路网的基本单位。因此, 我们从最小环形区域入手, 讨论如何根据Directons对象的寻径能力获取最小环形区域边界路段的方法。

作为路网的一种特例, 在最小环形区域内任意选择的两点 (v1和v2) 也能够通过Directions对象的函数获得其定位于该区域边界上的两点 (v1’和v2’) 。通过大量的研究发现:分别以v1’为起点v2’为终点和以v2’为起点v1’为终点进行路径查找 (以下简称为“正反路径查找”) , 分别得到的路径记为p1和p2, 将p1、p2合并为p, 则p仅可能存在以下3种情况:p为一条单线 (如图4中的pw1’和pw2’) ;p为一条单线与一条环线的组合 (如图4中的pw4’) ;p为一条环线 (如图4中的pw3和ppw4) 。进一步分析发现:

①当p为一条单线时, 观察p1和p2, 必为两条相同的路径。即若p1= , 则p2=, 且v1’和v2’必位于同一条twoWay路段或者由多条twoWay路段相连接的路径上。

②当p为一条单线与一条环线的组合时, 可以判断v1和v2必分别位于一条twoWay路段和一条doubleWay路段上。

③当p为一条环线时, 环线的类型可能有两类:一类是包含该最小环形区域的内部最小边界 (如图3中由pw1-pw5所有边界路段位于PL1内部的边构成的环路, 定义为ppw) , 另一类是仅包含环形边界上某一条或多条doubleWay类型的路段构成的环路 (如图3中pw2) 。以上两种情况可以通过v1与p的位置关系 (通过射线法可以判断出v1位于p的内部还是外部) 加以区分。

④若p为ppw时, v1’和v2’所在路段可能存在以下三种情况:a. v1’或v2’位于一条oneWay路段上;b. v1’和v2’分别位于两条doubleWay类型路段上;c. v1’和v2’分别位于相邻或不相邻的twoWay和doubleWay两条路段上。

⑤若p为pw类型环形路段时, v1’和v2’必分别位于两条可能相邻也可能不相邻的doubleWay路段上。

需要说明的是, 通过地图上任意两点 (v1, v2) 所定位到最近道路上的两点 (v1’, v2’) 不一定是路网中的节点 (比较图4和图5, v1’-v10’都不是路网中的原有的节点) , 因此需要在考虑提取策略和算法的时候剔除掉这些点, 以保证提取的道路与实际道路一致。

3.4局部交通道路提取策略

在探讨寻径策略和算法之前, 为便于区分, 需要对一些概念进行定义如下:

定义1:局部交通区域是指由用户所选定的一个需要提取交通道路网的区域, 记为PL, 如图3中所有交通路径所涵盖的区域。

定义2:最小环形区域是PL的最小组成单位, 即构成一个区域的封闭环形路径的所有边界所覆盖的区域, 记为PLn (n=1, 2, 3……) , 且PL=∪PLn (n=1, 2, ……) 。如图3中, PL1= pw1∪pw2∪pw3∪pw4∪pw5。

定义3:区域内部最小环路是指经过最小环形区域所有路段内部边界的一条环形路径, 记为ppw, 如图5中PL1所对应的ppw=。每个最小环形区域有且仅有一个ppw。

定义4:最小路段是组成最小环形区域边界的不能再分的完整路段, 记为pwn (n=1, 2, 3, …) 。最小路段有oneWay、twoWay和doubleWay三种类型, 如图3中pw1-pw10。

第3节的分析为我们获取局部交通地图提供了一条基本思路:①标记用户选择的矩形区域为PL;②将PL划分成最小环形区域 (即区域内不再包含更小的环形区域, 分别记为PL1, PL2, PL3……, 且PL=PL1∪PL2∪……。如图5中PL3=pw1∪pw2∪pw3∪pw4∪pw5) ;③通过限制搜索起点和搜索方向, 找到包含每个最小环形区域的内部最小环路ppw (如图5中对应于PL1的ppw= ) ;④依次通过对ppw中的相邻两点进行正反向路径查找, 获取各最小环形区域边界上的所有最小路段 (如PL3的边界路段由pw1、pw2、pw3、pw4、pw5构成, PL3=pw1∪pw2∪pw3∪pw4∪pw5) ;⑤尽管各最小环形区域内及各最小环形区域之间的最小路段可能会有重叠, 但消除重叠路段后, 可以获得包含在整个选择区域内的公路网。

4局部交通道路提取算法分析与设计

4.1提取算法分析

从以上分析可以看出, 实现局部交通道路提取策略的关键是如何通过任意定位在某个环形区域中的一点, 找到该环形区域中包含该环形区域内部边界的路径ppw。

但是, 在最小环形区域中选择的点 (记为v1) 的位置是随意的, 由于在自动处理的情况下事先并不知道环路的具体情况 (包括形状、边界路段类型等) , 因此需要从v1定位在边界路段上的点v1’开始, 通过适当的算法, 借助于该最小环形区域内的其它点, 逐步探索出边界路径, 最终合并为其ppw。

归纳3.3节的讨论结果, 可以得出如下结论:当且仅当v1’位于一条oneWay路段时, 可以通过在该区域内任意选择另一点v2, 通过正反向路径查找并合并找到的两条路径为p, 可以立即推断其获得的路径为ppw;在其他情况下, 通过另一点所在 (twoWay或doubleWay类型) 路段的类型无法直接推断所获取的路径必为ppw。

从整个Google Maps地图来看, 最小环形区域的形状和边界路段类型是多样的, 我们不能保证每一个最小环形区域都包含着oneWay类型路径, 因此无法通过仅在最小环形区域中选择两个特殊的点就一定能获取其ppw。通过反复实践, 我们发现单纯由twoWay类型路段构成的最小环形区域边界路径的提取是最复杂的, 其他类型的复合边界都可以通过在探索这种特殊情况的过程中获得其ppw。因此, 下面以图6为例, 讨论最小环形区域ppw的提取算法。

4.2提取算法设计

对仅由单线构成的最小环形区域而言, 其边界即为区域内部最小环路ppw。由于该环路仅由twoWay类型最小路段构成, 为无向环, 其寻径总是依据两点间距离最近的路径为基础的, 因此, 只需要在该ppw上找到任意两点 (记为v1’和v2’) , 找到其路径p0=;通过探索法让v1’或v2’沿边界路径延伸, 记录下新增长的路径p (如让v2’沿曲线v1’v2’的方向延伸, 延伸到新的点v3’, 得到p=) ;根据v1’和v3’获得新的路径p’, 判断p’是否包含原路径p0, 若包含, 继续延伸;否则, 一定能够找到不包含p0的一条新的路径p’, 该路径与p0和p、p’一起即可合并为ppw。

对不全由单线构成边界的ppw, 若环形区域边界上包含oneWay类型最小路段, 且在v1’或v2’延伸的过程中定位到oneWay类型路段上, 由4.1节的结论可知, 通过正反向查找即可获得该环路的ppw。

若v1’或v2’在延伸的过程中定位到doubleWay类型最小路段上, 由3.3节的结论可知, 通过正反向查找可能找到两种类型路径:ppw或者包含至少一条pw类型的环形边界, 通过射线法可以区分出这两类路径, 若不为ppw, 可以通过v1’和v2’的继续延伸最终获得ppw。

结合以上分析, 我们以图6所示的最小环形区域为例对提取算法说明如下:

①设环形区域内部任意确定一点v1, 通过Directions对象获取v1定位到最近边界路段上的点v1’, 求得的直线距离d及线段v1v1’的斜率k (若点v1和v1’具有相同的经度或纬度, 需要单独考虑, 以下只讨论一般情况) 。

②在以v1为起点, 沿v1’v1的垂线方向 (斜率为-1/k) 找到一个距离为d的点v2, 可以判断v2必在该环形区域内部 (因v1’是v1定位的最小距离点, 以v1为半径的圆上任一点都必在该环形区域内) 。

③以为两个顶点分别调用Directions对象进行正反向路径查找。设v2定位到边界路径上的位置为v2’, 查找获得的两条路径分别为p1’和p2’, 根据3.3节的结论, 可以判断出以v1’和v2’为端点的路径p的类型, 记n=0, p0<-p。若p为一条环路, v3’ <-v1’跳转至⑧。

④记n=n+1, 在v1’v2的正向寻找距v1距离为n×d的点, 并标记该点为v3。以 (v2, v3) 进行正反向路径查找, 合并所获得的两条路径p1’、p2’为新的路径p’, 记v3在边界上的定位点为v3’。

⑤在p’中找到距离v3’最近的点v4’, 求取直线v2v3与直线v4’v3’的交点v, 判断v2与v3是否在交点v的同侧, 即可判断出v3是否超出区域边界。若v3超出边界, 跳转至⑦;若p’为一条环路, p<-p’, 跳转至⑧;否则, 令p<-p’, 顺序执行。

⑥以 (v1, v3) 进行正反向查找, 合并所获得的两条新的路径p1’、p2’为新的路径p’, 若p’为一条环路, 令p<-p’, v2<-v1, v2’ <-v1’, 跳转至⑧;判断p’与p0、p的关系, 若p0与p、p’合并构成的路径p’’为不同路径, 则p0与p’’必构成一条环路, 该环路必为所要寻找的ppw, 记p0<-p0∪p’’, 算法终止;否则, p0<-p’’。置v1’ <-v2’, v2<-v3, v2’ <-v3’, 跳转执行④。

⑦求得的直线距离d1, 以v2为起点, 沿v1v2的垂线且与v2’v2不同侧方向, 寻找距离为d1的新点v, v也必在该环形区域内部。v1<-v2, v1’ <-v2, v2<-v, 跳转执行③。

⑧ (射线法) 在p中找到除端点v2’、v3’外且相邻的两点 (剔除非路网节点) , 选取线段v4’v5’的中点v6’。以为端点作一直线v1v6’, 求得v1v6’与p的另一交点v7’ (可能有多个交点, 只需要找到一个即可) 。判断v7’与v1是否在v6’的同侧, 即可判断v1在该环路p的内部还是外部。若v1在p的内部, 则该环路p必为ppw, 算法终止;否则该环路p必为一条或多条边界路径构成的边界环路 (如图7所示) , 跳转执行④。

5结束语

通过人工或者区域扩展的方法, 在局部区域内的各个最小环路上分别选取一点, 可以提取出各个最小环路, 进而通过合并、消除重复路径可以提取出该局部区域内的整个交通路径。限于篇幅, 本文没有给出如何判断重复路径和如何将各个最小环形区域合并为整个局部交通路径的细节。按照以上策略和算法, 对绵阳市高新区道路的提取结果显示, 该策略能够正确提取较小范围的交通路径, 并能够进一步获取各个交通路口、立交桥等特殊交通设施的交通寻径方式, 进而通过在Google Maps上叠加图层应用于具有特殊需求环境。但由于算法基于jsp编程, 程序通过Java虚拟机执行, 因此, 执行效率和存储空间受到较大的限制, 在提取较大范围交通路径时存在执行速度较慢的问题。

摘要：Google Maps API的开放策略将基于WebGIS的应用开发推向了一个新的高潮。针对具有特殊需求的路径规划、行车诱导等应用, 本文在深入分析路网特点及寻径方式的基础上, 探讨了如何根据Google Maps API提供的对象和方法提取局部交通道路数据的策略, 并针对一些特殊道路类型设计了提取算法。

关键词：Google Maps,交通道路数据,探索法,射线法

参考文献

[1]符海月, 赵军, 李春满.从GoogleMaps看我国全球化地理信息服务面临的挑战和对策[J].地理与地理信息科学, 2006, 22 (2)

[2]孙晓茹, 赵军.GoogleMapsAPI在WEBGIS中的应用[J].微计算机信息, 2006, 22 (19)

[3]Google地图API参考[EB/OL], http://code.google.com/intl/zh-CN/apis/maps/documentation/.

局部特征提取 篇5

随着电力系统的不断发展,电压等级越来越高,对大型电力设备的绝缘性能要求也越来越高。局部放电会逐渐腐蚀、损坏绝缘材料,使放电区域不断扩大,最终导致整个绝缘体击穿。电力设备内部早期故障产生的局部放电信号很微弱,而且现场存在着强烈的电磁干扰,局部放电信号往往被淹没在噪声中,对局部放电信号的检测和提取产生直接影响。因此,研究强噪声背景下的局部放电信号提取,对大型电力设备早期绝缘故障诊断有很重要的意义[1,2,3]。

局部放电的干扰信号按照其波形可分为窄带周期型干扰、平稳随机型干扰和随机脉冲干扰[4]。电力系统的载波通信、高频保护信号和无线电广播对监测的干扰属于窄带周期型干扰;各种设备以及监测线路中的白噪声属于平稳随机型干扰;输电线路和邻近电力设备的电晕放电和内部放电干扰信号属于脉冲型干扰。本文主要讨论平稳随机型干扰和窄带周期型干扰在局部放电信号提取中的影响。

目前研究局部放电信号提取的方法有小波阈值去噪法[5,6,7]、HHT方法[8,9,10]等,它们对局部放电信号的处理都能达到一定的效果,但是也都存在各自的问题。如小波的多分辨率特性能将信号在不同尺度下进行多分辨分解,但是小波变换的基函数的选取和最佳小波分层数的确定都存在问题,所以不能保证最优的分解效果。HHT方法可以分析非线性、非平稳的信号,具有良好的局部适应性,但这种方法会产生虚假分量,存在端点效应和模态混叠,且各分量的物理意义不是很明确。

本文提出将谱峭度SK(Spectral Kurtosis)[11,12]的方法应用到局部放电信号的提取中。谱峭度属于高阶统计量的范畴,是一个4阶累积量,理论上能够完全抑制高斯噪声,表征信号中的非平稳和非高斯信号,并且能够确定其在频带上的位置。谱峭度方法凭借其良好的统计特性已经在机械振动系统的故障诊断,如齿轮故障诊断[13]、轴承早期故障诊断[14,15]中得到广泛应用,并取得了很好的成果。本文通过求取谱峭度,根据其设计自适应最优带通滤波器,对局部放电信号进行滤波,从强噪声中提取出局部放电信号,并进行2~3层小波平滑去噪,以进一步提高局部放电信号特征的提取精度。通过仿真实验分析验证了本文方法的可行性和有效性。

1 谱峭度

1.1 谱峭度的定义

峭度是随机变量的4阶累计量,作为一个全局性指标不能反映特定信号分量的变化情况,因此,不适合强噪声环境下的状态检测问题。为了克服峭度在工程应用中的不足之处,Dwyer提出了谱峭度方法[16]。

考虑非平稳信号的Wold-Cramer分解,定义Y(t)为由信号X(t)激励的系统响应,则Y(t)可以表示为[17]:

其中,H(t,f)是系统的时变传递函数,可解释为信号Y(t)在频率f处的复包络。在实际系统中,H(t,f)是随机的,可表示为H(t,f,w),w表示滤波器时变性的随机变量。

S2n Y(t,f)为2n阶瞬时矩,是复包络能量的度量,定义为:

因此,当n分别取1、2时,谱峭度可定义为:

1.2 谱峭度的性质

谱峭度有很多性质[17],本文主要给出以下2个性质。

性质1一个纯的平稳过程Y(t)的谱峭度为:

特别地,一个纯的平稳高斯过程Y(t)的谱峭度为0。

性质2一个条件非平稳随机过程Z(t)=Y(t)+N(t),N(t)是添加的噪声信号,它与Y(t)相互独立。N(t)的谱峭度为:

其中,ρ(f)=S2N(f)/S2Y(f),为噪信比。根据式(5)可以得出,ρ(f)越大,谱峭度值越接近于0,而ρ(f)越小,谱峭度值近似等于KY(f),因此在整个频域进行搜索,可以找到谱峭度最大的频带。

2 基于谱峭度的局部放电去噪算法

本文根据局部放电信号的特征,提出一种基于谱峭度的局部放电信号提取新算法。该算法核心思想是首先求取出含噪局部放电信号的谱峭度,然后根据谱峭度设计自适应带通Wiener滤波器,再进行小波平滑去噪,提取出局部放电信号。

2.1 估计谱峭度

根据谱峭度的定义,本文采用短时傅里叶变换(STFT)的方法实现峭度谱的估计[18],如式(6)所示。

其中,Y(m)为一个随机过程;w(m)为窗函数,其长度为Nw。

定义Yw(u,f)的2n阶谱矩为:

其中,〈·〉l表示l阶时平均。因此,当n分别取1和2时,基于STFT的谱峭度的估计值为:

2.2 自适应最优带通滤波器设计

Wiener滤波器设计简单,性能优越,但是其参数是固定的。要设计这种滤波器,必须对信号和噪声的统计特性有先验知识,但是在实际中,常常无法预先知道这些统计特性,或者它们是随时间变化的,从而用Wiener滤波方法不能实现最优滤波。Wiener滤波器可表示为[18]:

根据1.2节中的谱峭度性质2,Wiener滤波器可以通过求取相应信号的谱峭度来实现,即:

其中,k为未知参数,其获取可以先假设一个比较小的值代入式(10),构造Wiener滤波器,再通过多次迭代,找到使滤波器滤波后相应输出信号的峭度值(可用MATLAB自带的峭度程序求取)的最大值,从而确定参数k[18]。

局部放电信号的窄带周期干扰可用一系列频率不同的正弦信号表示,平稳随机干扰可看作是白噪声。根据1.2节中的谱峭度性质1,理论上窄带周期干扰的谱峭度恒为-1,平稳随机干扰的谱峭度为0。但是基于STFT方法估计出的谱峭度值不是理论上的-1或者0,而是在它们附近波动。只根据式(10)设计出来的滤波器会包含其他频带的噪声,因此本文设置阈值σ来限定谱峭度的值。根据实际应用,σ可设置为:

其中,为信号谱峭度的最大值。谱峭度小于阈值σ可以认为是窄带周期干扰和白噪声信号,在设计滤波器时将其直接滤除。根据该方法设计出的滤波器是一个完全由数据驱动的自适应滤波器,不必知道信号的先验知识就能进行最优滤波。

2.3 小波平滑去噪

经过滤波过后的信号仍然存在与局部放电频率相同或相近的噪声信号,只需进一步通过2~3层的小波阈值平滑去噪,便可以得到更为精确的局部放电特征信号。

本文选用与局部放电信号波形匹配较好的db8小波进行分解[20],对经过自适应最优Wiener带通滤波器处理后的局放信号进行2~3层小波分解,采用无偏似然估计软阈值方法就可将局部放电信号提取出来。

综上所述,本文提出的基于谱峭度的局部放电信号提取算法的主要步骤如图1所示。

3 仿真验证

实际电力设备运行中,采集到的局部放电信号经常表现为指数衰减振荡形式,因此在理论分析中,通常可以采用以下2种数学模型来表示局部放电信号的形式[19]。

单指数衰减振荡型(图2(a)):

双指数衰减振荡型(图2(b)):

其中,A1、A2为信号幅值;τ1、τ2、τ3为衰减系数;fc为振荡频率。

在仿真中fc都取为1 MHz,τ1、τ2、τ3分别为2μs、2μs、4μs,A1、A2都为1 m V,采样频率为10 MHz。

3.1 抑制平稳随机噪声

平稳随机噪声可用白噪声模拟。理想局部放电信号见图3(a),加入白噪声(信噪比-2.029 0 d B)后的信号如图3(b)所示。

本文基于STFT方法求取含噪信号的谱峭度。大量的试验表明在窗函数长度确定的情况下,窗函数的类型(如汉明窗、海宁窗等)对谱峭度的估计值影响不大;但窗函数的类型确定时,窗函数的长度对求取出来的谱峭度影响比较大。这是由STFT本身决定的,窗函数太长不能保证时间分辨率,窗函数太短不能保证频率分辨率,只能在时间分辨率和频率分辨率之间折中处理。所以本文考虑到局部放电信号本身的频率特征,选择的窗函数的类型为汉明窗,窗函数每次移动一个数据点,讨论窗函数长度为分别为25、45、85、125,求取含噪局部放电信号的谱峭度的优劣性,如图4所示。

窗函数长度为25时,谱峭度频率分辨率不高;窗函数长度为125时,谱峭度频率分辨率达到要求,但是其中包含很多噪声成分;窗函数长度为45和85时,频率分辨率达到要求,其中的噪声成分比较少。从而对含噪局部放电信号,选取窗函数长度为45~85。本文选择窗函数的类型为汉明窗,每次移动1个数据点,窗函数长度为55,求出的谱峭度如图5所示。

由局部放电信号的谱峭度可以确定局部放电信号的频率,而且加入的白噪声的谱峭度几乎为0,从而就可以通过局部放电信号和噪声的谱峭度大小来设计滤波器。

根据式(10)、(11)可以设计Wiener最优带通滤波器,经过试验k的值使原信号谱峭度最大时,输出信号的峭度值最大,从而可得出最优带通滤波器,如图6所示,滤波后的信号如图7所示。

由图7可以看出,通过Wiener滤波器滤去了大部分噪声,能够清晰地看出局部放电信号的特征,但是还存在一些和局部放电信号频率相同或者相近的噪声信号。因此再进行小波平滑去噪处理,选用与本算例局部放电信号波形匹配较好的db8小波进行3层分解,就可以提取出光滑的局部放电信号。

分别用本文的方法和小波方法对局部放电信号进行去噪的结果比较如图8所示,小波方法的小波基选取的是与本算例局部放电信号波形匹配较好的db8小波,分解层数为8层。

原始信噪比为-2.029 0 d B的含噪局部放电信号去噪后的性能指标[20]为消噪后的信噪比、均方根误差和相关系数。在表1中列出本文方法进行3层分解和db8小波进行3层、8层分解的去噪性能指标。

从图8和表1中可见,在原始信噪比为-2.0290 d B时,本文方法对平稳随机型干扰有很好的消噪效果,不仅大幅提高了信号的信噪比,而且失真度小,与理想局部放电信号的波形相似度最高。通过消噪后的信噪比、均方根误差和相关系数这几个指标的比较,明显看出本文方法比db8小波的各方面性能都要好。

为了体现本文提出方法在低信噪比条件下的处理性能,在信噪比为-7.926 1 d B的情况下与传统小波去噪进行了对比分析,结果如图9所示,各个性能指标对比在表2中列出。

通过在不同信噪比条件下,对本文提出的方法和db8小波阈值去噪性能的比较,可以得出db8小波只有在信噪比比较高且分解层数较大的情况下,对局部放电信号的去噪效果较好,然而分解层数越多,去噪后将丢失更多的原始信息,随着信噪比降低,db8小波的去噪性能明显下降,不能准确地提取出局部放电信号;而本文提出的方法先通过谱峭度设计出来的自适应带通滤波器,提高了信噪比,然后进行3层小波平滑去噪,分解层数少,保留了更多的原始信息,在信噪比较小时也能很好地对局部放电信号进行去噪,对局部放电信号的波形特征提取效果更好。

3.2 抑制混有白噪声的窄带周期干扰

窄带周期干扰可以通过一组频率不同的正弦信号表示,其频率fc1、fc2、fc3、fc4分别设为150 k Hz、500 k Hz、1.75 MHz、2.5 MHz,干扰信号的数学表达式为:

理想局部放电信号中加入窄带周期干扰和白噪声(白噪声是由MATLAB中自带常用白噪声功能函数产生),原始信噪比为-12.7138d B,如图10所示。

窗函数的类型为汉明窗,每次移动一个数据点,窗函数长度为55,基于STFT谱峭度方法求出的谱峭度如图11(a)所示。

根据图11(a)可以得出窄带周期干扰的谱峭度值接近-1,而白噪声的谱峭度值在0附近。利用本文提出自适应带通滤波器的设计方法,可设计出最优带通Wiener滤波器,如图11(b)所示。

对含噪局部放电信号进行滤波,再进行小波平滑去噪,得到去噪后的局部放电信号,如图12(a)所示;图12(b)为直接用小波去噪后的局部放电信号。

原始信噪比为-12.7138 d B的含噪局部放电信号去噪后的各个性能指标在表3中列出。

从图11和表3中可以看出,本文方法能对窄带周期型干扰有很好的抑制效果,不仅大幅提高了信噪比,而且失真度小,与理想局部放电信号的波形相似度最大。通过消噪后的信噪比、均方根误差和相关系数这几个指标的比较,明显看出本文方法比db8小波方法的各方面性能都要好。

4 结论

局部特征提取 篇6

齿轮传动是一种常见的运动和动力传递方式,广泛地应用于机械设备中。齿轮故障会导致机械设备的整体性能下降,甚至引起严重设备事故从而造成重大经济损失[1,2],因此,研究齿轮故障诊断方法具有实际意义。而齿轮的振动信号蕴含了齿轮的运行状态信息,所以,通常通过振动信号分析进行齿轮故障诊断。

齿轮故障的振动信号一般具有非线性、非平稳的特性,为了准确地将信号分解,进而得到分量的局部特征,许多学者开展了时频分析技术的研究[3,4,5]。其中,经验模态分解(empirical mode decomposition,EMD)是一个研究的热点。作为一种有效的自适应时频分析方法,EMD与Hilbert变换(Hilbert transform,HT)相结合,被广泛地应用于齿轮故障诊断领域[6,7,8]。但EMD在使用过程中存在一些问题,如过包络、欠包络、端点效应、模态混叠等[6,7,8]。因此,国内外学者一直在寻求EMD的改进、替代方法。

程军圣等[9]在本征时间尺度分解(intrinsic time-scale decomposition,ITD)[10]的基础上创造性地提出了局部特征尺度分解(local characteristic-scale decomposition,LCD)方法,该方法在保证分解所得分量物理意义清晰的基础上,具有计算时间短、端点效应不明显等优势,开拓了自适应时频分析方法的新思路[11,12]。

齿轮故障振动信号可以看成多个调幅-调频信号的叠加,通过LCD可将信号分解成单个调幅-调频分量之和,对每个分量进行解调再进行故障诊断是一种有效的方法。常用的解调方法有Hilbert解调和经验调幅调频分解(empirical AM-FM decomposition,EAD)。Hilbert法的不足是:包络本身光滑性较差,且端点效应明显,求取信号瞬时频率时会出现负频率现象[13]。EAD法可避免Hilbert变换引起的端点效应及负频率,但EAD法只是一种经验方法,缺乏严密的理论推导[14]。因此,研究新的解调方法具有一定的工程意义。

为了准确分析齿轮故障振动信号蕴含的故障特征,本文提出了一种新的解调方法———局部均值解调法,并结合LCD进行故障诊断。

1 基本理论

1.1 LCD方法

LCD方法将一个信号x(t)分解成若干个内禀尺度分量(intrinsic scale components,ISC)之和[10?12]。每个ISC分量必须满足2个条件:①任意2个相邻的极值点符号互异;②考虑所有的极值点Xk及对应的时刻τk(k=1,2,…,M,其中M为极值点的个数)。

取两个相距最近且符号互异的极值点(τk,Xk)和(τk+2,Xk+2),按照下式定义一个τk+1时刻的函数值:

则Ak+1与极值点Xk+1的比值关系不变,即

式中,α∈(0,1)为一常量,典型地,α=0.5。

为了筛选出ISC分量,LCD方法构造了信号的均值曲线,将均值曲线不断地从原信号中分离,其迭代过程与EMD方法类似,不同之处在于均值曲线的构造方法。LCD方法用均值点构造均值曲线,均值点的定义为

式中,Lk+1为均值点,在迭代结束时数值为0。

标准LCD算法采用分段线性方法由均值点Lk计算均值曲线。为了得到更光滑的ISC分量,本文采用三次样条方法替代分段线性方法。

1.2 局部均值解调法

文献[14]提出了EAD法,指出该方法大多数情况下较Hilbert法有优势,笔者通过理论分析及数值仿真发现,在ISC分量的解调中,有必要对EAD法进行改进。

EAD法的基本思路是:将信号的所有极值点取绝对值,用三次样条函数插值形成包络函数,用原始信号除以包络函数得到标准化信号;若标准化信号为调频分量,则分解结束,否则将标准化信号作为原始信号,重复分解过程。通常情况下,迭代2~3次,标准化分解就会结束,信号可以表示为调幅分量(包络)和调频分量的乘积。

在上述过程中,有两点需要注意:

(1)EAD法认为信号是零均值对称的,所以将所有极值点取绝对值,进而求包络函数。但在实际应用中,LCD分解得到的分量往往不能满足EAD法的假设,只能“局部”满足“零均值对称”,这一点从均值点的定义及均值曲线的构造过程可以看出。为了从数值上说明这个问题,图1给出了调幅-调频信号x(t)的均值曲线。x(t)为

由图1可以看出,LCD方法定义的均值曲线是零均值的近似,在局部范围内与零值存在较大的偏差,按照均值曲线迭代得到的ISC分量并非是零均值对称的,而EAD法按照全局零均值处理求包络理论上存在误差,利用局部信息则能减小其影响。因此,在构造包络函数时,考虑用“局部均值曲线”代替“全局均值曲线”。

(2)EAD法采用三次样条插值方法求包络函数,存在“过包络”、“欠包络”的问题,从而导致解调精度下降。为了避免这种情况,考虑用滑动平均法代替三次样条插值。

基于以上两点考虑,参考局域均值分解(local mean decomposition,LMD)[15]的思想,提出了一种局部均值解调(local mean envelope,LME)法,将信号的调幅部分和调频部分分离。与EAD法相比,该方法在包络函数的构造方法上作了改进,其余步骤相同。包络函数的构造方法为:计算信号的局部均值点,采用滑动平均法求局部均值曲线,用原信号减去局部均值曲线,得到零均值对称曲线,用滑动平均法计算包络函数。计算过程如下。

(1)确定第j个ISC分量ISCj(t)的所有局部极值点ni及其时刻tn i,计算均值点mi、局部幅值ai:

(2)用滑动平均法得到均值函数m11(t)、包络估计函数a11(t)。

(3)将局部均值函数m11(t)从原始信号ISCj(t)中分离出来,得到

(4)用h11(t)除以包络估计函数对h11(t)进行解调,得到

(5)重复步骤(1)~(4)k次,直至得到一个纯调频信号s1k(t)(瞬时幅值均小于1)及k个包络信号a1k(t),按照下式计算包络信号:

从LME解调法的流程可以看出,其实质就是将一个ISC分量视为一个乘积函数(product function,PF)分量和一个剩余信号之和,然后按LMD方法迭代的第一步求出一个纯调频信号和一个包络信号,剩余信号视为计算误差舍去。ISC分量与PF分量的定义内涵类似,剩余信号很小,可作为误差舍去。

1.3 瞬时频率计算

对纯调频信号s1k(t)利用反正切函数计算相位:

将φ(t)展开并求导,可得到瞬时频率。

求导计算会导致局部极值点附近的瞬时频率出现畸点,可对瞬时频率作平滑或滤波处理。本文采用滑动平均法,对任何一个瞬时频率的计算结果取连续3个采样点的平均值。

1.4 基于瞬时频率频谱的齿轮故障诊断

当齿轮存在故障时,其振动信号会出现调幅、调频现象。忽略传递函数的影响,齿轮故障振动信号可表示为[8]:

式中,fr为轴的转频;z为齿轮的齿数;Xm为第m阶啮合频率谐波分量的幅值;φm为第m阶啮合频率谐波分量的初相位;dm和bm(t)为第m阶啮合频率谐波分量的幅值和相位调制函数,两者都是以fr及其倍频为重复频率的周期函数。

根据瞬时频率fm(t)的定义:

可以看出fm(t)可反映bm(t)的信息。

应用LCD方法对y(t)进行分解,对ISC分量求出瞬时频率fm(t)。fm(t)是一个以齿轮转轴转频及其倍频为中心的频率分量。对fm(t)进行频谱分析进而得到瞬时频率谱,由瞬时频率谱可以直观地判断fr及其倍频是否存在,从而进行故障诊断。

2 仿真数据分析

考察下式所示的仿真信号:

其中,x(t)为仿真信号,由调幅调频信号x1(t)、正弦信号x2(t)合成;采样频率为fs=1024Hz,仿真时间t∈[0,1]。对仿真信号进行LCD,得到2个ISC分量和一个剩余信号r(t),仿真信号及分解结果见图2。

由图2可以看出,LCD能将调频调幅信号、正弦信号分解出来,剩余信号幅值很小,具有良好的分解能力。

为了比较Hilbert法、EAD法和LME法的效果,分别采用三种方法对两个ISC分量进行包络解调,将瞬时幅值直接给出。对EAD法和LME法采用反正切法计算瞬时频率,三种方法的瞬时频率计算结果均采用三点平滑处理。

图3、图4所示是两个ISC分量瞬时幅值的计算结果,可以看出三种方法解调结果均能反映原始信号的变化趋势,其中,Hilbert法计算结果波动最大且端点处的计算值严重偏离理论值;EAD法效果优于Hilbert法;LME法求得的瞬时幅值最贴近理论值,端点效应最小,更符合原始信号的实际特征。

图5、图6所示是两个ISC分量瞬时频率的计算结果,可以看出不同方法计算结果差别较大:①在ISC1(t)的计算结果中,Hilbert法端点效应最明显,两端的计算误差已超过实际频率;EAD法的端点效应较Hilbert法有明显改善,但中间段的精确度有所下降;LME法计算结果最好;②在ISC2(t)的计算结果中,Hilbert法端点效应最小,EAD法端点效应最大;而中间段的计算效果,LME法最好,Hilbert法与EAD法相当。

总的来说,在瞬时幅值的计算方面,三种方法计算误差均不大,LME法优于EAD法,Hilbert法效果最差;在瞬时频率的计算方面,计算误差均较大,从端点和中间段的综合效果来看,LME法最优,Hilbert法最差。

3 试验数据分析

齿轮故障试验数据源于QPZZ-Ⅱ旋转机械振动分析及故障诊断试验平台系统,其齿轮箱传动结构如图7所示。变频调速电机通过联轴节驱动小齿轮,大齿轮与小齿轮直接啮合。大小齿轮均为圆柱齿轮,大齿轮齿数为75,小齿轮齿数为55。

人为地将小齿轮的一个齿尖切割约5mm,模拟小齿轮断齿中度故障。试验中,大齿轮输出轴负载为零。设置电机轴转速为880r/min,实测转速878r/min,则小齿轮的转频f1=14.6Hz,大齿轮的转频f2=10.7Hz,齿轮啮合频率为804.8Hz。采用加速度传感器采集振动信号,传感器在输出轴电机侧轴承处垂直于齿轮箱上表面安装,信号采样频率为fs=5120Hz,计算数据点N=1024。

图8给出了原始信号的时域波形图及包络谱。从时域波形图上可以看出原始信号具有调幅-调频特性,并且存在周期性冲击。为了验证LCD本身的能力,本文没有采用中值滤波、SVD等降噪技术,包络谱直接采用Hilbert法求出,为了观察特征频段的特征,只给出了低频段(0~1000Hz)的波形。从包络谱上可以看出20Hz、45Hz处有明显谱峰。20Hz可近似认为是大齿轮转频的2倍频,45Hz可近似认为是小齿轮转频的3倍频,由此可以判定齿轮箱出现了故障,但不能判定哪个齿轮出现了故障。

对原始采样信号进行LCD,得到4个ISC分量和1个剩余信号。图9是4个ISC分量的时域波形。可以看出:①LCD类似于自适应滤波器,4个ISC分量所包含的频率段逐渐降低;②每个ISC分量的幅值依次减小;③每个ISC分量均可近似认为是调幅-调频分量。因此,LCD是一种有效的信号自适应分解方法。

对4个ISC分量分别应用三种方法解调,并求瞬时频率谱,将低频段(0~1000Hz)的结果列于图10~图12。可以看出:①在分量1的瞬时频率谱中,小齿轮转频3倍频处存在相对明显的谱线;②在分量2、3的瞬时频率谱中,小齿轮转频2倍频处存在相对明显的谱线;③分量4的瞬时频率谱能清晰反映小齿轮转频及其2倍频。因此,LCD结合三种解调方法求出的瞬时频率谱均能得出小齿轮存在故障的正确结论。

但是Hilbert法求得的瞬时频率谱中存在明显的高频能量泄漏现象,其他频段的谱线杂乱,干扰较大;在EAD法求得的瞬时频率谱中,高频段的干扰得到明显的抑制;而LME法求得的瞬时频率谱高频干扰最小,信噪比最高。因此,LME法在齿轮故障诊断中具有一定的优越性。

4 结论

基于局部特征融合的人脸识别研究 篇7

人脸识别特指利用分析比较人脸视觉特征信息进行身份鉴别的计算机技术。目前对于人脸识别的方法, 用的较为成熟和广泛的方法是特征脸方法。特征脸方法是基于KL变换的人脸识别方法, KL变换是图像压缩的一种最优正交变换, 完全是基于图像灰度的统计特性的。

2 基于信息融合的人脸识别研究现状分析

过去传统的人脸识别是基于全局统计特征而开发的算法, 虽然具有一定的识别率, 然而基于全局统计特征的方法对于人脸表情、遮挡等局部变化是却很容易受到较大的影响。后来人们尝试将人脸分解表示, 提出基于局部特征的人脸表示方法, 例如局部特征分析算法、弹性图匹配算法、局部二值模式人脸识别方法、基于SIFT算子的人脸识别方法等等。

在上述的基于局部特征实现人脸识别的算法中, 基于SIFT算子的人脸识别方法获得了巨大成功, 继SIFT算子的成功应用, 人们关注到, SIFT在前期做了大量的工作如在尺度空间进行特征检测, 确定关键点 (Keypoints) , 并确定关键点的主方向, 最后才用关键点临近区域的统计梯度方向直方图来表示该关键点。为了简化计算, 各种优化算法开始得到大量研究与应用, 其中以融合多维局部特征信息的HOG特征描述子识别率最高, 2008年, Monzo等人借鉴基于Gabor小波的EBGM (Gabor-EBGM) 在人脸识别中的成功, 他们将EBGM算法中用Gabor小波提取局部特征的方法使用HOG特征代替, 并称为HOG-EBGM算法, 他们的实验指出HOG-EBGM相比Gabor-EBGM有更好的性能表现, 从而人脸识别算法研究进入到了新的多维信息融合的阶段。

3 基于局部特征融合的人脸识别研究

3.1 基于信息融合的Gabor特征描述子基本原理

由于Gabor特征对光照、姿态具有一定的鲁棒性, 因此Gabor特征在人脸识别领域的应用非常广泛, 目前已经成为了最为主流的人脸特征抽取方法。Gabor特征具有良好的空间局部性和方向选择性, 可以很好的描述图像的纹理信息, 因此许多人脸识别算法都采用Gabor特征作为识别特征。

提取单幅人脸图像的Gabor特征的基本步骤如下:

(1) 划分Gabor特征的采样点

一般而言, 对于被测图像, Gabor特征采样点划分的越多越好, 采样点划分越多, 表明有可能捕捉到的特征数越多, 但是对于人脸识别图像而言, 由于人脸的特征具有相似性, 因此划分过多的Gabor特征采样点反而会增加识别的计算量, 因此必须要在划分Gabor特征采样点后对采样点进行降维处理。

(2) 提取人脸图像的Gabor特征

对每一个经过降维之后的采样点, 提取在采样点上提取各个滤波器输出的值, 联合起来组成列向量, 从而得到每幅图像的Gabor特征向量, 并以此作为该人脸图像的局部特征向量。

与原图灰度值随位置的变换相比, 滤波器输出的幅值随位置的变化要更小更不敏感;Gabor相位信息随着空间位置呈周期性变化, 而幅值的变化相对平滑而稳定。因此, 即使采样点稍有偏移也不会使提取出的特征值发生太大变化;也就是说, Gabor特征具有良好的空间局部性和方向选择性, 而且对光照、姿态具有一定的鲁棒性, 可以容忍更大的视觉定位和对齐误差。

3.2 基于局部特征融合的人脸识别优化算法

(1) 采样点划分

这里以一幅40×30的点阵图像为例, 划分1200个采样点。

(2) 采样点降维

一个M×N的人脸图像经Gabor变换后就得到幅值特征和相位特征的维数都是40×M×N, 计算量很大, 而且Gabor特征在相邻像素间是高度相关的, 所以必须对Gabor特征进行降维。本论文所采用的降维方法就是对Gabor变换系数进行简单的下采样, 即对Gabor特征进行固定行固定列均匀分布采样, 从而实现对采样点的降维处理。

(3) 提取局部特征

假定全部的特征量都采用三角形隶属度函数提取特征, Xlj∈[0, 1], 权值∈[-1, 1]。采用的误差信号为:

提取特征的目的是使E最小, 但由于视觉识别系统的特殊性, 即使所有Yi=Ti, E也不会为零, 因此需要特殊的停止迭代的规则。这种规则可由如下方法得到。

设T1[0]=[t11, t12], 如果对于所有i都有Yi=Tii, 那么指定E[0]=[?λ, λ], λ由下式决定

给定一个可接受的偏差ε, 如果E落入区间[?λ?ε, λ+ε], 则停止迭代。

(4) 构建基于Gabor描述子的多维局部特征融合算子

二维Gabor小波的定义为:

式中:σ是与小波频率带宽有关的常数;Z= (X, Y) 为空间位置坐标;k确定了Gabor内核的方向和尺度。在采用8个方向和5个尺度的采样时, 方向和尺度上的k可以写为:

Kmax为最大频率, f是频域中的内核间隔因子。

令参数、、, 可以获得较好的小波表征和辨别效果。

下面构建基于Gabor描述子的多维特征融合算子, 对Gabor变换与Gabor内核做卷积, 可以得到:

在提取人脸图像的Gabor特征时, 通常采用多个在不同尺度和方向上的Gabor滤波器组成滤波器组, 并根据图像的特点和神经生理学的结论来选择参数。通常研究采用共包括8个方向 (n=8;u=0, 1, ……7) 和5个尺度 (Kmax=π/2;f=2;v=0, 1, 2, 3, 4) 的Gabor滤波器组, 并令σ=π, 使滤波器的带宽约为1倍频程。

结语

信息融合技术是一门集传感器技术、计算机技术、通信技术、信息处理技术和模式识别技术等于一起的多学科交叉技术, 它需要很丰富的理论知识和工程先验知识。目前, 我国的多传感器信息融合技术基本还处于理论研究阶段, 只研制出少量具有初步融合功能的多传感器系统, 尤其是应用在生物识别系统方面, 可应用的成果十分少。本论文结合局部特征融合探讨了人脸识别的算法, 并对具体的实现给出了具体的执行算法。

参考文献

[1]杨杰.基于数据融合和多重智能模型的目标识别和跟踪[J].高技术通讯, 1999, (4) :10-14.