function [ztrata u y b P]=ps(regulator_parameters,system,apriori,seed)

horizont=regulator_parameters.horizont;
prvni_faze=1;
druha_faze=horizont-prvni_faze;

y=zeros(1,horizont);
b=zeros(1,prvni_faze+1);
P=zeros(1,prvni_faze+1);
u=zeros(1,horizont-1);
ztrata=0;

y(1)=apriori(1);
b(1)=apriori(2);
P(1)=apriori(3);

for i=1:prvni_faze
 %   u(i)=10*system.sigma/sqrt(P(i));
     u(i)=sqrt(abs(100-1/P(i)))*system.sigma;
    [y(i+1) b(i+1) P(i+1)]=get_next([y(i) b(i) P(i)],u(i),system,seed(i));
    ztrata=ztrata+(y(i+1)-system.yr(i+1))^2;
end
b_1=b(end);

for i=1:druha_faze-1
    u(i+prvni_faze)=(system.yr(i+prvni_faze)-y(i+prvni_faze))/b_1;
    y(i+1+prvni_faze)=y(i+prvni_faze)+system.b*u(i+prvni_faze)+system.sigma*seed(i+prvni_faze);
    ztrata=ztrata+(y(i+1+prvni_faze)-system.yr(i+1+prvni_faze))^2;
end

end