王小君，鞏 超，和敬涵，李長城，張金國，吉 荔

（北京交通大學(xué) 國家能源主動配電網(wǎng)技術(shù)研發(fā)中心，北京 100044）

0 引言

隨著我國城市軌道交通的迅速發(fā)展，其對配電網(wǎng)電能質(zhì)量的影響日益受到人們的關(guān)注，其中牽引供電系統(tǒng)中廣泛應(yīng)用的多脈波整流電路產(chǎn)生的諧波電流和牽引電機(jī)頻繁起動、制動引起的電壓波動在諸多不利影響中尤為主要。本文充分考慮車輛各運行階段的動態(tài)過程對配電網(wǎng)電能質(zhì)量的影響，將牽引過程中的諧波電流時變特性和電壓波動細(xì)節(jié)分析作為研究重點。

城軌電流諧波的傳統(tǒng)評估方法往往采用線性全局的快速傅里葉變換（FFT）分析手段，選取測量時段內(nèi)各相實測值的95%概率值最大值作為判斷諧波是否超過允許值的依據(jù)。文獻(xiàn)[1]重點分析了城市軌道交通牽引供電系統(tǒng)諧波產(chǎn)生原因，詳細(xì)描述了諧波類型及其特性。文獻(xiàn)[2]基于實測數(shù)據(jù)，給出了城軌110 kV交流進(jìn)線諧波電流的計算結(jié)果，但未給出具體細(xì)節(jié)。文獻(xiàn)[3]基于PSCAD平臺仿真研究了不同運行工況對公共電網(wǎng)諧波含量的影響。但現(xiàn)有研究沒有充分考慮車輛運行動態(tài)過程對諧波含量的影響，即缺乏對諧波時變特性的研究。近年來，小波變換、Hilbert-Huang變換（HHT）等時頻分析方法的發(fā)展，為電流諧波時變特性的精細(xì)化描述提供了研究手段。文獻(xiàn)[4]將小波包變換（WPT）應(yīng)用于電網(wǎng)諧波信號的提取，表明其具有精確劃分頻帶的能力，但無法直接計算諧波的瞬時幅值。文獻(xiàn)[5]將HHT方法應(yīng)用于電力系統(tǒng)諧波的分析，對影響HHT方法分析精度的邊界問題和插值函數(shù)進(jìn)行改進(jìn)，最后得到各次諧波的瞬時頻率和瞬時幅值?？紤]到城軌牽引供電系統(tǒng)電流信號基波能量很大、諧波能量較小的特點，直接使用該算法會出現(xiàn)嚴(yán)重的模態(tài)混疊[6]，從而導(dǎo)致檢測結(jié)果失真。文獻(xiàn)[7]采用基于傅里葉變換的HHT方法較好地解決了上述問題，但濾波器的延時效應(yīng)使得瞬時參數(shù)檢測結(jié)果在起始處產(chǎn)生波動。

通過分析車輛的牽引制動曲線可以看出，其運行過程存在有功、無功的快速隨機(jī)波動，從而引起網(wǎng)側(cè)電壓的波動。傳統(tǒng)的波動程度評估指標(biāo)往往采用短時閃變值，文獻(xiàn)[2-3]基于實測數(shù)據(jù)，給出了列車運行過程中網(wǎng)側(cè)電壓波動和閃變的計算結(jié)果，但未給出具體細(xì)節(jié)?？紤]到城軌車輛牽引過程時間較短，因此在分析電壓波動時，基于統(tǒng)計的閃變指標(biāo)不再適用。文獻(xiàn)[8]指出電壓閃變包絡(luò)線對于電壓波動細(xì)節(jié)描述的重要意義，并提出基于瞬時無功功率理論求取電壓閃變參數(shù)的方法。文獻(xiàn)[9-10]通過仿真將瞬時電壓d-q分解算法與平方檢測法等算法進(jìn)行了對比，結(jié)果表明前者具有計算量小和檢測精度高的優(yōu)勢。因此，利用瞬時電壓d-q分解算法提取電壓閃變包絡(luò)線可以對車輛牽引過程引起的電壓波動情況進(jìn)行精細(xì)化描述。但實測信號頻率存在一定程度的波動，因此有必要對該方法在頻率波動情況下的工程適應(yīng)性進(jìn)行分析，而目前沒有相關(guān)文獻(xiàn)涉及這一問題。

本文充分考慮城軌電流諧波的特殊性，融合WPT精確劃分頻帶的優(yōu)勢和HHT準(zhǔn)確提取瞬時參數(shù)的能力，提出了基于WPT預(yù)處理的 HHT（WPTHHT）方法，分析諧波電流的時變特性；考慮到城軌電壓信號頻率存在一定的波動，首先分析了頻率波動對瞬時電壓d-q分解法計算精度的影響，給出了頻率波動情況下的誤差公式，然后將該方法應(yīng)用于車輛牽引過程中網(wǎng)側(cè)電壓閃變包絡(luò)線的提取。最后結(jié)合北京地鐵一號線某牽引變電站的實測錄波數(shù)據(jù)，研究了城軌供電系統(tǒng)雙邊供電方式下車輛牽引過程對配電網(wǎng)電能質(zhì)量的影響。

1 基于WPT-HHT的諧波電流分析

1.1 基于HHT的諧波檢測原理

HHT方法首先通過經(jīng)驗?zāi)B(tài)分解（EMD）將信號分解為若干個本征模態(tài)函數(shù)（IMF），其中IMF必須滿足2個條件：整個數(shù)據(jù)段內(nèi)，極值點個數(shù)和過零點個數(shù)相等或相差1個；任意時刻，由局部極大值點和極小值點形成的上、下包絡(luò)線平均值為0。然后對IMF進(jìn)行Hilbert變換得到信號的Hilbert譜。

對于電力系統(tǒng)諧波信號，EMD得到的各IMF對應(yīng)著不同的諧波分量[5]，此時原始信號 s（t）表示為：

其中，ci為第i個IMF分量，表示某次諧波成分；rn為殘余項，代表了信號的趨勢；n為IMF分量個數(shù)。

對滿足單分量信號要求的IMF分量ci（t）進(jìn)行Hilbert變換得：

構(gòu)造解析信號：

最后根據(jù)式（4）、（5）求得具有物理意義的瞬時頻率 fi（t）和瞬時幅值 ai（t），實現(xiàn)諧波參數(shù)辨識。

1.2 HHT存在的問題及處理方法

HHT方法能夠根據(jù)信號自身特性將隱藏在其中的各個模態(tài)依次剝離，實現(xiàn)電力系統(tǒng)諧波的自適應(yīng)提取，但在應(yīng)用過程中存在端點效應(yīng)、模態(tài)混疊等問題。

端點效應(yīng)表現(xiàn)在2個方面：在應(yīng)用EMD方法對非平穩(wěn)信號進(jìn)行分解時，數(shù)據(jù)的兩端會產(chǎn)生發(fā)散現(xiàn)象，并且這種發(fā)散的結(jié)果會逐漸向內(nèi)“污染”整個數(shù)據(jù)序列從而使結(jié)果嚴(yán)重失真；在對各IMF分量進(jìn)行Hilbert變換時，由于Hilbert變換算法不可避免的加窗效應(yīng)，信號的端點也會出現(xiàn)嚴(yán)重的端點效應(yīng)[11]。目前對于端點效應(yīng)抑制的研究主要包括基于神經(jīng)網(wǎng)絡(luò)的數(shù)據(jù)延拓方法、基于AR模型的時間序列線性預(yù)測方法、邊界波形匹配預(yù)測法、極值點延拓法等。文獻(xiàn)[12]提出基于支持向量回歸機(jī)（SVR）的端點延拓方法，通過與神經(jīng)網(wǎng)絡(luò)方法的比較，驗證了該算法的有效性。文獻(xiàn)[13]引入粒子群優(yōu)化（PSO）算法對SVR參數(shù)進(jìn)行優(yōu)化，提高了預(yù)測精度，較好地抑制了端點效應(yīng)。

利用EMD方法進(jìn)行篩分過程中，容易發(fā)生混疊的信號主要有3種：包含頻率很鄰近的分量；包含的2個分量頻率雖不是很鄰近，但是其中一個分量的能量遠(yuǎn)超過另一個分量；信號擾動有跳變發(fā)生[11]。實測電網(wǎng)信號基波能量很大，各次諧波能量相對較小，因此直接利用EMD方法進(jìn)行分解會產(chǎn)生模態(tài)混疊的問題。目前對于模態(tài)混疊改善的研究主要有噪聲輔助數(shù)據(jù)分析法、改進(jìn)掩膜信號法、差分運算和累計求和、歸一化幅值系數(shù)、基于傅里葉變換的EMD方法等。文獻(xiàn)[14-15]首先通過WPT將待分析信號分解為若干窄帶信號，然后進(jìn)行HHT分析，實現(xiàn)了模態(tài)混疊信號的有效篩分。

1.3 基于WPT-HHT的諧波電流時變特性分析方法

牽引變電站10 kV進(jìn)線電流信號基波能量很大，其他各次諧波能量較小，因此直接對其進(jìn)行EMD會出現(xiàn)模態(tài)混疊現(xiàn)象。WPT具有良好的分頻特性，而且分解精度不受諧波含量的影響，利用其對采集信號進(jìn)行預(yù)處理，可以準(zhǔn)確地提取任意頻率的諧波信號。在此基礎(chǔ)上對各次諧波信號進(jìn)行EMD，可避免模態(tài)混疊問題?；谝陨峡紤]，本文提出WPT-HHT方法，詳述如下：首先對實測電流信號進(jìn)行基于傅里葉變換的諧波特性分析，確定需要分析的諧波電流頻率；綜合考慮采樣頻率以及分析精度的要求，選取適當(dāng)?shù)男〔ò瘮?shù)和分解層數(shù)進(jìn)行小波包分解；分析各層所含諧波的成分進(jìn)而將待分析諧波成分所在層進(jìn)行重構(gòu)，采用基于PSO算法的SVR對重構(gòu)信號兩端進(jìn)行延拓；對延拓后的窄帶信號進(jìn)行EMD，通過相關(guān)系數(shù)法剔除虛假分量后得到真實的IMF分量，最后通過Hilbert變換計算各次諧波幅值。算法流程如圖1所示。

圖1 諧波電流分析流程Fig.1 Flowchart of harmonic current analysis

1.4 小波包基函數(shù)及分解層數(shù)的選取原則

在上述WPT-HHT方法中，WPT的預(yù)處理效果直接影響檢測精度，因此有必要對小波包基函數(shù)和分解層數(shù)選擇的關(guān)鍵問題進(jìn)行研究。

基于WPT的信號預(yù)處理實質(zhì)上就是利用小波的分頻特性刻畫信號在不同頻帶下各個時刻的信號特征，為避免劃分頻帶時的頻率混疊現(xiàn)象，需要優(yōu)先選擇頻域緊支撐正交小波[16]。Daubechies系列小波具備正交性、緊支撐性和可離散性等特征，滿足電力系統(tǒng)諧波分析的要求，而且隨著階次的增大，頻域消失矩階數(shù)增大，頻帶劃分效果越好，頻域分辨能力越強(qiáng)。但根據(jù)測不準(zhǔn)原理，此時時域緊支撐性減弱，而且計算量會相應(yīng)增加，實時性變差[17]。經(jīng)過多次實驗和綜合分析，本文最終選取40階的Daubechies小波函數(shù)（db40）作為分解和重構(gòu)算法的母小波。

小波包的分解層數(shù)與頻帶寬度滿足關(guān)系式：

其中，Δf為頻帶寬度；fs為采樣頻率；j為分解層數(shù)。

由式（6）可知：分解層數(shù)越大，頻帶劃分越精細(xì)。但隨著分解層數(shù)增多，級間濾波器的數(shù)目相應(yīng)增多，會造成信號位移增大；而且小波分解、合成都需要邊界延拓，層數(shù)越多則引起的邊界失真越嚴(yán)重。因此在保證頻率分辨率的前提下，分解層數(shù)越少則諧波分析精度越高[18]。另一方面，由于小波包分解在頻帶邊緣處的能量混疊情況最嚴(yán)重，因此在確定分解層數(shù)時，還要避免所要提取頻帶的能量出現(xiàn)在頻域劃分邊緣[19]。本文采用錄波數(shù)據(jù)的采樣頻率為5000 Hz，綜合考慮頻率分辨率和分析精度要求，分解層數(shù)確定為6層，形成了頻帶寬度為39.0625 Hz的64個均勻頻帶。

1.5 仿真信號分析

為驗證本文所提算法的有效性，以式（7）所示確定性的時域信號為分析對象，設(shè)定采樣頻率為5000 Hz，基于MATLAB平臺進(jìn)行了仿真實驗。

該仿真信號的2個頻率分量能量差別較大，不能滿足文獻(xiàn)[20]的完全分解條件，無法根據(jù)特征尺度有效分離不同的模態(tài)成分，對于模態(tài)混疊問題的分析具有典型代表性。采用傳統(tǒng)EMD算法對其進(jìn)行分解只能產(chǎn)生1個IMF分量，然后進(jìn)行Hilbert變換，得到瞬時參數(shù)如圖2所示，由于出現(xiàn)了模態(tài)混疊現(xiàn)象，分析結(jié)果不能正確地反映仿真信號的諧波信息。

為了消除EMD過程中的模態(tài)混疊現(xiàn)象，本文按照圖1所示的信號處理流程對仿真信號進(jìn)行了諧波分析，具體步驟如下。

a.選用db40小波作為WPT函數(shù)，進(jìn)行6層小波包分解與重構(gòu)。

圖2 基于傳統(tǒng)HHT方法的仿真信號諧波分析Fig.2 Harmonic analysis of simulation signal based on traditional HHT method

b.參照文獻(xiàn)[13]，選取徑向基函數(shù)為SVR的核函數(shù)，采用真實值向量與預(yù)測值向量的歐氏距離作為PSO算法的目標(biāo)函數(shù)對其進(jìn)行參數(shù)優(yōu)化。設(shè)定PSO算法的慣性權(quán)值為0.6，認(rèn)知學(xué)習(xí)因子為1.5，社會學(xué)習(xí)因子為1.7，迭代50次得到核函數(shù)最優(yōu)懲罰參數(shù)為10，核參數(shù)為0.037。利用上述參數(shù)建立SVR模型對各重構(gòu)信號兩端進(jìn)行延拓并進(jìn)行EMD。

c.計算各IMF與原始信號的相關(guān)系數(shù)，發(fā)現(xiàn)第1個IMF分量與原始信號相關(guān)性最高，剔除與原始信號相關(guān)性差的分量，從而消除了檢測儀器等外界因素產(chǎn)生的虛假頻率成分。

d.最后進(jìn)行Hilbert變換計算得到各次諧波瞬時頻率和瞬時幅值。

同時，為了說明利用PSO-SVR算法進(jìn)行端點延拓的重要性，將本文HHT分析結(jié)果與未進(jìn)行端點延拓的HHT分析結(jié)果進(jìn)行了對比，如圖3所示?？梢?，WPT-HHT方法有效地解決了EMD過程中不同頻率分量能量差別過大引起的模態(tài)混疊，且PSO-SVR算法具有良好的端點效應(yīng)抑制效果。

圖3 基于WPT-HHT方法的仿真信號諧波分析Fig.3 Harmonic analysis of simulation signal based on WPT-HHT method

2 基于瞬時電壓d-q分解法的電壓波動分析

瞬時電壓d-q分解法實質(zhì)是通過Park變換將abc坐標(biāo)系下的三相電壓轉(zhuǎn)換到dq0坐標(biāo)系下，進(jìn)而提取基波正序分量。文獻(xiàn)[21]對該方法的基本原理進(jìn)行了詳細(xì)闡述，并重點討論了諧波對于計算結(jié)果的影響。本文則重點研究頻率波動對該方法檢測精度的影響，討論其處理實測信號的工程適用性。

2.1 頻率波動引起的計算誤差分析

以a相電壓為分析對象，設(shè)其時域表達(dá)式為：

其中，U為電壓有效值；ω1為角頻率；初相角為0°。

將a相電壓延時120°，同時考慮到電壓瞬時頻率的變化，設(shè)由頻率波動造成的虛擬b相電壓相位偏移為 θ（t），得到 ub的時域表達(dá)式：

最后通過uc=-ua-ub，可得虛擬c相電壓。則abc坐標(biāo)系下電壓向量表達(dá)式為：

設(shè)d-q變換矩陣C的頻率為ω2，C可表示為：

則有：

化簡得：

其中，ur為頻率波動造成的誤差項系數(shù)。

由式（13）可以得出如下結(jié)論：

a.變換矩陣的頻率ω2不會對u2d+u2q結(jié)果產(chǎn)生影響，因此在利用瞬時電壓d-q分解法求取電壓閃變包絡(luò)線時可以直接設(shè)定ω2為100πrad／s；

b.電壓頻率波動導(dǎo)致的相位偏移量θ（t）會產(chǎn)生u2d+u2q的計算誤差項urU2，誤差項由urdU2和uraU2構(gòu)成，其中uraU2可以通過中心頻率為ω1／π的帶阻濾波器濾除。

考慮到uraU2為可控誤差，下文重點分析urdU2的影響。設(shè)測量時段內(nèi)頻率波動最大值為Δfmax，則最大相位偏移量θmax和urdmax滿足：

定義有效值計算的相對誤差為：

根據(jù)現(xiàn)行國家標(biāo)準(zhǔn)GB／T15945—2008《電能質(zhì)量 電力系統(tǒng)頻率偏差》規(guī)定：正常運行條件下頻率偏差限值為±0.2 Hz，系統(tǒng)容量較小時可放寬到±0.5Hz。本文選取國標(biāo)規(guī)定頻率偏差范圍內(nèi)頻率波動的不同場景，計算相位偏移最大值和相對誤差最大值如表1所示。表1中的結(jié)果表明，利用瞬時電壓d-q分解法提取電壓閃變包絡(luò)線可以滿足工程計算的要求。

表1 不同頻率波動條件下的相位偏移與相對誤差Table 1 Phase shift and relative error in different frequency fluctuation conditions

2.2 電壓閃變包絡(luò)線提取方法

為了研究車輛牽引過程中牽引變電站10 kV母線電壓波動特征，本文將瞬時電壓d-q分解法應(yīng)用于電壓閃變包絡(luò)線的提取，具體步驟如下：首先對實測電壓信號進(jìn)行帶通濾波以保留基波信息，利用式（12）計算u2d+u2q，對其進(jìn)行帶阻濾波以消除uraU2項，然后根據(jù)

得到電壓閃變包絡(luò)曲線。算法流程如圖4所示。

3 基于現(xiàn)場錄波數(shù)據(jù)的算例分析

為了研究車輛牽引過程中網(wǎng)側(cè)諧波電流時變特性，掌握網(wǎng)側(cè)電壓波動規(guī)律，本文以北京地鐵一號線典型牽引變電站為分析對象，采用第1、2節(jié)所述算法進(jìn)行了算例分析。

圖4 電壓波動分析流程Fig.4 Flowchart of voltage fluctuation analysis

3.1 數(shù)據(jù)來源

實驗使用博電MR1200便攜式多通道波形監(jiān)視記錄儀對牽引變電站內(nèi)部各中壓評估點的電流、電壓進(jìn)行了24 h錄波。圖5所示為牽引變電站主接線。通過對全天錄波數(shù)據(jù)的分析，本文最終選擇單車牽引典型工況下35 s時間長度的牽引站進(jìn)線開關(guān)201電流數(shù)據(jù)和PT49電壓數(shù)據(jù)作為分析樣本。測試時，北京電力公司10kV變電站電源處于檢修狀態(tài)，開關(guān)202退出運行，牽引站內(nèi)電能來源為開閉所。

圖5 牽引變電站電氣主接線圖Fig.5 Main electric connection diagram of traction substation

3.2 諧波電流和電壓波動分析

首先利用MATLAB中FFT函數(shù)計算了進(jìn)線電流的傅里葉幅值譜，結(jié)果如圖6所示。由圖6分析可知其主要諧波成分為 3、5、7、23、25 次諧波，其中 23、25次特征諧波含量最高。

圖6 進(jìn)線電流信號傅里葉幅值譜Fig.6 FFT amplitude spectrum of inlet current

本文以含量較高的7次非特征諧波和23次典型特征諧波為例，根據(jù)第1.4節(jié)所述小波包基函數(shù)及分解層數(shù)的確定原則，選用db40小波作為WPT函數(shù)，進(jìn)行6層小波包分解與重構(gòu)，保證了諧波檢測的嚴(yán)格頻帶分割和頻域緊支撐正交條件。當(dāng)采用6層小波包分解時，7次諧波信號剛好位于頻帶（8Δf，9Δf）和（9Δf，10Δf）邊緣，為了盡量減少能量泄漏，將2個頻帶信號重構(gòu)后進(jìn)行了疊加處理，作為7次諧波的提取結(jié)果；23 次諧波位于頻帶（29Δf，30Δf）中間，具有較好的頻域能量聚集性。根據(jù)文獻(xiàn)[22]對WPT頻帶劃分特性的分析，濾波器路徑（二進(jìn)制數(shù)字）與通頻帶（二進(jìn)制數(shù)字）的對應(yīng)關(guān)系如表2所示。按照第1.3節(jié)分析流程，計算各諧波分量小波包重構(gòu)信號如圖7所示，諧波電流瞬時幅值如圖8所示（圖7、8中由上至下的波形分別對應(yīng)基頻帶、7次諧波帶和23次諧波帶），諧波含有率時變曲線如圖9所示。

分析圖8、圖9可得以下結(jié)論。

a.本算例中，牽引過程7次諧波電流瞬時幅值最大值為0.89 A，其含有率最大值為6.3%；23次諧波電流瞬時幅值最大值為2.15 A，其含有率最大值為1.9%。這符合現(xiàn)行國家標(biāo)準(zhǔn)GB／T14549—1993《電能質(zhì)量 公用電網(wǎng)諧波》的相關(guān)規(guī)定。

表2 6層小波包分解樹通頻帶與濾波器路徑對應(yīng)關(guān)系Table 2 Relationship between passband and filter path of six-layer wavelet packet decomposition tree

圖7 諧波電流小波包重構(gòu)信號Fig.7 Signal reconstructed based on harmonic current wavelet packets

圖8 諧波電流瞬時幅值Fig.8 Instantaneous amplitude of harmonic current

圖9 電流諧波含有率時變曲線Fig.9 Time-varying harmonic ratio curve

b.在城軌車輛牽引過程的恒轉(zhuǎn)矩、恒功率、自然特性的不同階段，基波電流均存在較大波動。

c.10 kV進(jìn)線電流特征諧波幅值與基波電流存在正相關(guān)性，牽引過程諧波含有率基本不變；非特征諧波幅值與基波電流相關(guān)性較弱，牽引過程諧波含有率變化較大。

根據(jù)第2.2節(jié)電壓分析流程，本文計算的電壓閃變包絡(luò)曲線如圖10所示。由圖10可知：城軌車輛牽引過程中10 kV母線電壓波動能力有限，低于現(xiàn)行國家標(biāo)準(zhǔn)GB／T12326—2008《電能質(zhì)量 電壓波動和閃變》規(guī)定的限值3%，并且由于變頻調(diào)速技術(shù)的使用，電壓暫降現(xiàn)象已經(jīng)基本消除。

圖10 電壓閃變包絡(luò)曲線Fig.10 Envelope of voltage flicker

4 結(jié)論

本文提出基于WPT-HHT的諧波電流時變特性分析方法，研究了WPT基函數(shù)及分解層數(shù)的選取原則，并采用PSO-SVR算法進(jìn)行HHT端點延拓。該方法克服了傳統(tǒng)HHT方法的端點效應(yīng)和模態(tài)混疊問題，為諧波時變特性的研究提供了一種新思路。分析了檢測信號頻率波動對瞬時電壓d-q分解法計算精度的影響，給出頻率波動情況下的誤差公式，并將其應(yīng)用于電壓閃變包絡(luò)線的提取，對于分析沖擊性負(fù)荷接入配電網(wǎng)引起的電壓波動具有重要意義。最后，基于北京地鐵一號線某牽引變電站的實測錄波數(shù)據(jù)，利用所提算法分析了車輛牽引過程對配電網(wǎng)電能質(zhì)量的影響。算例分析結(jié)果表明，該算法是研究電流諧波時變特性和電壓波動規(guī)律的有效手段。