查 翔倪世宏 張 鵬

(空軍工程大學航空航天工程學院 西安 710038)

一類非線性信號去噪的奇異值分解有效迭代方法

查 翔*倪世宏 張 鵬

(空軍工程大學航空航天工程學院 西安 710038)

對于一類非線性信號的去噪問題，該文提出一種基于奇異值分解(Singular Value Decomposition, SVD)的有效迭代方法。對現(xiàn)有奇異值差分譜方法在兩類不同非線性信號上的去噪效果進行了對比，指出在信號不具有明顯特征頻率、非周期性變化時這一方法并不適用，并分析了現(xiàn)象產生的原因；然后針對該類信號的特點重新定義了Hankel矩陣結構，給出有效奇異值的確定方式，并通過SVD多次迭代過程實現(xiàn)對該類信號的有效去噪。對實際飛行數(shù)據(jù)去噪的實驗結果表明，該方法對提出的一類信號對象不僅去噪效果良好，而且可提高運算效率。

信號處理；去噪；奇異值分解；特征頻率； Hankel矩陣

1 引言

在工程測試中受系統(tǒng)本身和環(huán)境干擾的影響，現(xiàn)場采集的機械設備狀態(tài)信號往往含有各種噪聲，若能及時有效地將噪聲從有用信號中消除，還原設備運行的真實狀態(tài)，對于進一步的信號分析處理以及設備狀態(tài)監(jiān)控、故障診斷等具有重要意義。

當前業(yè)界流行的去噪方法主要涉及非線性濾波[1]、小波閾值法[2?4]以及奇異值分解(Singular Value Decomposition, SVD)[5]等。非線性濾波去噪效果的優(yōu)劣取決于濾波器結構是否設計合理，而且過程存在一定的時間延遲。小波閾值法過度依賴小波閾值以及分解層數(shù)的選擇，不同選擇方式下的去噪效果存在較大差異。SVD去噪的實質是線性加權式的正交化分解，得到的去噪信號不存在相位偏移和時間延遲，在1維振動信號[6]以及2維圖像去噪[7]、時間序列預測[8]中得到廣泛應用。

基于SVD的去噪方法主要解決兩個方面的問題：Hankel矩陣結構的確定以及有效奇異值的選擇。對Hankel矩陣結構的確定問題常被忽略，由于無有效理論依據(jù)，對其奇異值分解前的結構一般憑經(jīng)驗事先確定[9]。另一個關鍵問題是有效奇異值的選擇，SVD的基本思想就是通過有效奇異值對矩陣重構，從而實現(xiàn)信號去噪的目的，通常用信噪比經(jīng)驗值的估計選擇有效奇異值的閾值[10,11]，但在工程中難以實現(xiàn)。文獻[12]提出一種非常有代表性的方法：基于奇異值差分譜理論的SVD去噪，根據(jù)差分譜的最大峰值確定有效奇異值的個數(shù)，從而實現(xiàn)有效奇異值的自適應選擇。然而，對于一類無明顯特征頻率，不具有明顯周期性特點的信號而言，采用以上矩陣結構方式以及有效奇異值的選擇方法均不能有效實現(xiàn)信號去噪。主要表現(xiàn)在去噪后的信號只能描述出原始信號的大尺度變化趨勢，細節(jié)信息丟失過多，以致引起波形失真。文獻[13]所針對的具有明顯變化趨勢的信號對象亦可劃為這一類信號的范疇。

本文提出一種針對一類非線性信號去噪的奇異值分解有效迭代方法，處理后的信號與原始目標波形非常接近，細節(jié)特征基本被完好保留，能夠實現(xiàn)信號的有效去噪。

2 基于SVD的去噪原理

設1維離散信號x={x(i)}, i=1,2,…,N。為利用SVD對該信號進行去噪處理，首先要對其進行相空間重構，采用滑動窗口依次截取的方式獲得多個行向量，從而構造出Hankel矩陣：其中n為滑動窗口的寬度，且1＜n＜N。若令m= N?n+1，則A∈Rm×n。對A作奇異值分解，即必存在正交矩陣U∈Rm×m和正交矩陣V∈Rn×n，使得

其中對角陣Σ∈Rm×n, Σ中的各對角線元素σ1, σ2,…,σq稱為A的奇異值，q=min(m,n)≥2，并滿足σ1≥σ2≥…≥σq＞0，實際中常取m≥n。x對應的Hankel矩陣中任意兩行完全不相關，為滿秩矩陣，并且奇異值序列具有的特點為

式(3)中的“?”表示前k個奇異值遠大于之后的q?k個奇異值，在第k個與第k+1個奇異值之間存在較明顯的跳變。由于每個奇異值對應著一個分量信號，前k個奇異值代表著有用信號，而后q?k個奇異值則對應著噪聲。因此，可以選擇前k個奇異值作為有效奇異值，而將其它剩余奇異值置0，就可以依據(jù)Frobenius范數(shù)意義下的矩陣最佳逼近定理實現(xiàn)逆過程運算，以得到重構后的去噪信號。

信號去噪效果的優(yōu)劣與有效奇異值的選取有關?，F(xiàn)有基于奇異值差分譜的選取方法可較好描述奇異值序列的突變規(guī)律[12]，差分譜的定義為

其中k=1,2,…,q ?1，則序列c=(c1,c2,…,cq?1)構成了奇異值差分譜，若相鄰兩奇異值差值較大，則必會在差分譜序列中產生一個最大峰值，其對應式(3)中的突變位置k，從而實現(xiàn)有效奇異值的自適應選擇。

Hankel矩陣的不同結構直接影響到q的取值以及去噪效果。文獻[14]指出Hankel矩陣在最為接近方陣即q值最大時，濾波效果最好，具體結構為

其中rem(?,?)表示取余運算，這也是基于奇異值差分譜的SVD方法常采用的形式。

3 基于SVD的一類非線性信號去噪特點描述

以式(5)所示的Hankel矩陣結構為前提，基于奇異值差分譜的SVD去噪方法對于一類具有特殊屬性的非線性信號而言并不適用。為更好地說明問題，首先設計一個仿真信號x(t)=sin(4t)+sin(8t)，在t∈[0,2π]區(qū)間內共采樣512個點，同時加入信噪比約為2.24的白噪聲，將奇異值序列以及奇異值差分譜同時置于圖1(a)中(為能清楚顯示只繪出了前50個點)?？梢钥吹?，奇異值序列的前4個值相對較大，在第5個位置處序列出現(xiàn)陡降；奇異值差分譜的最大峰值位于第4個坐標，k應取為4。按照差分譜選擇奇異值的原則，取前4個奇異值重構后的信號如圖1(b)所示。對比去噪前后的信號，奇異值差分譜方法有效去除了噪聲，波形足夠平滑，并且能反映出x(t)的構造特點，即由兩個不同特征頻率的周期正弦信號疊加組成，且仍為周期信號。圖2給出了x(t)經(jīng)快速傅里葉變換得到的幅值譜，雖然存在噪聲頻率的干擾，但仍可發(fā)現(xiàn)大約在0.79 Hz f=以及f=1.57 Hz 處存在明顯的兩個尖峰，這反映了信號的兩個主特征頻率，也是信號能量的主要分布頻率。

然而，實際中的信號并不一定具有特征頻率或周期特性，此時奇異值差分譜方法的去噪效果可能很差。以典型的Block信號為例(如圖3(a)所示)，按式(5)對其構造Hankel矩陣并進行奇異值分解，得到奇異值與奇異值差分譜序列，取前50個位置點的結果如圖3(b)所示?？梢钥吹?，奇異值差分譜的最大峰值位于第1個位置，說明原始信號包含較強的直流成分，此時應以第2大峰值位置確定k，易知k=3。選擇前3個奇異值進行重構，恢復出的信號如圖3(c)所示。不難看出，按奇異值差分譜方法得到的去噪信號丟失了大部分的波形細節(jié)，信號被過度平滑且出現(xiàn)波形失真，僅反映出大時間尺度上的整體趨勢。

圖1 仿真信號波形及SVD處理結果

圖2 仿真信號的幅值譜

圖3 Block信號波形及SVD處理結果

觀察仿真信號與Block信號的時域波形，雖然混有一定程度的噪聲，但仍能看出仿真信號帶有明顯的周期性，而Block信號卻不具有這一特征，主體變化不存在明顯規(guī)律，帶有較強的隨意性。為使結論更可靠，從頻域角度分析兩種信號的頻率分布情況，給出Block信號的幅值譜如圖4所示。根據(jù)噪聲能量的分布特點，Block信號中的噪聲在整個頻域內分布較為均勻，而有用成分主要集中在低頻段，說明Block信號除噪聲外主要為低頻的直流成分，但與文獻[14]討論的直流成分有明顯區(qū)別。后者研究對象為直流、交流和噪聲3種成分疊加在一起的的混合信號，其所指的直流成分是相對于交流而言的，實際為帶有偏置特征的常數(shù)。從這個意義上來說，之前所設計的仿真信號可看作是由交流和噪聲成分組成的混合信號。為與文獻[14]中常數(shù)直流對象以示區(qū)別，本文將Block信號中這種低頻直流成分稱為廣義直流成分，其顯著特點是變化相對平緩，帶有一定的隨機性，并且無明顯周期規(guī)律。

4 針對一類非線性信號的SVD有效迭代去噪

4.1 方法引入

通過對比兩類信號的去噪效果可知，實際中一類無代表性特征頻率的含噪信號主要包含了廣義直流成分與噪聲，對于這類信號，基于奇異值差分譜的SVD方法無法實現(xiàn)有效去噪。具體原因可從以下兩個方面分析：

(1)Hankel矩陣結構的合理性 一般認為，奇異值數(shù)量q越大，奇異值序列的分布規(guī)律會更加明顯，如此才會在保留基本波形的前提下?lián)碛凶畲蟮脑肼暼コ?，這一結論在信號僅含有交流和噪聲成分，或同時含有直流、交流和噪聲3種成分時的去噪效果是十分顯著的，而對僅含有廣義直流和噪聲成分的信號則不然。可以推斷，不論原始信號中是否含有直流成分，交流成分與噪聲成分之間的分布差異需要較多的奇異值來體現(xiàn)，以實現(xiàn)兩種成分的分離；而對于僅含有廣義直流和噪聲成分的無代表性特征頻率的信號而言，這種對奇異值較多數(shù)量的要求未必適用。

文獻[12]推導出直流分量的能量僅由第1個奇異值所代表；若信號中同時存在交流分量，其能量則由第2個奇異值起的有限個奇異值所代表。雖然本文定義的廣義直流分量不一定為常數(shù)，但與噪聲能量的分布特性相比，能量仍然較為集中。由于本文的研究對象不存在明顯的交流成分，可以推斷，過大的q將導致廣義直流分量能量的分散，即有一部分能量可能散失至噪聲對應的奇異值中，在進行有效奇異值篩選的同時也將這部分能量進行了去除，從而丟失了部分細節(jié)。為使這種能量損失降到最低，并能同時完成去噪的目標， q應盡可能取小，同時考慮到Hankel矩陣本身的定義要求，n只能設置為2，且由m與n之間的約束關系，m應取為N?1。

(2)有效奇異值數(shù)量 k的選擇 奇異值差分譜方法雖然在有效奇異值與其它奇異值之間進行了隔斷，但不能保證有用分量與噪聲之間的精確分離，有效奇異值的過多或過少均會影響去噪效果。根據(jù)圖1(a)，第1大峰值與第2大峰值之間的差距很大，使得噪聲能被有效去除。作為對比，Block信號奇異值差分譜的兩最大峰值分別位于第3和第5個坐標(見圖3(b))，兩峰值之間的差距遠不如仿真信號的明顯，并且之后還連續(xù)有一些差距不大的峰值出現(xiàn)，而這些峰值攜帶了有用信號的某些細節(jié)信息。通過一系列實驗可以發(fā)現(xiàn)，若以其它峰值的位置為基準對Block信號繼續(xù)增加奇異值，發(fā)現(xiàn)去噪后的波形細節(jié)雖有一定程度的增加，但同時也會引入部分噪聲。根據(jù)這一分析，基于奇異值差分譜的SVD方法在面對一類無代表性特征頻率的含噪信號時，難以取得理想的去噪效果。

考慮到矩陣結構被重新定義，k也需重新選擇，但必須受q的約束，即k＜q。根據(jù)對Hankel矩陣結構合理性的討論以及q=min(m,n)可知，k僅能取1。

以上分析過程完成了對SVD去噪方法中相關參數(shù)的重新定義，仍以Block信號為例，觀察進行1次SVD去噪的效果(見圖5(a))。可以判定，相對于圖3(a)中的原始信號，進行1次去噪處理后的信號去除了部分噪聲，同時保留了大部分有用細節(jié)，說明重新定義的相關參數(shù)對提升Block信號的去噪效果具有一定作用。還需注意到，僅進行此處理過程還無法達到去噪所要求的理想波形，但不可否認的是，這一處理方式是有效的，關鍵在于如何進一步處理，以去除剩余噪聲。為此，自然聯(lián)想到繼續(xù)迭代處理的方式，即對前一次處理得到的信號再次按照重新定義的相關參數(shù)進行處理。為驗證這一思路的合理性，圖5(b)給出了重復迭代3次處理后的信號波形，與圖5(a)的波形相比，顯然更多的噪聲被去除，波形曲線也更為平滑，同時也保證了細節(jié)信息的保留。

4.2 方法的流程描述

基于以上分析，本文提出一種SVD迭代處理的去噪方法，具體步驟為：

步驟 1 對長度為N的目標信號x構造Hankel矩陣，其中設置矩陣列數(shù)n=2，則行數(shù)m=N?1；

步驟2 確定迭代次數(shù)T，并令當前迭代次數(shù)i=1；

步驟 3 對Hankel矩陣進行奇異值分解，選擇第1個奇異值為有效奇異值(因為k=1)，同時將其它奇異值置0，然后按平均法進行信號重構[15]，得到初次處理的去噪信號x'；

步驟 4 若i＜T，說明當前結果信號未達到結束條件，i=i+1，并將x'作為下次處理的目標信號，x'→x，返回步驟3，否則輸出x'作為最終的去噪信號。

一個需要解決的問題是如何確定迭代次數(shù)T，對不同信號進行迭代處理的結果表明，只需較少的迭代次數(shù)就會達到令人滿意的效果，此時噪聲去除量已趨近飽和，若繼續(xù)迭代，信號波形的變化將十分微小。從這一意義上來說，迭代次數(shù)的選擇具有較強的靈活性，可根據(jù)實際去噪效果進行合理選取。

4.3 時間復雜度分析

步驟1需要對信號序列以滑動窗口的方式依次截取m次，時間復雜度為O(m)；對于步驟2～步驟3，若矩陣維數(shù)為m×n，則分解與重構過程的時間復雜度分別為O(4m2n+8mn2+9n3), O(mn)；對于步驟4，每次迭代時信號長度均為N, Hankel矩陣維數(shù)也一致，故每次迭代所用時間復雜度相同。假若共進行T次迭代，則整個過程所用時間復雜度為(O(m) +O(4m2n+8mn2+9n3)+O(mn))×T≈O(Tm2n)。按本文方法n被設置為2, m=n或僅相差1。一般而言，信號長度N通常會很大，考慮到m+n=N+1為常數(shù)，根據(jù)數(shù)學常識，當m與n差距很小時兩者之積要遠大于差距很大時的兩者之積，并且由于T非常小，本文方法的復雜度要遠小奇異值差分譜方法。

5 實例驗證與分析

反映軍用飛機發(fā)動機狀態(tài)的飛行數(shù)據(jù)是一種典型的非平穩(wěn)、非線性信號，受各種內外部因素的影響，飛行數(shù)據(jù)中將不可避免地會混入各種噪聲。圖6(a)給出了某型飛機在處于飛行狀態(tài)時采集到的一段真實振動值Vib信號，從波形上看不具有明顯的周期性。求解它的幅值譜觀察整個頻帶區(qū)間內的頻率分布(見圖6(b))，發(fā)現(xiàn)能量主要集中在低頻段，不存在明顯尖峰，故符合本文研究對象的特點。

Vib原始采樣信號的長度N=1200，按照本文所提出的方法，首先構造Hankel矩陣，其列數(shù)n=2，行數(shù)m=1199；取不同的試驗迭代次數(shù)T=5和T=50。對構造好的Hankel矩陣作奇異值分解，并選擇第1個奇異值為有效奇異值，同時將其它奇異值置0，然后重構信號，不同迭代次數(shù)下的去噪結果如圖7。

從圖7可以看出，本文方法對于不規(guī)律變化的Vib信號去噪效果較好，在經(jīng)過5次迭代后的去噪波形已十分理想，信號較為純凈，不帶有明顯的噪聲；同時波形細節(jié)也非常豐富，在整個時間范圍內可明顯分辨信號的變化趨勢，能夠滿足振動值狀態(tài)監(jiān)控的實際需求。另外，當?shù)?0次后，相對于圖7(a), Vib信號的波形改變量極其微小，只表現(xiàn)在某些細節(jié)處，因而繼續(xù)迭代已無必要。在實際應用中，迭代次數(shù)的選擇不宜過多，應該根據(jù)實際信號的變化特點合理選擇迭代次數(shù)，提高處理效率。

圖4 Block信號的幅值譜

圖5 不同SVD處理次數(shù)的結果對比

圖6 Vib原始采樣信號及幅值譜

圖7 不同迭代次數(shù)下Vib信號的去噪效果對比

圖8 基于奇異值差分譜方法的Vib 序列處理結果及去噪波形

為體現(xiàn)本文方法在去噪上的優(yōu)勢，對Vib信號采用基于奇異值差分譜的SVD方法進行去噪實驗，得到奇異值序列以及奇異值差分譜如圖8(a)所示。由于第1個奇異值較大，為清晰顯示其它較小的奇異值，圖中并未將第1個奇異值繪出?？梢钥吹?，差分譜序列的第1大峰值位于第1個坐標處，按照差分譜方法的選擇原理，應以第2大峰值位置為準選擇有效奇異值，即k=13。選擇前13個奇異值作為有效奇異值進行信號重構，得到的結果如圖8(b)所示，與圖7相比，雖然奇異值差分譜方法繪出了Vib信號大致的變化趨勢，但仍存在有用細節(jié)丟失的現(xiàn)象。例如在150～230以及950～1050采樣點區(qū)間內，信號被過度平滑，從而引起波形失真。

計算時間方面的結果如表1所示?？梢钥闯?，本文方法在取得較好去噪效果的同時(即T=5)，所耗費時間要遠小于基于奇異值差分譜方法，而且當T=50時，所耗費時間才與基于奇異值差分譜方法相近，這一結論也印證了4.3節(jié)對時間復雜度分析的正確性。這同時也從另一個側面說明，過多的迭代次數(shù)并不能保證去噪質量的大幅度提升，反而會降低本文方法在時間效率上的優(yōu)勢。通過以上實驗結果表明，無論從去噪效果還是計算效率來看，本文方法均要優(yōu)于基于奇異值差分譜的方法。

6 結束語

對于一類非線性含噪信號，本文提出了一種迭代處理的有效去噪方法。首先為Hankel矩陣設置有效維數(shù)，對其進行奇異值分解后確定有效奇異值數(shù)目，然后完成信號的重構。若未達到迭代次數(shù)，對獲得的新信號重復上述過程，最終獲得理想的去噪信號。在真實飛行數(shù)據(jù)上的實驗結果顯示出本文方法只需較少的迭代次數(shù)即可達到理想去噪效果，運算時間也得到有效縮減，能夠滿足應用需求。

表1 不同方法的運行時間比較

[1] Massari C, Brocca L, Ciabatta L, et al.. A Wiener-wavelet based filter for de-noising satellite soil moisture retrievals[C]. Proceedings of the EGU General Assembly, Vienna, Austria, 2014: 1123-1148.

[2] Donoho D L. De-noising by soft-thresholding[J]. IEEE Transactions on Information Theory, 1995, 41(2): 613-618.

[3] He Wang-peng, Zi Yan-yang, Chen Bin-qiang, et al.. Tunable Q-factor wavelet transform denoising with neighboring coefficients and its application to rotating machinery fault diagnosis[J]. SCIENCE CHINA Technological Sciences, 2013, 56(8): 1956-1965.

[4] Maryam A and Rodrigo Q Q. Automatic denoising of single-trial evoked potentials[J]. NeuroImage, 2013, 66(10): 672-680.

[5] Pascal P M, Christian B, and Florence B. Denoising NMR time-domain signal by singular-value decomposition accelerated by graphics processing unites[J]. Solid State Nuclear Magnetic Resonance, 2014, 61(5): 28-34.

[6] Li Zhen-xing and Dai Wei-xiao. Local mean decomposition combined with SVD and application in telemetry vibration signal processing[J]. Applied Mechanics and Materials, 2013,347(2): 854-858.

[7] Ajit R, Anand R, and Arunava B. Image denoising using the higher order singular value decomposition[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2013, 35(4): 849-862.

[8] Cheng D H, Xu Q J, and Yao W F. Regional wind energy resource forecasting based on SVD and support vector machine[J]. Advanced Materials Research, 2014, 247(12): 1070-1072.

[9] 趙學智, 葉邦彥, 陳統(tǒng)堅. 奇異值差分譜理論及其在車床主軸箱故障診斷中的應用[J]. 機械工程學報, 2010, 46(1): 100-108.

Zhao Xue-zhi, Ye Bang-yan, and Chen Tong-jian. Difference spectrum theory of singular value and its application to the fault diagnosis of headstock of lathe[J]. Journal of Mechanical Engineering, 2010, 46(1): 100-108.

[10] 呂永樂, 郎榮玲. 基于奇異值分解的飛行數(shù)據(jù)降噪方法[J]. 計算機工程, 2010, 36(3): 260-262.

Lü Yong-le and Lang Rong-ling. Noise reduction method for flight data based on singular value decomposition[J]. Computer Engineering, 2010, 36(3): 260-262.

[11] Zhao H M, Shen H, Fu Y, et al.. Using singular value decomposition and high order spectrum for bearings fault diagnosis[C]. Proceedings of the IEEE Transportation Electrification Conference and Expo Asia-Pacific(ITEC Asia-Pacific), Beijing, China, 2014: 1-4.

[12] Zhao Xue-zhi and Ye Bang-yan. Selection of effective singular values using difference spectrum and its application to fault diagnosis of headstock[J]. Mechanical Systems and Signal Processing, 2011, 25(5): 1617-1631.

[13] 雷達, 鐘詩勝. 基于奇異值分解和經(jīng)驗模態(tài)分解的航空發(fā)動機健康信號降噪[J]. 吉林大學學報, 2013, 43(3): 764-770.

Lei Da and Zhong Shi-sheng. Aircraft engine health signal denoising based on singular value decomposition and empirical mode decomposition methods[J]. Journal of Jilin University, 2013, 43(3): 764-770.

[14] 趙學智, 葉邦彥, 陳統(tǒng)堅. 基于小波-奇異值分解差分譜的弱故障特征提取方法[J]. 機械工程學報, 2012, 48(7): 37-48.

Zhao Xue-zhi, Ye Bang-yan, and Chen Tong-jian. Extraction method of faint fault feature based on wavelet-SVD difference spectrum[J]. Journal of Mechanical Engineering, 2012, 48(7): 37-48.

[15] Brand M. Incremental singular value decomposition of uncertain data with missing values[C]. Proceeding of the 2002 European Conference on Computer Vision, Copenhagen, Denmark, 2002: 1-12.

查 翔： 男，1988年生，博士生，研究方向為飛機狀態(tài)監(jiān)控與故障診斷、人工智能及其應用等.

倪世宏： 男，1963年生，教授，博士生導師，研究方向為飛行數(shù)據(jù)智能處理、飛行器狀態(tài)監(jiān)控與健康管理等.

張 鵬： 男，1982年生，博士，講師，研究方向為航空發(fā)動機故障診斷、故障預測與健康管理等.

Effective Iteration Method of a Class of Nonlinear Signal Denoising Based on Singular Value Decomposition

Zha Xiang Ni Shi-hong Zhang Peng
(College of Aeronautics and Astronautics Engineering, Air Force Engineering University, Xi'an 710038, China)

To solve a class of nonlinear signal denoising, an effective iteration method based on the Singular Value Decomposition (SVD) is proposed. When the signals have no obvious characteristic frequency and non-periodic change, the current difference spectrum method is not applicable by comparing the results on the two class of nonlinear signal, and then the corresponding reason is analyzed. According to the signal feature, the structure of the Hankel matrix is defined again and the valid singular values are determined. The effective denoising is realized by the repeated iteration which is based on the SVD. The results of the flight data demonstrate that the proposed method can effectively reduce the noise and improve the computing efficiency as well.

Signal processing; Denoising; Singular Value Decomposition (SVD); Characteristic frequency; Hankel matrix

TN911.7

： A

：1009-5896(2015)06-1330-06

10.11999/JEIT141605

2014-12-15收到，2015-03-05改回

國家自然科學基金(61372167, 61379104)資助課題

*通信作者：查翔 zha_xiang@126.com