function [x,y]=modelmag(sys)
% Determines the logarithm of the spectral density 
% of the model sys
%
% plot(x,y) draws a plot of the logarithm of the 
% spectral density.
%
% Author: Per Enqvist, April 2007

sys=tf(sys);

s=sys.num{1};
a=sys.den{1};
n=length(a)-1;

x=(1:180)*pi/180;
cos_matrix=cos((0:n)'*x);
a2=conv(a,fliplr(a)); 			% Determines a(z)*z^n a(1/z)
a2=[a2(n+1) 2*a2(n+2:2*n+1)]; 		% let a2=a2(1)+a2(2)(z+1/z)/2+...
y1=a2*cos((0:n)'*x); 		        % let  y=a2(1)+a2(2)cos(x)+...
s2=conv(s,fliplr(s)); 			% Determines s(z)*z^n s(1/z)
s2=[s2(n+1) 2*s2(n+2:2*n+1)]; 		% let s2=s2(1)+s2(2)(z+1/z)/2+...
y2=s2*cos((0:n)'*x); 		        % let  y=s2(1)+s2(2)cos(x)+...
y=10*log10(y2)-10*log10(y1);
x=x/pi;

