fft of impulse response does not equal bode of transfer function

조회 수: 17 (최근 30일)
Lukas
Lukas 2018년 9월 11일
답변: Lukas 2018년 11월 2일
Hello
I am struggling with the problem, that the fourier transform of the impulse response of my system does not equal the bode plot of my transfer function. It would be great if somebody could help me to find the misstake.
clear all, clc, close all
%.........
% System
%.........
s=tf('s');
stoerung_G = 75 / (s + 25);
figure; bode(stoerung_G); grid;
stoerund_dt = 0.00001;
stoerung_t = [0 : stoerund_dt : 0.4]';
[stoerung_impulse] = impulse(stoerung_G, stoerung_t);
%.........
% Fourieranalyse
%.........
stoerung_y = stoerung_impulse;
stoerung_Fs = 1/stoerund_dt; % Hz
stoerung_L = size(stoerung_y,1);
stoerung_Y_ft = fft(stoerung_y);
P2 = abs(stoerung_Y_ft/stoerung_L);
P1 = P2(1:stoerung_L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
P1_log = 20 * log10(P1);
Phase2 = angle(stoerung_Y_ft);
Phase1 = Phase2(1:stoerung_L/2+1);
Phase1(2:end-1) = 2*Phase1(2:end-1);
Phase1 = Phase1*180/pi();
stoerung_f = stoerung_Fs*(0:(stoerung_L/2))/stoerung_L;
stoerung_w = 2*pi()*stoerung_f;
figure
subplot(211)
semilogx(stoerung_w,P1_log), grid, xlim([1 3000]);
title('Single-Sided Amplitude Spectrum of X(t)')
xlabel('w (rad/s)')
ylabel('|Mag| (dB)')
subplot(212)
semilogx(stoerung_w,Phase1), grid, xlim([1 3000]);
xlabel('w (rad/s)')
ylabel('Phase (degree)')

채택된 답변

Lukas
Lukas 2018년 11월 2일
Hello
Following link helped me to solve my issue: https://www.12000.org/my_notes/on_scaling_factor_for_ftt_in_matlab/
I have divided the fft by the sampling time. Then I have also deleted the factor 2 at the magnitude and at the phase calculation. It is still not reasonable to me why, because I am plotting just one half (as I understood: frequencies from 0 to fs/2 and not -fs/2 to fs/2) of the fft. In order to keep the "signal energy" constant I should multiply by 2. Anyways, it works now without the "2".
clear all, clc, close all
%.........
% System
%.........
s=tf('s');
stoerung_G = 75 / (s + 25);
figure; bode(stoerung_G); grid;
stoerung_dt = 0.00001;
stoerung_t = [0 : stoerung_dt : 0.4]';
[stoerung_impulse] = impulse(stoerung_G, stoerung_t);
%.........
% Fourieranalyse
%.........
stoerung_y = stoerung_impulse;
stoerung_Fs = 1/stoerung_dt; % Hz
stoerung_L = size(stoerung_y,1);
stoerung_Y_ft = fft(stoerung_y)*stoerung_dt;
P2 = abs(stoerung_Y_ft);
P1 = P2(1:stoerung_L/2+1);
P1_log = 20 * log10(P1);
Phase2 = angle(stoerung_Y_ft);
Phase1 = Phase2(1:stoerung_L/2+1);
Phase1 = Phase1*180/pi();
stoerung_f = stoerung_Fs*(0:(stoerung_L/2))/stoerung_L;
stoerung_w = 2*pi()*stoerung_f;
figure
subplot(211)
semilogx(stoerung_w,P1_log), grid, xlim([1 3000]);
title('Single-Sided Amplitude Spectrum of X(t)')
xlabel('w (rad/s)')
ylabel('|Mag| (dB)')
subplot(212)
semilogx(stoerung_w,Phase1), grid, xlim([1 3000]);
xlabel('w (rad/s)')
ylabel('Phase (degree)')

추가 답변 (1개)

Dimitris Kalogiros
Dimitris Kalogiros 2018년 9월 11일
You have to change the time interval that you used to define stoerung_t
Try a larger one:
stoerung_t = [0 : stoerund_dt : 40-stoerund_dt]';
Now, two bode diagrams will become identical
  댓글 수: 2
Lukas
Lukas 2018년 9월 18일
Unfortunately this is not solving my problem. It is also strange, that the fft leads to an offset in the magnitude plot of -26dB and that the phase shift is -180° instead of -90° for a PT1 system.
Dimitris Kalogiros
Dimitris Kalogiros 2018년 9월 18일
These artifacts stems from the fact that you compare a continous time system, discribed by laplace transform (stoerung_G), and a discrete time system, described by z-transform (stoerung_Y_ft).
For more details, have a look here: bilinear transformation

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Spectral Measurements에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by