# Estimation of transfer function

조회 수: 2 (최근 30일)
Kamil 2011년 6월 5일
Hello,
I would like to estimate a transfer function based on the experimental data (magnitude, phase and frequency. I was trying to use pem function but it does not work well. For example in my system I have a derivative and my transfer function should be something like: (as^2+bs)/(a1s^3+b1s^2+c1s+d1). I cannot get it in the estimation model. It gives me something like as^2+bs+c. I do not want "c" parameter. There is any estimation function that allows me to settle my searching function not only number of orders? I mean user0searching parameter? If is there any how I can do it?
Thank you for any response.
Kamil

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

### 답변 (5개)

Arnaud Miege 2011년 6월 6일
For frequency-domain data, there are three types of id models supported by the System Identification Toolbox(see documentation):
For your application, I would recommend the latter. Choose the correct order for your system and once you have your estimated state-space matrices, you can use the tf command from the Control System Toolbox to convert it to a transfer function.
HTH,
Arnaud
##### 댓글 수: 0이전 댓글 -2개 표시이전 댓글 -2개 숨기기

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

Rajiv Singh 2011년 6월 6일
There is also the "grey box" approach that lets you parametrize your model any way you like. For this, write a MATLAB function that takes in a, b, a1, b1, c1 and d1 as input arguments and returns a state-space form of the model as output arguments. Use this file to create an "idgrey" model object whose parameters (your a, b, a1,..) can be estimated to fit data. You will need to derive the expressions for SS matrices to begin with. A controls textbook (e.g., Linear Analysis by Kailath) would show you some common forms.
##### 댓글 수: 0이전 댓글 -2개 표시이전 댓글 -2개 숨기기

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

Kamil 2011년 6월 26일
Hello,
Thank you both of you for answers. I had much work and I could not answer at all. I will try to do it now. The links are to Matlab 2011 I have Matlab 2007b. I hope it will work as well. I will give you answer if I get what I want.
Kamil
##### 댓글 수: 0이전 댓글 -2개 표시이전 댓글 -2개 숨기기

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

Kamil 2011년 6월 27일
Hello again,
I did test of the estimation but I cannot get the right answer. Well, I can but it strictly depends on the initial points. As I know the function I know what to put. Could you check if I am going in the right direction?
G=tf([1 1],[1 1 1]);
[num,den] = tfdata(G,'v');
[a,b,c,d]=tf2ss(num,den);
[m,ph,f]=bode(G);
zfr = m(:).*exp(i*ph(:)*pi/180);
Ts = 0;
data = idfrd(zfr,f,Ts);
A = [1 1; 1 1];
B = [1; 1];
C = [1 1];
D = [0];
m = idss(A,B,C,D,'Ts',0);
m.As = [NaN NaN; NaN NaN];
m.Bs = [NaN; NaN];
m.Cs = [NaN NaN];
m.Ds = [NaN];
res = pem(data,m);
A=res.a;
B=res.b;
C=res.c;
D=res.d;
[num1, den1] = ss2tf(A, B, C, D);
tf(num1,den1)
Thank you.
Kamil
##### 댓글 수: 0이전 댓글 -2개 표시이전 댓글 -2개 숨기기

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

Kamil 2011년 6월 28일
Hello again,
I changed a file to my own measured data but unfortunately Matlab has problem to find something similar. Anybody has any idea how I can fix it?
function estimation_data
format short e
f=[6.2832e-001 1.2566e+000 1.8850e+000 2.5133e+000 3.1416e+000 3.7699e+000 4.3982e+000 5.0265e+000 5.6549e+000 6.2832e+000 6.9115e+000...
7.5398e+000 8.1681e+000 8.7965e+000 9.4248e+000 1.0053e+001 1.0681e+001 1.1310e+001 1.1938e+001 1.2566e+001 1.3195e+001 1.3823e+001...
1.4451e+001 1.5080e+001 1.5708e+001 1.6336e+001 1.6965e+001 1.7593e+001 1.8221e+001 1.8850e+001 1.9478e+001 2.0106e+001 2.0735e+001...
2.1363e+001 2.1991e+001 2.2619e+001 2.3248e+001 2.3876e+001 2.4504e+001 2.5133e+001 2.5761e+001 2.6389e+001 2.7018e+001 2.7646e+001...
2.8274e+001 2.8903e+001 2.9531e+001 3.0159e+001 3.0788e+001 3.1416e+001 3.2044e+001 3.2673e+001 3.3301e+001 3.3929e+001 3.4558e+001...
3.5186e+001 3.5814e+001 3.6442e+001 3.7071e+001 3.7699e+001 3.8327e+001 3.8956e+001 3.9584e+001 4.0212e+001 4.0841e+001 4.1469e+001...
4.2097e+001 4.2726e+001 4.3354e+001 4.3982e+001 4.4611e+001]; % [rad/s]
m2=[ -2.4555e+001 -1.8863e+001 -1.5273e+001 -1.3022e+001 -1.1368e+001 -1.0106e+001 -9.0412e+000 -8.1258e+000 -7.5362e+000 -6.8305e+000 -6.3075e+000...
-5.8365e+000 -5.3956e+000 -4.9887e+000 -4.6133e+000 -4.2957e+000 -3.9825e+000 -3.7159e+000 -3.5035e+000 -3.2543e+000 -3.0982e+000 -2.8406e+000...
-2.5244e+000 -2.1561e+000 -1.9043e+000 -1.5786e+000 -1.3299e+000 -1.0794e+000 -1.0280e+000 -1.0049e+000 -1.0973e+000 -1.1484e+000 -1.3472e+000...
-1.6581e+000 -1.7930e+000 -2.1920e+000 -2.0484e+000 -2.6549e+000 -2.3225e+000 -2.4648e+000 -2.8349e+000 -2.5328e+000 -2.5987e+000 -2.5688e+000...
-3.6216e+000 -3.6520e+000 -3.9749e+000 -3.6494e+000 -3.8424e+000 -3.3791e+000 -3.0330e+000 -3.3493e+000 -3.0616e+000 -3.6999e+000 -3.2204e+000...
-3.6160e+000 -3.4716e+000 -3.2962e+000 -3.5312e+000 -3.2992e+000 -3.5701e+000 -3.3806e+000 -3.7542e+000 -3.7192e+000 -3.5343e+000 -3.8899e+000...
-3.9252e+000 -4.1289e+000 -4.7280e+000 -4.2376e+000 -4.8412e+000]; % [dB]
ph=[7.9373e+001 7.3436e+001 6.8009e+001 6.2013e+001 5.6244e+001 5.0607e+001 4.6156e+001 4.1887e+001 3.7741e+001 3.3938e+001 3.0157e+001...
2.6500e+001 2.2878e+001 1.9852e+001 1.6389e+001 1.2625e+001 9.3336e+000 6.4858e+000 3.4874e+000 8.1324e-001 -2.2675e+000 -5.3157e+000...
-8.9050e+000 -1.3065e+001 -1.5734e+001 -1.9750e+001 -2.2638e+001 -2.6510e+001 -2.9596e+001 -3.2418e+001 -3.5896e+001 -4.0109e+001 -4.3411e+001...
-4.7108e+001 -5.0974e+001 -5.4830e+001 -5.9614e+001 -6.3986e+001 -6.6701e+001 -7.0808e+001 -7.3210e+001 -8.2986e+001 -8.3098e+001 -8.6777e+001...
-9.1118e+001 -9.3465e+001 -9.4334e+001 -9.8264e+001 -9.8835e+001 -9.9147e+001 -1.0102e+002 -1.0385e+002 -1.0682e+002 -1.0833e+002 -1.0766e+002...
-1.1003e+002 -1.0825e+002 -1.1121e+002 -1.0969e+002 -1.1306e+002 -1.1880e+002 -1.1594e+002 -1.2291e+002 -1.2404e+002 -1.2812e+002 -1.2867e+002...
-1.3082e+002 -1.3201e+002 -1.3627e+002 -1.3741e+002 -1.4003e+002];
zfr = m2.*exp(i*ph*pi/180);
Ts = 0;
data = idfrd(zfr,f,Ts);
A = [1 1 1 1 1 -1; 1 0 0 0 0 0; 0 1 0 0 0 0; 0 0 1 0 0 0; 0 0 0 1 0 0; 0 0 0 0 1 0];
B = [1; 0; 0; 0; 0; 0];
C = [1 1 1 1 1 0];
D = [0];
m = idss(A,B,C,D,'Ts',0);
m.As = [NaN NaN NaN NaN NaN -1; 1 0 0 0 0 0; 0 1 0 0 0 0; 0 0 1 0 0 0; 0 0 0 1 0 0; 0 0 0 0 1 0]
m.Bs = [1; 0; 0; 0; 0; 0];
m.Cs = [NaN NaN NaN NaN NaN 0];
m.Ds = [0];
res = pem(data,m);
A=res.a;
B=res.b;
C=res.c;
D=res.d;
[num1, den1] = ss2tf(A, B, C, D);
[m1,ph1,f1] = bode(num1,den1);
tf(num1, den1)
figure(1)
plot(f/(2*pi),m2,'o','LineWidth',2)
hold on
grid on
plot(f1/(2*pi),20*log10(m1(:)),'r','LineWidth',2)
axis([0 7.5 -26 50])
xlabel('Frequency [2*pi*f]','VerticalAlignment','top','fontsize',17)
ylabel('Magnitiude [dB]','VerticalAlignment','bottom','fontsize',17)
legend('Measured data','New',1)
figure(2)
plot(f/(2*pi),ph(:),'o','LineWidth',2)
hold on
plot(f1/(2*pi),ph1(:),'r','LineWidth',2)
grid on
axis([0 7.5 -150 100])
xlabel('Frequency [2*pi*f]','VerticalAlignment','top','fontsize',17)
ylabel('Phase [deg]','VerticalAlignment','bottom','fontsize',17)
legend('Measured data','New',1)
##### 댓글 수: 1이전 댓글 -1개 표시이전 댓글 -1개 숨기기
Kamil 2012년 3월 3일
Any solutions???? Anybody can help me?

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

### 카테고리

Help CenterFile Exchange에서 Linear Model Identification에 대해 자세히 알아보기

### Community Treasure Hunt

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

Start Hunting!

Translated by