FFT和小波包变换程序

%本程序为实现论文“基于FFT和小波包变换的电力系统谐波检测方法”
%%% 用小波包重构时，不能很好的分辨出50HZ、100HZ、200HZ的信号。重构时，信号出现在第5层1、2、4节点处。

clear all;
N=4800;   %采样点为480个
fs=4800;  %采样频率为4800Hz
n=0:N-1;  %序列
t=n/fs;  %时间序列
x=220*sin（2*pi*50*t）+80*sin（2*pi*100*t）+50*sin（2*pi*200*t）;  %信号表达式
figure（1）;
subplot（411）;plot（tx）;
title（‘原始信号‘）;
xlabel（‘时间‘）;ylabel（‘幅值‘）;
grid on;
%**************************************************对信号进行小波包分析

T=wpdec（x5‘db44‘）;  %对信号进行小波包分解，得到分解树
ca1=wprcoef（T[50]）;  %重构信号
subplot（412）;plot（tca1）;
ylabel（‘50Hz‘）;
ca2=wprcoef（T[51]）;
subplot（413）;plot（tca2）;
ylabel（‘100Hz‘）;
ca3=wprcoef（T[53]）;
subplot（414）;plot（tca3）;
ylabel（‘200Hz‘）;


%**********************************************************进行傅里叶频谱分析
y1=fft（ca1N）;  %对重构的50、100、200HZ的信号进行傅里叶频谱分析。
y2=fft（ca2N）;
y3=fft（ca3N）;
f=n*fs/N;
mag1=abs（y1）;  %取绝对值。傅里叶频谱分析的结果显示：重构的信号频率与预期信号相同，但幅值发生了衰减。
mag2=abs（y2）;
mag3=abs（y3）;
figure（2）;
subplot（311）;plot（f（1:N/10）mag1（1:N/10）*2/N）;
grid on;
title（‘50HZ的傅里叶变换频谱图‘）;
%xlim（[0 400]）;限制x轴的取值范围 ylim（[0 400]）限制y轴的取值范围，axis（[xmin xmax ymin
%ymax]）同时限制x、y的取值范围。
subplot（312）;plot（f（1:N/10）mag2（1:N/10）*2/N）;
grid on;
title（‘100HZ的傅里叶变换频谱图‘）;
subplot（313）;plot（f（1:N/10