Transfer Function to State Space Using Prediction Error Method

Hi Everyone,

I have Experiemental Data (Bode Plot), using MATLAB I found the Transfer Function of my plant. I have attached my E7i_CSV file.

T2 = readtable('E7i_CSV.csv');   
                                                                       % Fill With Actual Sampling Frequency
FHz = T2.F;
Ts = 1/(2*(max(FHz)+10000))
for k = 1:size(T2,1)-1
    if FHz(k+1) == FHz(k)
        FHz(k+1) = FHz(k+1)+0.5;                                                % 'Brute Force' Interpolation 
    end
end
Mag = T2.G;
PhDeg = T2.P;
Response = Mag.*exp(1j*deg2rad(PhDeg));                                         % Complex Vector
sysfr = idfrd(Response, FHz, Ts, 'FrequencyUnit','Hz')
%bode(sysfr)
%tfsys = tfest(sysfr,4,3) This is much better
tfsys = tfest(sysfr,4,3) %82.69%
%tfsys = tfest(sysfr,6,3)
figure
compare(sysfr, tfsys)
is_stable = isstable(tfsys);
disp(is_stable);

I want to convert this transfer function to statespace equations, which will be used for Model Predictive Control Design. Simple tf2ss command give me TF but it doesn't look very accrurate.

[A, B, C, D] = tf2ss(num,den)

So, I was reading a paper which i will follow, they have mentioned that they had FRF (Frequency Response Function) data and they find state space model/equations using Prediction Error Method (PEM). But I am unable to understand how I can used Bode Plot data and convert it state space using PEM. I will appreciate your help.

Secondy I verified the paper StateSpace Model by first converting it to Transfer Function and then using tf2ss command. which shows inaccurate results similiar to my results. I am looking for your help. Following is the verification code for paper i am following

%Rana Paper State Space Equations Coversion to
                            %Transfer Function for Validation Reason
% Define state-space matrices
Ax = [-49.3277, -980.7648, 618.2198;
      567.3557, -120.2273, 651.7342;
      -45.8885, -415.4194, -115.7522];

Bx = [-1.8122; 3.901370; 4.3011];

Cx = [14.9368, 16.9208, -23.6827];

Dx = 0;

% Create a state-space model
sys = ss(Ax, Bx, Cx, Dx);

% Display the state-space model
disp('State-Space Model:');
disp(sys);
% Convert the state-space model to a transfer function
tf_sys = tf(sys);

% Display the transfer function
disp('Transfer Function:');
disp(tf_sys);
                        %TF to StateSpace Rana Model
num_rana = [ -62.92 3.625*10^4 -1.065*10^8];
den_rana = [ 1 285.3 8.811*10^5 1.982*10^8];
tff_rana = tf(num_rana,den_rana);
tffss_rana = tf2ss(num_rana,den_rana);
[Arana,Brana,Crana,Drana] = tf2

Answer 1

[Mathieu NOE]

hello

a tried some coding around your data

I don't know why you have "accuracy" problem when you go from TF to SS models

here I use a low order model , you can see TF and SS are matching well

I also opted for a different sampling rate (Ts is one fourth of the nyquist rate)

T2 = readtable('E7i_CSV.csv');   
                                                                       % Fill With Actual Sampling Frequency
FHz = T2.F;
Mag = T2.G;
PhDeg = T2.P;
Response = Mag.*exp(1j*pi/180*PhDeg);                                         % Complex Vector

Ts = 0.25*1/(2*(max(FHz)));
% Ts = 1/(2*(max(FHz)+10000));

Fs = 1/Ts;

% TF  design
NA = 6;
NB = NA-2; 
W = 2*pi*FHz;
[num,den] = invfreqs(Response,W,NB,NA,[],100);
hverif = freqs(num,den,W);
% convert to SS
[A,B,C,D]=tf2ss(num,den);
[Ad,Bd,Cd,Dd] = c2dm(A,B,C,D,Ts,'tustin');
SYS = ss(Ad,Bd,Cd,Dd,Ts);
[M,P] = bode(SYS,W);
 
M = squeeze(M);
P = squeeze(P);
P = wrapTo180(P);

figure(1)
subplot(2,1,1),semilogx(FHz,20*log10(Mag),'b',FHz,20*log10(abs(hverif)),'r',FHz,20*log10(M),'c');
legend('data','TF ','SS '); 
xlabel('Frequency (Hz)');
ylabel('modulus (dB)');
subplot(2,1,2),semilogx(FHz,PhDeg,'b',FHz,180/pi*(angle(hverif)),'r',FHz,(P),'c');
legend('data','TF ','SS '); 
xlabel('Frequency (Hz)');
ylabel('modulus (dB)');

[Mathieu NOE]

here a higher order model - if it's what you're after :

T2 = readtable('E7i_CSV.csv');   
                                                                       % Fill With Actual Sampling Frequency
FHz = T2.F;
Mag = T2.G;
PhDeg = T2.P;
Response = Mag.*exp(1j*pi/180*PhDeg);                                         % Complex Vector

Ts = 0.25*1/(2*(max(FHz)));
% Ts = 1/(2*(max(FHz)+10000));

Fs = 1/Ts;

% TF  design
NA = 22;
NB = NA-2; 
W = 2*pi*FHz;
[num,den] = invfreqs(Response,W,NB,NA,[],100);
hverif = freqs(num,den,W);
% convert to SS
[A,B,C,D]=tf2ss(num,den);
[Ad,Bd,Cd,Dd] = c2dm(A,B,C,D,Ts,'tustin');
SYS = ss(Ad,Bd,Cd,Dd,Ts);
[M,P] = bode(SYS,W);
 
M = squeeze(M);
P = squeeze(P);
P = wrapTo180(P);

figure(1)
subplot(2,1,1),semilogx(FHz,20*log10(Mag),'b',FHz,20*log10(abs(hverif)),'r',FHz,20*log10(M),'c');
legend('data','TF ','SS '); 
xlabel('Frequency (Hz)');
ylabel('modulus (dB)');
subplot(2,1,2),semilogx(FHz,PhDeg,'b',FHz,180/pi*(angle(hverif)),'r',FHz,(P),'c');
legend('data','TF ','SS '); 
xlabel('Frequency (Hz)');
ylabel('modulus (dB)');

[Muhammad]

Dear @Mathieu NOE, Thanks for your answer. Sorry for reply late, I wrote reply to your answer but I check out today it was not submitted maybe or deleted.

I have only one question.

[Ad,Bd,Cd,Dd] = c2dm(A,B,C,D,Ts,'tustin');

In the above line of code from contibeous time to discrete you used 'tustin', I implemented it on another experimental data, actually my data has been changed because of some adittional components.

So, when I used c2dm command without tustin i.e. only [Ad,Bd,Cd,Dd] = c2dm(A,B,C,D,Ts] my ss model follows 100% of transfer function but when I used c2dm(A,B,C,D,Ts,'tustin') the ss follows only 80% of that model.

So, I wanted to know is it okay to use this command only without tustin? [Ad,Bd,Cd,Dd] = c2dm(A,B,C,D,Ts]

Why Tustin is used?

[Mathieu NOE]

hello

have a look here to see what method is best suited for your case :

Continuous-Discrete Conversion Methods - MATLAB & Simulink - MathWorks France

hope it helps !

[Mathieu NOE]

hello @Muhammad

problem solved ?

[Muhammad]

Hi @Mathieu NOE

I read the document and try to understand all the six discretization methods, for frequency domain they used Tustin Approximation but I was not able to understand it fully, and still I have confusion that why c2dm(A,B,C,D,Ts,'tustin'); doesn't follow my TF completely. For now I am using c2dm without tustin approximation.

Thanks for follow up. Please help me in this regard for choosing right one