function k = dcgain(a,b,c,d)
%DCGAIN D.C. gain of continuous system.
%	K = DCGAIN(A,B,C,D) computes the steady state (D.C. or low 
%	frequency) gain of the continuous state-space system (A,B,C,D).
%
%	K = DCGAIN(NUM,DEN) computes the steady state gain of the 
%	continuous polynomial transfer function system G(s)=NUM(s)/DEN(s)
%	where NUM and DEN contain the polynomial coefficients in 
%	descending powers of s.
%
%	See also: DDCGAIN.

%	Clay M. Thompson  7-6-90
%	Copyright (c) 1986-93 by the Mathworks, Inc.

error(nargchk(2,4,nargin));

if (nargin==2), % Transfer function description
  [num,den] = tfchk(a,b);
  n = length(den);
  if (n == 0 | length(num) == 0)
      k = [];
  else
      k = num(:,n)/den(n);
  end

elseif nargin==4, % State space description
  error(abcdchk(a,b,c,d));
  if any(abs(eig(a))<eps), % System is singular.
    % Use LU decomposition to solve system.  Put a very small number in place
    % of the singularity.  Set to INF if necessary after calculation.
    [l,u] = lu(a);
    ndx = find(abs(diag(u))<eps);
    if length(ndx)>0, 
       % Make a very small number.
       u(ndx,ndx) = -1.e-300*ones(length(ndx)); 
    end 
    k = -c*(u\(l\b)) + d;
	% Note there is caveat here. If the inverse has all infinites in it
	% and the C matrix has some zero elements, Nans may propogate when the 
	% true solution should have inifinites in it.
    out = find(abs(k)>1.e+200);
    if length(out)>0, k(out) = inf*sign(k(out)); end
  else
    k = -c/a*b + d;
  end
else
  error('Wrong number of input arguments.');
end