Nyquist plot rescaling is confusing

Hi - not sure I fully understand how the indentation around an imaginary axis pole is being done for the nyquist plot.

Looking at freqplot.py L 1449:1555 I see:

# Add points for half/quarter circle around pole frequency
# (these will get indented left or right below)
splane_contour = np.concatenate((
    splane_contour[0:first_point+1],
    (1j * np.linspace(
        start_freq, p.imag + indent_radius, indent_points)),
    splane_contour[last_point:]))

but this seems to make all the added points imaginary - what I would have expected to see would be a quarter/half circle with
real and imaginary parts -- something like:

  if indent_direction == 'right':
      direction = 1
  else:
      direction = -1
  angle_start = [0 if start_freq == 0 else -np.pi/2][0]
  angle_end = np.pi/2

  splane_contour = np.concatenate((
      splane_contour[0:first_point+1],
      p + direction*indent_radius*np.exp(direction*1j * np.linspace(
          angle_start, angle_end, indent_points)),
      splane_contour[last_point:]))

with the phase being from -pi/2:pi/2 (unless at origin) with a direction determined by the indent direction
and np.exp(...) creating the small complex valued quarter/half circle centered at pole p (on imag axis)

this change seems to break the encirclement math that follows in the code, so I cannot provide a working example

think I see my error - the code 1458-1482 is doing the indent - I thought it was correcting for something else

You can see exactly what the contour is using the contour attribute from the map returned by nyquist_response:

import matplotlib.pyplot as plt
import control as ct

sys = ct.zpk([], [0, +5j, -5j], 10, dt=0)  # dt required prior to v0.10.2
resp = ct.nyquist_response(sys)

plt.plot(resp.contour.real, resp.contour.imag)
plt.axis([-1e-4, 2e-4, 0, 2e-4])

Hi - yes, I see that now - the artifact that I was trying to track down seems to come from L1908
resp[rescale] *= max_curve_magnitude / abs(resp[rescale])

I think I understand the intent of this line, but compare the 2 nyquist plots you get if that line was changed to:
resp *= max_curve_magnitude / np.max(np.abs(resp)) # just a regular scaling of entire plot

I might be missing something, but I am not convinced the current "resp[rescale] *= max_curve_magnitude / abs(resp[rescale])"
is right

example case

from control import nyquist_plot, zpk, nyquist_response, tf
import matplotlib.pyplot as plt
import numpy as np

G2 = tf([1],[1, 1]) * tf([1],[1, 0])**2

fig, ax = plt.subplots(1,figsize=(4, 4),dpi=150)
resp = nyquist_response(G2,indent_points=100,indent_direction='left',indent_radius=0.1)
resp.plot(title='Nyquist Plot - left Indentation')
plt.xlim(-30,0)
plt.ylim(0,30)

I'll look into this further. I don't think we want to rescale the entire plot, since that could change the number of encirclements that show up on the plot.

agreed - I was just trying to find the easiest way to show the weird bit that appears in the bottom plot above (around -20,3)

I've been looking into what is going on in the rescaling and how to fix it.

First, what is going on is that the way in which the nonlinear rescaling is done can distort the shape of the Nyquist curve. Consider the system from above, with the Nyquist indentation around the poles set to be on the left:

import control as ct

G = ct.tf([1], [1, 1]) * ct.tf([1], [1, 0])**2
ct.nyquist_plot(G, indent_dir='left')

This is a system with two poles at the origin and a pole at -1. The default Nyquist plot for the system looks like this:

From the plot, it is not that clear what is happening at $-\infty$. In particular, the apparent reversal in the direction on the imaginary axis seems odd, since the magnitude and phase of the frequency response are both monotonic:

ct.bode(G)

We can get a better view if we turn off the rescaling and adjust the size of the indentation around the poles of the origin:

ct.nyquist_plot(G, indent_direction='left', indent_radius=0.05, maximum_curve_magnitude=1000)

In this view we can no longer see the -1 point, but we see that there is a clear gap between the Nyquist curve and the negative imaginary axis. This (quite large) gap is not at all event in the default plot. Note also that the inversion in the direction of the curve is now gone.

The reason for the distortion is the way in which the rescaling is done. Currently any point that has magnitude greater than max_curve_magnitude is linearly rescaled in magnitude to map to a circle of radius max_curve_magnitude. If we were to apply that rescaling to the plot with an indentation radius of 0.05, the point at roughly (-400, 18) maps to something like (-20, 0.9). However, the unscaled portion of the curve has a point at approximately (-19, 4) and so the rescaled point moves down in the rescaled plot, even though in the original plot the imaginary component of the curve is larger (18 vs 4).

An easy fix is to change the indentation radius so that the magnitude of the curve near the pole does not exceed the maximum curve magnitude. For a double pole at the origin, if we want the magnitude to be less than 20 when we evaluation at the point $\epsilon$, we should set $\epsilon = 1/\sqrt{20} \approx 0.225$. Using that value gives a much nicer looking plot:

ct.nyquist_plot(G, indent_direction='left', indent_radius=0.225)

We see that comparing this plot to the "unscaled" plot gives a quite nice representation with no major distortions in the shape.

Now that we understand the problem, the question is how to best fix it. I can see a couple of paths:

Set the indentation radius based on max_curve_magnitude and the location of any poles that are near the imaginary axis.
Independently rescale the real and imaginary axes when their values are bigger than max_curve_magnitude.
Figure out some sort of nonlinear rescaling of the outer portions of the Nyquist plot so that you can see the -1 point at the same time as not exceeding max_curve_magnitude.

Before talking about these options, one other important point to remember is that the generation of the Nyquist response is actually separate from the plotting of that response. In the commands above, I have used nyquist_plot for everything, but underneath the hood what is happening looks more like this:

resp = ct.nyquist_response(G, indent_direction='left', indent_radius=0.05)
cplt = resp.plot(max_curve_magnitude=1000)

The reason this matters is indicated in the use of the parameters above: the choice of the indentation direction and radius needs to happen at the time of generating the response, while the choice of the maximum magnitude that you plot happens at the time of creating the plot.

Now to the options:

I see two issues with this option:
- If the indentation radius is large, there might be situations where you have a closed loop pole "inside" the indentation circle and this would get missed (remembering that the encirclement count tells you about poles inside the Nyquist "D" curve => if you have a large indentation then you won't take poles in the indented region into account).
- This is a bit more minor, but this option requires knowledge of plotting parameters (max_curve_magnitude) at the time of response generation.
This approach works in principle, but the plots it generates are quite ugly. For the system above, the plot that I get from doing the independent scaling looks like this:
I like this idea the best, but I haven't yet been able to figure out how to implement it in a way that gives a good representation the curve using only the information that is in the response.

My current plan is to implement a two-pronged solution: allow automatic computation of indent_radius from max_curve_magnitude in the nyquist_response command (perhaps with a warning if there are closed loop poles inside the indentation circle) and redo the scaling in nyquist_plot so that it does something that doesn't distort the shape of the curve but generates a circle of radius max_curve_magnitude for points "at infinity".

I also think drawing negative frequencies in Nyquist plots is unnecessary if the system is not complex valued adds no value but confuses the directions when there are axis crossings. As far as I can look into it historically (which is not much), it's just out of habit which is often more powerful than any logical argument.

It also pairs better with one-sided Bode plots.

I also think drawing negative frequencies in Nyquist plots is unnecessary if the system is not complex valued adds no value but confuses the directions when there are axis crossings. As far as I can look into it historically (which is not much), it's just out of habit which is often more powerful than any logical argument.

It also pairs better with one-sided Bode plots.

ct.set_defaults('nyquist', mirror_style=False)
ct.nyquist_plot(G, indent_direction='left', indent_radius=0.225)

I think having the negative freq part in the plot is helpful for those learning to visualize encirclements - so I would at least leave it as an option to display both parts. What I had been doing recently was the equivalent of RM solution #1 and creating 2 graphs, one just a zoomed in version (around -1) of the larger. Suboptimal to have 2 plots, but when I tried to show the zoomed in part as an inset, I noticed another behavior I didn't understand regarding axes. For example, if you do this:

import control as ct
import numpy as np
import matplotlib.pyplot as plt

G = ct.tf([1], [1, 1]) * ct.tf([1], [1, 0])**2

fig, ax = plt.subplots(1,2,figsize=(8, 4),dpi=150)

resp = ct.nyquist_response(G, indent_direction='left', indent_radius=0.05)
cplt = resp.plot(max_curve_magnitude=1000,ax=ax[0])

resp = ct.nyquist_response(G, indent_direction='left', indent_radius=0.05)
cplt = resp.plot(max_curve_magnitude=100,ax=ax[1])

both plots appear on ax[1] as opposed to what I want, which is 1 each on ax[0] and ax[1].
resp.plot seems to accept ax as an input - not sure what it does with it though

Navigation Menu

Repository navigation

Nyquist plot rescaling is confusing #1143

example case

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Nyquist plot rescaling is confusing #1143

Description

Activity

JonHowMIT commented on Mar 24, 2025

murrayrm commented on Mar 24, 2025

JonHowMIT commented on Mar 24, 2025

example case

murrayrm commented on Mar 24, 2025

JonHowMIT commented on Mar 24, 2025

murrayrm commented on Apr 5, 2025

ilayn commented on Apr 5, 2025

murrayrm commented on Apr 5, 2025

JonHowMIT commented on Apr 5, 2025

Metadata

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Issue actions