function table = spec_table(Z,P)
% table = spec_table(Z,P);
% Constructs a specification table for 
% generalized digital Butterworth filter design.
% Z : total number of zeros
% P : total number of poles
% Example:
%       Z = 10; P = 2;
%       table = spec_table(Z,P);
%
%  required subprograms : choose.m

if Z <= P
	table = [Z 0 P 0 1];
	break
end
Fo = (1/2)^2;   % value of mag squared at half-mag freq.

%%%%%%% begin with all zeros at z = -1 %%%%%%%
L = Z;
M = 0;
N = P;
if N > 0
	%%%%%%% construct polynomial for checking boundary %%%%%%
        i = N:-1:0;
        k = L:-1:N+1;
        Y = [choose(L,k).*(-1).^k, (1-Fo)*choose(L,i).*(-1).^i];
        if rem(N,2) == 1 
		Y(L+1-N) = Y(L+1-N) - choose(L-1,N)*Fo;
        end
	%%% extract appropriate root %%%
	if rem(N,2) == 1
                % since it is known that in this case Y(x) has a root
                % at x=1, (see table ? in text) we choose to remove it directly.
		[Y,r] = deconv(Y,[1 -1]);
	end
	r = roots(Y);
	r = r(imag(r)==0);
	[temp,i] = min(abs(r-0.5));
	wo_max = acos(1-2*r(i));
	table = [Z, 0, P, 0, wo_max/pi];
end
%%%%%%%%% increment number of passband %%%%%%
%%%%%%%%% zeros by one in a loop %%%%%%%%%%%%
for M = 1:Z-P-1
	L = Z-M;
	k = M+L:-1:0;
        GR = choose(M+N-k-1,M-1).*choose(M+L-1,k).*(-1).^k;
        GT = choose(M+N-k-1,M-1).*choose(M+L-1,k-1).*(-1).^(k-1);
	AGR = GR(M+L+1-N:M+L+1);
	AGT = GT(M+L+1-N:M+L+1);
	%%% ------- wo_max ---------------------------------- %%%
        if rem(N,2) == 0
                c = (L-N)/(M+N);  % N even
        else
                c = (L-N)/N;      % N odd
        end
        Y = [zeros(1,M+L-N), Fo*AGR+c*Fo*AGT]-GR-c*GT;
	%%%%%%% extract appropriate root %%%
        if rem(N,2) == 1
                % since it is known that in this case Y(x) has a root
                % at x=1, (see table ? in text) we choose to remove it directly.
                [Y,r] = deconv(Y,[1 -1]);
        end
        r = roots(Y);
	r = r(imag(r)==0);
	[temp,i] = min(abs(r-0.5));
        wo_max = acos(1-2*r(i));        
        %%% ------- wo_min ---------------------------------- %%%
	if N > 0
		wo_min = table(M,5)*pi;
	elseif M > 1
		wo_min = table(M-1,5)*pi;
	else
                % for FIR filters, wo_min on the first row of the spec
                % table is not 0, so we need to compute it.
                Y = [zeros(1,M+L-N), Fo*(AGR-AGT)]-GR+GT;
        	r = roots(Y);
        	r = r(imag(r)==0);
        	[temp,i] = min(abs(r-0.5));
        	wo_min = acos(1-2*r(i));
	end
        table = [table; [L,M,N, wo_min/pi, wo_max/pi]];
end

if N > 0
	table = [table; [P,Z-P,N, table(Z-P,5), 1]];
end