Model estimation¶
The estimation of the measurand in the analysis of dynamic measurements typically corresponds to a deconvolution problem. Therefore, a digital filter can be designed whose input is the measured system output signal and whose output is an estimate of the measurand. The package Model estimation implements methods for the design of such filters given an array of frequency response values or the reciprocal of frequency response values with associated uncertainties for the measurement system.
The package Model estimation also contains a function for the identification of transfer function models.
The package consists of the following modules:
PyDynamic.model_estimation.fit_filter
: least-squares fit to a given complex frequency response or its reciprocalPyDynamic.model_estimation.fit_transfer
: identification of transfer function models
Fitting filters to frequency response or reciprocal¶
The module PyDynamic.model_estimation.fit_filter
contains several functions to
carry out a least-squares fit to a given complex frequency response and the design of
digital deconvolution filters by least-squares fitting to the reciprocal of a given
frequency response each with associated uncertainties.
This module contains the following functions:
LSIIR()
: Least-squares IIR filter fit to a given frequency responseLSFIR()
: Least-squares fit of a digital FIR filter to a given frequency responseinvLSFIR()
: Least-squares fit of a digital FIR filter to the reciprocal of a given frequency response.invLSFIR_unc()
: Design of FIR filter as fit to reciprocal of frequency response values with uncertaintyinvLSFIR_uncMC()
: Design of FIR filter as fit to reciprocal of frequency response values with uncertainty via Monte CarloinvLSIIR()
: Design of a stable IIR filter as fit to reciprocal of frequency response valuesinvLSIIR_unc()
: Design of a stable IIR filter as fit to reciprocal of frequency response values with uncertainty
- PyDynamic.model_estimation.fit_filter.LSFIR(H, N, tau, f, Fs, Wt=None)[source]¶
Least-squares fit of a digital FIR filter to a given frequency response.
- Parameters
H ((complex) frequency response values of shape (M,)) –
N (FIR filter order) –
tau (delay of filter) –
f (frequencies of shape (M,)) –
Fs (sampling frequency of digital filter) –
Wt ((optional) vector of weights of shape (M,) or shape (M,M)) –
- Returns
- Return type
filter coefficients bFIR (ndarray) of shape (N+1,)
- PyDynamic.model_estimation.fit_filter.LSIIR(Hvals, Nb, Na, f, Fs, tau=0, justFit=False)[source]¶
Least-squares IIR filter fit to a given frequency response.
This method uses Gauss-Newton non-linear optimization and pole mapping for filter stabilization
- Parameters
Hvals (numpy array of (complex) frequency response values of shape (M,)) –
Nb (integer numerator polynomial order) –
Na (integer denominator polynomial order) –
f (numpy array of frequencies at which Hvals is given of shape) –
(M –
) –
Fs (sampling frequency) –
tau (integer initial estimate of time delay) –
justFit (boolean, when true then no stabilization is carried out) –
- Returns
b,a (IIR filter coefficients as numpy arrays)
tau (filter time delay in samples)
References
Eichstädt et al. 2010 [Eichst2010]
Vuerinckx et al. 1996 [Vuer1996]
- PyDynamic.model_estimation.fit_filter.invLSFIR(H, N, tau, f, Fs, Wt=None)[source]¶
Least-squares fit of a digital FIR filter to the reciprocal of a given frequency response.
- Parameters
- Returns
bFIR – filter coefficients
- Return type
np.ndarray of shape (N,)
References
Elster and Link [Elster2008]
- PyDynamic.model_estimation.fit_filter.invLSFIR_unc(H, UH, N, tau, f, Fs, wt=None, verbose=True, trunc_svd_tol=None)[source]¶
Design of FIR filter as fit to reciprocal of frequency response values with uncertainty
Least-squares fit of a digital FIR filter to the reciprocal of a frequency response for which associated uncertainties are given for its real and imaginary part. Uncertainties are propagated using a truncated svd and linear matrix propagation.
- Parameters
H (np.ndarray of shape (M,)) – frequency response values
UH (np.ndarray of shape (2M,2M)) – uncertainties associated with the real and imaginary part
N (int) – FIR filter order
tau (float) – delay of filter
f (np.ndarray of shape (M,)) – frequencies
Fs (float) – sampling frequency of digital filter
wt (np.ndarray of shape (2M,) - optional) – array of weights for a weighted least-squares method
verbose (bool, optional) – whether to print statements to the command line
trunc_svd_tol (float) – lower bound for singular values to be considered for pseudo-inverse
- Returns
b (np.ndarray of shape (N+1,)) – filter coefficients of shape
Ub (np.ndarray of shape (N+1,N+1)) – uncertainties associated with b
References
Elster and Link [Elster2008]
- PyDynamic.model_estimation.fit_filter.invLSFIR_uncMC(H, UH, N, tau, f, Fs, verbose=True)[source]¶
Design of FIR filter as fit to reciprocal of frequency response values with uncertainty
Least-squares fit of a FIR filter to the reciprocal of a frequency response for which associated uncertainties are given for its real and imaginary parts. Uncertainties are propagated using a Monte Carlo method. This method may help in cases where the weighting matrix or the Jacobian are ill-conditioned, resulting in false uncertainties associated with the filter coefficients.
- Parameters
H (np.ndarray of shape (M,) and dtype complex) – frequency response values
UH (np.ndarray of shape (2M,2M)) – uncertainties associated with the real and imaginary part of H
N (int) – FIR filter order
tau (int) – time delay of filter in samples
f (np.ndarray of shape (M,)) – frequencies corresponding to H
Fs (float) – sampling frequency of digital filter
verbose (bool, optional) – whether to print statements to the command line
- Returns
b (np.ndarray of shape (N+1,)) – filter coefficients of shape
Ub (np.ndarray of shape (N+1, N+1)) – uncertainties associated with b
References
Elster and Link [Elster2008]
- PyDynamic.model_estimation.fit_filter.invLSIIR(Hvals, Nb, Na, f, Fs, tau, justFit=False, verbose=True)[source]¶
Design of a stable IIR filter as fit to reciprocal of frequency response values
Least-squares fit of a digital IIR filter to the reciprocal of a given set of frequency response values using the equation-error method and stabilization by pole mapping and introduction of a time delay.
- Parameters
Hvals (np.ndarray of shape (M,) and dtype complex) – frequency response values.
Nb (int) – order of IIR numerator polynomial.
Na (int) – order of IIR denominator polynomial.
f (np.ndarray of shape (M,)) – frequencies corresponding to Hvals
Fs (float) – sampling frequency for digital IIR filter.
tau (float) – initial estimate of time delay for filter stabilization.
justFit (bool) – if True then no stabilization is carried out.
verbose (bool) – If True print some more detail about input parameters.
- Returns
b (np.ndarray) – The IIR filter numerator coefficient vector in a 1-D sequence
a (np.ndarray) – The IIR filter denominator coefficient vector in a 1-D sequence
tau (int) – time delay (in samples)
References
Eichstädt, Elster, Esward, Hessling [Eichst2010]
- PyDynamic.model_estimation.fit_filter.invLSIIR_unc(H, UH, Nb, Na, f, Fs, tau=0)[source]¶
Design of stabel IIR filter as fit to reciprocal of given frequency response with uncertainty
Least-squares fit of a digital IIR filter to the reciprocal of a given set of frequency response values with given associated uncertainty. Propagation of uncertainties is carried out using the Monte Carlo method.
- Parameters
H (np.ndarray of shape (M,) and dtype complex) – frequency response values.
UH (np.ndarray of shape (2M,2M)) – uncertainties associated with real and imaginary part of H
Nb (int) – order of IIR numerator polynomial.
Na (int) – order of IIR denominator polynomial.
f (np.ndarray of shape (M,)) – frequencies corresponding to H
Fs (float) – sampling frequency for digital IIR filter.
tau (float) – initial estimate of time delay for filter stabilization.
- Returns
b,a (np.ndarray) – IIR filter coefficients
tau (int) – time delay (in samples)
Uba (np.ndarray of shape (Nb+Na+1, Nb+Na+1)) – uncertainties associated with [a[1:],b]
References
Eichstädt, Elster, Esward and Hessling [Eichst2010]
Identification of transfer function models¶
This module contains a function for the identification of transfer function models:
fit_som()
: Fit second-order model to complex-valued frequency response
- PyDynamic.model_estimation.fit_transfer.fit_som(f: numpy.ndarray, H: numpy.ndarray, UH: Optional[numpy.ndarray] = None, weighting: Optional[numpy.ndarray] = None, MCruns: Optional[int] = 10000, scaling: Optional[float] = 0.001, verbose: Optional[bool] = False)[source]¶
Fit second-order model to complex-valued frequency response
Fit second-order model (spring-damper model) with parameters \(S_0, delta\) and \(f_0\) to complex-valued frequency response with uncertainty associated with real and imaginary parts.
For a transformation of an uncertainty associated with amplitude and phase to one associated with real and imaginary parts, see
PyDynamic.uncertainty.propagate_DFT.AmpPhase2DFT
.- Parameters
f ((M,) np.ndarray) – vector of frequencies
H ((2M,) np.ndarray) – real and imaginary parts of measured frequency response values at frequencies f
UH ((2M,) or (2M,2M) np.ndarray, optional) – uncertainties associated with real and imaginary parts When UH is one-dimensional, it is assumed to contain standard uncertainties; otherwise it is taken as covariance matrix. When UH is not specified no uncertainties assoc. with the fit are calculated, which is the default behaviour.
weighting (str or (2M,) np.ndarray, optional) – Type of weighting associated with frequency responses, can be (‘diag’, ‘cov’) if UH is given, or Numpy array of weights, defaults to None, which means all values are considered equally important
MCruns (int, optional) – Number of Monte Carlo trials for propagation of uncertainties, defaults to 10000. When MCruns is set to ‘None’, matrix multiplication is used for the propagation of uncertainties. However, in some cases this can cause trouble.
scaling (float, optional) – scaling of least-squares design matrix for improved fit quality, defaults to 1e-3
verbose (bool, optional) – if True a progressbar will be printed to console during the Monte Carlo simulations, if False nothing will be printed out, defaults to False
- Returns
p (np.ndarray) – vector of estimated model parameters [S0, delta, f0]
Up (np.ndarray) – covariance associated with parameter estimate