一類半波對(duì)稱FFT改進(jìn)算法與電網(wǎng)諧波分析

湯天浩，鄭 慧

（上海海事大學(xué)電力傳動(dòng)與控制研究所，上海200135）

引言

由于電力電子技術(shù)的飛速發(fā)展，各種電力電子裝置普遍使用，諧波所造成電力系統(tǒng)污染等危害也日趨嚴(yán)重。因此，近年來(lái)，電能質(zhì)量控制研究已被人們逐漸重視，成為國(guó)際電氣工程領(lǐng)域的研究熱點(diǎn)之一。電能質(zhì)量控制的問題涉及面很廣，而諧波測(cè)量和諧波分析是研究電能質(zhì)量問題中的一個(gè)重要分支，也是解決諧波問題和治理電網(wǎng)污染的出發(fā)點(diǎn)和主要依據(jù)。

長(zhǎng)期以來(lái)，傅里葉變換（FT）一直是研究分析諧波的經(jīng)典方法和有效工具，也廣泛用于電力系統(tǒng)的諧波測(cè)量分析中。隨后又發(fā)展了用離散傅里葉變換（DFT）或快速傅里葉變換（FFT）進(jìn)行頻譜分析進(jìn)而得到各次諧波幅值和相位。在實(shí)際應(yīng)用當(dāng)中，電網(wǎng)諧波分析就是通過對(duì)信號(hào)進(jìn)行同步采樣，將周期信號(hào)展開為傅里葉級(jí)數(shù)，得到各次諧波系數(shù)進(jìn)而進(jìn)行各項(xiàng)諧波指標(biāo)的計(jì)算，如THD，HRn等。

然而，F(xiàn)T不具有時(shí)域分析能力，對(duì)含有短時(shí)高頻分量和長(zhǎng)時(shí)間低頻分量的電能質(zhì)量信號(hào)進(jìn)行分析時(shí)具有局限性。為了改善FT的局限性，1994年出現(xiàn)了窗口FT即后來(lái)的短時(shí)FT算法，比如：矩形窗、三角形窗、漢寧窗、海明窗、布萊克曼窗等。但是短時(shí)FT的時(shí)頻窗口是固定的，不適于分析信號(hào)的突變過程，且其離散形式?jīng)]有正交展開，難以實(shí)現(xiàn)高效算法。

本文首先通過對(duì)DFT與周期信號(hào)傅里葉級(jí)數(shù)以及非周期信號(hào)FT之間關(guān)系的研究，得到了用FFT對(duì)非周期信號(hào)進(jìn)行精確諧波測(cè)量的公式。進(jìn)而，針對(duì)一類電網(wǎng)波形具有半波對(duì)稱性的特點(diǎn)，提出了一種改進(jìn)的FFT算法。改進(jìn)的FFT算法對(duì)于只含奇次諧波分量的一類半波對(duì)稱性波形，在算法流程中省去了偶次諧波分量部分。因此，改進(jìn)的FFT算法的計(jì)算量是傳統(tǒng)FFT的一半，大大節(jié)省了運(yùn)算時(shí)間。

1 快速傅里葉變換算法

1.1 DFT基本算法[1]

對(duì)于對(duì)N點(diǎn)序列x(n)，其DFT變換定義為：

顯然，求出其N點(diǎn)X(k)需要N2次復(fù)數(shù)乘法及N(N－1)次復(fù)數(shù)加法。而實(shí)現(xiàn)一次復(fù)數(shù)乘法需要四次實(shí)數(shù)乘兩次實(shí)數(shù)加，實(shí)現(xiàn)一次復(fù)數(shù)加則需要兩次實(shí)數(shù)加。當(dāng)N很大時(shí)，其計(jì)算量是相當(dāng)可觀的。

由于上式中W因子的取值具有如下特點(diǎn)：

（1） W0=1，WN/2=-1；

（2） WNN+r=WrN，WN/2+r=-Wr。

因此，在DFT中包含了大量重復(fù)運(yùn)算，如果能夠巧妙地利用W因子的周期性及對(duì)稱性，將大大減少運(yùn)算量。

1.2 基于頻率抽取得的FFT算法

FFT算法最早由Cooley和Tukey于1965年提出，該算法使N點(diǎn)DFT的乘法計(jì)算量由N2次降為N/2log2N。自Cooley-Tukey的算法提出之后，F(xiàn)FT的發(fā)展方向有兩個(gè)：一是針對(duì)N等于2的整數(shù)次冪的算法，如基2算法、基4算法等；二是以Winograd為代表的一類N不等于2的整數(shù)次冪的算法。近年來(lái)，又不斷涌現(xiàn)許多新的算法，比如：由A.M.Grigoryan等人提出的流圖簡(jiǎn)化算法[2]。

如果采用時(shí)間抽?。―IT）基2 FFT算法，將頻域 X(k)的序號(hào) k按奇、偶分開，對(duì)式（1）的 DFT，先將x(n)按序號(hào)分成上、下兩部分，得：

式中：WNNk/2=（-1）k；k＝0，1，…，N－1。

分別令 k＝2r，k＝2r＋1，可以推導(dǎo)出奇、偶兩個(gè)計(jì)算序列：

這樣就將一個(gè)N點(diǎn)DFT分成了兩個(gè)N/2點(diǎn)的DFT，分的辦法就是將X(k)按序號(hào)k的奇、偶分開。

2 基于半波對(duì)稱的DIF基2FFT改進(jìn)算法

2.1 FFT算法及特點(diǎn)

由上節(jié)對(duì)頻率抽?。―IF）基2 FFT算法的分析可知，該算法是將每一組計(jì)算得到的輸出序列按奇、偶分開而得到的一種快速傅里葉變換算法。而且，通常FFT算法都是全波分析和計(jì)算，即同時(shí)計(jì)算各次諧波，包括奇次諧波和偶次諧波。

2.2 基于半波對(duì)稱FFT算法的提出

為研究電力系統(tǒng)諧波，需先分析電力系統(tǒng)諧波的產(chǎn)生根源及其特點(diǎn)。電力系統(tǒng)的諧波主要是由于電力電子裝置等非線性負(fù)載的大量使用造成的，且具有如下特性[3]：

（1） 奇對(duì)稱性的特點(diǎn)是 f(－t)＝－f(t)，展開為傅里葉級(jí)數(shù)時(shí)沒有余弦相。

（2） 偶對(duì)稱性的特點(diǎn)是 f(－t)＝f(t)，展開為傅里葉級(jí)數(shù)時(shí)沒有正弦相而只有余弦相。

（3） 半波對(duì)稱的特點(diǎn)是 f(t±T/2)＝－f(t)，沒有直流分量且偶次諧波被抵消。

因?yàn)殡娏ο到y(tǒng)往往是由雙向?qū)ΨQ的元件組成的，這些元件產(chǎn)生的電壓和電流具有半波對(duì)稱性[4]，這類半波對(duì)稱特點(diǎn)使我們?cè)谥C波分析時(shí)可以忽略電力系統(tǒng)中的偶次諧波。

如果考慮一類電力系統(tǒng)諧波的半波對(duì)稱特性，并利用DIF的基2 FFT算法的奇、偶序列的可分性，對(duì)頻率抽?。―IF）基2 FFT算法進(jìn)行改進(jìn)，可得到一個(gè)只計(jì)算奇次諧波的用于電網(wǎng)諧波分析的快速傅里葉變換算法。

2.3 一類半波對(duì)稱的DIF基2FFT計(jì)算方法

考慮到一類電力系統(tǒng)的電壓和電流具有半波對(duì)稱性，其偶次諧波被抵消。因此僅需在上節(jié)給出的DIF基2FFT算法中，僅需進(jìn)行奇序列計(jì)算。但X（2r+1）仍然是高復(fù)合數(shù)（N/2）的 DFT，對(duì) X（2r+1）仍然按照式（2）進(jìn)行分解，即：

分別令 r＝2l，r＝2l+1，l＝0，1，…，N/4－1。 有：

式中：l＝0，1，…，N/4－1。

例如：N＝8，則有 l＝0，1，這時(shí) f(n)和 i(n)都是兩點(diǎn)的 DFT，無(wú)需再分。 將 l＝0，1 代入式（10）和（11）有：

這樣，我們就得到了8點(diǎn)DFT的奇次諧波。若N＝16，32或2的更高的冪，可按上述方法繼續(xù)分下去，直到兩點(diǎn)的DFT為止。

2.4 算法的實(shí)現(xiàn)

2.4.1 算法的實(shí)現(xiàn)

改進(jìn)的算法仍然是基于頻率抽取 （DIF）基2 FFT算法的，即：在每一級(jí)都對(duì)頻域序號(hào)k按照奇、偶分開，直到兩點(diǎn)的DFT為止。其抽取過程可用一個(gè)樹形圖來(lái)表示，如圖1所示，由于改進(jìn)算法是不計(jì)算偶次諧波的，因此圖2中只對(duì)奇次諧波進(jìn)行抽取。

圖1 改進(jìn)算法的樹形圖

整個(gè)算法的關(guān)鍵是計(jì)算每一級(jí)所對(duì)應(yīng)的g(n)，h(n)，f(n)，i(n)，…，對(duì)照?qǐng)D 1，可以畫出每一級(jí)所對(duì)應(yīng)的 g(n)，h(n)，f(n)，i(n)，…，的樹形圖，如圖 2 所示。 由于只需計(jì)算奇次諧波部分，因此可以直接從X(2k+1)作為第一級(jí)開始往下分，每一級(jí)所對(duì)應(yīng)的h(n)，f(n)和 i(n)，…，均用 h(n)，h1(n)，h2(n)，…，來(lái)表示。 圖中m 表示第 m 級(jí)，m＝1，2，log2N－2。 l表示第 m 級(jí)所需要計(jì)算的 h(n)的個(gè)數(shù)，即 h1(n)，h2(n)，…，hl(n)，l＝1，2，…，2m。

2.4.2 碼位倒置

圖2 改進(jìn)算法的樹形圖

由頻率抽取信號(hào)流圖，其變換后輸出序列X(k)的次序不再是原來(lái)的自然順序，如果將輸出X(k)的序號(hào)寫成二進(jìn)制，并將二進(jìn)制數(shù)碼翻轉(zhuǎn)，它們對(duì)應(yīng)的十進(jìn)制序號(hào)這正是我們需要的自然順序輸出。在改進(jìn)算法中，只計(jì)算了奇次諧波分量，對(duì)于直流分量和偶次諧波分量可以直接令其為零，這是因?yàn)樵诎氩▽?duì)稱波形中不含直流分量和偶次諧波分量。這樣在編程的時(shí)候，前N/2點(diǎn)數(shù)據(jù)直接賦值為0，而后面N/2點(diǎn)就等于前面計(jì)算出來(lái)的各奇次諧波分量，對(duì)于這樣的N個(gè)數(shù)據(jù)，我們就可以利用碼位倒置的規(guī)律，最后得到X(k)的自然輸出序列。圖3給出了只有奇次諧波的16點(diǎn)DIF基2 FFT改進(jìn)算法流程圖。

圖3 輸出只有奇次諧波的16點(diǎn)DIF基2 FFT改進(jìn)算法

2.5 算法比較

在MATLAB下編寫了該改進(jìn)算法的程序——half_dif.m，并通過對(duì)方波信號(hào)的仿真分析，證明了該程序的正確性。表1給出了方波信號(hào)諧波的理論值、用一般FFT分析的諧波值以及用改進(jìn)算法進(jìn)行分析的諧波值，這里，為方便比較，只給出了其中直流分量、基波分量和2到15次諧波的幅值。

表1 方波信號(hào)（0—15次）諧波幅值表

通過比較可以發(fā)現(xiàn)，利用改進(jìn)算法進(jìn)行諧波分析，其奇次諧波跟一般FFT分析的結(jié)果是一致的，而偶次諧波則跟理論值是一致的，與前面的理論分析是相吻合的，從而證明了該程序的正確性。

3 電力系統(tǒng)諧波分析實(shí)驗(yàn)

本文利用matlab里面的GUI向?qū)гO(shè)計(jì)器設(shè)計(jì)了一個(gè)用于電網(wǎng)諧波測(cè)量分析的人機(jī)交互界面，可以實(shí)現(xiàn)如下功能：

（1）可以由用戶確定濾波器截止頻率、階數(shù)以及濾波器型號(hào)，其中包含了三種基本濾波器：巴特沃思濾波器、切比雪夫I型濾波器和橢圓濾波器。

（2）由用戶輸入一周期采樣點(diǎn)數(shù)N。

（3）實(shí)現(xiàn)對(duì)一些基本信號(hào)的諧波測(cè)量分析，

（1），（2）所確定的濾波器以及采樣點(diǎn)數(shù)都是針對(duì)這些基本信號(hào)的?；拘盘?hào)包括：方波信號(hào)、矩形信號(hào)和正弦加3次、7次諧波信號(hào)，其中正弦信號(hào)幅值為1，相位為0，而3次、7次諧波幅值均為0.5，3次諧波相位為π/2，7次諧波相位為0，如圖4所示。（1），（2）都是用于基本信號(hào)的諧波測(cè)量分析。

圖4 三種基本信號(hào)波形選擇界面圖

實(shí)現(xiàn)對(duì)實(shí)時(shí)信號(hào)的諧波測(cè)量分析，如圖5所示。其中圖5（a）為對(duì)逆變器一相輸出電壓的測(cè)量波形，圖5（b）為諧波分析結(jié)果。

圖5（a） 逆變器輸出電壓時(shí)域波形顯示界面圖

圖5(b) 各次諧波含有率波形顯示界面圖（正弦加3，7次諧波信號(hào)）

4 結(jié)論

本文根據(jù)一類電網(wǎng)信號(hào)具有半波對(duì)稱性的特點(diǎn)，提出了一種用于電網(wǎng)諧波測(cè)量改進(jìn)的快速傅里葉變換算法，該算法的計(jì)算量是傳統(tǒng)FFT算法的一半，因此大大減少了運(yùn)行時(shí)間；在對(duì)改進(jìn)算法進(jìn)行詳細(xì)推導(dǎo)的基礎(chǔ)上，詳細(xì)介紹了算法的實(shí)現(xiàn)過程；并完成了在matlab下對(duì)改進(jìn)算法的程序的編寫，該程序在matlab里面可以像其工具箱里面的FFT函數(shù)一樣直接被調(diào)用；進(jìn)而通過電網(wǎng)諧波分析實(shí)驗(yàn)，驗(yàn)證了所提方法的可行性。

[1] 胡廣書.數(shù)字信號(hào)處理――理論、算法與實(shí)現(xiàn)[M].北京：清華大學(xué)出版社，2003.

[2] Artyom M Grigoryan，Veerabhadra S Bhamidipati.Method of Flow Graph Simplification for the 16-Point Discrete Fourier Transform[J].IEEE Transaction on Signal Process,2005,53(1):384-359.

[3] George J Wakileh（奧地利）著，徐政譯.電力系統(tǒng)諧波——基本原理、分析方法和濾波器設(shè)計(jì)[M].北京：機(jī)械工業(yè)出版社，2003:4-21.

[4]Marek ROGó z˙.Model of the Harmonic Analyzer for the Needs of Power Quality Assessment[A].In:Proceedings of the 8thInternational Conference on Electrical Power Quality and Utilisation[C].Cracow,Poland,2005:21-23.

[5]Hui Zheng,Zhao Chen,Dongkai Peng and Tianhao Tang.A Variable N Synchronous Sampling Method for Grid Aperiodic Harmonic Detection and Analysis[A].In:Proceedings of ICFPEE2010[C].Shenzhen,China,2010.