Repository navigation
Expand file tree
/
Copy pathtimeresp.py
More file actions
2273 lines (1854 loc) · 86.3 KB
/
Copy pathtimeresp.py
File metadata and controls
2273 lines (1854 loc) · 86.3 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
"""
timeresp.py - time-domain simulation routines.
The :mod:`~control.timeresp` module contains a collection of
functions that are used to compute time-domain simulations of LTI
systems.
Arguments to time-domain simulations include a time vector, an input
vector (when needed), and an initial condition vector. The most
general function for simulating LTI systems the
:func:`forced_response` function, which has the form::
t, y = forced_response(sys, T, U, X0)
where `T` is a vector of times at which the response should be
evaluated, `U` is a vector of inputs (one for each time point) and
`X0` is the initial condition for the system.
See :ref:`time-series-convention` for more information on how time
series data are represented.
Copyright (c) 2011 by California Institute of Technology
All rights reserved.
Copyright (c) 2011 by Eike Welk
Copyright (c) 2010 by SciPy Developers
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions
are met:
1. Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright
notice, this list of conditions and the following disclaimer in the
documentation and/or other materials provided with the distribution.
3. Neither the name of the California Institute of Technology nor
the names of its contributors may be used to endorse or promote
products derived from this software without specific prior
written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH
OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
SUCH DAMAGE.
Initial Author: Eike Welk
Date: 12 May 2011
Modified: Sawyer B. Fuller (minster@uw.edu) to add discrete-time
capability and better automatic time vector creation
Date: June 2020
Modified by Ilhan Polat to improve automatic time vector creation
Date: August 17, 2020
Modified by Richard Murray to add TimeResponseData class
Date: August 2021
$Id$
"""
import warnings
from copy import copy
import numpy as np
import scipy as sp
from numpy import einsum, maximum, minimum
from scipy.linalg import eig, eigvals, matrix_balance, norm
from . import config
from .exception import pandas_check
from .iosys import NamedSignal, isctime, isdtime
from .timeplot import time_response_plot
__all__ = ['forced_response', 'step_response', 'step_info',
'initial_response', 'impulse_response', 'TimeResponseData',
'TimeResponseList']
class TimeResponseData:
"""A class for returning time responses.
This class maintains and manipulates the data corresponding to the
temporal response of an input/output system. It is used as the return
type for time domain simulations (step response, input/output response,
etc).
A time response consists of a time vector, an output vector, and
optionally an input vector and/or state vector. Inputs and outputs can
be 1D (scalar input/output) or 2D (vector input/output).
A time response can be stored for multiple input signals (called traces),
with the output and state indexed by the trace number. This allows for
input/output response matrices, which is mainly useful for impulse and
step responses for linear systems. For multi-trace responses, the same
time vector must be used for all traces.
Time responses are accessed through either the raw data, stored as
:attr:`t`, :attr:`y`, :attr:`x`, :attr:`u`, or using a set of properties
:attr:`time`, :attr:`outputs`, :attr:`states`, :attr:`inputs`. When
accessing time responses via their properties, squeeze processing is
applied so that (by default) single-input, single-output systems will have
the output and input indices suppressed. This behavior is set using the
``squeeze`` keyword.
Attributes
----------
t : 1D array
Time values of the input/output response(s). This attribute is
normally accessed via the :attr:`time` property.
y : 2D or 3D array
Output response data, indexed either by output index and time (for
single trace responses) or output, trace, and time (for multi-trace
responses). These data are normally accessed via the :attr:`outputs`
property, which performs squeeze processing.
x : 2D or 3D array, or None
State space data, indexed either by output number and time (for single
trace responses) or output, trace, and time (for multi-trace
responses). If no state data are present, value is ``None``. These
data are normally accessed via the :attr:`states` property, which
performs squeeze processing.
u : 2D or 3D array, or None
Input signal data, indexed either by input index and time (for single
trace responses) or input, trace, and time (for multi-trace
responses). If no input data are present, value is ``None``. These
data are normally accessed via the :attr:`inputs` property, which
performs squeeze processing.
squeeze : bool, optional
By default, if a system is single-input, single-output (SISO)
then the outputs (and inputs) are returned as a 1D array
(indexed by time) and if a system is multi-input or
multi-output, then the outputs are returned as a 2D array
(indexed by output and time) or a 3D array (indexed by output,
trace, and time). If ``squeeze=True``, access to the output
response will remove single-dimensional entries from the shape
of the inputs and outputs even if the system is not SISO. If
``squeeze=False``, the output is returned as a 2D or 3D array
(indexed by the output [if multi-input], trace [if multi-trace]
and time) even if the system is SISO. The default value can be
set using config.defaults['control.squeeze_time_response'].
transpose : bool, optional
If True, transpose all input and output arrays (for backward
compatibility with MATLAB and :func:`scipy.signal.lsim`). Default
value is False.
issiso : bool, optional
Set to ``True`` if the system generating the data is single-input,
single-output. If passed as ``None`` (default), the input data
will be used to set the value.
ninputs, noutputs, nstates : int
Number of inputs, outputs, and states of the underlying system.
input_labels, output_labels, state_labels : array of str
Names for the input, output, and state variables.
success : bool, optional
If ``False``, result may not be valid (see
:func:`~control.input_output_response`). Defaults to ``True``.
message : str, optional
Informational message if ``success`` is ``False``.
sysname : str, optional
Name of the system that created the data.
params : dict, optional
If system is a nonlinear I/O system, set parameter values.
plot_inputs : bool, optional
Whether or not to plot the inputs by default (can be overridden in
the plot() method)
ntraces : int, optional
Number of independent traces represented in the input/output
response. If ntraces is 0 (default) then the data represents a
single trace with the trace index surpressed in the data.
title : str, optional
Set the title to use when plotting.
trace_labels : array of string, optional
Labels to use for traces (set to sysname it ntraces is 0)
trace_types : array of string, optional
Type of trace. Currently only 'step' is supported, which controls
the way in which the signal is plotted.
Notes
-----
The responses for individual elements of the time response can be
accessed using integers, slices, or lists of signal offsets or the
names of the appropriate signals::
sys = ct.rss(4, 2, 1)
resp = ct.initial_response(sys, X0=[1, 1, 1, 1])
plt.plot(resp.time, resp.outputs['y[0]'])
In the case of multi-trace data, the responses should be indexed using
the output signal name (or offset) and the input signal name (or
offset)::
sys = ct.rss(4, 2, 2, strictly_proper=True)
resp = ct.step_response(sys)
plt.plot(resp.time, resp.outputs[['y[0]', 'y[1]'], 'u[0]'].T)
For backward compatibility with earlier versions of python-control,
this class has an ``__iter__`` method that allows it to be assigned to
a tuple with a variable number of elements. This allows the following
patterns to work::
t, y = step_response(sys)
t, y, x = step_response(sys, return_x=True)
When using this (legacy) interface, the state vector is not affected
by the `squeeze` parameter.
For backward compatibility with earlier version of python-control, this
class has ``__getitem__`` and ``__len__`` methods that allow the return
value to be indexed:
response[0]: returns the time vector
response[1]: returns the output vector
response[2]: returns the state vector
When using this (legacy) interface, the state vector is not affected
by the `squeeze` parameter.
The default settings for ``return_x``, ``squeeze`` and ``transpose``
can be changed by calling the class instance and passing new values::
response(tranpose=True).input
See :meth:`TimeResponseData.__call__` for more information.
"""
def __init__(
self, time, outputs, states=None, inputs=None, issiso=None,
output_labels=None, state_labels=None, input_labels=None,
title=None, transpose=False, return_x=False, squeeze=None,
multi_trace=False, trace_labels=None, trace_types=None,
plot_inputs=True, sysname=None, params=None, success=True,
message=None
):
"""Create an input/output time response object.
Parameters
----------
time : 1D array
Time values of the output. Ignored if None.
outputs : ndarray
Output response of the system. This can either be a 1D array
indexed by time (for SISO systems or MISO systems with a specified
input), a 2D array indexed by output and time (for MIMO systems
with no input indexing, such as initial_response or forced
response) or trace and time (for SISO systems with multiple
traces), or a 3D array indexed by output, trace, and time (for
multi-trace input/output responses).
states : array, optional
Individual response of each state variable. This should be a 2D
array indexed by the state index and time (for single trace
systems) or a 3D array indexed by state, trace, and time.
inputs : array, optional
Inputs used to generate the output. This can either be a 1D
array indexed by time (for SISO systems or MISO/MIMO systems
with a specified input), a 2D array indexed either by input and
time (for a multi-input system) or trace and time (for a
single-input, multi-trace response), or a 3D array indexed by
input, trace, and time.
title : str, optonal
Title of the data set (used as figure title in plotting).
squeeze : bool, optional
By default, if a system is single-input, single-output (SISO)
then the inputs and outputs are returned as a 1D array (indexed
by time) and if a system is multi-input or multi-output, then
the inputs are returned as a 2D array (indexed by input and
time) and the outputs are returned as either a 2D array (indexed
by output and time) or a 3D array (indexed by output, trace, and
time). If squeeze=True, access to the output response will
remove single-dimensional entries from the shape of the inputs
and outputs even if the system is not SISO. If squeeze=False,
keep the input as a 2D or 3D array (indexed by the input (if
multi-input), trace (if single input) and time) and the output
as a 3D array (indexed by the output, trace, and time) even if
the system is SISO. The default value can be set using
config.defaults['control.squeeze_time_response'].
Other parameters
----------------
input_labels, output_labels, state_labels: array of str, optional
Optional labels for the inputs, outputs, and states, given as a
list of strings matching the appropriate signal dimension.
sysname : str, optional
Name of the system that created the data.
transpose : bool, optional
If True, transpose all input and output arrays (for backward
compatibility with MATLAB and :func:`scipy.signal.lsim`).
Default value is False.
return_x : bool, optional
If True, return the state vector when enumerating result by
assigning to a tuple (default = False).
plot_inputs : bool, optional
Whether or not to plot the inputs by default (can be overridden
in the plot() method)
multi_trace : bool, optional
If ``True``, then 2D input array represents multiple traces. For
a MIMO system, the ``input`` attribute should then be set to
indicate which trace is being specified. Default is ``False``.
success : bool, optional
If ``False``, result may not be valid (see
:func:`~control.input_output_response`).
message : str, optional
Informational message if ``success`` is ``False``.
"""
#
# Process and store the basic input/output elements
#
# Time vector
self.t = np.atleast_1d(time)
if self.t.ndim != 1:
raise ValueError("Time vector must be 1D array")
self.title = title
self.sysname = sysname
self.params = params
#
# Output vector (and number of traces)
#
self.y = np.array(outputs)
if self.y.ndim == 3:
multi_trace = True
self.noutputs = self.y.shape[0]
self.ntraces = self.y.shape[1]
elif multi_trace and self.y.ndim == 2:
self.noutputs = 1
self.ntraces = self.y.shape[0]
elif not multi_trace and self.y.ndim == 2:
self.noutputs = self.y.shape[0]
self.ntraces = 0
elif not multi_trace and self.y.ndim == 1:
self.noutputs = 1
self.ntraces = 0
# Reshape the data to be 2D for consistency
self.y = self.y.reshape(self.noutputs, -1)
else:
raise ValueError("Output vector is the wrong shape")
# Check and store labels, if present
self.output_labels = _process_labels(
output_labels, "output", self.noutputs)
# Make sure time dimension of output is the right length
if self.t.shape[-1] != self.y.shape[-1]:
raise ValueError("Output vector does not match time vector")
#
# State vector (optional)
#
# If present, the shape of the state vector should be consistent
# with the multi-trace nature of the data.
#
if states is None:
self.x = None
self.nstates = 0
else:
self.x = np.array(states)
self.nstates = self.x.shape[0]
# Make sure the shape is OK
if multi_trace and \
(self.x.ndim != 3 or self.x.shape[1] != self.ntraces) or \
not multi_trace and self.x.ndim != 2 :
raise ValueError("State vector is the wrong shape")
# Make sure time dimension of state is the right length
if self.t.shape[-1] != self.x.shape[-1]:
raise ValueError("State vector does not match time vector")
# Check and store labels, if present
self.state_labels = _process_labels(
state_labels, "state", self.nstates)
#
# Input vector (optional)
#
# If present, the shape and dimensions of the input vector should be
# consistent with the trace count computed above.
#
if inputs is None:
self.u = None
self.ninputs = 0
self.plot_inputs = False
else:
self.u = np.array(inputs)
self.plot_inputs = plot_inputs
# Make sure the shape is OK and figure out the nuumber of inputs
if multi_trace and self.u.ndim == 3 and \
self.u.shape[1] == self.ntraces:
self.ninputs = self.u.shape[0]
elif multi_trace and self.u.ndim == 2 and \
self.u.shape[0] == self.ntraces:
self.ninputs = 1
elif not multi_trace and self.u.ndim == 2 and \
self.ntraces == 0:
self.ninputs = self.u.shape[0]
elif not multi_trace and self.u.ndim == 1:
self.ninputs = 1
# Reshape the data to be 2D for consistency
self.u = self.u.reshape(self.ninputs, -1)
else:
raise ValueError("Input vector is the wrong shape")
# Make sure time dimension of output is the right length
if self.t.shape[-1] != self.u.shape[-1]:
raise ValueError("Input vector does not match time vector")
# Check and store labels, if present
self.input_labels = _process_labels(
input_labels, "input", self.ninputs)
# Check and store trace labels, if present
self.trace_labels = _process_labels(
trace_labels, "trace", self.ntraces)
self.trace_types = trace_types
# Figure out if the system is SISO
if issiso is None:
# Figure out based on the data
if self.ninputs == 1:
issiso = (self.noutputs == 1)
elif self.ninputs > 1:
issiso = False
else:
# Missing input data => can't resolve
raise ValueError("Can't determine if system is SISO")
elif issiso is True and (self.ninputs > 1 or self.noutputs > 1):
raise ValueError("Keyword `issiso` does not match data")
# Set the value to be used for future processing
self.issiso = issiso
# Keep track of whether to squeeze inputs, outputs, and states
if not (squeeze is True or squeeze is None or squeeze is False):
raise ValueError("Unknown squeeze value")
self.squeeze = squeeze
# Keep track of whether to transpose for MATLAB/scipy.signal
self.transpose = transpose
# Store legacy keyword values (only needed for legacy interface)
self.return_x = return_x
# Information on the whether the simulation result may be incorrect
self.success = success
self.message = message
def __call__(self, **kwargs):
"""Change value of processing keywords.
Calling the time response object will create a copy of the object and
change the values of the keywords used to control the ``outputs``,
``states``, and ``inputs`` properties.
Parameters
----------
squeeze : bool, optional
If squeeze=True, access to the output response will remove
single-dimensional entries from the shape of the inputs, outputs,
and states even if the system is not SISO. If squeeze=False, keep
the input as a 2D or 3D array (indexed by the input (if
multi-input), trace (if single input) and time) and the output and
states as a 3D array (indexed by the output/state, trace, and
time) even if the system is SISO.
transpose : bool, optional
If True, transpose all input and output arrays (for backward
compatibility with MATLAB and :func:`scipy.signal.lsim`).
Default value is False.
return_x : bool, optional
If True, return the state vector when enumerating result by
assigning to a tuple (default = False).
input_labels, output_labels, state_labels: array of str
Labels for the inputs, outputs, and states, given as a
list of strings matching the appropriate signal dimension.
"""
# Make a copy of the object
response = copy(self)
# Update any keywords that we were passed
response.transpose = kwargs.pop('transpose', self.transpose)
response.squeeze = kwargs.pop('squeeze', self.squeeze)
response.return_x = kwargs.pop('return_x', self.return_x)
# Check for new labels
input_labels = kwargs.pop('input_labels', None)
if input_labels is not None:
response.input_labels = _process_labels(
input_labels, "input", response.ninputs)
output_labels = kwargs.pop('output_labels', None)
if output_labels is not None:
response.output_labels = _process_labels(
output_labels, "output", response.noutputs)
state_labels = kwargs.pop('state_labels', None)
if state_labels is not None:
response.state_labels = _process_labels(
state_labels, "state", response.nstates)
# Make sure there were no extraneous keywords
if kwargs:
raise TypeError("unrecognized keywords: ", str(kwargs))
return response
@property
def time(self):
"""Time vector.
Time values of the input/output response(s).
:type: 1D array"""
return self.t
# Getter for output (implements squeeze processing)
@property
def outputs(self):
"""Time response output vector.
Output response of the system, indexed by either the output and time
(if only a single input is given) or the output, trace, and time
(for multiple traces). See :attr:`TimeResponseData.squeeze` for a
description of how this can be modified using the `squeeze` keyword.
:type: 1D, 2D, or 3D array
"""
# TODO: move to __init__ to avoid recomputing each time?
y = _process_time_response(
self.y, issiso=self.issiso,
transpose=self.transpose, squeeze=self.squeeze)
return NamedSignal(y, self.output_labels, self.input_labels)
# Getter for states (implements squeeze processing)
@property
def states(self):
"""Time response state vector.
Time evolution of the state vector, indexed indexed by either the
state and time (if only a single trace is given) or the state, trace,
and time (for multiple traces). See :attr:`TimeResponseData.squeeze`
for a description of how this can be modified using the `squeeze`
keyword.
:type: 2D or 3D array
"""
# TODO: move to __init__ to avoid recomputing each time?
if self.x is None:
return None
elif self.squeeze is True:
x = self.x.squeeze()
elif self.ninputs == 1 and self.noutputs == 1 and \
self.ntraces == 1 and self.x.ndim == 3 and \
self.squeeze is not False:
# Single-input, single-output system with single trace
x = self.x[:, 0, :]
else:
# Return the full set of data
x = self.x
# Transpose processing
if self.transpose:
x = np.transpose(x, np.roll(range(x.ndim), 1))
return NamedSignal(x, self.state_labels, self.input_labels)
# Getter for inputs (implements squeeze processing)
@property
def inputs(self):
"""Time response input vector.
Input(s) to the system, indexed by input (optiona), trace (optional),
and time. If a 1D vector is passed, the input corresponds to a
scalar-valued input. If a 2D vector is passed, then it can either
represent multiple single-input traces or a single multi-input trace.
The optional ``multi_trace`` keyword should be used to disambiguate
the two. If a 3D vector is passed, then it represents a multi-trace,
multi-input signal, indexed by input, trace, and time.
See :attr:`TimeResponseData.squeeze` for a description of how the
dimensions of the input vector can be modified using the `squeeze`
keyword.
:type: 1D or 2D array
"""
# TODO: move to __init__ to avoid recomputing each time?
if self.u is None:
return None
u = _process_time_response(
self.u, issiso=self.issiso,
transpose=self.transpose, squeeze=self.squeeze)
return NamedSignal(u, self.input_labels, self.input_labels)
# Getter for legacy state (implements non-standard squeeze processing)
# TODO: remove when no longer needed
@property
def _legacy_states(self):
"""Time response state vector (legacy version).
Time evolution of the state vector, indexed indexed by either the
state and time (if only a single trace is given) or the state,
trace, and time (for multiple traces).
The `legacy_states` property is not affected by the `squeeze` keyword
and hence it will always have these dimensions.
:type: 2D or 3D array
"""
if self.x is None:
return None
elif self.ninputs == 1 and self.noutputs == 1 and \
self.ntraces == 1 and self.x.ndim == 3:
# Single-input, single-output system with single trace
x = self.x[:, 0, :]
else:
# Return the full set of data
x = self.x
# Transpose processing
if self.transpose:
x = np.transpose(x, np.roll(range(x.ndim), 1))
return x
# Implement iter to allow assigning to a tuple
def __iter__(self):
if not self.return_x:
return iter((self.time, self.outputs))
return iter((self.time, self.outputs, self._legacy_states))
# Implement (thin) getitem to allow access via legacy indexing
def __getitem__(self, index):
# See if we were passed a slice
if isinstance(index, slice):
if (index.start is None or index.start == 0) and index.stop == 2:
return (self.time, self.outputs)
# Otherwise assume we were passed a single index
if index == 0:
return self.time
if index == 1:
return self.outputs
if index == 2:
return self._legacy_states
raise IndexError
# Implement (thin) len to emulate legacy testing interface
def __len__(self):
return 3 if self.return_x else 2
# Convert to pandas
def to_pandas(self):
"""Convert response data to pandas data frame.
Creates a pandas data frame using the input, output, and state
labels for the time response.
"""
if not pandas_check():
raise ImportError("pandas not installed")
import pandas
# Create a dict for setting up the data frame
data = {'time': self.time}
if self.ninputs > 0:
data.update(
{name: self.u[i] for i, name in enumerate(self.input_labels)})
if self.noutputs > 0:
data.update(
{name: self.y[i] for i, name in enumerate(self.output_labels)})
if self.nstates > 0:
data.update(
{name: self.x[i] for i, name in enumerate(self.state_labels)})
return pandas.DataFrame(data)
# Plot data
def plot(self, *args, **kwargs):
"""Plot the time response data objects.
This method calls :func:`time_response_plot`, passing all arguments
and keywords.
"""
return time_response_plot(self, *args, **kwargs)
#
# Time response data list class
#
# This class is a subclass of list that adds a plot() method, enabling
# direct plotting from routines returning a list of TimeResponseData
# objects.
#
class TimeResponseList(list):
"""This class consist of a list of :class:`TimeResponseData` objects.
It is a subclass of the Python `list` class, with a `plot` method that
plots the individual :class:`TimeResponseData` objects.
"""
def plot(self, *args, **kwargs):
from .ctrlplot import ControlPlot
lines = None
label = kwargs.pop('label', [None] * len(self))
for i, response in enumerate(self):
cplt = TimeResponseData.plot(
response, *args, label=label[i], **kwargs)
if lines is None:
lines = cplt.lines
else:
# Append the lines in the new plot to previous lines
for row in range(cplt.lines.shape[0]):
for col in range(cplt.lines.shape[1]):
lines[row, col] += cplt.lines[row, col]
return ControlPlot(lines, cplt.axes, cplt.figure)
# Process signal labels
def _process_labels(labels, signal, length):
"""Process time response signal labels.
Parameters
----------
labels : list of str or dict
Description of the labels for the signal. This can be a list of
strings or a dict giving the index of each signal (used in iosys).
signal : str
Name of the signal being processed (for error messages).
length : int
Number of labels required.
Returns
-------
labels : list of str
List of labels.
"""
if labels is None or len(labels) == 0:
return None
# See if we got passed a dictionary (from iosys)
if isinstance(labels, dict):
# Form inverse dictionary
ivd = {v: k for k, v in labels.items()}
try:
# Turn into a list
labels = [ivd[n] for n in range(len(labels))]
except KeyError:
raise ValueError("Name dictionary for %s is incomplete" % signal)
# Convert labels to a list
if isinstance(labels, str):
labels = [labels]
else:
labels = list(labels)
# Make sure the signal list is the right length and type
if len(labels) != length:
raise ValueError("List of %s labels is the wrong length" % signal)
elif not all([isinstance(label, str) for label in labels]):
raise ValueError("List of %s labels must all be strings" % signal)
return labels
# Helper function for checking array-like parameters
def _check_convert_array(in_obj, legal_shapes, err_msg_start, squeeze=False,
transpose=False):
"""Helper function for checking array_like parameters.
* Check type and shape of ``in_obj``.
* Convert ``in_obj`` to an array if necessary.
* Change shape of ``in_obj`` according to parameter ``squeeze``.
* If ``in_obj`` is a scalar (number) it is converted to an array with
a legal shape, that is filled with the scalar value.
The function raises an exception when it detects an error.
Parameters
----------
in_obj : array like object
The array or matrix which is checked.
legal_shapes : list of tuple
A list of shapes that in_obj can legally have.
The special value "any" means that there can be any
number of elements in a certain dimension.
* ``(2, 3)`` describes an array with 2 rows and 3 columns
* ``(2, "any")`` describes an array with 2 rows and any number of
columns
err_msg_start : str
String that is prepended to the error messages, when this function
raises an exception. It should be used to identify the argument which
is currently checked.
squeeze : bool
If True, all dimensions with only one element are removed from the
array. If False the array's shape is unmodified.
For example:
``array([[1,2,3]])`` is converted to ``array([1, 2, 3])``
transpose : bool, optional
If True, assume that 2D input arrays are transposed from the standard
format. Used to convert MATLAB-style inputs to our format.
Returns
-------
out_array : array
The checked and converted contents of ``in_obj``.
"""
# convert nearly everything to an array.
out_array = np.asarray(in_obj)
if (transpose):
out_array = np.transpose(out_array)
# Test element data type, elements must be numbers
legal_kinds = set(("i", "f", "c")) # integer, float, complex
if out_array.dtype.kind not in legal_kinds:
err_msg = "Wrong element data type: '{d}'. Array elements " \
"must be numbers.".format(d=str(out_array.dtype))
raise TypeError(err_msg_start + err_msg)
# If array is zero dimensional (in_obj is scalar):
# create array with legal shape filled with the original value.
if out_array.ndim == 0:
for s_legal in legal_shapes:
# search for shape that does not contain the special symbol any.
if "any" in s_legal:
continue
the_val = out_array[()]
out_array = np.empty(s_legal, 'd')
out_array.fill(the_val)
break
# Test shape
def shape_matches(s_legal, s_actual):
"""Test if two shape tuples match"""
# Array must have required number of dimensions
if len(s_legal) != len(s_actual):
return False
# All dimensions must contain required number of elements. Joker: "all"
for n_legal, n_actual in zip(s_legal, s_actual):
if n_legal == "any":
continue
if n_legal != n_actual:
return False
return True
# Iterate over legal shapes, and see if any matches out_array's shape.
for s_legal in legal_shapes:
if shape_matches(s_legal, out_array.shape):
break
else:
legal_shape_str = " or ".join([str(s) for s in legal_shapes])
err_msg = "Wrong shape (rows, columns): {a}. Expected: {e}." \
.format(e=legal_shape_str, a=str(out_array.shape))
raise ValueError(err_msg_start + err_msg)
# Convert shape
if squeeze:
out_array = np.squeeze(out_array)
# We don't want zero dimensional arrays
if out_array.shape == tuple():
out_array = out_array.reshape((1,))
return out_array
# Forced response of a linear system
def forced_response(sysdata, T=None, U=0., X0=0., transpose=False, params=None,
interpolate=False, return_x=None, squeeze=None):
"""Compute the output of a linear system given the input.
As a convenience for parameters `U`, `X0`:
Numbers (scalars) are converted to constant arrays with the correct shape.
The correct shape is inferred from arguments `sys` and `T`.
For information on the **shape** of parameters `U`, `T`, `X0` and
return values `T`, `yout`, `xout`, see :ref:`time-series-convention`.
Parameters
----------
sysdata : I/O system or list of I/O systems
I/O system(s) for which forced response is computed.
T : array_like, optional for discrete LTI `sys`
Time steps at which the input is defined; values must be evenly spaced.
If None, `U` must be given and `len(U)` time steps of sys.dt are
simulated. If sys.dt is None or True (undetermined time step), a time
step of 1.0 is assumed.
U : array_like or float, optional
Input array giving input at each time `T`.
If `U` is None or 0, `T` must be given, even for discrete
time systems. In this case, for continuous time systems, a direct
calculation of the matrix exponential is used, which is faster than the
general interpolating algorithm used otherwise.
X0 : array_like or float, default=0.
Initial condition.
params : dict, optional
If system is a nonlinear I/O system, set parameter values.
transpose : bool, default=False
If True, transpose all input and output arrays (for backward
compatibility with MATLAB and :func:`scipy.signal.lsim`).
interpolate : bool, default=False
If True and system is a discrete time system, the input will
be interpolated between the given time steps and the output
will be given at system sampling rate. Otherwise, only return
the output at the times given in `T`. No effect on continuous
time simulations.
return_x : bool, default=None
Used if the time response data is assigned to a tuple:
* If False, return only the time and output vectors.
* If True, also return the the state vector.