MATLAB Answers

Help with creating my own Nyquist plotting function

4 views (last 30 days)
Hi everyone,
I am trying to make a function that takes a transfer function L(s) in the form of an array of numerator and denominator coefficients and spits out the nyquist plot, avoiding imaginary poles by drawing a semi-circle on the right half plane with radius epsilon around those imaginary poles and ending with an s-plane contour that is a semi-circle on the right half plane with radius R (see attached image)
So far, this is my code:
function n = nyquill(N,D,R,epsilon)
%N is numerator coefficients, D is denominator coefffs
L = tf(N,D);
t = pole(L)
r = NaN(3,1);
d = zeros(3,1);
for n = 1:size(t)
if real(t(n)) == 0
t(n) = r(n);
d(n) = NaN;
else
t(n) = d(n);
r(n) = NaN;
end
end
if r(1)==NaN & r(2)==NaN & r(3)==NaN
g1 = [-R:0.1:R]*i;
g2 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2];
elseif r(1)==NaN & r(2)==NaN
g1 = [-R;0.1:imag(r(3))-epsilon]*i;
g2 = imag(r(3)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g3 = [imag(r(3)) + epsilon:.01:R]*i;
g4 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2 g3 g4];
elseif r(1)==NaN
g1 = [-R;0.1:imag(r(2))-epsilon]*i;
g2 = imag(r(2)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g3 = [imag(r(2)) + epsilon:.01:imag(r(3))-epsilon]*i;
g4 = imag(r(3)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g5 = [imag(r(3)) + epsilon:.01:R]*i;
g6 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2 g3 g4 g5 g6];
elseif r(2)==NaN
g1 = [-R;0.1:imag(r(1))-epsilon]*i;
g2 = imag(r(1)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g3 = [imag(r(1)) + epsilon:.01:imag(r(3))-epsilon]*i;
g4 = imag(r(3)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g5 = [imag(r(3)) + epsilon:.01:R]*i;
g6 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2 g3 g4 g5 g6];
elseif r(3)==NaN
g1 = [-R;0.1:imag(r(1))-epsilon]*i;
g2 = imag(r(1)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g3 = [imag(r(1)) + epsilon:.01:imag(r(2))-epsilon]*i;
g4 = imag(r(2)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g5 = [imag(r(2)) + epsilon:.01:R]*i;
g6 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2 g3 g4 g5 g6];
elseif r(2)==NaN & r(3)==NaN
g1 = [-R;0.1:imag(r(1))-epsilon]*i;
g2 = imag(r(1)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g3 = [imag(r(1)) + epsilon:.01:R]*i;
g4 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2 g3 g4];
elseif r(1)==NaN & r(3)==NaN
g1 = [-R;0.1:imag(r(2))-epsilon]*i;
g2 = imag(r(2)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g3 = [imag(r(2)) + epsilon:.01:R]*i;
g4 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2 g3 g4];
elseif r(1)~=NaN & r(2)~=NaN & r(3)~=NaN
g1 = [-R;0.1:imag(r(1))-epsilon]*i;
g2 = imag(r(1)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g3 = [imag(r(1)) + epsilon:.01:imag(r(2))-epsilon]*i;
g4 = imag(r(2)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g5 = [imag(r(2)) + epsilon:.01:imag(r(3))-epsilon]*i;
g6 = imag(r(3)) + epsilon*exp(i*[-pi/2:pi/90:pi/2]);
g7 = [imag(r(3)) + epsilon:.01:R]*i;
g8 = R*exp(i*[pi/2:-pi/90:-pi/2]);
g = [g1 g2 g3 g4 g5 g6 g7 g8];
end
c = poly2sym(N);
d = poly2sym(D);
ln = c/d;
x = g;
subs(ln);
plot(real(ln),imag(ln));grid on;axis('equal')
hold on
plot(-1,0, '-o')
end
And I end up receiving the error:
Error using plot
Data must be numeric, datetime, duration or an array convertible to double.
Error in nyquill (line 83)
plot(real(ln),imag(ln));grid on;axis('equal')
I don't know what this error means or exactly how to fix it, any help would be appreciated!

  2 Comments

Daniel Beeson
Daniel Beeson on 24 May 2020
The inputs N and D are the numerator and denominator coefficients
Daniel Beeson
Daniel Beeson on 24 May 2020
The idea is to avoid the poles that lie on the imaginary axis by going around them with a small semi circle of radius epsilon that only goes into the right half plane

Sign in to comment.

Accepted Answer

Daniel Beeson
Daniel Beeson on 24 May 2020
Hey everyone, I figured it out. My conditions for NaN were incorrect. instead of using ==NaN, I changed it to isnan(r(...))==1 for if the value of r was NaN or ==0 for if it is not. Thank you everyone for the help!

  0 Comments

Sign in to comment.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!