# Source code for PyDynamic.uncertainty.propagate_convolution

```"""This module assists in uncertainty propagation for the convolution operation

The convolution operation is a common operation in signal and data
processing. Convolving signals is mathematically similar to a filter
application.

This module contains the following function:

* :func:`convolve_unc`: Convolution with uncertainty propagation based on FIR-filter
"""

__all__ = ["convolve_unc"]

import numpy as np

from .propagate_filter import _fir_filter
from ..misc.tools import is_vector

[docs]def convolve_unc(x1, U1, x2, U2, mode="full"):
"""Discrete convolution of two signals with uncertainty propagation

This function supports the convolution modes of :func:`numpy.convolve` and
:func:`scipy.ndimage.convolve1d`.

.. note::
The option to provide the uncertainties as 1D-arrays of standard uncertainties
is given for convenience only. It does not result in any performance benefits,
as they are internally just converted into a diagonal covariance matrix.
Moreover, the output will always be a full covariance matrix (and will almost
always have off-diagonal entries in practical scenarios).

Parameters
----------
x1 : np.ndarray, (N,)
first input signal
U1 : np.ndarray, (N,) or (N, N)
- 1D-array: standard uncertainties associated with x1 (corresponding to uncorrelated entries of x1)
- 2D-array: full 2D-covariance matrix associated with x1
- None: corresponds to a fully certain signal x1, results in more efficient calculation (compared to using np.zeros(...))
x2 : np.ndarray, (M,)
second input signal
U2 : np.ndarray, (M,) or (M, M)
- 1D-array: standard uncertainties associated with x2 (corresponding to uncorrelated entries of x2)
- 2D-array: full 2D-covariance matrix associated with x2
- None: corresponds to a fully certain signal x2, results in more efficient calculation (compared to using np.zeros(...))
mode : str, optional
:func:`numpy.convolve`-modes:

- full:  len(y) == N+M-1 (default)
- valid: len(y) == max(M, N) - min(M, N) + 1
- same:  len(y) == max(M, N) (value+covariance are padded with zeros)

:func:`scipy.ndimage.convolve1d`-modes:

- nearest: len(y) == N (value+covariance are padded with by stationary assumption)
- reflect:  len(y) == N
- mirror:   len(y) == N

Returns
-------
conv : np.ndarray
convoluted output signal
Uconv : np.ndarray
full 2D-covariance matrix of y

References
----------
.. seealso::
:func:`numpy.convolve`
:func:`scipy.ndimage.convolve1d`
"""

# if a numpy-mode is chosen, x1 is expected to be the longer signal
# remember that pure convolution is commutative
if len(x1) < len(x2) and mode in ["valid", "full", "same"]:
x1, x2 = x2, x1
U1, U2 = U2, U1

# convert 1d array of standard uncertainties to covariance matrix (if necessary)
U1, U2 = _ensure_cov_matrix(U1), _ensure_cov_matrix(U2)

# actual computation
if mode == "valid":
# apply _fir_filter directly
y, Uy = _fir_filter(x=x1, theta=x2, Ux=U1, Utheta=U2, initial_conditions="zero")

# compensate boundary adjustments from _fir_filter
conv = y[len(x2) - 1 :]
Uconv = Uy[len(x2) - 1 :, len(x2) - 1 :]

elif mode == "full":
# append zeros to adapt to _fir_filter mechanism
if isinstance(U1, np.ndarray):
)
else:
U1_mod = None

# apply _fir_filter
y, Uy = _fir_filter(
x=x1_mod, theta=x2, Ux=U1_mod, Utheta=U2, initial_conditions="zero"
)

# use output directly
conv = y
Uconv = Uy

elif mode == "same":
# append zeros to adapt to _fir_filter mechanism
pad_len = (len(x2) - 1) // 2
if isinstance(U1, np.ndarray):
)
else:
U1_mod = None

# apply _fir_filter
y, Uy = _fir_filter(
x=x1_mod, theta=x2, Ux=U1_mod, Utheta=U2, initial_conditions="zero"
)

# compensate boundary adjustments from _fir_filter

elif mode in ["nearest", "reflect", "mirror"]:

# scipy.ndimage.convolve1d and numpy.pad use different (but overlapping) terminology
mode_translation = {
"nearest": "edge",
"reflect": "symmetric",
"mirror": "reflect",
}

# prepend and append to x1 and U1 to get requested boundary effect
n1 = len(x1)
n2 = len(x2)
pad_len = (n2 + 1) // 2
# we assume that U1 is an array or None
if isinstance(U1, np.ndarray):
else:
U1_mod = None

# apply _fir_filter
y, Uy = _fir_filter(
x=x1_mod, theta=x2, Ux=U1_mod, Utheta=U2, initial_conditions="zero"
)

# compensate boundary adjustments from _fir_filter
conv = y[n2 : n2 + n1]
Uconv = Uy[n2 : n2 + n1, n2 : n2 + n1]

else:
raise ValueError(f'convolve_unc: Mode "{mode}" is not supported.')

return conv, Uconv

def _ensure_cov_matrix(unc_array):
"""
Converts 1D-arrays of standard uncertainties into the corresponding
(diagonal) covariance matrix by *square*+diag.

Does not modify inputs which are 2D-arrays or None.
"""

# squeeze() ensures proper execution on arrays with only one dimension of non-zero length
if unc_array is not None and is_vector(unc_array):
unc_array = np.diag(np.square(unc_array))

return unc_array
```