演習3解答(12 月 1 日まで) Q3.m -----------------------------------------------------------------------------% 孤立三角パルスのフーリエ変換 clear all close all df=0.02; fs=10; ts=1/fs; t=-5:ts:5; % % % % 周波数分解能 [Hz] 標本化周波数 [Hz] 標本化間隔 [sec] 時間 [sec] % 三角パルス x1(t) x1=zeros(size(t)); x1(41:51)=(0:1:10)/10; x1(51:60)=(10:-1:1)/10; % 三角パルス x2(t) x2=zeros(size(t)); x2(51:61)=(0:10)/10; x2(61:70)=(10:-1:1)/10; figure() plot(t,x1,'-',t,x2,'--') axis([-5 5 -1 2]) legend('x1(t)','x2(t)') xlabel('sec') % 三角パルス x1(t),x2(t)をプロット % x軸とy軸の範囲を指定 % 凡例を指定 [X1,x11,df1]=fourier_analysis(x1,ts,df); [X2,x21,df2]=fourier_analysis(x2,ts,df); X11=fftshift(X1)/fs; X21=fftshift(X2)/fs; f=[0:df1:df1*(length(x11)-1)]-fs/2; % % % % % 三角パルス x1(t)をMATLAB DFT(FFT) 三角パルス x2(t)をMATLAB DFT(FFT) MATLAB DFT(FFT)とフーリエ変換の関係 MATLAB DFT(FFT)とフーリエ変換の関係 スペクトルに対応した周波数軸を作成 figure() plot(f,abs(X11),'k-',f,abs(X21),'k--') % 振幅スペクトル(実線と破線が重なる) legend('x1(t)の振幅スペクトル','x2(t)の振幅スペクトル') xlabel('Hz') figure() % 位相スペクトル [位相シフトの確認] plot(f(250:264),angle(X11(250:264)),'k-',f(250:264),angle(X21(250:264)),'k--') legend('x1(t)の位相スペクトル','x2(t)の位相スペクトル') axis([-0.15 0.15 -4 4]) xlabel('Hz') fourier_analysis.m -----------------------------------------------------------------function [M,m,df]=fourier_analysis(m,ts,df) % m: 入力信号/出力信号 % ts: 標本時間間隔 % df: 必要周波数分解能/実際の周波数分解能 % M: 周波数スペクトル fs=1/ts; % 入力信号の標本周波数 n1=fs/df; % 必要周波数分解能を得るためのサンプル数 n2=length(m); % 入力信号の長さ n=2^max(nextpow2(n1),nextpow2(n2)); % 必要周波数分解能を得るためのfftサンプル数 M=fft(m,n); % 周波数スペクトル m=[m,zeros(1,n-n2)]; % 0をつけたデータ df=fs/n; % 実際の周波数分解能
© Copyright 2024 ExpyDoc