供电系统中谐波的产生，危害及治理(4)

for k=1 : m

x=x+Am(k)*cos(2*pi*f0*k*t+PH(k)/180*pi);

end; figure(1); plot(x);

图5.3.1是用Matlab仿真出的谐波电压信号图：

图5.3.1 谐波电压信号图

（2）当取采样频率为3000HZ，数据长度为1024采样点时，利用FFT编写Matlab语言得到信号频谱的程序如下： %dierge.m N=length(x); f=50; m=13; n=1:N; Fn=zeros(9,3);

An=zeros(9,2); Pn=zeros(9,2); for i=1 : m

Fn(i,1) = i; Fn(i,2) = f * i; An(i,1) = i; Pn(i,1) = i; end fsN=fs/N; f0=50; X=fft(x); X=X(1:N); Xabs=abs(X); for i= 1 : m

[Amax,index]=TriFind(Xabs,floor((i*f0-15)/fsN),ceil((i*f0+15)/fsN)); if(index==-1) Fn(i,3) = 0; An(i,2) = 0; Pn(i,2) = 0; else

if(Xabs(index-1) > Xabs(index+1))

a1 = Xabs(index-1) / Xabs(index); r1 = 1/(1+a1); k01 = index -1; else

a1 = Xabs(index) / Xabs(index+1); r1 = 1/(1+a1); k01 = index; end

Fn(i,3) = (k01+r1-1)*fs/N;

An(i,2) = 2*pi*r1*Xabs(k01)/(N*sin(r1*pi)); Pn(i,2) = phase(X(k01))-pi*r1; Pn(i,2) = Pn(i,2)*180/pi; Pn(i,2) = mod(Pn(i,2),180); end end Fn An Pn

dlmwrite('An.txt',An,'delimiter', '\\t','precision', '%.10f','newline', 'pc'); dlmwrite('Fn.txt',Fn,'delimiter', '\\t','precision', '%.10f','newline', 'pc');

dlmwrite('Pn.txt',Pn,'delimiter', '\\t','precision', '%.10f','newline', 'pc');

经过FFT变换之后，得到的谐波频率（Fn）、幅值（An）、相位（Pn）如下：

将仿真后得到的结果和表1中的数据进行比对，发现通过快速傅里叶变换得到的谐波频率（Fn）、幅值（An）、相位（Pn）都与原测量值十分接近，这说明，快速傅里叶变换法的运算精度是非常高的，它也是一种非常实用的谐波分析的方法。