% SCHEBYSHEV Symbolic counterpart of chebyshev.m
%
function [ab,normsq]=schebyshev(dig,N,mom,abm)
syms sig
digits(dig);
if N<=0, error('N out of range'), end
if N>size(mom,2)/2, N=size(mom,2)/2; end
if nargin<4, abm=zeros(2*N-1,2); end
if N>(size(abm,1)+1)/2, N=(size(abm,1)+1)/2; end
smom=vpa(mom,dig); sabm=vpa(abm,dig);
ab(1,1)=sabm(1,1)+smom(2)/smom(1); ab(1,2)=smom(1);
if N==1, normsq(1)=smom(1); return, end
sig(1,1:2*N)=0; sig(2,:)=smom(1:2*N);
for n=3:N+1
  for m=n-1:2*N-n+2
    sig(n,m)=sig(n-1,m+1)-(ab(n-2,1)-sabm(m,1))*sig(n-1,m) ...
      -ab(n-2,2)*sig(n-2,m)+sabm(m,2)*sig(n-1,m-1);
  end
  ab(n-1,1)=sabm(n-1,1)+sig(n,n)/sig(n,n-1)-sig(n-1,n-1)/ ...
    sig(n-1,n-2);
  ab(n-1,2)=sig(n,n-1)/sig(n-1,n-2);
end
for n=1:N, normsq(n)=sig(n+1,n); end; normsq=normsq';