How can I find the Transfer Function having Magnitude(dB), Phase(degree) and Frequency(Hz)?

36 visualizaciones (últimos 30 días)
Liang Kar Yan
Liang Kar Yan el 10 de Dic. de 2021
Comentada: Liang Kar Yan el 14 de Dic. de 2021
I have 3 individual files which are Magnitude(dB), Phase(degrees) and Frequency(Hz) in excel. I need to find the transfer function with these data. After getting the transfer function, I need to plot back the graph (magnitude and phase) from transfer function to compare with my data.
I wish to get something like the picture attached below.
  5 comentarios
Mostrar 3 comentarios más antiguos
Mathieu NOE
Mathieu NOE el 14 de Dic. de 2021
hello @Star Strider
I believe Magn is given in dB (like in the plot)
so first think is to convert back to linear magnitude
Magn = 10.^(Magn/20) , and then
cplxv = Magn .* exp(1j*deg2rad(Phse))

Iniciar sesión para comentar.

Iniciar sesión para responder a esta pregunta.

Respuesta aceptada

Mathieu NOE
Mathieu NOE el 14 de Dic. de 2021
hello
I tried a few options , IIR or FIR filters fit.
As the phase plot shows , there is a quite significant phase roll rate, idicating the presence of a huge delay in the system.
I assumed a sampling rate of Fs = 1000 hz and found out that more or less we can fit either a FIR or a IIr high pass filter in series with almost 200 samples of delay (hudge !!)
maybe there are more powerful tools then the simple invfreqz (not giving any good results here) or the manual fit I am doing here
FIR fit plot :
IIR fit plot :
clc
clearvars
Freq = readmatrix('Frequency.xlsx');
Magn_dB = readmatrix('Magnitude.xlsx');
Phse = readmatrix('Phase.xlsx');
figure(1)
subplot(211),plot(Freq,Magn_dB);
subplot(212),plot(Freq,Phse);
Magn = 10.^(Magn_dB/20) ;
%% high pass filter model (IIR)
Fs = 1000; % ? to be confirmed
N = 8; % filter order
dc_gain = Magn(end); % asymptotic value
[val,ind] = min(abs(Magn_dB - Magn_dB(end) + 5)); % - 5dB (vs dc_gain) cut off frequency index search
fc = Freq(ind); % - 5dB (vs dc_gain) cut off frequency
[b,a] = butter(N,2*fc/Fs,'high');
b = b*dc_gain; % apply dc gain on numerator
[g,p] = dbode(b,a,1/Fs,2*pi*Freq);
% adding delay due to sampling
nd = 200; % delay (samples)
rpd = -360*nd*Freq/Fs;
p = p+rpd; % adding filter phase to samples delay phase
p = mod(p,360);
p = p -180; % polarity correction
figure(1)
subplot(211),plot(Freq,Magn_dB,Freq,20*log10(g));
title('IIR model fit')
ylabel('Modulus (dB)');
subplot(212),plot(Freq,Phse,Freq,p);
xlabel('Frequency (Hz)');
ylabel('Phase (°)');% return
%% high pass filter model (FIR)
Fs = 1000; % ? to be confirmed
N = 20;
dc_gain = Magn(end); % asymptotic value
[val,ind] = min(abs(Magn_dB - Magn_dB(end) + 3)); % - 3dB (vs dc_gain) cut off frequency index search
fc = Freq(ind); % - 3dB (vs dc_gain) cut off frequency
[b,a] = fir1(N,2*fc/Fs,'high');
b = b*dc_gain; % apply dc gain on numerator
[g,p] = dbode(b,a,1/Fs,2*pi*Freq);
% adding delay due to sampling
Fs = 1000; % ? to be confirmed
nd = 200; % delay (samples)
rpd = -360*nd*Freq/Fs;
p = p+rpd; % adding filter phase to samples delay phase
p = mod(p,360);
p = p -180; % polarity correction
figure(2)
subplot(211),plot(Freq,Magn_dB,Freq,20*log10(g));
title('FIR model fit')
ylabel('Modulus (dB)');
subplot(212),plot(Freq,Phse,Freq,p);
xlabel('Frequency (Hz)');
ylabel('Phase (°)');
%% Id with invfreqz (FIR)
h = Magn .* exp(1j*180/pi*(Phse));
nb = 40+nd;
na = 1;
iter = 1000;
[bb,aa] = invfreqz(h,pi*Freq/Fs,nb,na,[],iter); % stable approximation to system
[g,p] = dbode(bb,aa,1/Fs,2*pi*Freq);
p = mod(p,360);
figure(3)
subplot(211),plot(Freq,Magn_dB,Freq,20*log10(g));
subplot(212),plot(Freq,Phse,Freq,p);
  10 comentarios
Mostrar 8 comentarios más antiguos

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Digital Filter Analysis en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by