function eta=thinit(eta0,R,par,sp);
%THINIT gives (random) initial values to the parameters, with guaranteed 
%	predictor and/or system stability
%
%	TH = thinit(TH0)
%or
%	TH = thinit(TH0,R,PAR,SP)
%
%	TH0:Original theta corresponding to state-space model (See help th(ss)) 
%	TH: Modified TH0 with new initial parameter values
%	R: variance of random initial parameters 
%	   (row vector with dim = dimension of parameter vector)
%	PAR: mean of initial parameters (PAR =[] gives default)
%	AP: SP='p' :stability of predictor assured
%	    SP='s' :stability of system is assured
%	    SP='b' :both predictor and system will be stable
%
%	Defaults: SP='p', PAR=the parameter values of TH0, R=[1 ... 1]

%	L. Ljung 10-2-1990
%	Copyright (c) 1990 by the MathWorks, Inc.
%	All Rights Reserved

if ~isthss(eta0),error('This routine is for state-space model structures only!'),end
[nr,nc]=size(eta0);
sspm=getmfth(eta0);
T=eta0(1,2); Tmod=abs(T);
nd=eta0(1,5);
ms=getargth(eta0);
if nargin<2;R=[];end
if nargin<3;par=[];end
if nargin<4;sp=[];end

if isempty(R),R=ones(nd,1);end,if isempty(sp),sp='p';end
if isempty(par),par=eta0(3,1:nd);end

if length(R)~=nd,error('The length of the vector R must equal the number of free parameters in th!'),end
if length(par)~=nd,error('The length of the parameter vector must equal the number of free parameters in th!'),end
test=2;nr=1;if norm(R)==0,noR=1;else noR=0;end
while test>1 & nr<50
	parval=rand(1,nd)*diag(sqrt(R))+par;
	if any(eta0(2,8)==[2 3]),Tmod=-1;end
	[a,b,c,d,k]=feval(sspm,parval,Tmod,ms);

	if sp=='s',test=max(abs(eig(a)));
	elseif sp=='p',test=max(abs(eig(a-k*c)));
	else test=max([max(abs(eig(a))), max(abs(eig(a-k*c)))]);
	end
nr=nr+1; if nr==50,error('50 samples drawn without satisfying the stability conditions'),end
if noR,nr=50;if test>1,disp('WARNING: System/predictor unstable'),end,end
end
eta=eta0;
eta(3,1:nd)=parval;
eta(2,7)=27;
if eta(2,8)==3,eta(2,7)=37;end
