Libraries

pyfda contains the following libraries:

pyfda_lib

Library with various general functions and variables needed by the pyfda routines

pyfda.libs.pyfda_lib.H_mag(num, den, z, H_max, H_min=None, log=False, div_by_0='ignore')[source]

Calculate |H(z)| at the complex frequency(ies) z (scalar or array-like). The function H(z) is given in polynomial form with numerator and denominator. When log == True, \(20 \log_{10} (|H(z)|)\) is returned.

The result is clipped at H_min, H_max; clipping can be disabled by passing None as the argument.

Parameters:
  • num (float or array-like) – The numerator polynome of H(z).

  • den (float or array-like) – The denominator polynome of H(z).

  • z (float or array-like) – The complex frequency(ies) where H(z) is to be evaluated

  • H_max (float) – The maximum value to which the result is clipped

  • H_min (float, optional) – The minimum value to which the result is clipped (default: None)

  • log (boolean, optional) – When true, return 20 * log10 (|H(z)|). The clipping limits have to be given as dB in this case.

  • div_by_0 (string, optional) – What to do when division by zero occurs during calculation (default: ‘ignore’). As the denomintor of H(z) becomes 0 at each pole, warnings are suppressed by default. This parameter is passed to numpy.seterr(), hence other valid options are ‘warn’, ‘raise’ and ‘print’.

Returns:

H_mag – The magnitude |H(z)| for each value of z.

Return type:

float or ndarray

pyfda.libs.pyfda_lib.ceil_even(x) int[source]

Return the smallest even integer not less than x. x can be integer or float.

pyfda.libs.pyfda_lib.ceil_odd(x) int[source]

Return the smallest odd integer not less than x. x can be integer or float.

pyfda.libs.pyfda_lib.clean_ascii(arg)[source]

Remove non-ASCII-characters (outside range 0 … x7F) from arg when it is a str. Otherwise, return arg unchanged.

Parameters:

arg (str) – This is a unicode string under Python 3

Returns:

arg – Input string, cleaned from non-ASCII characters when arg is a string

or

Unchanged parameter arg when not a string

Return type:

str

pyfda.libs.pyfda_lib.cmp_version(mod: str, version: str) int[source]

Compare version number of installed module mod against value version (str) and return 1, 0 or -1 if the installed version is greater, equal or less than the number in version. If mod is not installed, return -2.

Parameters:
  • mod (str) – name of the module to be compared

  • version (str) – version number in the form e.g. “0.1.6”

Returns:

result

one of the following error codes:

-3:

version number could not be determined

-2:

module is not installed

-1:

version of installed module is lower than the specified version

0:

version of installed module is equal to specied version

1:

version of installed module is higher than specified version

Return type:

int

pyfda.libs.pyfda_lib.cmplx_sort(p)[source]

sort roots based on magnitude.

pyfda.libs.pyfda_lib.cround(x, n_dig=0)[source]

Round complex number to n_dig digits. If n_dig == 0, don’t round at all, just convert complex numbers with an imaginary part very close to zero to real.

pyfda.libs.pyfda_lib.dB(lin: float, power: bool = False) float[source]

Calculate dB from linear value. If power = True, calculate 10 log …, else calculate 20 log …

pyfda.libs.pyfda_lib.expand_lim(ax, eps_x: float, eps_y: float = None) None[source]

Expand the xlim and ylim-values of passed axis by eps

Parameters:
  • ax (axes object)

  • eps_x (float) – factor by which x-axis limits are expanded

  • eps_y (float) – factor by which y-axis limits are expanded. If eps_y is None, eps_x is used for eps_y as well.

Return type:

None

pyfda.libs.pyfda_lib.fil_convert(fil_dict: dict, format_in) None[source]

Convert between poles / zeros / gain, filter coefficients (polynomes) and second-order sections and store all formats not generated by the filter design routine in the passed dictionary fil_dict.

Parameters:
  • fil_dict (dict) – filter dictionary containing a.o. all formats to be read and written.

  • format_in (string or set of strings) –

    format(s) generated by the filter design routine. Must be one of

    ’sos’:

    a list of second order sections - all other formats can easily be derived from this format

    ’zpk’:

    [z,p,k] where z is the array of zeros, p the array of poles and k is a scalar with the gain - the polynomial form can be derived from this format quite accurately

    ’ba’:

    [b, a] where b and a are the polynomial coefficients - finding the roots of the a and b polynomes may fail for higher orders

Returns:

  • None

  • Exceptions

  • ———-

  • ValueError for Nan / Inf elements or other unsuitable parameters

pyfda.libs.pyfda_lib.fil_save(fil_dict: dict, arg, format_in: str, sender: str, convert: bool = True) None[source]

Save filter design arg given in the format specified as format_in in the dictionary fil_dict. The format can be either poles / zeros / gain, filter coefficients (polynomes) or second-order sections.

Convert the filter design to the other formats if convert is True.

Parameters:
  • fil_dict (dict) – The dictionary where the filter design is saved to.

  • arg (various formats) – The actual filter design

  • format_in (string) –

    Specifies how the filter design in ‘arg’ is passed:

    ’ba’:

    Coefficient form: Filter coefficients in FIR format (b, one dimensional) are automatically converted to IIR format (b, a).

    ’zpk’:

    Zero / pole / gain format: When only zeroes are specified, poles and gain are added automatically.

    ’sos’:

    Second-order sections

  • sender (string) – The name of the method that calculated the filter. This name is stored in fil_dict together with format_in.

  • convert (boolean) – When convert = True, convert arg to the other formats.

Return type:

None

pyfda.libs.pyfda_lib.floor_even(x) int[source]

Return the largest even integer not larger than x. x can be integer or float.

pyfda.libs.pyfda_lib.floor_odd(x) int[source]

Return the largest odd integer not larger than x. x can be integer or float.

pyfda.libs.pyfda_lib.format_ticks(ax, xy: str, scale: float = 1.0, format: str = '%.1f') None[source]

Reformat numbers at x or y - axis. The scale can be changed to display e.g. MHz instead of Hz. The number format can be changed as well.

Parameters:
  • ax (axes object)

  • xy (string, either 'x', 'y' or 'xy') – select corresponding axis (axes) for reformatting

  • scale (float (default: 1.)) – rescaling factor for the axes

  • format (string (default: %.1f)) – define C-style number formats

Return type:

None

Examples

Scale all numbers of x-Axis by 1000, e.g. for displaying ms instead of s.

>>> format_ticks('x',1000.)

Two decimal places for numbers on x- and y-axis

>>> format_ticks('xy',1., format = "%.2f")
pyfda.libs.pyfda_lib.lin2unit(lin_value: float, filt_type: str, amp_label: str, unit: str = 'dB') float[source]

Convert linear amplitude specification to dB or W, depending on filter type (‘FIR’ or ‘IIR’) and whether the specifications belong to passband or stopband. This is determined by checking whether amp_label contains the strings ‘PB’ or ‘SB’ :

  • Passband:
    \[ \begin{align}\begin{aligned}\begin{split}\\text{IIR:}\quad A_{dB} &= -20 \log_{10}(1 - lin\_value)\end{split}\\\begin{split}\\text{FIR:}\quad A_{dB} &= 20 \log_{10}\\frac{1 + lin\_value}{1 - lin\_value}\end{split}\end{aligned}\end{align} \]
  • Stopband:
    \[A_{dB} = -20 \log_{10}(lin\_value)\]

Returns the result as a float.

pyfda.libs.pyfda_lib.mod_version(mod: str = '') str[source]

Return the version of the module ‘mod’. If the module is not found, return empty string. When no module is specified, return a string with all modules and their versions sorted alphabetically.

pyfda.libs.pyfda_lib.qstr(text)[source]

Convert text (QVariant, QString, string) or numeric object to plain string.

In Python 3, python Qt objects are automatically converted to QVariant when stored as “data” (itemData) e.g. in a QComboBox and converted back when retrieving to QString. In Python 2, QVariant is returned when itemData is retrieved. This is first converted from the QVariant container format to a QString, next to a “normal” non-unicode string.

Parameters:

text (QVariant, QString, string or numeric data type that can be converted) – to string

Return type:

The current text data as a unicode (utf8) string

pyfda.libs.pyfda_lib.round_even(x) int[source]

Return the nearest even integer from x. x can be integer or float.

pyfda.libs.pyfda_lib.round_odd(x) int[source]

Return the nearest odd integer from x. x can be integer or float.

pyfda.libs.pyfda_lib.safe_eval(expr, alt_expr=0, return_type: str = 'float', sign: str = None)[source]

Try … except wrapper around numexpr to catch various errors When evaluation fails or returns None, try evaluating alt_expr. When this also fails, return 0 to avoid errors further downstream.

Parameters:
  • expr (str or scalar) – Expression to be evaluated, is cast to a string

  • alt_expr (str or scalar) – Expression to be evaluated when evaluation of first string fails, is cast to a string.

  • return_type (str) – Type of returned variable [‘float’ (default) / ‘cmplx’ / ‘int’ / ‘’ or ‘auto’]

  • sign (str) –

    enforce positive / negative sign of result [‘pos’, ‘poszero’ / None (default)

    ’negzero’ / ‘neg’]

Returns:

the evaluated result or 0 when both arguments fail

Return type:

float (default) / complex / int

Function attribute err contains number of errors that have occurred during evaluation (0 / 1 / 2)

pyfda.libs.pyfda_lib.set_dict_defaults(d: dict, default_dict: dict) None[source]

Add the key:value pairs of default_dict to dictionary d in-place for all missing keys.

pyfda.libs.pyfda_lib.sos2zpk(sos)[source]

Taken from scipy/signal/filter_design.py - edit to eliminate first order section

Return zeros, poles, and gain of a series of second-order sections

Parameters:

sos (array_like) – Array of second-order filter coefficients, must have shape (n_sections, 6). See sosfilt for the SOS filter format specification.

Returns:

  • z (ndarray) – Zeros of the transfer function.

  • p (ndarray) – Poles of the transfer function.

  • k (float) – System gain.

Notes

Added in version 0.16.0.

pyfda.libs.pyfda_lib.to_html(text: str, frmt: str = None) str[source]
Convert text to HTML format:
  • pretty-print logger messages

  • convert “n” to “<br />

  • convert “< “ and “> “ to “&lt;” and “&gt;”

  • format strings with italic and / or bold HTML tags, depending on parameter frmt. When frmt=None, put the returned string between <span> tags to enforce HTML rendering downstream

  • replace ‘_’ by HTML subscript tags. Numbers 0 … 9 are never set to italic format

Parameters:
  • text (str) – Text to be converted

  • frmt (str) –

    define text style

    • ’b’ : bold text

    • ’i’ : italic text

    • ’bi’ or ‘ib’ : bold and italic text

Returns:

HTML - formatted text

Return type:

str

Examples

>>> to_html("F_SB", frmt='bi')
"<b><i>F<sub>SB</sub></i></b>"
>>> to_html("F_1", frmt='i')
"<i>F</i><sub>1</sub>"
pyfda.libs.pyfda_lib.unique_roots(p, tol: float = 0.001, magsort: bool = False, rtype: str = 'min', rdist: str = 'euclidian')[source]

Determine unique roots and their multiplicities from a list of roots.

Parameters:
  • p (array_like) – The list of roots.

  • tol (float, default tol = 1e-3) – The tolerance for two roots to be considered equal. Default is 1e-3.

  • magsort (Boolean, default = False) – When magsort = True, use the root magnitude as a sorting criterium (as in the version used in numpy < 1.8.2). This yields false results for roots with similar magniudes (e.g. on the unit circle) but is signficantly faster for a large number of roots (factor 20 for 500 double roots.)

  • rtype ({'max', 'min, 'avg'}, optional) – How to determine the returned root if multiple roots are within tol of each other. - ‘max’ or ‘maximum’: pick the maximum of those roots (magnitude ?). - ‘min’ or ‘minimum’: pick the minimum of those roots (magnitude?). - ‘avg’ or ‘mean’ : take the average of those roots. - ‘median’ : take the median of those roots

  • dist ({'manhattan', 'euclid'}, optional) – How to measure the distance between roots: ‘euclid’ is the euclidian distance. ‘manhattan’ is less common, giving the sum of the differences of real and imaginary parts.

Returns:

  • pout (list) – The list of unique roots, sorted from low to high (only for real roots).

  • mult (list) – The multiplicity of each root.

Notes

This utility function is not specific to roots but can be used for any sequence of values for which uniqueness and multiplicity has to be determined. For a more general routine, see numpy.unique.

Examples

>>> vals = [0, 1.3, 1.31, 2.8, 1.25, 2.2, 10.3]
>>> uniq, mult = unique_roots(vals, tol=2e-2, rtype='avg')

Check which roots have multiplicity larger than 1:

>>> uniq[mult > 1]
array([ 1.305])

Find multiples of complex roots on the unit circle: >>> vals = np.roots(1,2,3,2,1) uniq, mult = unique_roots(vals, rtype=’avg’)

pyfda.libs.pyfda_lib.unit2lin(unit_value: float, filt_type: str, amp_label: str, unit: str = 'dB') float[source]

Convert amplitude specification in dB or W to linear specs:

  • Passband:
    \[ \begin{align}\begin{aligned}\begin{split}\\text{IIR:}\quad A_{PB,lin} &= 1 - 10 ^{-unit\_value/20}\end{split}\\\begin{split}\\text{FIR:}\quad A_{PB,lin} &= \\frac{10 ^ {unit\_value/20} - 1}{10 ^ {unit\_value/20} + 1}\end{split}\end{aligned}\end{align} \]
  • Stopband:
    \[A_{SB,lin} = -10 ^ {-unit\_value/20}\]

Returns the result as a float.

pyfda_sig_lib

Library with various signal processing related functions

pyfda.libs.pyfda_sig_lib.angle_zero(X, n_eps=1000.0, mode='auto', wrapped='auto')[source]

Calculate angle of argument X when abs(X) > n_eps * machine resolution. Otherwise, zero is returned.

pyfda.libs.pyfda_sig_lib.div_safe(num, den, n_eps: float = 1.0, i_scale: float = 1.0, verbose: bool = False)[source]

Perform elementwise array division after treating singularities, meaning:

  • check whether denominator (den) coefficients approach zero

  • check whether numerator (num) or denominator coefficients are non-finite, i.e. one of nan, ìnf or ninf.

At each singularity, replace denominator coefficient by 1 and numerator coefficient by 0

Parameters:
  • num (array_like) – numerator coefficients

  • den (array_like) – denominator coefficients

  • n_eps (float) – n_eps * machine resolution is the limit for the denominator below which the ratio is set to zero. The machine resolution in numpy is given by np.spacing(1), the distance to the nearest number which is equivalent to matlab’s “eps”.

  • i_scale (float) – The scale for the index i for num, den for printing the index of the singularities.

  • verbose (bool, optional) – whether to print The default is False.

Returns:

ratio – The ratio of num and den (zero at singularities)

Return type:

array_like

pyfda.libs.pyfda_sig_lib.group_delay(b, a=1, nfft=512, whole=False, analog=False, verbose=True, fs=6.283185307179586, sos=False, alg='scipy', n_eps=100)[source]

Calculate group delay of a discrete time filter, specified by numerator coefficients b and denominator coefficients a of the system function H ( z).

When only b is given, the group delay of the transversal (FIR) filter specified by b is calculated.

Parameters:
  • b (array_like) – Numerator coefficients (transversal part of filter)

  • a (array_like (optional, default = 1 for FIR-filter)) – Denominator coefficients (recursive part of filter)

  • whole (boolean (optional, default : False)) – Only when True calculate group delay around the complete unit circle (0 … 2 pi)

  • verbose (boolean (optional, default : True)) – Print warnings about frequency points with undefined group delay (amplitude = 0) and the time used for calculating the group delay

  • nfft (integer (optional, default: 512)) – Number of FFT-points

  • fs (float (optional, default: fs = 2*pi)) – Sampling frequency.

  • alg (str (default: "scipy")) –

    The algorithm for calculating the group delay:
    • ”scipy” The algorithm used by scipy’s grpdelay,

    • ”jos”: The original J.O.Smith algorithm; same as in “scipy” except that the frequency response is calculated with the FFT instead of polyval

    • ”diff”: Group delay is calculated by differentiating the phase

    • ”Shpakh”: Group delay is calculated from second-order sections

  • n_eps (integer (optional, default : 100)) – Minimum value in the calculation of intermediate values before tau_g is set to zero.

Returns:

  • tau_g (ndarray) – group delay

  • w (ndarray) – angular frequency points where group delay was computed

Notes

The following explanations follow [JOS].

Definition and direct calculation (‘diff’)

The group delay \(\tau_g(\omega)\) of discrete time (DT) and continuous time (CT) systems is the rate of change of phase with respect to angular frequency. In the following, derivative is always meant w.r.t. \(\omega\):

\[\tau_g(\omega) = -\frac{\partial }{\partial \omega}\angle H( \omega) = -\frac{\partial \phi(\omega)}{\partial \omega} = - \phi'(\omega)\]

With numpy / scipy, the group delay can be calculated directly with

w, H = sig.freqz(b, a, worN=nfft, whole=whole)
tau_g = -np.diff(np.unwrap(np.angle(H)))/np.diff(w)

The derivative can create numerical problems for e.g. phase jumps at zeros of frequency response or when the complex frequency response becomes very small e.g. in the stop band.

This can be avoided by calculating the group delay from the derivative of the logarithmic frequency response in polar form (amplitude response and phase):

\[ \begin{align}\begin{aligned}\ln ( H( \omega)) = \ln \left({H_A( \omega)} e^{j \phi(\omega)} \right) = \ln \left({H_A( \omega)} \right) + j \phi(\omega)\\ \Rightarrow \; \frac{\partial }{\partial \omega} \ln ( H( \omega)) = \frac{H_A'( \omega)}{H_A( \omega)} + j \phi'(\omega)\end{aligned}\end{align} \]

where \(H_A(\omega)\) is the amplitude response. \(H_A(\omega)\) and its derivative \(H_A'(\omega)\) are real-valued, therefore, the group delay can be calculated by separating real and imginary components (and discarding the real part):

\[\begin{align} \Re \left\{\frac{\partial }{\partial \omega} \ln ( H( \omega))\right\} &= \frac{H_A'( \omega)}{H_A( \omega)} \ \Im \left\{\frac{\partial }{\partial \omega} \ln ( H( \omega))\right\} &= \phi'(\omega) \end{align}\]

and hence

\[\tau_g(\omega) = -\phi'(\omega) = -\Im \left\{ \frac{\partial }{\partial \omega} \ln ( H( \omega)) \right\} =-\Im \left\{ \frac{H'(\omega)}{H(\omega)} \right\}\]

Note: The last term contains the complex response \(H(\omega)\), not the amplitude response \(H_A(\omega)\)!

In the following, it will be shown that the derivative of birational functions (like DT and CT filters) can be calculated very efficiently and from this the group delay.

J.O. Smith’s basic algorithm for FIR filters (‘scipy’)

An efficient form of calculating the group delay of FIR filters based on the derivative of the logarithmic frequency response has been described in [JOS] and [Lyons08] for discrete time systems.

A FIR filter is defined via its polyome \(H(z) = \sum_k b_k z^{-k}\) and has the following derivative:

\[\frac{\partial }{\partial \omega} H(z = e^{j \omega T}) = \frac{\partial }{\partial \omega} \sum_{k = 0}^N b_k e^{-j k \omega T} = -jT \sum_{k = 0}^{N} k b_{k} e^{-j k \omega T} = -jT H_R(e^{j \omega T})\]

where \(H_R\) is the “ramped” polynome, i.e. polynome \(H\) multiplied with a ramp \(k\), yielding

\[\tau_g(e^{j \omega T}) = -\Im \left\{ \frac{H'(e^{j \omega T})} {H(e^{j \omega T})} \right\} = -\Im \left\{ -j T \frac{H_R(e^{j \omega T})} {H(e^{j \omega T})} \right\} = T \, \Re \left\{\frac{H_R(e^{j \omega T})} {H(e^{j \omega T})} \right\}\]

scipy’s grpdelay directly calculates the complex frequency response \(H(e^{j\omega T})\) and its ramped function at the frequency points using the polyval function.

When zeros of the frequency response are on or near the data points of the DFT, this algorithm runs into numerical problems. Hence, it is neccessary to check whether the magnitude of the denominator is less than e.g. 100 times the machine eps. In this case, \(\tau_g\) is set to zero.

J.O. Smith’s basic algorithm for IIR filters (‘scipy’)

IIR filters are defined by

\[H(z) = \frac {B(z)}{A(z)} = \frac {\sum b_k z^k}{\sum a_k z^k},\]

their group delay can be calculated numerically via the logarithmic frequency response as well.

The derivative of \(H(z)\) w.r.t. \(\omega\) is calculated using the quotient rule and by replacing the derivatives of numerator and denominator polynomes with their ramp functions:

\[\begin{split}\begin{align} \frac{H'(e^{j \omega T})}{H(e^{j \omega T})} &= \frac{\left(B(e^{j \omega T})/A(e^{j \omega T})\right)'}{B(e^{j \omega T})/A(e^{j \omega T})} = \frac{B'(e^{j \omega T}) A(e^{j \omega T}) - A'(e^{j \omega T})B(e^{j \omega T})} { A(e^{j \omega T}) B(e^{j \omega T})} \\ &= \frac {B'(e^{j \omega T})} { B(e^{j \omega T})} - \frac { A'(e^{j \omega T})} { A(e^{j \omega T})} = -j T \left(\frac { B_R(e^{j \omega T})} {B(e^{j \omega T})} - \frac { A_R(e^{j \omega T})} {A(e^{j \omega T})}\right) \end{align}\end{split}\]

This result is substituted once more into the log. derivative from above:

\[\begin{split}\begin{align} \tau_g(e^{j \omega T}) =-\Im \left\{ \frac{H'(e^{j \omega T})}{H(e^{j \omega T})} \right\} &=-\Im \left\{ -j T \left(\frac { B_R(e^{j \omega T})} {B(e^{j \omega T})} - \frac { A_R(e^{j \omega T})} {A(e^{j \omega T})}\right) \right\} \\ &= T \Re \left\{\frac { B_R(e^{j \omega T})} {B(e^{j \omega T})} - \frac { A_R(e^{j \omega T})} {A(e^{j \omega T})} \right\} \end{align}\end{split}\]

If the denominator of the computation becomes too small, the group delay is set to zero. (The group delay approaches infinity when there are poles or zeros very close to the unit circle in the z plane.)

J.O. Smith’s algorithm for CT filters

The same process can be applied for CT systems as well: The derivative of a CT polynome \(P(s)\) w.r.t. \(\omega\) is calculated by:

\[\frac{\partial }{\partial \omega} P(s = j \omega) = \frac{\partial }{\partial \omega} \sum_{k = 0}^N c_k (j \omega)^k = j \sum_{k = 0}^{N-1} (k+1) c_{k+1} (j \omega)^{k} = j P_R(s = j \omega)\]

where \(P_R\) is the “ramped” polynome, i.e. its k th coefficient is multiplied by the ramp k + 1, yielding the same form as for DT systems (but the ramped polynome has to be calculated differently).

\[\tau_g(\omega) = -\Im \left\{ \frac{H'(\omega)}{H(\omega)} \right\} = -\Im \left\{j \frac{H_R(\omega)}{H(\omega)} \right\} = -\Re \left\{\frac{H_R(\omega)}{H(\omega)} \right\}\]

J.O. Smith’s improved algorithm for IIR filters (‘jos’)

J.O. Smith gives the following speed and accuracy optimizations for the basic algorithm:

  • convert the filter to a FIR filter with identical phase and group delay (but with different magnitude response)

  • use FFT instead of polyval to calculate the frequency response

The group delay of an IIR filter \(H(z) = B(z)/A(z)\) can also be calculated from an equivalent FIR filter \(C(z)\) with the same phase response (and hence group delay) as the original filter. This filter is obtained by the following steps:

  • The zeros of \(A(z)\) are the poles of \(1/A(z)\), its phase response is \(\angle A(z) = - \angle 1/A(z)\).

  • Transforming \(z \rightarrow 1/z\) mirrors the zeros at the unit circle, correcting the negative phase response. This can be performed numerically by “flipping” the order of the coefficients and multiplying by \(z^{-N}\) where \(N\) is the order of \(A(z)\). This operation also conjugates the coefficients (?) which mirrors the zeros at the real axis. This effect has to be compensated, yielding the polynome \(\tilde{A}(z)\). It is the “flip-conjugate” or “Hermitian conjugate” of \(A(z)\).

    Frequently (e.g. in the scipy and until recently in the Matlab implementation) the conjugate operation is omitted which gives wrong results for complex coefficients.

  • Finally, \(C(z) = B(z) \tilde{A}(z)\):

\[C(z) = B(z)\left[ z^{-N}{A}^{*}(1/z)\right] = B(z)\tilde{A}(z)\]

where

\[\begin{split}\begin{align} \tilde{A}(z) &= z^{-N}{A}^{*}(1/z) = {a}^{*}_N + {a}^{*}_{N-1}z^{-1} + \ldots + {a}^{*}_1 z^{-(N-1)}+z^{-N}\\ \Rightarrow \tilde{A}(e^{j\omega T}) &= e^{-jN \omega T}{A}^{*}(e^{-j\omega T}) \\ \Rightarrow \angle\tilde{A}(e^{j\omega T}) &= -\angle A(e^{j\omega T}) - N\omega T \end{align}\end{split}\]

In Python, the coefficients of \(C(z)\) are calculated efficiently by convolving the coefficients of \(B(z)\) and \(\tilde{A}(z)\):

c = np.convolve(b, np.conj(a[::-1]))

where \(b\) and \(a\) are the coefficient vectors of the original numerator and denominator polynomes. The actual group delay is then calculated from the equivalent FIR filter as described above.

Calculating the frequency response with the np.polyval(p,z) function at the NFFT frequency points along the unit circle, \(z = \exp(-j \omega)\), seems to be numerically less robust than using the FFT for the same task, it is also much slower.

This measure fixes already most of the problems described for narrowband IIR filters in scipy issues [Scipy_9310] and [Scipy_1175]. In my experience, these problems occur for all narrowband IIR response types.

Shpak algorithm for IIR filters

The algorithm described above is numerically efficient but not robust for narrowband IIR filters. Especially for filters defined by second-order sections, it is recommended to calculate the group delay using the D. J. Shpak’s algorithm.

Code is available at [Endolith_5828333] (GPL licensed) or at [SPA] (MIT licensed).

This algorithm sums the group delays of the individual sections which is much more robust as only second-order functions are involved. However, converting (b,a) coefficients to SOS coefficients introduces inaccuracies.

Examples

>>> b = [1,2,3] # Coefficients of H(z) = 1 + 2 z^2 + 3 z^3
>>> tau_g, td = pyfda_lib.grpdelay(b)
pyfda.libs.pyfda_sig_lib.group_delayz(b, a, w, plot=None, fs=6.283185307179586)[source]

Compute the group delay of digital filter.

Parameters:
  • b (array_like) – Numerator of a linear filter.

  • a (array_like) – Denominator of a linear filter.

  • w (array_like) – Frequencies in the same units as fs.

  • plot (callable) – A callable that takes two arguments. If given, the return parameters w and gd are passed to plot.

  • fs (float, optional) – The angular sampling frequency of the digital system.

Returns:

  • w (ndarray) – The frequencies at which gd was computed, in the same units as fs.

  • gd (ndarray) – The group delay in seconds.

pyfda.libs.pyfda_sig_lib.impz(b, a=1, FS=1, N=0, step=False)[source]

Calculate impulse response of a discrete time filter, specified by numerator coefficients b and denominator coefficients a of the system function H(z).

When only b is given, the impulse response of the transversal (FIR) filter specified by b is calculated.

Parameters:
  • b (array_like) – Numerator coefficients (transversal part of filter)

  • a (array_like (optional, default = 1 for FIR-filter)) – Denominator coefficients (recursive part of filter)

  • FS (float (optional, default: FS = 1)) – Sampling frequency.

  • N (float (optional)) – Number of calculated points. Default: N = len(b) for FIR filters, N = 100 for IIR filters

  • step (bool (optional, default False)) – return the step response instead of the impulse response

Returns:

  • hn (ndarray) – impulse or step response with length N (see above)

  • td (ndarray) – contains the time steps with same length as hn

Examples

>>> b = [1,2,3] # Coefficients of H(z) = 1 + 2 z^2 + 3 z^3
>>> h, n = dsp_lib.impz(b)
pyfda.libs.pyfda_sig_lib.impz_len(system, zpk: bool = False, level: float = -40) int[source]

Calculate length of impulse response for FIR and IIR filters.

Parameters:
  • system (array-like) – The discrete-time system, either specified by its coefficients or its poles and zeros.

  • zpk (bool) – When False (default), the system is specified by its coefficients.

  • level (float) –

    The relative level in dB to which the impulse response has decayed (relevant

    IIR filters only).

Returns:

N – The length of the impulse response.

Return type:

int

Notes

For FIR filters, the lenngth of the impulse response is equal to the number of filter taps (= filter order + 1).

For IIR filters, it is only possible to approximate the number of steps until the relative level of the impulse response is below the specified level.

The most simple approximation for IIR filters is to only regard the pole \(p_{max}\) nearest to the unit circle. For many real-world filters and systems, this pole dominates the transient response (“dominant pole approximation”) [JOS_time_constant] with a time constant

\[\tau\approx \frac{T_S}{1-\abs{p}_{max}}\]

When calculating the number of samples instead of an absolute time, T_S = 1.

The number of samples required to reach level in dBs is calculated with

\[N \approx -level/20 * \ln(10) \tau\]

When poles are specified (zpk == True), it is of course easy to find the dominant pole. When coefficients of the system are specified, poles can be calculated as the roots of the system

Alternatively, partial fraction expansion (PFE) yields poles and residues. These can be used for a more general solution, see [dsp_stackexchange_2021], [dsp_stackexchange_2022].

pyfda.libs.pyfda_sig_lib.quadfilt_group_delayz(b, w, fs=6.283185307179586)[source]

Compute group delay of 2nd-order digital filter.

Parameters:
  • b (array_like) – Coefficients of a 2nd-order digital filter.

  • w (array_like) – Frequencies in the same units as fs.

  • fs (float, optional) – The sampling frequency of the digital system.

Returns:

  • w (ndarray) – The frequencies at which gd was computed.

  • gd (ndarray) – The group delay in seconds.

pyfda.libs.pyfda_sig_lib.sos_group_delayz(sos, w, plot=None, fs=6.283185307179586)[source]

Compute group delay of digital filter in SOS format.

Parameters:
  • sos (array_like) – Array of second-order filter coefficients, must have shape (n_sections, 6). Each row corresponds to a second-order section, with the first three columns providing the numerator coefficients and the last three providing the denominator coefficients.

  • w (array_like) – Frequencies in the same units as fs.

  • plot (callable, optional) – A callable that takes two arguments. If given, the return parameters w and gd are passed to plot.

  • fs (float, optional) – The sampling frequency of the digital system.

Returns:

  • w (ndarray) – The frequencies at which gd was computed.

  • gd (ndarray) – The group delay in seconds.

pyfda.libs.pyfda_sig_lib.validate_sos(sos)[source]

Helper to validate a SOS input

Copied from scipy.signal._filter_design._validate_sos()

pyfda.libs.pyfda_sig_lib.zeros_with_val(N: int, val: float = 1.0, pos: int = 0)[source]

Create a 1D array of N zeros where the element at position pos has the value val.

Parameters:
  • N (int) – number of elements

  • val (scalar) – value to be inserted at position pos (default: 1)

  • pos (int) – Position of val to be inserted (default: 0)

Returns:

arr – Array with zeros except for element at position pos

Return type:

np.ndarray

pyfda.libs.pyfda_sig_lib.zorp_group_delayz(zorp, w, fs=1)[source]

Compute group delay of digital filter with a single zero/pole.

Parameters:
  • zorp (complex) – Zero or pole of a 1st-order linear filter

  • w (array_like) – Frequencies in the same units as fs.

  • fs (float, optional) – The sampling frequency of the digital system.

Returns:

  • w (ndarray) – The frequencies at which gd was computed.

  • gd (ndarray) – The group delay in seconds.

pyfda.libs.pyfda_sig_lib.zpk2array(zpk)[source]

Test whether Z = zpk[0] and P = zpk[1] have the same length, if not, equalize the lengths by adding zeros.

Test whether gain = zpk[2] is a scalar or a vector and whether it (or the first element of the vector) is != 0. If the gain is 0, set gain = 1.

Finally, convert the gain into an vector with the same length as P and Z and return the the three vectors as one array.

Parameters:

zpk (list, tuple or ndarray) – Zeros, poles and gain of the system

Return type:

zpk as an array or an error string

pyfda.libs.pyfda_sig_lib.zpk_group_delay(z, p, k, w, plot=None, fs=6.283185307179586)[source]

Compute group delay of digital filter in zpk format.

Parameters:
  • z (array_like) – Zeroes of a linear filter

  • p (array_like) – Poles of a linear filter

  • k (scalar) – Gain of a linear filter

  • w (array_like) – Frequencies in the same units as fs.

  • plot (callable, optional) – A callable that takes two arguments. If given, the return parameters w and gd are passed to plot.

  • fs (float, optional) – The sampling frequency of the digital system.

Returns:

  • w (ndarray) – The frequencies at which gd was computed.

  • gd (ndarray) – The group delay in seconds.

pyfda_qt_lib

Library with various helper functions for Qt widgets

class pyfda.libs.pyfda_qt_lib.EventTypes[source]

https://stackoverflow.com/questions/62196835/how-to-get-string-name-for-qevent-in-pyqt5 Events in Qt5: https://doc.qt.io/qt-5/qevent.html

Stores a string name for each event type.

With PySide2 str() on the event type gives a nice string name, but with PyQt5 it does not. So this method works with both systems.

Example usage (simultaneous initialization and method call / translation) > event_str = EventTypes().as_string(QEvent.UpdateRequest) > assert event_str == “UpdateRequest”

Example usage, separate initialization and method call > event types = EventTypes() > event_str = event_types.as_string(event.type())

as_string(event: Type) str[source]

Return the string name for this event.

class pyfda.libs.pyfda_qt_lib.PushButton(parent=None, text: str = '', icon: QIcon = None, N_x: int = 8, checkable: bool = True, checked: bool = False, objectName='', **kwargs)[source]

Create a QPushButton with a width fitting the label with bold font as well

Parameters:
  • text (str) – Text for button (optional)

  • rtf (bool) – Render text as rich text

  • icon (QIcon) – Icon for button. Either text or icon must be defined.

  • N_x (int) – Width in number of “x”

  • margin (int) – margin for ?

  • checkable (bool) – Whether button is checkable

  • checked (bool) – Whether initial state is checked

eventFilter(self, a0: QObject | None, a1: QEvent | None) bool[source]
minimumSizeHint(self) QSize[source]
setCheckable(self, a0: bool)[source]
setChecked(self, a0: bool)[source]
sizeHint(self) QSize[source]
class pyfda.libs.pyfda_qt_lib.PushButtonRT(parent=None, text: str = '', N_x: int = 8, margin: int = 10, checkable: bool = True, checked: bool = False, objectName='', **kwargs)[source]

Subclass QPushButton using QLabel to render rich text

Parameters:
  • text (str) – Text for button (optional)

  • rtf (bool) – Render text as rich text

  • N_x (int) – Width in number of “x”

  • margin (int) – margin for ?

  • checkable (bool) – Whether button is checkable

  • checked (bool) – Whether initial state is checked

eventFilter(self, a0: QObject | None, a1: QEvent | None) bool[source]
minimumSizeHint(self) QSize[source]
setCheckable(self, a0: bool)[source]
setChecked(self, a0: bool)[source]
setText(self, text: str | None)[source]
sizeHint(self) QSize[source]
class pyfda.libs.pyfda_qt_lib.QHLine(width=1)[source]

Create a thin horizontal line utilizing the HLine property of QFrames Usage:

> myline = QHLine() > mylayout.addWidget(myline)

class pyfda.libs.pyfda_qt_lib.QLabelVert(text, orientation='west', forceWidth=True)[source]

Create a vertical label

Adapted from https://pyqtgraph.readthedocs.io/en/latest/_modules/pyqtgraph/widgets/VerticalLabel.html

https://stackoverflow.com/questions/34080798/pyqt-draw-a-vertical-label

check https://stackoverflow.com/questions/29892203/draw-rich-text-with-qpainter

minimumSizeHint(self) QSize[source]
paintEvent(self, a0: QPaintEvent | None)[source]
sizeHint(self) QSize[source]
class pyfda.libs.pyfda_qt_lib.QVLine(width=2)[source]

Create a thin vertical line utilizing the HLine property of QFrames Usage:

> myline = QVLine() > mylayout.addWidget(myline)

class pyfda.libs.pyfda_qt_lib.RotatedButton[source]

Create a rotated QPushButton

Taken from

https://forum.qt.io/topic/9279/moved-how-to-rotate-qpushbutton-63/7

minimumSizeHint(self) QSize[source]
paintEvent(self, a0: QPaintEvent | None)[source]
sizeHint(self) QSize[source]
pyfda.libs.pyfda_qt_lib.emit(self, dict_sig: dict = {}, sig_name: str = 'sig_tx') None[source]

Emit a signal self.<sig_name> (defined as a class attribute) with a dict dict_sig using Qt’s emit().

  • Add the keys ‘id’ and ‘class’ with id resp. class name of the calling instance if not contained in the dict

  • If key ‘ttl’ is in the dict and its value is less than one, terminate the signal. Otherwise, reduce the value by one.

  • If the sender has passed an objectName, add it with the key “sender_name” to the dict.

pyfda.libs.pyfda_qt_lib.popup_warning(self, param1: int = 0, param2: str = '', message: str = '') bool[source]

Pop-up a warning box and require a user prompt. When message == “”, warn of very large filter orders, otherwise display the passed message

pyfda.libs.pyfda_qt_lib.qcmb_box_add_item(cmb_box, item_list, data=True, fireSignals=False, caseSensitive=False)[source]

Add an entry in combobox with text / data / tooltipp from item_list. When the item is already in combobox (searching for data or text item, depending data), do nothing. Signals are blocked during the update of the combobox unless fireSignals is set True. By default, the search is case insensitive, this can be changed by passing caseSensitive=False.

Parameters:
  • item_list (list) – List with [“new_data”, “new_text”, “new_tooltip”] to be added.

  • data (bool (default: False)) – Whether the string refers to the data or text fields of the combo box

  • fireSignals (bool (default: False)) – When True, fire a signal if the index is changed (useful for GUI testing)

  • caseInsensitive (bool (default: False)) – When true, perform case sensitive search.

Returns:

  • The index of the found item with string / data. When not found in the

  • combo box, return index -1.

pyfda.libs.pyfda_qt_lib.qcmb_box_add_items(cmb_box: QComboBox, items_list: list) None[source]

Add items to combo box cmb_box with text, data and tooltip from the list items_list.

Text and tooltip are prepared for translation via self.tr()

Parameters:
  • cmb_box (instance of QComboBox) – Combobox to be populated

  • items_list (list) –

    List of combobox entries, in the format

    (“data 1st item”, “text 1st item”, “tooltip for 1st item” # [optional]), (“data 2nd item”, “text 2nd item”, “tooltip for 2nd item”)]

Return type:

None

pyfda.libs.pyfda_qt_lib.qcmb_box_del_item(cmb_box: QComboBox, string: str, data: bool = False, fireSignals: bool = False, caseSensitive: bool = False) int[source]

Try to find the entry in combobox corresponding to string in a text field (data = False) or in a data field (data=True) and delete the item. When string is not found,do nothing. Signals are blocked during the update of the combobox unless fireSignals is set True. By default, the search is case insensitive, this can be changed by passing caseSensitive=False.

Parameters:
  • string (str) – The label in the text or data field to be deleted.

  • data (bool (default: False)) – Whether the string refers to the data or text fields of the combo box

  • fireSignals (bool (default: False)) – When True, fire a signal if the index is changed (useful for GUI testing)

  • caseInsensitive (bool (default: False)) – When true, perform case sensitive search.

Returns:

  • The index of the item with string / data. When not found in the combo box,

  • return index -1.

pyfda.libs.pyfda_qt_lib.qcmb_box_populate(cmb_box: QComboBox, items_list: list, item_init: str) int[source]

Clear and populate combo box cmb_box with text, data and tooltip from the list items_list with initial selection of init_item (data).

Text and tooltip are prepared for translation via self.tr()

Parameters:
  • cmb_box (instance of QComboBox) – Combobox to be populated

  • items_list (list) –

    List of combobox entries, in the format [“Tooltip for Combobox”, (“data 1st item”, “text 1st item”, “tooltip for 1st item”), (“data 2nd item”, “text 2nd item”, “tooltip for 2nd item”)]

    Tooltipps are optional.

  • item_init (str) – data for initial setting of combobox. When data is not found, set combobox to first item.

Returns:

ret – Index of item_init in combobox. If index == -1, item_init was not in items_list

Return type:

int

pyfda.libs.pyfda_qt_lib.qget_cmb_box(cmb_box: QComboBox, data: bool = True) str[source]

Get current itemData or Text of comboBox and convert it to string.

In Python 3, python Qt objects are automatically converted to QVariant when stored as “data” e.g. in a QComboBox and converted back when retrieving. In Python 2, QVariant is returned when itemData is retrieved. This is first converted from the QVariant container format to a QString, next to a “normal” non-unicode string.

Returns:

The current text or data of combobox as a string

pyfda.libs.pyfda_qt_lib.qget_selected(table, select_all=False, reverse=True)[source]

Get selected cells in table and return a dictionary with the following keys:

‘idx’: indices of selected cells as an unsorted list of tuples

‘sel’: list of lists of selected cells per column, by default sorted in reverse

‘cur’: current cell selection as a tuple

Parameters:
  • select_all (bool) – select all table items and create a list when True

  • reverse (bool) – return selected fields upside down when True

pyfda.libs.pyfda_qt_lib.qset_cmb_box(cmb_box: QComboBox, string: str, data: bool = False, fireSignals: bool = False, caseSensitive: bool = False) int[source]

Set combobox to the index corresponding to string in a text field (data = False) or in a data field (data=True). When string is not found in the combobox entries, select the first entry. Signals are blocked during the update of the combobox unless fireSignals is set True. By default, the search is case insensitive, this can be changed by passing caseSensitive=False.

Parameters:
  • string (str) – The label in the text or data field to be selected. When the string is not found, select the first entry of the combo box.

  • data (bool (default: False)) – Whether the string refers to the data or text fields of the combo box

  • fireSignals (bool (default: False)) – When True, fire a signal if the index is changed (useful for GUI testing)

  • caseSensitive (bool (default: False)) – When true, perform case sensitive search.

Returns:

  • The index of the string. When the string was not found in the combo box,

  • select first entry of combo box and return index -1.

pyfda.libs.pyfda_qt_lib.qstyle_widget(widget, state)[source]

Apply the “state” defined in pyfda_rc.py to the widget, e.g.: Color the >> DESIGN FILTER << button according to the filter design state.

This requires setting the property, “unpolishing” and “polishing” the widget and finally forcing an update.

  • ‘normal’ : default, no color styling

  • ‘ok’: green, filter has been designed, everything ok

  • ‘changed’: yellow, filter specs have been changed

  • ‘running’: orange, simulation is running

  • ‘error’ : red, an error has occurred during filter design

  • ‘failed’ : pink, filter fails to meet target specs (not used yet)

  • ‘u’ or ‘unused’ : grey text color

  • ‘d’ or ‘disabled’: background color darkgrey

  • ‘a’ or ‘active’ : no special style defined

pyfda.libs.pyfda_qt_lib.qtext_height(text: str = 'X', font=None) int[source]

Calculate size of text in points`.

The actual size of the string is calculated using fontMetrics and the default or the passed font

Parameters:

test (str) – string to calculate the height for (default: “X”)

Returns:

lineSpacing – The height of the text (line spacing) in points

Return type:

int

Notes

This is based on https://stackoverflow.com/questions/27433165/how-to-reimplement-sizehint-for-bold-text-in-a-delegate-qt

and

https://stackoverflow.com/questions/47285303/how-can-i-limit-text-box-width-of- # qlineedit-to-display-at-most-four-characters/47307180#47307180

https://stackoverflow.com/questions/56282199/fit-qtextedit-size-to-text-size-pyqt5

pyfda.libs.pyfda_qt_lib.qtext_width(text: str = '', N_x: int = 17, bold: bool = True, font=None) int[source]

Calculate width of text in points`. When text=`, calculate the width of number N_x of characters ‘x’.

The actual width of the string is calculated by creating a QTextDocument with the passed text and retrieving its idealWidth()

Parameters:
  • test (str) – string to calculate the width for

  • N_x (int) – When text == ‘’, calculate the width from N_x * width(‘x’)

  • bold (bool (default: True)) – When True, determine width based on bold font

Returns:

width – The width of the text in points

Return type:

int

Notes

This is based on https://stackoverflow.com/questions/27433165/how-to-reimplement-sizehint-for-bold-text-in-a-delegate-qt

and

https://stackoverflow.com/questions/47285303/how-can-i-limit-text-box-width-of- # qlineedit-to-display-at-most-four-characters/47307180#47307180

pyfda.libs.pyfda_qt_lib.qwindow_stay_on_top(win: QDialog, top: bool) None[source]

Set flags for a window such that it stays on top (True) or not

On Windows 7 the new window stays on top anyway. Additionally setting WindowStaysOnTopHint blocks the message window when trying to close pyfda.

On Windows 10 and Linux, WindowStaysOnTopHint needs to be set.

pyfda_io_lib

pyfda_fix_lib