Design of deconvolution filters¶
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. This module implements methods for the design of such filters given an array of frequency response values with associated uncertainties for the measurement system.
-
PyDynamic.deconvolution.fit_filter.
LSFIR
(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: - H (np.ndarray) – frequency response values
- N (int) – FIR filter order
- tau (float) – delay of filter
- f (np.ndarray) – frequencies
- Fs (float) – sampling frequency of digital filter
- Wt (np.ndarray, optional) – vector of weights
Returns: bFIR – filter coefficients
Return type: np.ndarray
References
- Elster and Link [Elster2008]
-
PyDynamic.deconvolution.fit_filter.
LSFIR_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) – frequency response values
- UH (np.ndarray) – uncertainties associated with the real and imaginary part
- N (int) – FIR filter order
- tau (float) – delay of filter
- f (np.ndarray) – frequencies
- Fs (float) – sampling frequency of digital filter
- wt (np.ndarray, 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) – filter coefficients of shape
- Ub (np.ndarray) – uncertainties associated with b
References
- Elster and Link [Elster2008]
-
PyDynamic.deconvolution.fit_filter.
LSIIR
(Hvals, Nb, Na, f, Fs, tau=1, 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) – frequency response values.
- Nb (int) – order of IIR numerator polynomial.
- Na (int) – order of IIR denominator polynomial.
- f (np.ndarray) – 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.
Returns: b,a – IIR filter coefficients, int tau – time delay (in samples)
Return type: np.ndarray
References
- Eichstädt, Elster, Esward, Hessling [Eichst2010]
-
PyDynamic.deconvolution.fit_filter.
LSIIR_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) – frequency response values.
- UH (np.ndarray) – 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) – 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) – uncertainties associated with [a[1:],b]
References
- Eichstädt, Elster, Esward and Hessling [Eichst2010]
-
PyDynamic.deconvolution.fit_filter.
LSFIR_uncMC
(H, UH, N, tau, f, Fs, wt=None, 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) – frequency response values
- UH (np.ndarray) – 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) – frequencies corresponding to H
- Fs (float) – sampling frequency of digital filter
- wt (np.ndarray, optional) – vector of weights
- verbose (bool, optional) – whether to print statements to the command line
Returns: - b (np.ndarray) – filter coefficients of shape
- Ub (np.ndarray) – uncertainties associated with b
References
- Elster and Link [Elster2008]