Miscellaneous¶
Tools for 2nd order systems¶
A collection of modules and methods that are used throughout the whole package. Methods specialized for second order dynamic systems, such as the ones used for high-class accelerometers.
This module contains the following functions:
- sos_FreqResp: Calculation of the system frequency response
- sos_phys2filter: Calculation of continuous filter coefficients from physical parameters
- sos_absphase: Propagation of uncertainty from physical parameters to real and imaginary part of system’s transfer function using GUM S2 Monte Carlo
- sos_realimag: Propagation of uncertainty from physical parameters to real and imaginary part of system’s transfer function using GUM S2 Monte Carlo
-
PyDynamic.misc.SecondOrderSystem.
sos_FreqResp
(S, d, f0, freqs)[source]¶ Calculation of the system frequency response
The frequency response is calculated from the continuous physical model of a second order system given by
\(H(f) = \frac{4S\pi^2f_0^2}{(2\pi f_0)^2 + 2jd(2\pi f_0)f - f^2}\)
If the provided system parameters are vectors then \(H(f)\) is calculated for each set of parameters. This is helpful for Monte Carlo simulations by using draws from the model parameters
Parameters: - S (float or ndarray shape (K,)) – static gain
- d (float or ndarray shape (K,)) – damping parameter
- f0 (float or ndarray shape (K,)) – resonance frequency
- freqs (ndarray shape (N,)) – frequencies at which to calculate the freq response
Returns: H – complex frequency response values
Return type: ndarray shape (N,) or ndarray shape (N,K)
-
PyDynamic.misc.SecondOrderSystem.
sos_phys2filter
(S, d, f0)[source]¶ Calculation of continuous filter coefficients from physical parameters.
If the provided system parameters are vectors then the filter coefficients are calculated for each set of parameters. This is helpful for Monte Carlo simulations by using draws from the model parameters
Parameters: - S (float) – static gain
- d (float) – damping parameter
- f0 (float) – resonance frequency
Returns: b,a – analogue filter coefficients
Return type: ndarray
-
PyDynamic.misc.SecondOrderSystem.
sos_absphase
(S, d, f0, uS, ud, uf0, f, runs=10000)[source]¶ Propagation of uncertainty from physical parameters to real and imaginary part of system’s transfer function using GUM S2 Monte Carlo.
Parameters: - S (float) – static gain
- d (float) – damping
- f0 (float) – resonance frequency
- uS (float) – uncertainty associated with static gain
- ud (float) – uncertainty associated with damping
- uf0 (float) – uncertainty associated with resonance frequency
- f (ndarray, shape (N,)) – frequency values at which to calculate amplitue and phase
Returns: - Hmean (ndarray, shape (N,)) – best estimate of complex frequency response values
- Hcov (ndarray, shape (2N,2N)) – covariance matrix [ [u(abs,abs), u(abs,phase)], [u(phase,abs), u(phase,phase)] ]
-
PyDynamic.misc.SecondOrderSystem.
sos_realimag
(S, d, f0, uS, ud, uf0, f, runs=10000)[source]¶ Propagation of uncertainty from physical parameters to real and imaginary part of system’s transfer function using GUM S2 Monte Carlo.
Parameters: - S (float) – static gain
- d (float) – damping
- f0 (float) – resonance frequency
- uS (float) – uncertainty associated with static gain
- ud (float) – uncertainty associated with damping
- uf0 (float) – uncertainty associated with resonance frequency
- f (ndarray, shape (N,)) – frequency values at which to calculate real and imaginary part
Returns: - Hmean (ndarray, shape (N,)) – best estimate of complex frequency response values
- Hcov (ndarray, shape (2N,2N)) – covariance matrix [ [u(real,real), u(real,imag)], [u(imag,real), u(imag,imag)] ]
Tools for digital filters¶
A collection of methods which are related to filter design.
This module contains the following functions:
- db: Calculation of decibel values \(20\log_{10}(x)\) for a vector of values
- ua: Shortcut for calculation of unwrapped angle of complex values
- grpdelay: Calculation of the group delay of a digital filter
- mapinside: Maps the roots of polynomial with coefficients \(a\) to the unit circle
- kaiser_lowpass: Design of a FIR lowpass filter using the window technique with a Kaiser window.
- isstable: Determine whether a given IIR filter is stable
- savitzky_golay: Smooth (and optionally differentiate) data with a Savitzky-Golay filter
-
PyDynamic.misc.filterstuff.
grpdelay
(b, a, Fs, nfft=512)[source]¶ Calculation of the group delay of a digital filter
Parameters: - b (ndarray) – IIR filter numerator coefficients
- a (ndarray) – IIR filter denominator coefficients
- Fs (float) – sampling frequency of the filter
- nfft (int) – number of FFT bins
Returns: - group_delay (np.ndarray) – group delay values
- frequencies (ndarray) – frequencies at which the group delay is calculated
References
- Smith, online book [Smith]
-
PyDynamic.misc.filterstuff.
mapinside
(a)[source]¶ Maps the roots of polynomial to the unit circle.
Maps the roots of polynomial with coefficients \(a\) to the unit circle.
Parameters: a (ndarray) – polynomial coefficients Returns: a – polynomial coefficients with all roots inside or on the unit circle Return type: ndarray
-
PyDynamic.misc.filterstuff.
kaiser_lowpass
(L, fcut, Fs, beta=8.0)[source]¶ Design of a FIR lowpass filter using the window technique with a Kaiser window.
This method uses a Kaiser window. Filters of that type are often used as FIR low-pass filters due to their linear phase.
Parameters: - L (int) – filter order (window length)
- fcut (float) – desired cut-off frequency
- Fs (float) – sampling frequency
- beta (float) – scaling parameter for the Kaiser window
Returns: - blow (ndarray) – FIR filter coefficients
- shift (int) – delay of the filter (in samples)
-
PyDynamic.misc.filterstuff.
isstable
(b, a, ftype='digital')[source]¶ Determine whether IIR filter (b,a) is stable
Parameters: - b (ndarray) – filter numerator coefficients
- a (ndarray) – filter denominator coefficients
- ftype (string) – type of filter (digital or analog)
Returns: stable – whether filter is stable or not
Return type: bool
-
PyDynamic.misc.filterstuff.
savitzky_golay
(y, window_size, order, deriv=0, rate=1)[source]¶ Smooth (and optionally differentiate) data with a Savitzky-Golay filter.
The Savitzky-Golay filter removes high frequency noise from data. It has the advantage of preserving the original shape and features of the signal better than other types of filtering approaches, such as moving averages techniques.
Source obtained from scipy cookbook (online), downloaded 2013-09-13
Parameters: - y (ndarray, shape (N,)) – the values of the time history of the signal
- window_size (int) – the length of the window. Must be an odd integer number
- order (int) – the order of the polynomial used in the filtering. Must be less then window_size - 1.
- deriv (int) – the order of the derivative to compute (default = 0 means only smoothing)
Returns: ys – the smoothed signal (or it’s n-th derivative).
Return type: ndarray, shape (N,)
Notes
The Savitzky-Golay is a type of low-pass filter, particularly suited for smoothing noisy data. The main idea behind this approach is to make for each point a least-square fit with a polynomial of high order over a odd-sized window centered at the point.
References
- Savitzky et al. [Savitzky]
- Numerical Recipes [NumRec]
Test signals¶
Collection of test signals which can be used to simulate dynamic measurements and test methods.
This module contains the following functions:
- shocklikeGaussian: signal that resembles a shock excitation as a Gaussian followed by a smaller Gaussian of opposite sign
- GaussianPulse: Generates a Gaussian pulse at \(t_0\) with height \(m_0\) and std \(sigma\)
- rect: Rectangular signal of given height and width \(t_1 - t_0\)
- squarepulse: Generates a series of rect functions to represent a square pulse signal
-
PyDynamic.misc.testsignals.
shocklikeGaussian
(time, t0, m0, sigma, noise=0.0)[source]¶ Generates a signal that resembles a shock excitation as a Gaussian followed by a smaller Gaussian of opposite sign.
Parameters: - time (np.ndarray of shape (N,)) – time instants (equidistant)
- t0 (float) – time instant of signal maximum
- m0 (float) – signal maximum
- sigma (float) – std of main pulse
- noise (float, optional) – std of simulated signal noise
Returns: x – signal amplitudes at time instants
Return type: np.ndarray of shape (N,)
-
PyDynamic.misc.testsignals.
GaussianPulse
(time, t0, m0, sigma, noise=0.0)[source]¶ Generates a Gaussian pulse at t0 with height m0 and std sigma
Parameters: - time (np.ndarray of shape (N,)) – time instants (equidistant)
- t0 (float) – time instant of signal maximum
- m0 (float) – signal maximum
- sigma (float) – std of pulse
- noise (float, optional) – std of simulated signal noise
Returns: x – signal amplitudes at time instants
Return type: np.ndarray of shape (N,)
-
PyDynamic.misc.testsignals.
rect
(time, t0, t1, height=1, noise=0.0)[source]¶ Rectangular signal of given height and width t1-t0
Parameters: - time (np.ndarray of shape (N,)) – time instants (equidistant)
- t0 (float) – time instant of rect lhs
- t1 (float) – time instant of rect rhs
- height (float) – signal maximum
- noise (float, optional) – std of simulated signal noise
Returns: x – signal amplitudes at time instants
Return type: np.ndarray of shape (N,)
-
PyDynamic.misc.testsignals.
squarepulse
(time, height, numpulse=4, noise=0.0)[source]¶ Generates a series of rect functions to represent a square pulse signal
Parameters: - time (np.ndarray of shape (N,)) – time instants
- height (float) – height of the rectangular pulses
- numpulse (int) – number of pulses
- noise (float, optional) – std of simulated signal noise
Returns: x – signal amplitude at time instants
Return type: np.ndarray of shape (N,)
Miscellaneous useful helper functions¶
Collection of miscellaneous helper functions.
This module contains the following functions:
- print_vec: Print vector (1D array) to the console or return as formatted string
- print_mat: Print matrix (2D array) to the console or return as formatted string
- make_semiposdef: Make quadratic matrix positive semi-definite
- FreqResp2RealImag: Calculate real and imaginary parts from frequency response
- make_equidistant: Interpolate non-equidistant time series to equidistant
-
PyDynamic.misc.tools.
print_mat
(matrix, prec=5, vertical=False, retS=False)[source]¶ Print matrix (2D array) to the console or return as formatted string
Parameters: - matrix ((M,N) array_like) –
- prec (int) – the precision of the output
- vertical (bool) – print out vertical or not
- retS (bool) – print or return string
Returns: s – if retS is True
Return type: str
-
PyDynamic.misc.tools.
print_vec
(vector, prec=5, retS=False, vertical=False)[source]¶ Print vector (1D array) to the console or return as formatted string
Parameters: - vector ((M,) array_like) –
- prec (int) – the precision of the output
- vertical (bool) – print out vertical or not
- retS (bool) – print or return string
Returns: s – if retS is True
Return type: str
-
PyDynamic.misc.tools.
make_semiposdef
(matrix, maxiter=10, tol=1e-12, verbose=False)[source]¶ Make quadratic matrix positive semi-definite by increasing its eigenvalues
Parameters: - matrix ((N,N) array_like) –
- maxiter (int) – the maximum number of iterations for increasing the eigenvalues
- tol (float) – tolerance for deciding if pos. semi-def.
- verbose (bool) – If True print some more detail about input parameters.
Returns: quadratic positive semi-definite matrix
Return type: (N,N) array_like
-
PyDynamic.misc.tools.
FreqResp2RealImag
(Abs, Phase, Unc, MCruns=10000.0)[source]¶ Calculate real and imaginary parts from frequency response
Calculate real and imaginary parts from amplitude and phase with associated uncertainties.
Parameters: - Abs ((N,) array_like) – amplitude values
- Phase ((N,) array_like) – phase values in rad
- Unc ((2N, 2N) or (2N,) array_like) – uncertainties
- MCruns (bool) – Iterations for Monte Carlo simulation
Returns: - Re, Im ((N,) array_like) – real and imaginary parts (best estimate)
- URI ((2N, 2N) array_like) – uncertainties assoc. with Re and Im
-
PyDynamic.misc.tools.
make_equidistant
(t, y, uy, dt=0.05, kind='linear')[source]¶ Interpolate non-equidistant time series to equidistant
Interpolate measurement values and propagate uncertainties accordingly.
Parameters: - t ((N,) array_like) – timestamps in ascending order
- y ((N,) array_like) – corresponding measurement values
- uy ((N,) array_like) – corresponding measurement values’ uncertainties
- dt (float, optional) – desired interval length in seconds
- kind (str, optional) – Specifies the kind of interpolation for the measurement values as a string (‘previous’, ‘next’, ‘nearest’ or ‘linear’).
Returns: - t_new ((N,) array_like) – timestamps
- y_new ((N,) array_like) – measurement values
- uy_new ((N,) array_like) – measurement values’ uncertainties
References
- White [White2017]