frep2tf

パラメータ

frq

周波数ベクトル (単位: Hz).

repf

周波数応答ベクトル

dg

線形システムの次数

dom

時間領域 ('c' または 'd' または dt)

tols

長さ3のベクトルで, 相対および絶対許容誤差, 最大反復回数を指定します(デフォルト値: rtol=1.e-2; atol=1.e-4, N=10).

weight

周波数の重みベクトル

h

SISO 伝達関数

err

誤差 (例えば dom='c'の場合, sum(abs(h(2i*pi*frq) - rep)^2)/size(frq,*))

説明

周波数応答を伝達関数に変換します. hの次数は,必須の引数 dgで事前に指定されます. 以下の線形システムが最小二乗解として解かれます.

weight(k)*(n( phi_k) - d(phi_k)*rep_k)=0, k=1,..,n

ただし,dom='c'の場合は phi_k= 2*%i*%pi*frq, そうでない場合は phi_k=exp(2*%i*%pi*dom*frq. weight ベクトルが指定されない場合, デフォルトの重み付けが使用されます(dom='c'の場合).

安定で最小位相のシステムがfactors関数により 得られます.

例

s=poly(0,'s');
h=syslin('c',(s-1)/(s^3+5*s+20))
frq=0:0.05:3;repf=repfreq(h,frq);
clean(frep2tf(frq,repf,3))

Sys=ssrand(1,1,10); 
frq=logspace(-3,2,200);
[frq,rep]=repfreq(Sys,frq);  //Sysの周波数応答
[Sys2,err]=frep2tf(frq,rep,10);Sys2=clean(Sys2)//Sysの周波数応答からSys2を得る
[frq,rep2]=repfreq(Sys2,frq); //Sys2の周波数応答
clf();bode(frq,[rep;rep2])   //Sys と Sys2の応答
[gsort(spec(Sys('A'))),gsort(roots(Sys2('den')))] //極

dom=1/1000; // サンプリング時間
z=poly(0,'z');
h=syslin(dom,(z^2+0.5)/(z^3+0.1*z^2-0.5*z+0.08))
frq=(0:0.01:0.5)/dom;repf=repfreq(h,frq);
[Sys2,err]=frep2tf(frq,repf,3,dom);
[frq,rep2]=repfreq(Sys2,frq); //Sys2の周波数応答
clf();plot2d1("onn",frq',abs([repf;rep2])');

参照

• imrep2ss — インパルス応答の状態空間実現
• arl2 — L2伝達関数近似による SISO モデル実現
• time_id — SISO 最小二乗同定
• armax — armax 同定
• frfit — 周波数応答フィッティング