pyepri.apodization

This module contains functions for computing frequency attenuation (a.k.a., apodization) profiles and performing frequency attenuation operations by circular convolution.

Functions

`rfftconvol`(u, rfft_filter[, backend, dim, notest])	1D Fourier convolution (a.k.a. Circular convolution) for real input signal.
`smoothstep`(t, m, sig[, backend, notest])	Apodization profile computed as the convolution between a rectangular step function and a Gaussian function.
`_check_nd_inputs_`(caller, backend[, u, rfft_filter, ...])	Factorized consistency checks for functions in the `pyepri.apodization` submodule.

Module Contents

pyepri.apodization.rfftconvol(u, rfft_filter, backend=None, dim=-1, notest=False)[source]

1D Fourier convolution (a.k.a. Circular convolution) for real input signal.

Parameters:

u (array_like (with type backends.cls)) – Input signal, can be multidimensional, the filtering operation will be applied along a single dimension (broadcasting).
rfft_filter (array_like (with type backend.cls)) –
Input convolution filter in Fourier domain, can be either mono-dimensional with length u.shape[dim]//2+1 OR multidimensional with shape satisfying
- rfft_filter.shape[k] = u.shape[k]//2+1 if k == dim;
- rfft_filter.shape[k] = u.shape[k] otherwise.
backend (<class 'pyepri.backends.Backend'> or None, optional) –
A numpy, cupy or torch backend (see pyepri.backends module).

When backend is None, a default backend is inferred from the input arrays (u, rfft_filter).
dim (int, optional) – Dimension along which the one-dimensional direct/inverse real FFTs (see how output is computed below) will be computed.
notest (bool, optional) – Set notest=True to disable consistency checks.

Returns:

out – Filtered signal computed using

v = backend.irfft(backend.rfft(u, dim=dim) * rfft_filter, n=u.shape[dim], dim=dim).

Return type:

array_like (with type backend.cls)

Example

>>> ##################
>>> # import modules #
>>> ##################
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> import pyepri.apodization as apodization
>>> import pyepri.backends as backends
>>>
>>> ######################################################
>>> # create a numpy backend (you can try other choices) #
>>> ######################################################
>>> backend = backends.create_numpy_backend()
>>>
>>> #################################
>>> # compute a random input signal #
>>> #################################
>>> u = backend.rand(100,)
>>>
>>> ###########################################################
>>> # compute filter in the frequency domain over half of the #
>>> # full frequency domain (here a smoothed step profile)    #
>>> ###########################################################
>>> t = backend.arange(len(u)//2+1, dtype='float32')
>>> m, sig = 40, 2.
>>> rfft_filter = apodization.smoothstep(t, m, sig, backend=backend)
>>>
>>> ###########################
>>> # compute filtered signal #
>>> ###########################
>>> v = apodization.rfftconvol(u, rfft_filter, backend=backend)
>>>
>>> ###################
>>> # display signals #
>>> ###################
>>> plt.figure(figsize=(10,5))
>>> plt.subplot(1,2,1)
>>> plt.plot(backend.to_numpy(u))
>>> plt.plot(backend.to_numpy(v))
>>> plt.title("input & filtered signals")
>>> plt.legend(["input signal", "filtered signal"])
>>> plt.subplot(1,2,2)
>>> plt.plot(backend.to_numpy(t),backend.to_numpy(rfft_filter))
>>> plt.xlabel("frequency indexes (half of the full domain)")
>>> plt.title("filter in Fourier domain (frequency attenuation profile)")
>>> plt.show()

pyepri.apodization.smoothstep(t, m, sig, backend=None, notest=False)[source]

Apodization profile computed as the convolution between a rectangular step function and a Gaussian function.

The computed profile roughly (but not exactly) satisfies:

h(t) = 1 for abs(t) < m/2 - 5*sig
h(t) = 0 for abs(t) > m/2 + 5*sig

leading to the following “transition intervals” (where the profile increases from 0 to 1 or decreases from 1 to 0):

-m/2 + [-5*sig, 5*sig] : the profile increases from roughly 0 to roughly 1;
m/2 + [-5*sig, 5*sig] : the profile decreases from roughly 1 to roughly 0.

Parameters:

t (array_like (with type backend.cls)) – Input sampling points for the apodization profile.
m (float) – Support length for the step function (function f in the above mathematical description).
sig (float) – Smoothing parameter of the apodization profile (= standard deviation of the Gaussian profile g in the above mathematical description).
backend (<class 'pyepri.backends.Backend'>) –
A numpy, cupy or torch backend (see pyepri.backends module).

When backend is None, a default backend is inferred from the input array t.
notest (bool, optional) – Set notest=True to disable consistency checks.

Returns:

h – Computed values of h(t) (same shape as input t).

Return type:

array_like (with type backend.cls)

Notes

This module computes

\[h(t) = (f * g) (t)\]

where \(*\) denotes the convolution product, \(f\) is the rectangular step function with support length m defined by

\[\begin{split}\forall x\in\mathbb{R}\,, \quad f(x) = \left\{\begin{array}{cl}x&\text{if } |x| \leq \frac{m}{2}\\0&\text{otherwise,}\end{array}\right.\end{split}\]

and \(g\) is the Gaussian function with standard deviation \(\sigma\) (corresponding the the input parameter sig) defined by

\[\forall x\in\mathbb{R}\,,\quad g(x) = \frac{1}{\sigma \sqrt{2\pi}} \cdot \exp{\left(-\frac{x^2}{2\sigma^2}\right)}\,.\]

The practical evaluation of the output smooth step apodization profile h is done using

\[\forall t\in \mathbb{R}\,,\quad h(t) = \frac{1}{2} \cdot \left(\mathrm{erfc}\left(\frac{t-m/2}{\sigma \sqrt{2}}\right) - \mathrm{erfc}\left(\frac{t+m/2}{\sigma \sqrt{2}}\right)\right)\]

where

\[\forall t \in \mathbb{R}\,,\quad \mathrm{erfc}(t) = \frac{2}{\sqrt{\pi}} \int_{t}^{+\infty} e^{-x^2} \, dx\]

is the complementary error function.

Example

>>> ##################
>>> # import modules #
>>> ##################
>>> import numpy as np
>>> import pyepri.backends as backends
>>> import pyepri.apodization as apodization
>>> import matplotlib.pyplot as plt
>>>
>>> ##########################
>>> # prepare figure display #
>>> ##########################
>>> plt.rcParams.update({'axes.xmargin' : 0, 'axes.ymargin' : 0})
>>> plt.figure(figsize=(9.5, 5))
>>> plt.xlabel("t")
>>> plt.ylabel("h(t): smooth step profile")
>>>
>>> ###############################
>>> # compute apodization profile #
>>> ###############################
>>> backend = backends.create_numpy_backend()
>>> t = np.linspace(-500, 500, 10000)
>>> m = 350
>>> sig = 10
>>> h = apodization.smoothstep(t, m, sig, backend=backend)
>>>
>>> #############################
>>> # draw the computed profile #
>>> #############################
>>> plt.plot(backend.to_numpy(t), backend.to_numpy(h))
>>>
>>> #########################################################
>>> # draw boundaries on the right-side transition interval #
>>> #########################################################
>>> plt.plot((m/2-5*sig)*np.array([1, 1]), np.array([-.1, 1.1]), 'k:')
>>> plt.plot((m/2+5*sig)*np.array([1, 1]), np.array([-.1, 1.1]), 'k:')
>>>
>>> #########################################################
>>> # draw boundaries on the right-side transition interval #
>>> #########################################################
>>> plt.plot((-m/2+5*sig)*np.array([1, 1]), np.array([-.1, 1.1]), 'k:')
>>> plt.plot((-m/2-5*sig)*np.array([1, 1]), np.array([-.1, 1.1]), 'k:')
>>> plt.legend(('smooth step profile', 'boundaries of the' + '\n' +
... 'transition intervals'), loc='upper right')
>>> plt.show()

pyepri.apodization._check_nd_inputs_(caller, backend, u=None, rfft_filter=None, dim=None, t=None, m=None, sig=None)[source]: Factorized consistency checks for functions in the pyepri.apodization submodule.