function [u,l]=cont(theta,q,y,u)
% explicit LQG controller with input penalty q>=0 for first order system
% yt=theta*[y(t-1),u(t),u(t-1)] with state s:=theta1*y+theta3*u
% see report Karny no 1012 p. 36 (where is an error in state definition)
 
  d1=theta(1,1)*theta(2,1)+theta(3,1);
  m0=q+theta(2,1)^2;
  m1=0.5*(m0+q*theta(1,1)^2+theta(3,1)^2);
  l1=theta(2,1)*theta(3,1)-theta(1,1)*q;
  if(abs(d1)<1e-6)
     l=theta(2,1)/m0;
  else
     l=(theta(1,1)+l1/(m1+sqrt(m1^2-l1^2)))/d1;
  end
  u=-l*(theta(1,1)*y+theta(3,1)*u);