I get very good match:
t=0:0.01:100;
u=6*cos(2*t)';
y=(-exp(-t)-3*exp(-3*t)+4*cos(2*t)+4*sin(2*t))';
data=iddata(y,u,t(2));
sys=(9*s+14)/(3*(s+1)*(s+3));
sys_id=tfest(data,2,1);
Checking:
>>sys
sys =
9 s + 14
----------------
3 s^2 + 12 s + 9
>> sys_id
sys_id =
From input "u1" to output "y1":
3.023 s + 4.606
---------------
s^2 + 4 s + 3
Continuous-time identified transfer function.
As you can see, sys_id gets the poles exactly, the zero and gain are slightly off.
If we do the bode and step, the differences are really small:
bode(sys);
hold on;
bode(sys_id,'r');
figure;
step(sys,sys_id,'r');
