Frequency response of a system operating at a fixed frequency
17 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Alberto Moretti
am 10 Jun. 2024
Kommentiert: Mathieu NOE
am 12 Jun. 2024
I am trying to simulate the behavior of a system in simulink. I am providing a sinusoidal input with a fixed frequency (e.g. 0.3 Hz) and I am exporting the results in matlab. I want to verify that at this operating frequency, the frequency response is flat. Once done that I'll increase the frequency and so on untill I reach a decrease in the frequency response of -3dB.
However I am facing dissiculties in evaluating this.
I tried doing something like this:
Y_measured = fft(measured);
Y_target = fft(target);
Y = Y_measured./Y_target;
amplitude = abs(Y);
phase = unwrap(angle(Y(1:N/2+1)));
I then plotted the amplitude and the phase over time.
The fact is that being the system evaluated at a fixed frequency, the frequency response diagram loses a lot of significance and I become only intrested in what would be a specific point of the diagram (I guess?).
So in conclusion:
Is this a correct way of evaluating the frequency response of my system?
What would be an efficien way of plotting and showing the results?
Thank you.
0 Kommentare
Akzeptierte Antwort
Mathieu NOE
am 10 Jun. 2024
Bearbeitet: Mathieu NOE
am 11 Jun. 2024
hello
we could use the general fft approach , but you could also use the stepped sine identification method. We only need to generate a cos/sin waveforms at the desired frequency, generate a test signal that must be fed to the process being identified then project (à la DFT) the output on the cos/sin signals and from there you get the real and imaginary parts of your process
here a little demo - the process is here a simple butterworth filter
EDIT : updated the code so the ID process is operated on a fixed an round number of periods (as it should be);
the fact that the previous version would give correct results is because the samples amount was relatively high so that the last non complete period would have minimal impact on the result.
But this new version is the one and only correct implementation
ThemeCopy
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% stepped sinus identification
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%% inputs %%%%%%%%%%%%%%%
Fs = 500;
freq = linspace(0,Fs/2.56,300);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%tranfer function (process model)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% process is simulated here with a bandpass butterworth filter
fl = 10; % lower corner freq
fh = 100; % upper corner freq
N = 2; % order
[nz,dz] = butter(N,[fl fh]/(Fs/2)); % "process"
[g,p] = dbode(nz,dz,1/Fs,2*pi*freq);
%% Stepped sine main loop
f_id = (10:5:150); % frequencies at which we want to identify our process
period_ramp_up = 10;
period_ID = 20;
amplitude_ID = 0.1; %amplitude of the test signal
dt = 1/Fs;
FRF_mod_max = 0;
for k = 1:numel(f_id)
fc = f_id(k); % select frequency
% compute samples needed
samples_ramp_up = round(period_ramp_up*Fs/fc); % ramp up portion of the test signal - no identification is done during this period to avoid transients
samples_ID = round(period_ID*Fs/fc); % samples used for the identification process itself
samples = samples_ramp_up + samples_ID; % total amount of samples
% cos and sin DFT tables (coefficients)
t = (0:samples-1)*dt;
tcos = cos(2*pi*fc*t); %
tsin = -sin(2*pi*fc*t); %
w = [linspace(0,1,samples_ramp_up) ones(1,samples_ID)]; % ramp up "window"
signal_ID = amplitude_ID*tcos.*w; % signal with ramp up
% uncomment the line below if you want to have a look at the ID signal
% plot(t(1:samples_ramp_up),signal_ID(1:samples_ramp_up),'r')
% hold on
% plot(t(1+samples_ramp_up:samples),signal_ID(1+samples_ramp_up:samples),'b')
% hold off
% generate / measure the process output
signal_out = filter(nz,dz,signal_ID);
% TF real and imag DFT computation (once ramp up is finished)
ind = (samples_ramp_up+1:samples);
FRF_real(k) = sum(signal_out(ind).*tcos(ind))*2/samples_ID/amplitude_ID;
FRF_imag(k) = sum(signal_out(ind).*tsin(ind))*2/samples_ID/amplitude_ID;
end
% convert identified TF real and imaginary part into module / phase
FRF_mod = sqrt(FRF_real.^2 + FRF_imag.^2);
FRF_phase = atan2(FRF_imag,FRF_real);
figure(1);
subplot(2,1,1),semilogx(freq,20*log10(g),'b',f_id(1:k),20*log10(FRF_mod),'dr');grid on
title(' Stepped sine Identification');
ylabel('dB ');
legend('model','ID');
subplot(2,1,2),semilogx(freq,wrapTo180(p),'b',f_id(1:k),wrapTo180(FRF_phase*180/pi),'dr');grid on
xlabel('Frequency (Hz)');
ylabel(' phase (°)')
legend('model','ID');
14 Kommentare
Mathieu NOE
am 12 Jun. 2024
maybe then you could downsample the data to a much more appropriate rate (use "rate converter" or a "ZOH" block to transition to a lower sampling rate) - otherwise you gonna handle very large time data and make a lot of unnecessary computations ...
Weitere Antworten (1)
Paul
am 10 Jun. 2024
Hi Alberto,
Conceptually, what you're doing is correct (though there are other alternatives to the input signal). Keep in mind that the only element of Y that you care about is the element that corresponds as closely as possible to the target frequency (0.3 Hz in this example). Or you could interpolate Y to that frequency, or you could zero-pad the FFTs to ensure you get a frequency sample exactly at the target frequency. Another thing to consider is the effect of the startup transient; you may want to only do the FFT processing on the data after the system reaches steady state. And you should make sure that you have enough steady state cycles of data. If the system is nonlinear, you may want to consider checking sensitivity of the results to the amplitude of the input. And there may be other complications is you're system is unstable.
Also may be interested in frestimate, for which there is is also a GUI interface (I think in the Model Linearizer).
Whatever path you take, I suggest first working with a model for which you know what the frequency response should be to help iron out the workflow.
Siehe auch
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!