1

2

3

4

5

6

7

8

9

10

11

12

%% FFT

clear;clc;close all

Fs=1000;    % 採集頻率

T=1/Fs;     % 採集時間間隔

N=2000;     % 採集訊號的長度--取樣點數

f1=33;      % 第一個餘弦訊號的頻率

f2=200;     % 第二個餘弦訊號的頻率

t=(0:1:N-1)*T;  % 定義整個採集時間點

t=t';           % 轉置成列向量

y=1.2+2.7*cos(2*pi*f1*t+pi/4)+5*cos(2*pi*f2*t+pi/6);  % 時域訊號

1

2

3

4

5

6

%% 繪製時域訊號

figure

plot(t,y)

xlabel('時間')

ylabel('訊號值')

title('時域訊號')

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

%% fft變換

Y=fft(y);    % Y為fft變換的結果，為複數向量

A=abs(Y);    % 複數的幅值（模）

RE=real(Y);  % 複數的實部

IM=imag(Y);  % 複數的虛部

%% 繪製fft變換結果（幅值，實部，虛部）

figure

subplot(3,1,1)

plot(0:1:N-1,A)

xlabel('序號 0 ~ N-1')

ylabel('幅值')

grid on

%% 頻域只讀取0:1:N/2

subplot(3,1,2)

plot(0:1:N-1,RE)

xlabel('序號 0 ~ N-1')

ylabel('實部')

grid on

subplot(3,1,3)

plot(0:1:N-1,IM)

xlabel('序號 0 ~ N-1')

ylabel('虛部')

grid on

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

%% fft變換

Y=fft(y);          % Y為fft變換結果，複數向量

Y=Y(1:N/2+1);      % 只看變換結果的一半即可

A=abs(Y);          % 複數的幅值（模）

f=(0:1:N/2)*Fs/N;  % 生成頻率範圍

f=f';              % 轉置成列向量

%% 幅值修正

A_adj=zeros(N/2+1,1);

A_adj(1)=A(1)/N;        % 頻率為0的位置

A_adj(end)=A(end)/N;    % 頻率為Fs/2的位置

A_adj(2:end-1)=2*A(2:end-1)/N;

%% 繪製頻率幅值圖

figure

subplot(2,1,1)

plot(f,A_adj)

xlabel('頻率 (Hz)')

ylabel('幅值 (修正後)')

title('FFT變換幅值圖')

grid on

%% 繪製頻譜相點陣圖

subplot(2,1,2)

phase_angle=angle(Y);   % angle函式的返回結果為弧度

phase_angle=rad2deg(phase_angle);

plot(f,phase_angle)

xlabel('頻率 (Hz)')

ylabel('相位角 (degree)')

title('FFT變換相點陣圖')

grid on

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

clear;clc;close all

load('example')

Fs=2500;      % 採集頻率

T=1/Fs;       % 採集時間間隔

N=length(y);  % 採集訊號的長度

t=(0:1:N-1)*T;   % 定義整個採集時間點

t=t';            % 轉置成列向量

% 繪製時域訊號

figure

plot(t,y)

xlabel('時間')

ylabel('訊號值')

title('時域訊號')

% fft變換

Y=fft(y);          % Y為fft變換結果，複數向量

Y=Y(1:N/2+1);      % 只看變換結果的一半即可

A=abs(Y);          % 複數的幅值（模）

f=(0:1:N/2)*Fs/N;  % 生成頻率範圍

f=f';              % 轉置成列向量

% 幅值修正

A_adj=zeros(N/2+1,1);

A_adj(1)=A(1)/N;       % 頻率為0的位置

A_adj(end)=A(end)/N;   % 頻率為Fs/2的位置

A_adj(2:end-1)=2*A(2:end-1)/N;

% 繪製頻率幅值圖

figure

subplot(2,1,1)

plot(f,A_adj)

xlabel('頻率 (Hz)')

ylabel('幅值 (修正後)')

title('FFT變換幅值圖')

grid on

% 繪製頻譜相點陣圖

subplot(2,1,2)

phase_angle=angle(Y);    % angle函式的返回結果為弧度

phase_angle=rad2deg(phase_angle);

plot(f,phase_angle)

xlabel('頻率 (Hz)')

ylabel('相位角 (degree)')

title('FFT變換相點陣圖')

grid on