极点零点法分析频响特性，哪位大神能指出其中的错误，急啊！

clear
Hz=sym('2/(z-1)+1/(z-3)+4/(z-5)');
[E,F]=numden(Hz); %求分子、分母多项式
Z0=solve(E); %求零点
P0=solve(F); %求极点
z=double(Z0); %零点z双精度
p=double(P0); %极点p双精度
M=length(z); %零点个数M
N=length(p); %极点个数N
C=sym2poly(E); %分母多项式系数
D=sym2poly(F); %分子多项式系数
A=C(1)/D(1); %求幅度A
syms w H t g f g1 f1 t1 y
for loop1=1:M
x(loop1)=abs(exp(i*w)-z(loop1)); %零点向量的模x
end
for loop2=1:N
y(loop2)=abs(exp(i*w)-p(loop2)); %极点向量的模y
end
for loop3=1:M
H=A.*x(loop3);
end
for loop4=1:N
H=H./y(loop4); %H为幅频特性
end
for loop5=1:M
g=sin(w)./(cos(w)-z(loop5));
f=atan(t);
a(loop5)=compose(f,g,t,w); %零点向量的相角a
end
for loop6=1:N
g1=sin(w)./(cos(w)-p(loop6));
f1=atan(t1);
b(loop6)=compose(f1,g1,t1,w); %极点向量的相角b
end
for loop7=1:M
y=(N-M)*w+a(loop7);
end
for loop8=1:N
y=y-b(loop8); %y为相频特性
end
w=0:pi/50:2*pi;
h=eval(H); %幅频特性数据化
o=eval(y); %相频特性数据化
subplot(2,2,1)
zplane(z,p);
subplot(2,2,2)
plot(w,h);
subplot(2,2,3)
plot(w,o);