%*************************************************************************%
%                              FFT實踐及頻譜分析                          %
%*************************************************************************%
%*************************************************************************%
%***************1.正弦波****************%
fs=100;%設定采樣頻率
N=128;
n=0:N-1;
t=n/fs;
f0=10;%設定正弦信號頻率
%生成正弦信號
x=sin(2*pi*f0*t);
figure(1);
subplot(231);
plot(t,x);%作正弦信號的時域波形
xlabel('t');
ylabel('y');
title('正弦信號y=2*pi*10t時域波形');
grid;
%進行FFT變換并做頻譜圖
y=fft(x,N);%進行fft變換
mag=abs(y);%求幅值
f=(0:length(y)-1)'*fs/length(y);%進行對應的頻率轉換
figure(1);
subplot(232);
plot(f,mag);%做頻譜圖
axis([0,100,0,80]);
xlabel('頻率(Hz)');
ylabel('幅值');
title('正弦信號y=2*pi*10t幅頻譜圖N=128');
grid;
%求均方根譜
sq=abs(y);
figure(1);
subplot(233);
plot(f,sq);
xlabel('頻率(Hz)');
ylabel('均方根譜');
title('正弦信號y=2*pi*10t均方根譜');
grid;
%求功率譜
power=sq.^2;
figure(1);
subplot(234);
plot(f,power);
xlabel('頻率(Hz)');
ylabel('功率譜');
title('正弦信號y=2*pi*10t功率譜');
grid;
%求對數譜
ln=log(sq);
figure(1);
subplot(235);
plot(f,ln);
xlabel('頻率(Hz)');
ylabel('對數譜');
title('正弦信號y=2*pi*10t對數譜');
grid;
%用IFFT恢復原始信號
xifft=ifft(y);
magx=real(xifft);
ti=[0:length(xifft)-1]/fs;
figure(1);
subplot(236);
plot(ti,magx);
xlabel('t');
ylabel('y');
title('通過IFFT轉換的正弦信號波形');
grid;
%****************2.矩形波****************%
fs=10;%設定采樣頻率
t=-5:0.1:5;
x=rectpuls(t,2);
x=x(1:99);
figure(2);
subplot(231);
plot(t(1:99),x);%作矩形波的時域波形
xlabel('t');
ylabel('y');
title('矩形波時域波形');
grid;
%進行FFT變換并做頻譜圖
y=fft(x);%進行fft變換
mag=abs(y);%求幅值
f=(0:length(y)-1)'*fs/length(y);%進行對應的頻率轉換
figure(2);
subplot(232);
plot(f,mag);%做頻譜圖
xlabel('頻率(Hz)');
ylabel('幅值');
title('矩形波幅頻譜圖');
grid;
%求均方根譜
sq=abs(y);
figure(2);
subplot(233);
plot(f,sq);
xlabel('頻率(Hz)');
ylabel('均方根譜');
title('矩形波均方根譜');
grid;
%求功率譜
power=sq.^2;
figure(2);
subplot(234);
plot(f,power);
xlabel('頻率(Hz)');
ylabel('功率譜');
title('矩形波功率譜');
grid;
%求對數譜
ln=log(sq);
figure(2);
subplot(235);
plot(f,ln);
xlabel('頻率(Hz)');
ylabel('對數譜');
title('矩形波對數譜');
grid;
%用IFFT恢復原始信號
xifft=ifft(y);
magx=real(xifft);
ti=[0:length(xifft)-1]/fs;
figure(2);
subplot(236);
plot(ti,magx);
xlabel('t');
ylabel('y');
title('通過IFFT轉換的矩形波波形');
grid;
%****************3.白噪聲****************%
fs=10;%設定采樣頻率
t=-5:0.1:5;
x=zeros(1,100);
x(50)=100000;
figure(3);
subplot(231);
plot(t(1:100),x);%作白噪聲的時域波形
xlabel('t');
ylabel('y');
title('白噪聲時域波形');
grid;
%進行FFT變換并做頻譜圖
y=fft(x);%進行fft變換
mag=abs(y);%求幅值
f=(0:length(y)-1)'*fs/length(y);%進行對應的頻率轉換
figure(3);
subplot(232);
plot(f,mag);%做頻譜圖
xlabel('頻率(Hz)');
ylabel('幅值');
title('白噪聲幅頻譜圖');
grid;
%求均方根譜
sq=abs(y);
figure(3);
subplot(233);
plot(f,sq);
xlabel('頻率(Hz)');
ylabel('均方根譜');
title('白噪聲均方根譜');
grid;
%求功率譜
power=sq.^2;
figure(3);
subplot(234);
plot(f,power);
xlabel('頻率(Hz)');
ylabel('功率譜');
title('白噪聲功率譜');
grid;
%求對數譜
ln=log(sq);
figure(3);
subplot(235);
plot(f,ln);
xlabel('頻率(Hz)');
ylabel('對數譜');
title('白噪聲對數譜');
grid;
%用IFFT恢復原始信號
xifft=ifft(y);
magx=real(xifft);
ti=[0:length(xifft)-1]/fs;
figure(3);
subplot(236);
plot(ti,magx);
xlabel('t');
ylabel('y');
title('通過IFFT轉換的白噪聲波形');
grid;