close all;
clc;
clear all;


T=[15 10 5];
tau=[40 20 10];
load('pomiary_3out.mat');
time=1:1:300;
figure(1)
plot(time,pomiary_3out);

y1=pomiary_3out(:,1)- pomiary_3out(1,1);
y2=pomiary_3out(:,2)- pomiary_3out(1,2);
y3=pomiary_3out(:,3)- pomiary_3out(1,3);

figure(2)
plot(time,y1);
figure(3)
plot(time,y2);
figure(4)
plot(time,y3);

k1=(y1(300,1)-y1(1,1))/1.0
k2=(y2(300,1)-y2(1,1))/1.0
k3=(y3(300,1)-y3(1,1))/1.0

[ld1,md1]=pade(tau(1),10);
[l1,m1]=series([k1],[T(1) 1],ld1,md1);
ym1=step(l1,m1,time);
MSE1=sum((y1-ym1).^2);

[ld2,md2]=pade(tau(2),10);
[l2,m2]=series([k2],[T(2) 1],ld2,md2);
ym2=step(l2,m2,time);
MSE2=sum((y2-ym2).^2);

[ld3,md3]=pade(tau(1),10);
[l3,m3]=series([k3],[T(1) 1],ld3,md3);
ym3=step(l3,m3,time);
MSE3=sum((y3-ym3).^2);

figure(5)
plot(ym1);
figure(6)
plot(ym2);
figure(7)
plot(ym3);