Skip to content

Commit fed6c8d

Browse files
authored
Merge pull request #1 from murrayrm/fix_scalar_tfops-06Dec2024
Fix scalar tfops 06 dec2024
2 parents 1cc84a7 + 9acbfbb commit fed6c8d

15 files changed

Lines changed: 666 additions & 196 deletions

‎control/exception.py‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -52,6 +52,10 @@ class ControlArgument(TypeError):
5252
"""Raised when arguments to a function are not correct"""
5353
pass
5454

55+
class ControlIndexError(IndexError):
56+
"""Raised when arguments to an indexed object are not correct"""
57+
pass
58+
5559
class ControlMIMONotImplemented(NotImplementedError):
5660
"""Function is not currently implemented for MIMO systems"""
5761
pass

‎control/frdata.py‎

Lines changed: 127 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,7 @@
1010
FRD data.
1111
"""
1212

13+
from collections.abc import Iterable
1314
from copy import copy
1415
from warnings import warn
1516

@@ -20,7 +21,8 @@
2021

2122
from . import config
2223
from .exception import pandas_check
23-
from .iosys import InputOutputSystem, _process_iosys_keywords, common_timebase
24+
from .iosys import InputOutputSystem, NamedSignal, _process_iosys_keywords, \
25+
_process_subsys_index, common_timebase
2426
from .lti import LTI, _process_frequency_response
2527

2628
__all__ = ['FrequencyResponseData', 'FRD', 'frd']
@@ -33,8 +35,8 @@ class FrequencyResponseData(LTI):
3335
3436
The FrequencyResponseData (FRD) class is used to represent systems in
3537
frequency response data form. It can be created manually using the
36-
class constructor, using the :func:~~control.frd` factory function
37-
(preferred), or via the :func:`~control.frequency_response` function.
38+
class constructor, using the :func:`~control.frd` factory function or
39+
via the :func:`~control.frequency_response` function.
3840
3941
Parameters
4042
----------
@@ -65,6 +67,28 @@ class constructor, using the :func:~~control.frd` factory function
6567
frequency point.
6668
dt : float, True, or None
6769
System timebase.
70+
squeeze : bool
71+
By default, if a system is single-input, single-output (SISO) then
72+
the outputs (and inputs) are returned as a 1D array (indexed by
73+
frequency) and if a system is multi-input or multi-output, then the
74+
outputs are returned as a 2D array (indexed by output and
75+
frequency) or a 3D array (indexed by output, trace, and frequency).
76+
If ``squeeze=True``, access to the output response will remove
77+
single-dimensional entries from the shape of the inputs and outputs
78+
even if the system is not SISO. If ``squeeze=False``, the output is
79+
returned as a 3D array (indexed by the output, input, and
80+
frequency) even if the system is SISO. The default value can be set
81+
using config.defaults['control.squeeze_frequency_response'].
82+
ninputs, noutputs, nstates : int
83+
Number of inputs, outputs, and states of the underlying system.
84+
input_labels, output_labels : array of str
85+
Names for the input and output variables.
86+
sysname : str, optional
87+
Name of the system. For data generated using
88+
:func:`~control.frequency_response`, stores the name of the system
89+
that created the data.
90+
title : str, optional
91+
Set the title to use when plotting.
6892
6993
See Also
7094
--------
@@ -89,6 +113,20 @@ class constructor, using the :func:~~control.frd` factory function
89113
the imaginary access). See :meth:`~control.FrequencyResponseData.__call__`
90114
for a more detailed description.
91115
116+
A state space system is callable and returns the value of the transfer
117+
function evaluated at a point in the complex plane. See
118+
:meth:`~control.StateSpace.__call__` for a more detailed description.
119+
120+
Subsystem response corresponding to selected input/output pairs can be
121+
created by indexing the frequency response data object::
122+
123+
subsys = sys[output_spec, input_spec]
124+
125+
The input and output specifications can be single integers, lists of
126+
integers, or slices. In addition, the strings representing the names
127+
of the signals can be used and will be replaced with the equivalent
128+
signal offsets.
129+
92130
"""
93131
#
94132
# Class attributes
@@ -223,14 +261,20 @@ def __init__(self, *args, **kwargs):
223261

224262
# create interpolation functions
225263
if smooth:
264+
# Set the order of the fit
265+
if self.omega.size < 2:
266+
raise ValueError("can't smooth with only 1 frequency")
267+
degree = 3 if self.omega.size > 3 else self.omega.size - 1
268+
226269
self.ifunc = empty((self.fresp.shape[0], self.fresp.shape[1]),
227270
dtype=tuple)
228271
for i in range(self.fresp.shape[0]):
229272
for j in range(self.fresp.shape[1]):
230273
self.ifunc[i, j], u = splprep(
231274
u=self.omega, x=[real(self.fresp[i, j, :]),
232275
imag(self.fresp[i, j, :])],
233-
w=1.0/(absolute(self.fresp[i, j, :]) + 0.001), s=0.0)
276+
w=1.0/(absolute(self.fresp[i, j, :]) + 0.001),
277+
s=0.0, k=degree)
234278
else:
235279
self.ifunc = None
236280

@@ -243,21 +287,72 @@ def __init__(self, *args, **kwargs):
243287

244288
@property
245289
def magnitude(self):
246-
return np.abs(self.fresp)
290+
"""Magnitude of the frequency response.
291+
292+
Magnitude of the frequency response, indexed by either the output
293+
and frequency (if only a single input is given) or the output,
294+
input, and frequency (for multi-input systems). See
295+
:attr:`FrequencyResponseData.squeeze` for a description of how this
296+
can be modified using the `squeeze` keyword.
297+
298+
Input and output signal names can be used to index the data in
299+
place of integer offsets.
300+
301+
:type: 1D, 2D, or 3D array
302+
303+
"""
304+
return NamedSignal(
305+
np.abs(self.fresp), self.output_labels, self.input_labels)
247306

248307
@property
249308
def phase(self):
250-
return np.angle(self.fresp)
309+
"""Phase of the frequency response.
310+
311+
Phase of the frequency response in radians/sec, indexed by either
312+
the output and frequency (if only a single input is given) or the
313+
output, input, and frequency (for multi-input systems). See
314+
:attr:`FrequencyResponseData.squeeze` for a description of how this
315+
can be modified using the `squeeze` keyword.
316+
317+
Input and output signal names can be used to index the data in
318+
place of integer offsets.
319+
320+
:type: 1D, 2D, or 3D array
321+
322+
"""
323+
return NamedSignal(
324+
np.angle(self.fresp), self.output_labels, self.input_labels)
251325

252326
@property
253327
def frequency(self):
328+
"""Frequencies at which the response is evaluated.
329+
330+
:type: 1D array
331+
332+
"""
254333
return self.omega
255334

256335
@property
257336
def response(self):
258-
return self.fresp
337+
"""Complex value of the frequency response.
338+
339+
Value of the frequency response as a complex number, indexed by
340+
either the output and frequency (if only a single input is given)
341+
or the output, input, and frequency (for multi-input systems). See
342+
:attr:`FrequencyResponseData.squeeze` for a description of how this
343+
can be modified using the `squeeze` keyword.
344+
345+
Input and output signal names can be used to index the data in
346+
place of integer offsets.
347+
348+
:type: 1D, 2D, or 3D array
349+
350+
"""
351+
return NamedSignal(
352+
self.fresp, self.output_labels, self.input_labels)
259353

260354
def __str__(self):
355+
261356
"""String representation of the transfer function."""
262357

263358
mimo = self.ninputs > 1 or self.noutputs > 1
@@ -303,7 +398,12 @@ def __add__(self, other):
303398

304399
# Convert the second argument to a frequency response function.
305400
# or re-base the frd to the current omega (if needed)
306-
other = _convert_to_frd(other, omega=self.omega)
401+
if isinstance(other, (int, float, complex, np.number)):
402+
other = _convert_to_frd(
403+
other, omega=self.omega,
404+
inputs=self.ninputs, outputs=self.noutputs)
405+
else:
406+
other = _convert_to_frd(other, omega=self.omega)
307407

308408
# Check that the input-output sizes are consistent.
309409
if self.ninputs != other.ninputs:
@@ -593,9 +693,25 @@ def __iter__(self):
593693
return iter((self.omega, fresp))
594694
return iter((np.abs(fresp), np.angle(fresp), self.omega))
595695

596-
# Implement (thin) getitem to allow access via legacy indexing
597-
def __getitem__(self, index):
598-
return list(self.__iter__())[index]
696+
def __getitem__(self, key):
697+
if not isinstance(key, Iterable) or len(key) != 2:
698+
# Implement (thin) getitem to allow access via legacy indexing
699+
return list(self.__iter__())[key]
700+
701+
# Convert signal names to integer offsets (via NamedSignal object)
702+
iomap = NamedSignal(
703+
self.fresp[:, :, 0], self.output_labels, self.input_labels)
704+
indices = iomap._parse_key(key, level=1) # ignore index checks
705+
outdx, outputs = _process_subsys_index(indices[0], self.output_labels)
706+
inpdx, inputs = _process_subsys_index(indices[1], self.input_labels)
707+
708+
# Create the system name
709+
sysname = config.defaults['iosys.indexed_system_name_prefix'] + \
710+
self.name + config.defaults['iosys.indexed_system_name_suffix']
711+
712+
return FrequencyResponseData(
713+
self.fresp[outdx, :][:, inpdx], self.omega, self.dt,
714+
inputs=inputs, outputs=outputs, name=sysname)
599715

600716
# Implement (thin) len to emulate legacy testing interface
601717
def __len__(self):

‎control/iosys.py‎

Lines changed: 94 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -13,9 +13,10 @@
1313
import numpy as np
1414

1515
from . import config
16+
from .exception import ControlIndexError
1617

17-
__all__ = ['InputOutputSystem', 'issiso', 'timebase', 'common_timebase',
18-
'isdtime', 'isctime']
18+
__all__ = ['InputOutputSystem', 'NamedSignal', 'issiso', 'timebase',
19+
'common_timebase', 'isdtime', 'isctime']
1920

2021
# Define module default parameter values
2122
_iosys_defaults = {
@@ -33,6 +34,69 @@
3334
}
3435

3536

37+
# Named signal class
38+
class NamedSignal(np.ndarray):
39+
def __new__(cls, input_array, signal_labels=None, trace_labels=None):
40+
# See https://numpy.org/doc/stable/user/basics.subclassing.html
41+
obj = np.asarray(input_array).view(cls) # Cast to our class type
42+
obj.signal_labels = signal_labels # Save signal labels
43+
obj.trace_labels = trace_labels # Save trace labels
44+
obj.data_shape = input_array.shape # Save data shape
45+
return obj # Return new object
46+
47+
def __array_finalize__(self, obj):
48+
# See https://numpy.org/doc/stable/user/basics.subclassing.html
49+
if obj is None:
50+
return
51+
self.signal_labels = getattr(obj, 'signal_labels', None)
52+
self.trace_labels = getattr(obj, 'trace_labels', None)
53+
self.data_shape = getattr(obj, 'data_shape', None)
54+
55+
def _parse_key(self, key, labels=None, level=0):
56+
if labels is None:
57+
labels = self.signal_labels
58+
try:
59+
if isinstance(key, str):
60+
key = labels.index(item := key)
61+
if level == 0 and len(self.data_shape) < 2:
62+
raise ControlIndexError
63+
elif isinstance(key, list):
64+
keylist = []
65+
for item in key: # use for loop to save item for error
66+
keylist.append(
67+
self._parse_key(item, labels=labels, level=level+1))
68+
if level == 0 and key != keylist and len(self.data_shape) < 2:
69+
raise ControlIndexError
70+
key = keylist
71+
elif isinstance(key, tuple) and len(key) > 0:
72+
keylist = []
73+
keylist.append(
74+
self._parse_key(
75+
item := key[0], labels=self.signal_labels,
76+
level=level+1))
77+
if len(key) > 1:
78+
keylist.append(
79+
self._parse_key(
80+
item := key[1], labels=self.trace_labels,
81+
level=level+1))
82+
if level == 0 and key[:len(keylist)] != tuple(keylist) \
83+
and len(keylist) > len(self.data_shape) - 1:
84+
raise ControlIndexError
85+
for i in range(2, len(key)):
86+
keylist.append(key[i]) # pass on remaining elements
87+
key = tuple(keylist)
88+
except ValueError:
89+
raise ValueError(f"unknown signal name '{item}'")
90+
except ControlIndexError:
91+
raise ControlIndexError(
92+
"signal name(s) not valid for squeezed data")
93+
94+
return key
95+
96+
def __getitem__(self, key):
97+
return super().__getitem__(self._parse_key(key))
98+
99+
36100
class InputOutputSystem(object):
37101
"""A class for representing input/output systems.
38102
@@ -965,3 +1029,31 @@ def _parse_spec(syslist, spec, signame, dictname=None):
9651029
ValueError(f"signal index '{index}' is out of range")
9661030

9671031
return system_index, signal_indices, gain
1032+
1033+
1034+
#
1035+
# Utility function for processing subsystem indices
1036+
#
1037+
# This function processes an index specification (int, list, or slice) and
1038+
# returns a index specification that can be used to create a subsystem
1039+
#
1040+
def _process_subsys_index(idx, sys_labels, slice_to_list=False):
1041+
if not isinstance(idx, (slice, list, int)):
1042+
raise TypeError("system indices must be integers, slices, or lists")
1043+
1044+
# Convert singleton lists to integers for proper slicing (below)
1045+
if isinstance(idx, (list, tuple)) and len(idx) == 1:
1046+
idx = idx[0]
1047+
1048+
# Convert int to slice so that numpy doesn't drop dimension
1049+
if isinstance(idx, int):
1050+
idx = slice(idx, idx+1, 1)
1051+
1052+
# Get label names (taking care of possibility that we were passed a list)
1053+
labels = [sys_labels[i] for i in idx] if isinstance(idx, list) \
1054+
else sys_labels[idx]
1055+
1056+
if slice_to_list and isinstance(idx, slice):
1057+
idx = range(len(sys_labels))[idx]
1058+
1059+
return idx, labels

‎control/nlsys.py‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -193,8 +193,8 @@ def __call__(sys, u, params=None, squeeze=None):
193193

194194
# Evaluate the function on the argument
195195
out = sys._out(0, np.array((0,)), np.asarray(u))
196-
_, out = _process_time_response(
197-
None, out, issiso=sys.issiso(), squeeze=squeeze)
196+
out = _process_time_response(
197+
out, issiso=sys.issiso(), squeeze=squeeze)
198198
return out
199199

200200
def __mul__(self, other):

0 commit comments

Comments
 (0)