clc;
clear all;
clf;
Fs=100000;%采样频率
n=10;%周期数。采样时间为0.1s时n=5;采样时间为0.2s时n=10。
t=[0:Fs/50*n]/Fs;
%x=(6*exp(-20*t).*sin(1000*pi*(t-0.02))+15*sin(100*pi*t)).*(t>=0.02);%输入信号
% x=5*cos(100*pi*t+pi/6)%+3*cos(300*pi*t+pi/2).*exp(-1200*t);
x=sqrt(2)*sin(100*pi*t+pi/6)+0.5*sin(300*pi*t+pi/3)+0.3*sin(475*pi*t)...
+0.2*sin(1100*pi*t+pi/4);
N=Fs*n/50;%采样点总数
Y = fft(x); %做FFT变换
Ayy = (abs(Y)); %取模
plot(Ayy(1:N)); %显示原始的FFT模值结果
title('FFT 模值');
figure
plot(t,x);
figure;
Ayy=Ayy/(N/2); %换算成实际的幅度
Ayy(1)=Ayy(1)/2;
F=([1:N]-1)*Fs/N; %换算成实际的频率值
plot(F(1:N/2),Ayy(1:N/2)); %显示换算后的FFT模值结果
title('幅度-频率曲线图');
% set(gca,'XLim',[0 200]);%Y轴的数据显示范围
% %傅里叶变换
% y = fft(x);
% m = abs(y)/length(y);
% m = m(1:floor(length(m)/2)+1);
% m(2:end) = m(2:end)*2;
% p = angle(y);
% p = p(1:floor(length(p)/2)+1);
% p=p*180/pi;
% f = [0:length(m)-1]/length(y)*Fs;
%
% subplot(2,3,4);
% % bar(f(1:15),m(1:15))
% maxf=600;
% index = find(f<=maxf);
% bar(f(index),p(index))
figure;
Pyy=[1:N/2];
for i=1:N/2
% Pyy(i)=phase(Y(i)); %计算相位
Pyy(i)=angle(Y(i)); %计算相位
Pyy(i)=Pyy(i)*180/pi; %换算为角度
end;
bar(F(1:N/2),Pyy(1:N/2)); %显示相位图
title('相位-频率曲线图');
grid on
% set(gca,'XLim',[0 200]);%Y轴的数据显示范围