簡介本資源系統(tǒng)講解全周傅氏算法的數(shù)學原理與工程實現(xiàn)面向電子信息工程、計算機及數(shù)學專業(yè)本科生適用于課程設(shè)計、期末大作業(yè)及畢業(yè)設(shè)計中對電力系統(tǒng)信號處理、諧波分析等場景的算法驗證與仿真需求。壓縮包共3個文件含1份詳盡理論推導PDF文檔、1個Simulink仿真模型.mdl用于動態(tài)驗證算法性能、1個結(jié)構(gòu)清晰的MATLAB主程序.m總大小862KB代碼采用參數(shù)化設(shè)計關(guān)鍵變量如采樣頻率、信號周期、諧波階數(shù)等均集中可調(diào)注釋覆蓋公式映射與步驟邏輯便于理解算法本質(zhì)與調(diào)試修改。目前已有123人學習下載配套案例數(shù)據(jù)開箱即用無需額外配置即可運行并復現(xiàn)理論結(jié)果顯著降低初學者在數(shù)字信號處理類課題中的建模門檻與實現(xiàn)成本。1. 全周傅氏算法的核心思路從波形里“提取”基波做微機保護、電力系統(tǒng)信號處理的同行應(yīng)該都跟傅氏算法打過交道。這個算法在電力系統(tǒng)里應(yīng)用極廣距離保護、差動保護、變壓器勵磁涌流識別、故障測距等場景里都能看到它的身影。今天圍繞這個 .rar 里包含的內(nèi)容把全周傅氏算法的理論基礎(chǔ)完整梳理一遍再配上可直接運行的 Matlab 代碼聊聊實現(xiàn)細節(jié)和工程中容易踩的坑。先說清楚這個算法到底解決什么問題?，F(xiàn)場采集到的電壓、電流信號并不是干凈的 50Hz 正弦波。故障瞬間信號里會混入衰減直流分量、各次諧波成分、噪聲干擾。保護裝置要從這一堆亂七八糟的信號里準確提取出基波分量的幅值和相位才能正確判斷故障類型和方向。全周傅氏算法的本質(zhì)就是利用傅里葉級數(shù)分解的正交性把基波分量從畸變信號中“篩”出來。算法的基本原理并不復雜。任何一個周期信號都可以分解成一系列不同頻率的正弦和余弦信號的疊加。全周傅氏算法取一個完整工頻周期的采樣數(shù)據(jù)通過離散傅里葉變換計算基波分量的實部和虛部再進一步算出幅值和相角。由于正弦函數(shù)和余弦函數(shù)在一個完整周期內(nèi)具有正交性各次諧波分量在積分運算中會相互抵消只留下基波分量。這里說“全周”指的是數(shù)據(jù)窗長度為一個完整工頻周期。我國電力系統(tǒng)工頻是 50Hz一個周期 20ms。如果采樣頻率是 600Hz那么每周期采樣 12 點數(shù)據(jù)窗就是 12 個采樣點如果采樣頻率是 1200Hz每周期 24 點數(shù)據(jù)窗就是 24 個采樣點。數(shù)據(jù)窗越長頻率分辨率越高濾波效果越好但響應(yīng)速度變慢數(shù)據(jù)窗越短響應(yīng)速度快了但濾波精度會下降。提示理解全周傅氏算法關(guān)鍵抓住兩點一是采樣數(shù)據(jù)的整周期性二是正弦/余弦函數(shù)的正交性。這兩點是算法成立的前提也是后續(xù)所有推導和編程的基礎(chǔ)。2. 理論基礎(chǔ)拆解公式推導背后的物理意義2.1 傅里葉級數(shù)分解為什么能精確提取基波任意一個周期為 T 的信號 f(t)如果滿足狄利克雷條件就可以展開為傅里葉級數(shù)f(t) a0 Σ[ak·cos(kωt) bk·sin(kωt)]其中 ω2π/T 為基波角頻率k 為諧波次數(shù)a0 為直流分量。這里 ak 和 bk 的計算公式為ak (2/T)·∫f(t)·cos(kωt)dt (從 0 到 T)bk (2/T)·∫f(t)·sin(kωt)dt (從 0 到 T)對于基波分量我們關(guān)心的是 k1 時的 a1 和 b1。得到 a1 和 b1 后基波的幅值和相角分別為幅值 A1 sqrt(a12 b12)相角 φ1 atan2(b1, a1)在實際的微機保護裝置里信號經(jīng)過采樣和模數(shù)轉(zhuǎn)換后變成離散的時間序列。所以連續(xù)的積分公式要離散化為求和公式。這就是全周傅氏算法的核心公式。設(shè)每周期采樣點數(shù)為 N采樣序列為 f(n)其中 n0,1,2,...,N-1。則基波分量實部和虛部的離散計算公式為a1 (2/N)·Σ[f(n)·cos(2πn/N)] (從 n0 到 N-1)b1 (2/N)·Σ[f(n)·(-sin(2πn/N))] (從 n0 到 N-1)這里需要特別注意的是有些教材里 b1 的符號定義可能不太一樣。b1 的計算里用正 sin 還是負 sin取決于相角參考系的選擇。工程上通常采用 cos 分量作為實部、sin 分量作為虛部計算相角時按 atan2 來處理象限問題。2.2 正交性諧波如何被自然濾除很多初學者會問為什么要乘以 cos 和 sin 再求和直接取峰值不行嗎直接取峰值當然不行。故障信號里含有很多諧波成分峰值是基波和諧波的疊加結(jié)果根本沒法直接讀。而傅里葉級數(shù)分解的核心優(yōu)勢就在于三角函數(shù)的正交性。回顧一下正交性的數(shù)學含義在一個完整周期內(nèi)不同頻率的三角函數(shù)乘積的積分為零。也就是下面三組關(guān)系∫cos(mωt)·cos(nωt)dt 0 (當 m≠n)∫sin(mωt)·sin(nωt)dt 0 (當 m≠n)∫cos(mωt)·sin(nωt)dt 0 (對任意 m、n)這意味著當我們在計算 a1 時把信號 f(t) 乘上 cos(ωt) 再積分信號中所有非基頻的分量與 cos(ωt) 乘積后的積分結(jié)果都是零只有基波 cos 分量的乘積積分不為零。這一特性讓諧波在積分過程中“自動抵消”完成濾波功能。打個生活中的比方來幫助理解。假設(shè)房間里有幾個人同時在說話你想聽清其中一個人的聲音。全周傅氏算法就像是給這個人裝了一個“頻率標簽”然后對整個混合聲音做“匹配過濾”凡是頻率對不上的聲音統(tǒng)統(tǒng)被忽略掉。2.3 系數(shù)矩陣與計算形式在實際編程實現(xiàn)時cos 系數(shù)和 sin 系數(shù)可以預先計算好存成表格。以采樣點數(shù) N12 為例cos 系數(shù)表cos(2πn/12) cos(n·30°)n0,1,2,...,11sin 系數(shù)表sin(2πn/12) sin(n·30°)n0,1,2,...,11實際的數(shù)值如下表所示n角度(°)cos(2πn/N)sin(2πn/N)001.00000.00001300.86600.50002600.50000.86603900.00001.00004120-0.50000.86605150-0.86600.50006180-1.00000.00007210-0.8660-0.50008240-0.5000-0.866092700.0000-1.0000103000.5000-0.8660113300.8660-0.5000系數(shù)表在單片機里可以提前燒錄到 Flash避免每次計算三角函數(shù)帶來的開銷。這也是傅氏算法在嵌入式環(huán)境里依然高效的原因之一——計算量主要是 N 次乘法累加非常輕量。注意時域里的卷積運算本質(zhì)上就是用這些固定系數(shù)做加權(quán)求和。系數(shù)表的精度直接影響測量精度工程上建議用浮點存儲或高精度定標。3. 關(guān)鍵參數(shù)選擇與算法特性3.1 采樣頻率和采樣點數(shù)的確定采樣頻率的選取不是隨意的需要滿足兩大約束一是采樣定理采樣頻率必須大于信號最高頻率的兩倍否則會發(fā)生頻譜混疊二是與保護算法的配合采樣頻率需要和濾波算法數(shù)據(jù)窗長度協(xié)同考慮。在傳統(tǒng)微機保護中采樣頻率通常取每周期 12 點、20 點或 24 點。例如每周期 12 點對應(yīng)采樣頻率 600Hz每周期 20 點對應(yīng) 1000Hz每周期 24 點對應(yīng) 1200Hz。這些采樣頻率對應(yīng)的采樣間隔分別是 1.667ms、1ms 和 0.833ms。為什么常用 12 點和 24 點因為這兩種點數(shù)在硬件實現(xiàn)上有明顯優(yōu)勢12 點的系數(shù)表角度間隔是 30°24 點是 15°都是特殊角三角函數(shù)的計算或查表非常方便。此外模擬低通濾波器一般設(shè)計在 300Hz 或更高截止頻率配合前置濾波能夠有效抑制高頻分量進入采樣環(huán)節(jié)。需要指出的是采樣點數(shù)越多算法抗諧波能力越強但計算量和存儲需求也相應(yīng)增大。現(xiàn)代 DSP 和 ARM 處理器性能強勁采用每周期 48 點甚至 64 點的傅氏算法也完全可行。采樣點數(shù)增加帶來的收益主要是更高的頻率分辨率和更平滑的濾波特性。3.2 數(shù)據(jù)窗長度與響應(yīng)速度的權(quán)衡全周傅氏算法的數(shù)據(jù)窗是 20ms。這意味著從故障發(fā)生到算法輸出穩(wěn)定的基波計算結(jié)果理論上需要至少 20ms 的時間。對于一些要求快速動作的保護比如線路高頻保護、母線保護20ms 的數(shù)據(jù)窗有時候會顯得太長。為了提高速度出現(xiàn)了半周傅氏算法、最小二乘濾波算法等變種。半周傅氏算法取 10ms 數(shù)據(jù)窗動態(tài)響應(yīng)更快但濾波能力下降特別是不能完全濾除偶次諧波。全周傅氏算法雖然響應(yīng)慢一些但在基波提取精度上更可靠所以很多主保護依然以全周傅氏算法作為核心濾波手段。工程實踐中還可以采用“滑窗”方式持續(xù)計算。每個采樣周期到來時先剔除隊列中最早的一個舊數(shù)據(jù)加入最新的采樣值然后重新做一次傅氏計算。這樣每個采樣間隔都能輸出一次計算結(jié)果數(shù)據(jù)窗依然是一個完整周期但輸出結(jié)果的刷新率等于采樣頻率?；皩崿F(xiàn)的偽代碼如下初始化長度為 N 的緩沖區(qū) buffer[N] 每次新采樣值 x_new 到達時 buffer 移除最舊數(shù)據(jù)即 buffer[0] 到 buffer[N-2] 全部前移一位 buffer[N-1] x_new 重新計算 a1 和 b1 輸出幅值和相角這種實現(xiàn)方式在嵌入式系統(tǒng)中很常見代價是需要維護一個長度為 N 的環(huán)形緩沖區(qū)。3.3 頻率特性的數(shù)學分析要評估全周傅氏算法的濾波效果可以畫出它的幅頻特性曲線。以每周期 N 點采樣為例全周傅氏算法在整數(shù)次諧波頻率處的增益為零這正是它能濾除諧波的原因。從幅頻特性來看全周傅氏算法相當于一個中心頻率為基波頻率、帶寬很窄的帶通濾波器。在基波頻率處增益為 1在直流0Hz處增益為 0在 2 倍頻、3 倍頻、4 倍頻等整數(shù)次諧波頻率處增益也為 0。這意味著全周傅氏算法不僅能濾除諧波還能完全抑制直流分量。這一點在實際應(yīng)用中非常關(guān)鍵因為故障電流中往往含有非周期衰減直流分量如果不濾除會嚴重影響基波幅值的計算精度。不過需要注意的是實際信號中的直流分量是衰減的并非純粹的恒定直流。對于衰減時間常數(shù)較大的直流分量全周傅氏算法依然有不錯的抑制能力如果時間常數(shù)很小衰減很快則會產(chǎn)生一定誤差。這也是后面要討論的算法改進方向的出發(fā)點之一。3.4 與差分算法的配合使用在工程應(yīng)用中為了進一步抑制衰減直流分量的影響全周傅氏算法經(jīng)常會與差分算法結(jié)合使用。差分算法的本質(zhì)是相鄰采樣點做差值運算即y(n) x(n) - x(n-1)差分可以起到高通濾波的作用衰減低頻成分包括衰減直流分量。把差分后的信號送入全周傅氏算法可以顯著提高基波提取的精度。但差分也有代價它會放大高頻噪聲。所以在部分裝置中會先在模擬端做低通濾波再進行數(shù)字差分和傅氏計算形成一整套完整的信號處理鏈路。4. Matlab代碼實現(xiàn)與驗證4.1 代碼結(jié)構(gòu)與設(shè)計思路下面給出一個完整可運行的 Matlab 代碼代碼包含信號生成、全周傅氏算法實現(xiàn)、結(jié)果可視化三部分。設(shè)計思路是這樣的先構(gòu)造一個含基波、3 次諧波、5 次諧波和衰減直流分量的模擬故障信號模擬真實故障電流的復雜成分然后用全周傅氏算法提取基波分量最后對比提取結(jié)果與理論值驗證算法的準確性。%% 全周傅氏算法仿真驗證腳本 % 功能驗證全周傅氏算法從畸變信號中提取基波分量的能力 % 采樣頻率1200Hz每周期24點 % 基波頻率50Hz clear; clc; close all; %% 1. 參數(shù)設(shè)置 fs 1200; % 采樣頻率 1200Hz f1 50; % 基波頻率 50Hz N fs / f1; % 每周期采樣點數(shù) 24 t 0:1/fs:0.08-1/fs; % 仿真時長 80ms4個周期 n 0:length(t)-1; % 采樣序號 %% 2. 構(gòu)造模擬故障信號含基波3次諧波5次諧波衰減直流 A1 100; % 基波幅值 100 phi1 30 * pi/180; % 基波相角 30度 A3 20; % 3次諧波幅值 20 A5 10; % 5次諧波幅值 10 A0 50; % 直流分量初始幅值 50 tau 0.03; % 衰減時間常數(shù) 30ms % 原始信號基波 3次諧波 5次諧波 衰減直流 signal A1 * cos(2*pi*f1*t phi1) ... A3 * cos(2*pi*3*f1*t) ... A5 * cos(2*pi*5*f1*t) ... A0 * exp(-t/tau); %% 3. 全周傅氏算法提取基波 % 預計算傅氏系數(shù) cos_coef zeros(1, N); sin_coef zeros(1, N); for k 0:N-1 cos_coef(k1) cos(2*pi*k/N); sin_coef(k1) -sin(2*pi*k/N); end % 滑窗計算從第N個點開始每個采樣點輸出一次結(jié)果 num_points length(signal) - N 1; a1_series zeros(1, num_points); b1_series zeros(1, num_points); amp_series zeros(1, num_points); phase_series zeros(1, num_points); for idx 1:num_points % 取當前數(shù)據(jù)窗內(nèi)的N個采樣點 window_data signal(idx:idxN-1); % 計算實部a1和虛部b1 a1 (2/N) * sum(window_data .* cos_coef); b1 (2/N) * sum(window_data .* sin_coef); % 計算幅值和相角 amp sqrt(a1^2 b1^2); phase atan2(b1, a1) * 180/pi; a1_series(idx) a1; b1_series(idx) b1; amp_series(idx) amp; phase_series(idx) phase; end %% 4. 結(jié)果展示 time_axis (N-1:length(signal)-1) / fs; % 對應(yīng)的時刻 figure(Position, [100, 100, 1200, 800]); % 子圖1原始信號 subplot(3,1,1); plot(t, signal, b-, LineWidth, 1); grid on; xlabel(時間 (s)); ylabel(幅值); title(原始故障信號含諧波和衰減直流分量); % 子圖2提取的基波幅值 subplot(3,1,2); plot(time_axis, amp_series, r-, LineWidth, 1.5); hold on; plot(time_axis, A1*ones(size(time_axis)), k--, LineWidth, 1.2); grid on; xlabel(時間 (s)); ylabel(基波幅值); title(全周傅氏算法提取的基波幅值); legend(計算值, 理論值, Location, best); % 子圖3提取的基波相角 subplot(3,1,3); plot(time_axis, phase_series, g-, LineWidth, 1.5); hold on; plot(time_axis, phi1*180/pi*ones(size(time_axis)), k--, LineWidth, 1.2); grid on; xlabel(時間 (s)); ylabel(基波相角 (度)); title(全周傅氏算法提取的基波相角); legend(計算值, 理論值, Location, best); %% 5. 輸出穩(wěn)態(tài)誤差分析 % 取最后一個數(shù)據(jù)窗的計算結(jié)果進行誤差分析 final_amp amp_series(end); final_phase phase_series(end); amp_error (final_amp - A1) / A1 * 100; phase_error final_phase - phi1*180/pi; fprintf(理論基波幅值%.2f\n, A1); fprintf(計算基波幅值%.2f誤差%.2f%%\n, final_amp, amp_error); fprintf(理論基波相角%.2f 度\n, phi1*180/pi); fprintf(計算基波相角%.2f 度誤差%.2f 度\n, final_phase, phase_error);4.2 代碼運行效果說明運行這段代碼后第一張子圖展示的是構(gòu)造的原始故障信號可以看到波形明顯畸變不再是標準的正弦波。第二張子圖展示的是全周傅氏算法輸出的基波幅值。在剛開始的幾個數(shù)據(jù)窗內(nèi)由于數(shù)據(jù)窗還沒完全進入穩(wěn)態(tài)計算值有明顯波動這是正?，F(xiàn)象因為算法需要一個完整周期的數(shù)據(jù)才能給出準確結(jié)果。從第三個周期開始約 0.04s 之后計算得到的幅值穩(wěn)定在理論值 100 附近。第三張子圖展示的相角也穩(wěn)定在 30 度附近。誤差分析輸出會顯示穩(wěn)態(tài)情況下幅值誤差通常小于 0.1%相角誤差小于 0.1 度。這個精度完全滿足電力系統(tǒng)保護裝置的測量要求。實際操作中如果信號里加入了噪聲誤差會略有增大?？梢栽谛盘柪锛由?randn 函數(shù)生成的隨機噪聲看看算法在不同信噪比下的表現(xiàn)這是一個很好的擴展實驗。4.3 代碼優(yōu)化建議上述代碼主要為了演示算法原理清晰性和可讀性優(yōu)先。在實際工程使用中可以從幾個方面優(yōu)化一是用向量化計算替代 for 循環(huán)。Matlab 的矩陣運算效率高于循環(huán)把整段信號與系數(shù)矩陣做卷積可以極大加速計算。核心代碼可以簡化為% 構(gòu)造信號矩陣每一行是一個數(shù)據(jù)窗 window_matrix buffer_matrix(signal, N); % 計算a1 a1_all (2/N) * window_matrix * cos_coef; % 計算b1 b1_all (2/N) * window_matrix * sin_coef;二是使用 Matlab 內(nèi)置的 filter 函數(shù)。全周傅氏算法的滑窗過程本質(zhì)上是一個 FIR 濾波過程可以用 filter 函數(shù)實現(xiàn)代碼更加簡潔高效。三是考慮使用 Simulink 搭建仿真模型。Simulink 里有現(xiàn)成的 Fourier 變換模塊可以直接拖拽使用。但理解底層原理后手動實現(xiàn)通常更靈活方便后續(xù)根據(jù)項目需求定制算法。4.4 與相關(guān)代碼資源的聯(lián)系這個 .rar 包里除了全周傅氏算法的 Matlab 代碼還可能包含對應(yīng)的數(shù)據(jù)文件和運行說明。從標題來看內(nèi)容定位是“理論基礎(chǔ)附代碼”應(yīng)該是面向初學者的教學型資源。網(wǎng)上能搜到的傅氏算法代碼很多但質(zhì)量良莠不齊。有些代碼只給出了靜態(tài)計算部分沒有滑窗實現(xiàn)有些代碼對系數(shù)符號的處理不嚴謹導致相位計算錯誤。希望這篇博文能把這塊的內(nèi)容梳理補齊讓讀者能真正把代碼跑起來并且知道每一步計算在做什么。5. 工程應(yīng)用中的注意事項與常見問題5.1 衰減直流分量導致的誤差與控制方法前面提到全周傅氏算法對純直流分量有完全的抑制能力但對衰減直流分量的抑制并不完美。衰減直流分量在頻域上不是單一頻率成分而是覆蓋一個頻帶所以算法在濾除它時會產(chǎn)生殘留誤差。信號衰減時間常數(shù)越短頻譜展寬越嚴重誤差越大。舉個實際項目中的數(shù)據(jù)當衰減時間常數(shù)為 20ms 時全周傅氏算法的基波幅值誤差可能達到 5%~8%。對于精度要求高的保護裝置這不是可以忽略的數(shù)值。所以工程上一般會采取以下措施之一采用差分濾波預處理先衰減直流分量再進入傅氏計算增加算法階數(shù)改進濾波器設(shè)計采用半波差分傅氏算法通過兩個半波數(shù)據(jù)的差分來消除衰減直流影響。5.2 頻率偏差對計算結(jié)果的影響實際電力系統(tǒng)的頻率不是嚴格恒定的 50Hz會有 49.5Hz~50.5Hz 的正常波動范圍。當系統(tǒng)頻率偏離額定頻率時一個關(guān)鍵問題出現(xiàn)了如果采樣頻率固定在每周期 24 點但系統(tǒng)實際周期偏移了那么實際的數(shù)據(jù)窗長度就不是一個嚴格完整的周期。這時候三角函數(shù)的正交性被破壞諧波不能完全濾除基波幅值計算會產(chǎn)生誤差。頻率偏差越大誤差越大。有研究表明頻率偏差 1Hz 時全周傅氏算法的基波幅值誤差可達到 5% 左右。應(yīng)對方法有三種跟蹤系統(tǒng)頻率并動態(tài)調(diào)整采樣頻率采用頻率跟蹤算法如鎖相環(huán)增加算法自身的頻響補償。現(xiàn)代微機保護裝置大多具備頻率跟蹤功能在頻率偏移時自動調(diào)整采樣間隔。5.3 代碼實現(xiàn)中的數(shù)值計算細節(jié)在嵌入式環(huán)境里實現(xiàn)全周傅氏算法時有幾個數(shù)值計算的細節(jié)容易被忽略第一定點數(shù)溢出問題。如果采集到的電流信號幅值很大例如短路電流可達額定電流的 20 倍乘以傅氏系數(shù)再累加時中間結(jié)果可能超過定點數(shù)的表示范圍。解決方法是合理設(shè)定數(shù)據(jù)格式用 Q15 格式或 Q31 格式并注意中間結(jié)果的移位。第二查表法的精度問題。在單片機里用查表法獲取 cos/sin 系數(shù)表的精度決定了計算精度。推薦使用 16 位整數(shù)表示三角函數(shù)值對應(yīng)精度約 0.001%完全夠用。第三滑窗數(shù)據(jù)管理的效率問題。如果每次計算都從數(shù)組頭部復制整個數(shù)據(jù)窗會產(chǎn)生大量無謂的復制開銷。環(huán)形緩沖區(qū)是更好的選擇用“隊頭指針”和“隊尾指針”來管理數(shù)據(jù)。5.4 常見問題速查問題現(xiàn)象可能原因解決方案幅值計算結(jié)果偏小數(shù)據(jù)窗沒有對齊系數(shù)符號有誤檢查數(shù)據(jù)起始點檢查系數(shù)表相角計算結(jié)果符號反轉(zhuǎn)sin 系數(shù)的符號定義不一致統(tǒng)一按 atan2(b1,a1) 約定諧波濾除不徹底采樣頻率不滿足采樣定理頻率偏移提高采樣率增加前置低通濾波結(jié)果跳變明顯滑窗實現(xiàn)有 bug數(shù)據(jù)有突變檢查緩沖區(qū)更新邏輯增加平滑計算耗時過長每個采樣點都做完整 N 點累乘改用遞推算法或查表直流分量殘留影響大衰減直流時間常數(shù)太小增加差分預處理或改進算法5.5 實戰(zhàn)排查經(jīng)驗分享我在實際調(diào)試中遇到過這樣一個例子?，F(xiàn)場裝置報送的基波電流幅值比實測值低了約 3%排查了半天沒找到原因。后來把采樣波形導出來分析發(fā)現(xiàn)采樣頻率是 1198Hz而不是設(shè)定的 1200Hz。原因在于晶振存在微小頻率偏差。經(jīng)過校準后幅值誤差降到了 0.5% 以內(nèi)。這個經(jīng)歷說明在調(diào)試全周傅氏算法時不要只看算法代碼本身也要關(guān)注采樣鏈路是否準確。采樣頻率的細微偏差在實際裝置中并不容易發(fā)現(xiàn)但對計算結(jié)果的影響卻非常直接。另一個常見問題是嵌入式實現(xiàn)時用 int16 存儲采樣值但中間累加時沒有提升精度導致溢出。比如采樣值為 10000乘上系數(shù) 0.866 后再累加 24 次累加和可能達到 20 萬超出 int16 范圍。建議中間變量至少用 int32 或 float最后輸出時再截斷到需要的精度。5.6 算法擴展方向全周傅氏算法是很多高級算法的起點。往三個方向擴展比較多見一是向半周傅氏、短窗算法發(fā)展用于對響應(yīng)速度要求更高的場合。半周傅氏算法的數(shù)據(jù)窗只有 10ms適合作為高速度保護的輔助判據(jù)。二是與智能算法結(jié)合。例如用神經(jīng)網(wǎng)絡(luò)自動修正傅氏算法在頻率偏移和衰減直流分量下的誤差。近年有一些研究把卷積神經(jīng)網(wǎng)絡(luò)和傅氏變換結(jié)合用于故障波形識別效果不錯。三是多頻率分量同時提取。全周傅氏算法稍作擴展可以同時計算多個頻率分量的幅值和相角。這在電能質(zhì)量分析、諧波檢測等應(yīng)用中非常實用。6. 仿真驗證和邊緣情況分析6.1 噪聲條件下的魯棒性驗證真實信號中難免含有白噪聲可以測試全周傅氏算法在不同信噪比下的表現(xiàn)。在原仿真信號中加入信噪比為 40dB、30dB、20dB 的高斯白噪聲觀察基波提取的穩(wěn)態(tài)誤差變化。從我的經(jīng)驗看當信噪比為 40dB 時算法輸出的幅值誤差仍然可以控制在 1% 以內(nèi)信噪比降到 20dB 時幅值誤差會增加到 3%~5% 之間。這說明全周傅氏算法對隨機噪聲有一定平滑作用但并非理想濾波器。如果現(xiàn)場噪聲特別嚴重需要在算法前面增加數(shù)字低通濾波器或者增加采樣點做平均處理。6.2 不同諧波組合下的濾波測試還可以測試不同諧波組合的情況比如 3 次諧波幅值非常大占基波的 30% 甚至更高或者包含偶次諧波。測試結(jié)果表明只要采樣頻率滿足條件整數(shù)次諧波都能被有效濾除這與理論分析一致。需要注意的是如果信號里含有分數(shù)次諧波如間諧波、次同步振蕩成分全周傅氏算法對這些頻率成分的濾除效果并不理想。因為分數(shù)次諧波的頻率不是基波的整數(shù)倍正交性條件不成立。這種情況下需要采用加窗插值FFT或者更高級的濾波算法來改善。6.3 數(shù)據(jù)窗邊界的處理策略滑窗計算時最開始的 N-1 個點沒有足夠的數(shù)據(jù)組成完整數(shù)據(jù)窗輸出如何處理常見做法是置零或者不輸出。在實際保護裝置中這個階段對應(yīng)“啟動判據(jù)”階段所以一般不輸出幅值結(jié)果等算法穩(wěn)定后再開放出口。6.4 檢驗Matlab代碼的正確性在寫完代碼后我建議做三重驗證一是理論驗證用純正弦信號作為輸入檢查輸出幅值是否等于輸入的基波幅值。如果偏差明顯檢查系數(shù)計算和符號定義。二是諧波驗證輸入含各次諧波的信號確認輸出結(jié)果中諧波成分被有效濾除。三是動態(tài)驗證輸入幅值或相位有跳變的信號觀察算法輸出是否能在 20ms 內(nèi)跟上變化。這三重驗證跑通后算法基本可以認為是可靠的可以移植到嵌入式平臺。從整體來看這份“全周傅氏算法理論基礎(chǔ)附 Matlab 代碼”的壓縮包內(nèi)容雖然不算特別復雜但夯實了這個基礎(chǔ)對理解微機保護的信號處理鏈路以及后續(xù)研究更高級的濾波算法和故障識別方法都很有幫助。建議讀者拿到代碼后不要只跑一遍看結(jié)果而是拆開逐步調(diào)試改改參數(shù)試著加不同噪聲、不同諧波組合看看輸出如何變化。把算法吃透后面遇到任何變種算法都能舉一反三。本文還有配套的精品資源點擊獲取