function nyq(a,b,c,r,p,q)

% Blood Model (z + a)(z + b) + c*e^(-r*z)
% p = x_max, q = y_max

h=p/100;
k=q/50;
x(1) = 0;
y(1) = q;
u1(1) = (x(1)+a)*(x(1)+b)-y(1)^2+c*exp(-r*x(1))*cos(r*y(1));
w1(1) = y(1)*(2*x(1)+a+b)-c*exp(-r*x(1))*sin(r*y(1));
for i=2:101
    x(i) = 0;
    y(i) = y(i-1)-k;
    u1(i) = (x(i)+a)*(x(i)+b)-y(i)^2+c*exp(-r*x(i))*cos(r*y(i));
    w1(i) = y(i)*(2*x(i)+a+b)-c*exp(-r*x(i))*sin(r*y(i));
end
for i=102:201
    x(i) = x(i-1)+h;
    y(i) = -q;
    u2(i-101) = (x(i)+a)*(x(i)+b)-y(i)^2+c*exp(-r*x(i))*cos(r*y(i));
    w2(i-101) = y(i)*(2*x(i)+a+b)-c*exp(-r*x(i))*sin(r*y(i));
end
for i=202:301
    x(i) = p;
    y(i) = y(i-1)+k;
    u3(i-201) = (x(i)+a)*(x(i)+b)-y(i)^2+c*exp(-r*x(i))*cos(r*y(i));
    w3(i-201) = y(i)*(2*x(i)+a+b)-c*exp(-r*x(i))*sin(r*y(i));
end
for i=302:401
    x(i) = x(i-1)-h;
    y(i) = q;
    u4(i-301) = (x(i)+a)*(x(i)+b)-y(i)^2+c*exp(-r*x(i))*cos(r*y(i));
    w4(i-301) = y(i)*(2*x(i)+a+b)-c*exp(-r*x(i))*sin(r*y(i));
end
plot(u1,w1,'b-',u2,w2,'r-',u3,w3,'g-',u4,w4,'m-');grid

% Normal
% nyq(2.8,0.001,0.17,6,2,1)
% Hemolytic
% nyq(6.65,0.08,4.73,4,3,3)
