data, omega=None, *fmt, grid=None, title=None, ax=None,

label=None, **kwargs):

"""Nichols plot for a system.

# Get parameter values

grid = config._get_param('nichols', 'grid', grid, True)

label = _process_line_labels(label)

rcParams = config._get_param('ctrlplot', 'rcParams', kwargs, pop=True)

# If argument was a singleton, turn it into a list

if not isinstance(data, (tuple, list)):

data = [data]

# If we were passed a list of systems, convert to data

if all([isinstance(

sys, (StateSpace, TransferFunction)) for sys in data]):

data = frequency_response(data, omega=omega)

# Make sure that all systems are SISO

if any([resp.ninputs > 1 or resp.noutputs > 1 for resp in data]):

raise NotImplementedError("MIMO Nichols plots not implemented")

fig, ax_nichols = _process_ax_keyword(ax, rcParams=rcParams, squeeze=True)

legend_loc, _, show_legend = _process_legend_keywords(

kwargs, None, 'upper left')

# Create a list of lines for the output

out = np.empty(len(data), dtype=object)

for idx, response in enumerate(data):

# Get the magnitude and phase of the system

mag = np.squeeze(response.magnitude)

phase = np.squeeze(response.phase)

omega = response.omega

# Convert to Nichols-plot format (phase in degrees,

# and magnitude in dB)

x = unwrap(np.degrees(phase), 360)

y = 20*np.log10(mag)

# Decide on the system name and label

sysname = response.sysname if response.sysname is not None \

else f"Unknown-sys_{idx}"

label_ = sysname if label is None else label[idx]

# Generate the plot

out[idx] = ax_nichols.plot(x, y, *fmt, label=label_, **kwargs)

# Label the plot axes

ax_nichols.set_xlabel('Phase [deg]')

ax_nichols.set_ylabel('Magnitude [dB]')

# Mark the -180 point

ax_nichols.plot([-180], [0], 'r+')

# Add grid

if grid:

nichols_grid(ax=ax_nichols)

# List of systems that are included in this plot

lines, labels = _get_line_labels(ax_nichols)

# Add legend if there is more than one system plotted

if show_legend == True or (show_legend != False and len(labels) > 1):

with plt.rc_context(rcParams):

legend = ax_nichols.legend(lines, labels, loc=legend_loc)

else:

legend = None

# Add the title

if ax is None:

if title is None:

title = "Nichols plot for " + ", ".join(labels)

_update_plot_title(

title, fig=fig, rcParams=rcParams, use_existing=False)

# intersection of data and view extents

# if intersection empty, return view extents

_inner = matplotlib.transforms.Bbox.intersection(ax.viewLim, ax.dataLim)

if _inner is None:

return ax.ViewLim.extents

else:

label_cl_phases=True):

"""Plot Nichols chart grid.

if ax is None:

ax = plt.gca()

# Default chart size

ol_phase_min = -359.99

ol_phase_max = 0.0

ol_mag_min = -40.0

ol_mag_max = default_ol_mag_max = 50.0

if ax.has_data():

# Find extent of intersection the current dataset or view

ol_phase_min, ol_mag_min, ol_phase_max, ol_mag_max = _inner_extents(ax)

# M-circle magnitudes.

if cl_mags is None:

# Default chart magnitudes

# The key set of magnitudes are always generated, since this

# guarantees a recognizable Nichols chart grid.

key_cl_mags = np.array([-40.0, -20.0, -12.0, -6.0, -3.0, -1.0, -0.5,

0.0, 0.25, 0.5, 1.0, 3.0, 6.0, 12.0])

# Extend the range of magnitudes if necessary. The extended arange

# will end up empty if no extension is required. Assumes that

# closed-loop magnitudes are approximately aligned with open-loop

# magnitudes beyond the value of np.min(key_cl_mags)

cl_mag_step = -20.0 # dB

extended_cl_mags = np.arange(np.min(key_cl_mags),

ol_mag_min + cl_mag_step, cl_mag_step)

cl_mags = np.concatenate((extended_cl_mags, key_cl_mags))

# a minimum 360deg extent containing the phases

phase_round_max = 360.0*np.ceil(ol_phase_max/360.0)

phase_round_min = min(phase_round_max-360,

360.0*np.floor(ol_phase_min/360.0))

# N-circle phases (should be in the range -360 to 0)

if cl_phases is None:

# aim for 9 lines, but always show (-360+eps, -180, -eps)

# smallest spacing is 45, biggest is 180

phase_span = phase_round_max - phase_round_min

spacing = np.clip(round(phase_span / 8 / 45) * 45, 45, 180)

key_cl_phases = np.array([-0.25, -359.75])

other_cl_phases = np.arange(-spacing, -360.0, -spacing)

cl_phases = np.unique(np.concatenate((key_cl_phases, other_cl_phases)))

elif not ((-360 < np.min(cl_phases)) and (np.max(cl_phases) < 0.0)):

raise ValueError('cl_phases must between -360 and 0, exclusive')

# Find the M-contours

m = m_circles(cl_mags, phase_min=np.min(cl_phases),

phase_max=np.max(cl_phases))

m_mag = 20*np.log10(np.abs(m))

m_phase = np.mod(np.degrees(np.angle(m)), -360.0) # Unwrap

# Find the N-contours

n = n_circles(cl_phases, mag_min=np.min(cl_mags), mag_max=np.max(cl_mags))

n_mag = 20*np.log10(np.abs(n))

n_phase = np.mod(np.degrees(np.angle(n)), -360.0) # Unwrap

# Plot the contours behind other plot elements.

# The "phase offset" is used to produce copies of the chart that cover

# the entire range of the plotted data, starting from a base chart computed

# over the range -360 < phase < 0. Given the range

# the base chart is computed over, the phase offset should be 0

# for -360 < ol_phase_min < 0.

phase_offsets = 360 + np.arange(phase_round_min, phase_round_max, 360.0)

cl_mag_lines = []

cl_phase_lines = []

cl_mag_labels = []

cl_phase_labels = []

for phase_offset in phase_offsets:

# Draw M and N contours

cl_mag_lines.extend(

ax.plot(m_phase + phase_offset, m_mag, color='lightgray',

linestyle=line_style, zorder=0))

cl_phase_lines.extend(

ax.plot(n_phase + phase_offset, n_mag, color='lightgray',

linestyle=line_style, zorder=0))

# Add magnitude labels

for x, y, m in zip(m_phase[:][-1] + phase_offset, m_mag[:][-1],

cl_mags):

align = 'right' if m < 0.0 else 'left'

cl_mag_labels.append(

ax.text(x, y, str(m) + ' dB', size='small', ha=align,

color='gray', clip_on=True))

# phase labels

if label_cl_phases:

for x, y, p in zip(n_phase[:][0] + phase_offset,

n_mag[:][0],

cl_phases):

if p > -175:

align = 'right'

elif p > -185:

align = 'center'

else:

align = 'left'

cl_phase_labels.append(

ax.text(x, y, f'{round(p)}\N{DEGREE SIGN}',

size='small',

ha=align,

va='bottom',

color='gray',

clip_on=True))

# Fit axes to generated chart

ax.axis([phase_round_min,

phase_round_max,

np.min(np.concatenate([cl_mags,[ol_mag_min]])),

np.max([ol_mag_max, default_ol_mag_max])])

"""Contours of the function Gcl = Gol/(1+Gol), where

# Compute the contours in Gcl-space. Since we're given closed-loop

# magnitudes and phases, this is just a case of converting them into

# a complex number.

Gcl = Gcl_mags*np.exp(1.j*Gcl_phases)

# Invert Gcl = Gol/(1+Gol) to map the contours into the open-loop space

"""Constant-magnitude contours of the function Gcl = Gol/(1+Gol), where

# Convert magnitudes and phase range into a grid suitable for

# building contours

phases = np.radians(np.linspace(phase_min, phase_max, 2000))

Gcl_mags, Gcl_phases = np.meshgrid(10.0**(mags/20.0), phases)

"""Constant-phase contours of the function Gcl = Gol/(1+Gol), where

# Convert phases and magnitude range into a grid suitable for

# building contours

mags = np.linspace(10**(mag_min/20.0), 10**(mag_max/20.0), 2000)

Gcl_phases, Gcl_mags = np.meshgrid(np.radians(phases), mags)