ltfatpy.gabor package - Gabor analysis¶
Basic Time/Frequency analysis¶
s0norm¶
Module of S0 norm calculation
Ported from ltfat_2.1.0/gabor/s0norm.m
-
ltfatpy.gabor.s0norm.
s0norm
(f, rel=False, dim=None)[source]¶ S0-norm of signal
Usage:
y = s0norm(f)
Input parameters:
Parameters: - f (numpy.ndarray) – is the signal
- rel (bool) – True if returning result should be relative to the norm (the energy) of the signal. False by default
- dim (int) – dimension along which norm is applied (first non-singleton dimension as default)
- Output parameters:
Returns: s0-norm Return type: float s0norm(f)
computes the -norm of a vector.If the input is a matrix or ND-array, the -norm is computed along the first (non-singleton) dimension, and a vector of values is returned.
Warning
The -norm is computed by computing a full short-time Fourier transform of a signal, which can be quite time-consuming. Use this function with care for long signals.
Gabor systems¶
dgt¶
Module of dgt calculation
Ported from ltfat_2.1.0/gabor/dgt.m
-
ltfatpy.gabor.dgt.
dgt
(f, g, a, M, L=None, pt='freqinv')[source]¶ Discrete Gabor Transform
Usage:
(c, Ls, g) = dgt(f, g, a, M)
(c, Ls, g) = dgt(f, g, a, M, L)
(c, Ls, g) = dgt(f, g, a, M, L, pt)
Input parameters:
Parameters: - f (numpy.ndarray) – Input data. f dtype has to be float64 or complex128.
- g (str, dict or numpy.ndarray) – Window function.
- a (int) – Length of time shift.
- M (int) – Number of channels.
- L (int) – Length of transform to do. Default is None.
- pt (str) – ‘freqinv’ or ‘timeinv’. Default is ‘freqinv’
- Output parameters:
Returns: (c, Ls, g)
Return type: Variables: - c (numpy.ndarray) – array of gabor transform coefficients. Its dtype is complex128.
- Ls (int) – length of input signal
- g (numpy.ndarray) – updated window function. Its dtype is float64 or complex128 depending on f dtype.
dgt(f, g, a, M)
computes the Gabor coefficients (also known as a windowed Fourier transform) of the input signal f with respect to the window g and parameters a and M. The output is a one or two dimensionalnumpy.ndarray
in a rectangular layout.The length of the transform will be the smallest multiple of a and M that is larger than the signal. f will be zero-extended to the length of the transform. If f is a 2d array, the transformation is applied to each column. The length of the transform done can be obtained by
L = c.shape[1] * a
The window g may be an array of numerical values, a text string or a dictionary. See the help of
gabwin()
for more details.dgt(f, g, a, M, L)
computes the Gabor coefficients as above, but does a transform of length L. f will be cut or zero-extended to length L before the transform is done.(c, Ls) = dgt(f, g, a, M)
or(c, Ls) = dgt(f, g, a, M, L)
returns the length of the input signal f. This is handy for reconstruction:- Examples:
In the following example we create a Hermite function, which is a complex-valued function with a circular spectrogram, and visualize the coefficients using both
imshow()
andplotdgt()
:>>> import numpy as np >>> import matplotlib.pyplot as plt >>> from ltfatpy import plotdgt >>> from ltfatpy import pherm >>> a = 10 >>> M = 40 >>> L = a * M >>> h, _ = pherm(L, 4) # 4th order hermite function. >>> c = dgt(h, 'gauss', a, M)[0] >>> # Simple plot: The squared modulus of the coefficients on >>> # a linear scale >>> _ = plt.imshow(np.abs(c)**2, interpolation='nearest', origin='upper') >>> plt.show() >>> # Better plot: zero-frequency is displayed in the middle, >>> # and the coefficients are show on a logarithmic scale. >>> _ = plotdgt(c, a, dynrange=50) >>> plt.show()
(c, Ls, g)=dgt(...)
outputs the window used in the transform. This is useful if the window was generated from a description in a string or dictionary.The Discrete Gabor Transform is defined as follows: Consider a window g and a one-dimensional signal f of length L and define . The output from
c = dgt(f, g, a, M)
is then given by:where , and are computed modulo L.
Additional parameters:
dgt
takes the following keyword at the end of the line of input arguments:- pt = ‘freqinv’
Compute a DGT using a frequency-invariant phase. This is the default convention described above.
- pt = ‘timeinv’
Compute a DGT using a time-invariant phase. This convention is typically used in FIR-filter algorithms.
See also
idgt()
,gabwin()
,dwilt()
,gabdual()
,phaselock()
- References:
- [FS98][Grochenig01]
idgt¶
Module of idgt calculation
Ported from ltfat_2.1.0/gabor/idgt.m
-
ltfatpy.gabor.idgt.
idgt
(coef, g, a, Ls=None, pt='freqinv')[source]¶ Inverse discrete Gabor transform
Usage:
(f, g) = idgt(c, g, a)
(f, g) = idgt(c, g, a, Ls)
(f, g) = idgt(c, g, a, Ls, pt)
Input parameters:
Parameters: - c (numpy.ndarray) – Array of coefficients
- g (str, dict or numpy.ndarray) – Window function
- a (int) – Length of time shift
- Ls (int) – Length of signal
- pt (str) – ‘freqinv’ or ‘timeinv’. Default is ‘freqinv’.
- Output parameters:
Returns: Signal (dtype = complex128) Return type: numpy.ndarray idgt(c, g, a)
computes the Gabor expansion of the input coefficients c with respect to the window g and time shift a. The number of channels is deduced from the size of the coefficients c.idgt(c, g, a, Ls)
does as above but cuts or extends f to length Ls.(f, g)=idgt(...)
additionally outputs the window used in the transform. This is useful if the window was generated from a description in a string or cell array.For perfect reconstruction, the window used must be a dual window of the one used to generate the coefficients.
The window g may be a vector of numerical values, a text string or a cell array. See the help of
gabwin()
for more details.If g is a row vector, then the output will also be a row vector. If c is 3-dimensional, then
idgt
will return a matrix consisting of one column vector for each of the TF-planes in c.Assume that
f=idgt(c, g, a, L)
for an array c of size . Then the following holds for :Additional parameters:
idgt
takes the following keyword at the end of the line of input arguments:- pt=’freqinv’
Compute a DGT using a frequency-invariant phase. This is the default convention described above.
- pt=’timeinv’
Compute a DGT using a time-invariant phase. This convention is typically used in FIR-filter algorithms.
Examples:
The following example demonstrates the basic principles for getting perfect reconstruction (short version):
>>> from ltfatpy import greasy >>> from ltfatpy import dgt >>> f = greasy()[0] # Input test signal >>> a = 32 # time shift >>> M = 64 # frequency shift >>> gs = {'name': 'blackman', 'M': 128} # synthesis window >>> # analysis window >>> ga = {'name' : ('dual', gs['name']), 'M' : 128} >>> (c, Ls) = dgt(f, ga, a, M)[0:2] # analysis >>> # ... do interesting stuff to c at this point ... >>> r = idgt(c, gs, a, Ls)[0] # synthesis >>> np.linalg.norm(f-r) < 1e-10 # test True
See also
dgt()
,gabwin()
,dwilt()
,gabtight()
dgtreal¶
Module of dgtreal calculation
Ported from ltfat_2.1.0/gabor/dgtreal.m
-
ltfatpy.gabor.dgtreal.
dgtreal
(f, g, a, M, L=None, pt='freqinv')[source]¶ Discrete Gabor transform for real signals
- Usage:
(c, Ls, g) = dgtreal(f, g, a, M)
(c, Ls, g) = dgtreal(f, g, a, M, L)
(c, Ls, g) = dgtreal(f, g, a, M, L, pt)
Input parameters:
Parameters: - Output parameters:
Returns: (c, Ls, g)
Return type: Variables: - c (numpy.ndarray) – array of complex128 coefficients.
- Ls (int) – length of input signal
- g (numpy.ndarray) – updated window function. dtype is float64.
dgtreal(f, g, a, M)
computes the Gabor coefficients (also known as a windowed Fourier transform) of the real-valued input signal f with respect to the real-valued window g and parameters a and M. The output is a vector/matrix in a rectangular layout.As opposed to
dgt()
only the coefficients of the positive frequencies of the output are returned.dgtreal
will refuse to work for complex valued input signals.The length of the transform will be the smallest multiple of a and M that is larger than the signal. f will be zero-extended to the length of the transform. If f is a matrix, the transformation is applied to each column. The length of the transform done can be obtained by
L = c.shape[1] * a
.The window g may be a vector of numerical values, a text string or a dictionary. See the help of
gabwin()
for more details.dgtreal(f, g, a, M, L)
computes the Gabor coefficients as above, but does a transform of length L. f will be cut or zero-extended to length L before the transform is done.The
dgtreal
function returns the length of the input signal f. This is handy for reconstruction:>>> (c, Ls, g) = dgtreal(f, g, a, M) >>> fr = idgtreal(c, gd, a, M, Ls)
will reconstruct the signal f no matter what the length of f is, provided that gd is a dual window of g.
It also outputs the window used in the transform. This is useful if the window was generated from a description in a string or dictionary.
See the help on
dgt()
for the definition of the discrete Gabor transform. This routine will return the coefficients for channel frequencies from 0 tofloor(M/2)
.dgtreal
optionnaly takes a pt argument:- ‘freqinv’
- Compute a
dgtreal
using a frequency-invariant phase. This is the default convention described in the help fordgt()
. - ‘timeinv’
- Compute a
dgtreal
using a time-invariant phase. This convention is typically used in filter bank algorithms.
dgtreal
can be used to manually compute a spectrogram, if you want full control over the parameters and want to capture the output.- Example:
>>> import matplotlib.pyplot as plt >>> from ltfatpy import greasy >>> (f,fs) = greasy() # Input test signal >>> a = 10 # Downsampling factor in time >>> M = 200 # Total number of channels, only 101 will be computed >>> # Compute the coefficients using a 20 ms long Hann window >>> c = dgtreal(f, {'name' : 'hann', 'M' : 0.02*fs}, a, M)[0] >>> # Visualize the coefficients as a spectrogram >>> dynrange = 90 # 90 dB dynamical range for the plotting >>> from ltfatpy import plotdgtreal >>> coef = plotdgtreal(c, a, M, fs=fs, dynrange=dynrange) >>> plt.show()
See also
dgt()
,idgtreal()
,gabwin()
,dwilt()
,gabtight()
,plotdgtreal()
- References:
- [FS98][Grochenig01]
idgtreal¶
Module of idgtreal calculation
Ported from ltfat_2.1.0/gabor/idgtreal.m
-
ltfatpy.gabor.idgtreal.
idgtreal
(coef, g, a, M, Ls=None, pt='freqinv')[source]¶ Inverse discrete Gabor transform for real-valued signals
- Usage:
(f, g) = idgtreal(c, g, a, M)
(f, g) = idgtreal(c, g, a, M, Ls)
(f, g) = idgtreal(c, g, a, M, Ls, pt)
Input parameters:
Parameters: - c (numpy.ndarray) – Array of coefficients
- g (str, dict or numpy.ndarray) – Window function
- a (int) – Length of time shift
- M (int) – Number of channels
- Ls (int) – Length of signal
- pt (str) – ‘freqinv’ or ‘timeinv’. Default is ‘freqinv’.
- Output parameters:
Returns: (f, g)
Return type: Variables: - f (numpy.ndarray) – signal
- g (numpy.ndarray) – window
idgtreal(c, g, a, M)
computes the Gabor expansion of the input coefficients c with respect to the real-valued window g, time shift a and number of channels M. c is assumed to be the positive frequencies of the Gabor expansion of a real-valued signal.It must hold that
c.shape[0] == np.floor(M/2)+1
. Note that since the correct number of channels cannot be deduced from the input,idgtreal
takes an additional parameter as opposed toidgt()
.The window g may be a vector of numerical values, a text string or a dictionary. See the help of
gabwin()
for more details.idgtreal(c, g, a, M, Ls)
does as above but cuts or extends f to length Ls.(f, g) = idgtreal(...)
outputs the window used in the transform. This is useful if the window was generated from a description in a string or dictionary.For perfect reconstruction, the window used must be a dual window of the one used to generate the coefficients.
If g is a row vector, then the output will also be a row vector. If c is 3-dimensional, then
idgtreal
will return a matrix consisting of one column vector for each of the TF-planes in c.See the help on
idgt()
for the precise definition of the inverse Gabor transform.- Additional parameters
idgtreal
optionnaly takes a pt arguments that can take the following values:‘freqinv’ Compute a idgtreal
using a frequency-invariant phase. This is the default convention described in the help fordgt()
.‘timeinv’ Compute a idgtreal
using a time-invariant phase. This convention is typically used in filter bank algorithms.Examples
The following example demonstrates the basic principles for getting perfect reconstruction (short version):
>>> from ltfatpy import greasy >>> from ltfatpy import dgtreal >>> f = greasy()[0] # Input test signal >>> a = 32 # time shift >>> M = 64 # frequency shift >>> gs = {'name': 'blackman', 'M': 128} # synthesis window >>> # analysis window >>> ga = {'name' : ('dual', gs['name']), 'M' : gs['M']} >>> (c, Ls) = dgtreal(f, ga, a, M)[0:2] # analysis >>> r = idgtreal(c, gs, a, M, Ls)[0] # synthesis >>> np.linalg.norm(f-r) < 1e-10 # test True
dgtlength¶
Module of dgtlength calculation
Ported from ltfat_2.1.0/gabor/dgtlength.m
-
ltfatpy.gabor.dgtlength.
dgtlength
(Ls, a, M)[source]¶ DGT length from signal
Usage:
L = dgtlength(Ls, a, M)
Input parameters:
Parameters: - Output parameters:
Returns: Corrected signal length Return type: int dgtlength(Ls, a, M)
returns the length of a Gabor system that is long enough to expand a signal of length Ls. Please see the help ondgt()
for an explanation of the parameters a and M.See also
gabwin¶
Module of Gabor window calculation
Ported from ltfat_2.1.0/gabor/gabwin.m
-
ltfatpy.gabor.gabwin.
gabwin
(g, a, M, L=None)[source]¶ Computes a Gabor window
Usage:
(g, info) = gabwin(g, a, M, L)
Input parameters:
Parameters: - Output parameters:
Returns: (g, info)
Return type: Variables: - g (numpy.ndarray) – the computed gabor window
- info (dict) – the information dictionary (see before)
With the transform length parameter L specified, it computes a window that fits well with the specified number of channels M, time shift a and transform length L. The window itself is specified by a text description or a dictionary containing additional parameters.
The window can be specified directly as a vector of numerical values. In this case,
gabwin
only checks assumptions about transform sizes etc.Without the transform length parameter L, it does the same, but the window must be a FIR window.
The window can be specified as one of the following text strings:
‘gauss’ Gaussian window fitted to the lattice, i.e. . ‘dualgauss’ Canonical dual of Gaussian window. ‘tight’ Tight window generated from a Gaussian. In these cases, a long window is generated with a length of L.
It is also possible to specify one of the window names from
firwin()
. In such a case,gabwin()
will generate the specified FIR window with a length of M.The window can also be specified as a dictionary. The possibilities are:
- {‘name’: ‘gauss’, ‘tfr’: 1.0, ‘fs’: 0.0, ‘width’: 0.0,
- ‘bw’: 0.0, ‘c_f’: 0.0, ‘centering’: ‘hp’, ‘delay’: 0.0, ‘norm’: 2}:
Additional parameters are passed to
pgauss()
.
- {‘name’: (‘dual’,…), …}:
- Canonical dual window of whatever follows. See the examples below.
- {‘name’: (‘tight’,…), …}:
- Canonical tight window of whatever follows.
It is also possible to specify one of the window names from
firwin()
as a string. In this case, the remaining entries of the cell array are passed directly tofirwin()
.Some examples:
To compute a Gaussian window of length L fitted for a system with time-shift a and M channels use:
g = gabwin('gauss', a, M, L)[0]
To compute Gaussian window with equal time and frequency support irrespective of a and M:
g = gabwin({'name': 'gauss', 'tfr': 1}, a, M, L)[0]
To compute the canonical dual of a Gaussian window fitted for a system with time-shift a and M channels:
gd = gabwin('gaussdual', a, M, L)[0]
To compute the canonical tight window of the Gaussian window fitted for the system:
gd = gabwin({'name': ('tight','gauss')}, a, M, L)[0]
To compute the dual of a Hann window of length 20:
g = gabwin({'name': ('dual', 'hann'), 'M': 20}, a, M, L)[0]
The returned info dictionary provides some information about the computed window:
keys Values ‘gauss’ True if the window is a Gaussian. ‘tfr’ Time/frequency support ratio of the window. Set whenever it makes sense. ‘wasrow’ True if input was a row window ‘isfir’ True if input is an FIR window ‘isdual’ True if output is the dual window of the auxiliary window. ‘istight’ True if output is known to be a tight window. ‘auxinfo’ Info about auxiliary window. ‘gl’ Length of window.
Reconstructing windows¶
gabdual¶
Module of Canonical dual window calculation
Ported from ltfat_2.1.0/gabor/gabdual.m
-
ltfatpy.gabor.gabdual.
gabdual
(g, a, M, L=None)[source]¶ Canonical dual window of Gabor frame
Usage:
gd = gabdual(g, a, M)
gd = gabdual(g, a, M, L)
Input parameters:
Parameters: - Output parameters:
Returns: the canonical dual window Return type: numpy.ndarray gabdual(g, a, M)
computes the canonical dual window of the discrete Gabor frame with window g and parameters a, M.The window g may be a vector of numerical values, a text string or a dictionary.
If the length of g is equal to M, then the input window is assumed to be an FIR window. In this case, the canonical dual window also has length of M. Otherwise the smallest possible transform length is chosen as the window length.
gabdual(g, a, M, L)
returns a window that is the dual window for a system of length L. Unless the dual window is a FIR window, the dual window will have length L.If then the dual window of the Gabor Riesz sequence with window g and parameters a and M will be calculated.
- Example:
The following example shows the canonical dual window of the Gaussian window.
>>> import matplotlib.pyplot as plt >>> from ltfatpy import pgauss, gabdual >>> a = 20 >>> M = 30 >>> L = 300 >>> g = pgauss(L, a*M/L)[0] >>> gd = gabdual(g, a, M) >>> # Plot in the time-domain >>> _ = plt.plot(gd) >>> plt.show()
See also
gabtight¶
Module of canonical tight windows calculation
Ported from ltfat_2.1.0/gabor/gabtight.m
-
ltfatpy.gabor.gabtight.
gabtight
(g, a, M, L=None)[source]¶ Canonical tight window of Gabor frame
Usage:
gt = gabtight(None, a, M, L)
gt = gabtight(g, a, M)
gt = gabtight(g, a, M, L)
Input parameters:
Parameters: - Output parameters:
Returns: Canonical tight window Return type: numpy.ndarray gabtight(None, a, M, L)
computes a nice tight window of length L for a lattice with parameters a, M. The window is not an FIR window, meaning that it will only generate a tight system if the system length is equal to L.gabtight(g, a, M)
computes the canonical tight window of the Gabor frame with window g and parameters a, M.The window g may be a vector of numerical values, a text string or a dictionary. See the help of
gabwin()
for more details.If the length of g is equal to M, then the input window is assumed to be a FIR window. In this case, the canonical dual window also has length of M. Otherwise the smallest possible transform length is chosen as the window length.
gabtight(g, a, M, L)
returns a window that is tight for a system of length L. Unless the input window g is a FIR window, the returned tight window will have length L.If
a > M
then an orthonormal window of the Gabor Riesz sequence with window g and parameters a and M will be calculated.Examples:
The following example shows the canonical tight window of the Gaussian window. This is calculated by default by
gabtight()
if no window is specified:>>> import matplotlib.pyplot as plt >>> from ltfatpy import gabtight >>> a = 20 >>> M = 30 >>> L = 300 >>> gt = gabtight(None, a, M, L) >>> # Plot in the time-domain >>> _ = plt.plot(gt) >>> plt.show()
See also
Conditions numbers¶
gabframediag¶
Module that computes Diagonal of Gabor frame operator
Ported from ltfat_2.1.0/gabor/gabframediag.m
-
ltfatpy.gabor.gabframediag.
gabframediag
(g, a, M, L)[source]¶ Diagonal of Gabor frame operator
Usage:
d = gabframediag(g, a, M, L)
Input parameters:
Parameters: - g (numpy.ndarray) – Window function
- a (int) – Length of time shift
- M (int) – Number of channels
- L (int) – Length of transform to do
- Output parameters:
Returns: Diagonal stored as a column vector Return type: numpy.ndarray gabframediag(g, a, M, L)
computes the diagonal of the Gabor frame operator with respect to the window g and parameters a and M. The diagonal is stored a as column vector of length L.The diagonal of the frame operator can for instance be used as a preconditioner.
See also
Phase gradient methods and reassignment¶
gabphasegrad¶
Module of phase gradient computation
Ported from ltfat_2.1.0/gabor/gabphasegrad.m
-
ltfatpy.gabor.gabphasegrad.
gabphasegrad
(method, *args, **kwargs)[source]¶ Phase gradient of the discrete Gabor transform
Usage:
(tgrad, fgrad, c) = gabphasegrad('dgt', f, g, a, M, L=None)
(tgrad, fgrad) = gabphasegrad('phase', cphase, a)
(tgrad, fgrad) = gabphasegrad('abs', s, g, a, difforder=4)
Input parameters:
Parameters: - method (str) – Method used to compute the phase gradient, see the possible values below
- f (numpy.ndarray) – (defined if
method='dgt'
) Input signal - cphase (numpy.ndarray) – (defined if
method='phase'
) Phase of adgt()
of the signal - s (numpy.ndarray) – (defined if
method='abs'
) Spectrogram of the signal - g (numpy.ndarray) – (defined if
method='dgt'
ormethod='phase'
) Window function - a (int) – (defined if
method='dgt'
ormethod='phase'
ormethod='abs'
) Length of time shift - M (int) – (defined if
method='dgt'
) Number of channels - L (int) – (defined if
method='dgt'
, optional) Length of transform to do - difforder (int) – (defined if
method='abs'
, optional) Order of the centered finite difference scheme used to perform the needed numerical differentiation
- Output parameters:
Returns: (tgrad, fgrad, c)
ifmethod='dgt'
, or(tgrad, fgrad)
ifmethod='phase'
ormethod='abs'
Return type: Variables: - tgrad (numpy.ndarray) – Instantaneous frequency
- fgrad (numpy.ndarray) – Local group delay
- c (numpy.ndarray) – Gabor coefficients
gabphasegrad
computes the time-frequency gradient of the phase of thedgt()
of a signal. The derivative in time tgrad is the instantaneous frequency while the frequency derivative fgrad is the local group delay.tgrad and fgrad measure the deviation from the current time and frequency, so a value of zero means that the instantaneous frequency is equal to the center frequency of the considered channel.
tgrad is scaled such that distances are measured in samples. Similarly, fgrad is scaled such that the Nyquist frequency (the highest possible frequency) corresponds to a value of
L/2
.The computation of tgrad and fgrad is inaccurate when the absolute value of the Gabor coefficients is low. This is due to the fact the the phase of complex numbers close to the machine precision is almost random. Therefore, tgrad and fgrad may attain very large random values when
abs(c)
is close to zero.The computation can be done using three different methods:
(tgrad, fgrad, c) = gabphasegrad('dgt', f, g, a, M)
computes the time-frequency gradient using adgt()
of the signal f. Thedgt()
is computed using the window g on the lattice specified by the time shift a and the number of channels M. The algorithm used to perform this calculation computes several DGTs, and therefore this routine takes the exact same input parameters asdgt()
.The window g may be specified as in
dgt()
. If the window used is'gauss'
, the computation will be done by a faster algorithm.(tgrad, fgrad, c) = gabphasegrad('dgt', f, g, a, M)
additionally returns the Gabor coefficientsc
, as they are always computed as a byproduct of the algorithm.(tgrad, fgrad) = gabphasegrad('phase', cphase, a)
computes the phase gradient from the phase cphase of adgt()
of the signal. The originaldgt()
from which the phase is obtained must have been computed using a time-shift of a.(tgrad, fgrad) = gabphasegrad('abs', s, g, a)
computes the phase gradient from the spectrogram s. The spectrogram must have been computed using the window g and time-shift a.(tgrad, fgrad) = gabphasegrad('abs', s, g, a, difforder=ord)
uses a centered finite difference scheme of orderord
to perform the needed numerical differentiation. Default is to use a 4th order scheme.Currently the ‘abs’ method only works if the window g is a Gaussian window specified as a string or cell array.
See also
resgram()
,gabreassign()
,dgt()
Phase reconstruction¶
Phase conversions¶
phaselock¶
Module of phaselocking
Ported from ltfat_2.1.0/gabor/phaselock.m
-
ltfatpy.gabor.phaselock.
phaselock
(c, a)[source]¶ Phaselock Gabor coefficients
Usage:
c_out = phaselock(c, a)
Input parameters:
Parameters: - c (numpy.ndarray) – non-phaselocked Gabor coefficients
- a (int) – Length of time shift
- Output parameters:
Returns: phaselocked Gabor coefficients Return type: numpy.ndarray phaselock(c, a)
phaselocks the Gabor coefficients c. The coefficients must have been obtained from adgt()
with parameter a.Phaselocking the coefficients modifies them so as if they were obtained from a time-invariant Gabor system. A filter bank produces phase locked coefficients.
Phaselocking of Gabor coefficients corresponds to the following transform:
Consider a signal
f
of lengthL
and defineN = L/a
.The output from
c = phaselock(dgt(f, g, a, M), a)
is given bywhere
m = 0,..., M-1
andn = 0,..., N-1
andl-a*n
are computed moduloL
.See also
dgt()
,phaseunlock()
,symphase()
- References:
- [Puc95]
phaseunlock¶
Module of phaseunlocking
Ported from ltfat_2.1.0/gabor/phaseunlock.m
-
ltfatpy.gabor.phaseunlock.
phaseunlock
(c, a)[source]¶ Undo phase lock of Gabor coefficients
Usage:
c_out = phaseunlock(c, a)
Input parameters:
Parameters: - c (numpy.ndarray) – phaselocked Gabor coefficients
- a (int) – Length of time shift
- Output parameters:
Returns: non-phaselocked Gabor coefficients Return type: numpy.ndarray phaseunlock(c, a)
removes phase locking from the Gabor coefficients c. The coefficient must have been obtained from adgt()
with parameter a.Phase locking the coefficients modifies them so as if they were obtained from a time-invariant Gabor system. A filter bank produces phase locked coefficients.
See also
dgt()
,phaselock()
,symphase()
- References:
- [Puc95]
Plots¶
gabimagepars¶
Module to find Gabor parameters to generate image
Ported from ltfat_2.1.0/gabor/gabimagepars.m
-
ltfatpy.gabor.gabimagepars.
gabimagepars
(Ls, x, y)[source]¶ Find Gabor parameters to generate image
Usage:
(a, M, L, N, Ngood) = gabimagepars(Ls, x, y)
Input parameters:
Parameters: - Output parameters:
Returns: (a, M, L, N, Ngood)
Return type: Variables: - a (int) – Length of time shift
- M (int) – Number of frequency channels
- L (int) – Length of transform to do
- N (int) – Total number of time steps
- Ngood (int) – Number of time steps (columns in the coefficients
matrix) that contain relevant information. The columns from
Ngood-1
untilN-1
only contain information from a zero-extension of the signal.
(a, M, L, N, Ngood) = gabimagepars(Ls, x, y)
will compute a reasonable set of parameters a, M and L to produce a nice Gabor ‘image’ of a signal of length Ls.If you use this function to calculate a grid size for analysis of a real-valued signal (using
dgtreal()
), please input twice of the desired size y. This is becausedgtreal()
only returns half as many coefficients in the frequency direction asdgt()
.For this function to work properly, the specified numbers for x and y must not be large prime numbers.
Example:
We wish to compute a Gabor image of a real valued signal
f
of length 7500. The image should have an approximate resolution of 600 x 800 pixels:>>> from matplotlib.pyplot import show >>> from ltfatpy import linus, gabimagepars, dgtreal, plotdgtreal >>> f, fs = linus() >>> f = f[4000:4000+7500] >>> a, M, L, N, Ngood = gabimagepars(7500, 800, 2*600) >>> c = dgtreal(f, 'gauss', a, M)[0] >>> _ = plotdgtreal(c, a, M, fs=fs, dynrange=90) >>> show()
The size of
c
is(M/2)+1 x N
equal to 601 x 700 pixels.
instfreqplot¶
Module of instantaneous frequency plot
Ported from ltfat_2.1.0/gabor/instfreqplot.m
-
ltfatpy.gabor.instfreqplot.
instfreqplot
(f, fs=None, tfr=1.0, wlen=None, nf=None, thr=None, fmax=None, xres=800, yres=600, method='dgt', climsym=None, **kwargs)[source]¶ Plot of instantaneous frequency
- Input parameters:
Parameters: - f (numpy.ndarray) – Analyzed signal
- fs (float) – Sampling rate in Hz of the analyzed signal
- tfr (float) – Set the ratio of frequency resolution to time resolution.
A value
tfr = 1.0
is the default. Settingtfr > 1.0
will give better frequency resolution at the expense of a worse time resolution. A value of0.0 < tfr < 1.0
will do the opposite. - wlen (int) – Window length. Specifies the length of the window
measured in samples. See help of
pgauss()
on the exact details of the window length (parameter width). - nf (bool) – If
True
, display negative frequencies, with the zero-frequency centered in the middle. For real signals, this will just mirror the upper half plane. This is standard for complex signals. - thr (float) – Keep the coefficients with a magnitude larger than thr times the largest magnitude. Set the instantaneous frequency of the rest of the coefficients to zero.
- fmax (float) – Display fmax as the highest frequency.
- xres (int) – Approximate number of pixels in the time direction
- yres (int) – Number of pixels in the frequency direction
- method (str) – Specify the method used for instantaneous frequency
computation. Possible values are ‘dgt’, ‘phase’ and ‘abs’. See the help
of
gabphasegrad()
for details. - climsym (float) – Use a colormap ranging from
-climsym
to+climsym
- **kwargs –
instfreqplot
supports all the optional parameters oftfplot()
. Please see the help oftfplot()
for an exhaustive list.
- Output parameters:
Returns: The instantaneous frequency data used in the plotting Return type: numpy.ndarray instfreqplot(f)
plots the instantaneous frequency of f using a Discrete Gabor Transform (dgt()
).The instantaneous frequency contains extreme spikes in regions where the spectrogram is close to zero. These points are usually uninteresting and destroy the visibility of the plot. Use the thr or clim or climsym options to remove these points.
An example:
>>> from ltfatpy import instfreqplot, greasy >>> from matplotlib.pyplot import show >>> _ = instfreqplot(greasy()[0], 16000., thr=.03, climsym=100.) >>> show()
will produce a nice instantaneous frequency plot of the
greasy()
signal.
See also
sgram()
,resgram()
,phaseplot()
phaseplot¶
Module of phase plot
Ported from ltfat_2.1.0/gabor/phaseplot.m
-
ltfatpy.gabor.phaseplot.
phaseplot
(f, fs=None, tfr=1.0, wlen=None, nf=None, thr=None, fmax=None, phase='timeinv', norm='2', **kwargs)[source]¶ Phase plot
- Input parameters:
Parameters: - f (numpy.ndarray) – Analyzed signal
- fs (float) – Sampling rate in Hz of the analyzed signal
- tfr (float) – Set the ratio of frequency resolution to time resolution.
A value
tfr = 1.0
is the default. Settingtfr > 1.0
will give better frequency resolution at the expense of a worse time resolution. A value of0.0 < tfr < 1.0
will do the opposite. - wlen (int) – Window length. Specifies the length of the window
measured in samples. See help of
pgauss()
on the exact details of the window length (parameter width). - nf (bool) – If
True
, display negative frequencies, with the zero-frequency centered in the middle. For real signals, this will just mirror the upper half plane. This is standard for complex signals. - thr (float) – Keep the coefficients with a magnitude larger than thr times the largest magnitude. Set the phase of the rest of the coefficients to zero. This is useful, because for small amplitude the phase values can be meaningless.
- fmax (float) – Display fmax as the highest frequency.
- phase (str) – Phase convention for the dgt computation. Use phase=’timeinv’ to display the phase as computed by a time-invariant dgt (this is the default) or phase=’freqinv’ to display the phase as computed by a frequency-invariant dgt.
- norm (string) – Window normalization. See parameter norm in the
help of
pgauss()
. - **kwargs –
phaseplot
supports all the optional parameters oftfplot()
. Please see the help oftfplot()
for an exhaustive list.
- Output parameters:
Returns: The phase data used in the plotting Return type: numpy.ndarray phaseplot(f)
plots the phase of f using a dgt.phaseplot
should only be used for short signals (shorter than the resolution of the screen), as there will otherwise be some visual aliasing, such that very fast changing areas will look very smooth.phaseplot
always calculates the phase of the full time/frequency plane (as opposed tosgram()
), and you therefore risk running out of memory for long signals.For the best result when using phaseplot, use a circulant color map, for instance hsv.
Examples:
The following code shows the phaseplot of a periodic, hyperbolic secant visualized using the hsv colormap:
>>> from matplotlib.pyplot import hsv, show >>> from ltfatpy import phaseplot, psech >>> hsv() >>> _ = phaseplot(psech(200)[0], tc=True, nf=True) >>> show()
The following phaseplot shows the phase of white, Gaussian noise:
>>> from numpy.random import randn >>> from matplotlib.pyplot import hsv, show >>> from ltfatpy import phaseplot >>> hsv() >>> _ = phaseplot(randn(200)) >>> show()
See also
- References:
- [CHTorresani98]
plotdgt¶
Module of dgt coefficients plotting
Ported from ltfat_2.1.0/gabor/plotdgt.m
-
ltfatpy.gabor.plotdgt.
plotdgt
(coef, a, **kwargs)[source]¶ Plot dgt coefficients
- Input parameters:
Parameters: - coef (numpy.ndarray) – Gabor coefficients
- a (int) – Length of time shift used when computing coef
- **kwargs –
plotdgt
supports all the optional parameters oftfplot()
. Please see the help oftfplot()
for an exhaustive list.
- Output parameters:
Returns: The processed image data used in the plotting. Inputting this data directly to matshow()
or similar functions will create the plot. This is useful for custom post-processing of the image data.Return type: numpy.ndarray plotdgt(coef, a)
plots the Gabor coefficients coef. The coefficients must have been produced with a time shift of a.The figure generated by this function places the zero-frequency in the center of the y-axis, with positive frequencies above and negative frequencies below.
See also
plotdgtreal¶
Module of dgtreal coefficients plotting
Ported from ltfat_2.1.0/gabor/plotdgtreal.m
-
ltfatpy.gabor.plotdgtreal.
plotdgtreal
(coef, a, M, **kwargs)[source]¶ Plot dgtreal coefficients
Parameters: - coef (numpy.ndarray) – Gabor coefficients
- a (int) – Length of time shift used when computing coef
- M (int) – Number of modulations used when computing coef
- **kwargs –
plotdgtreal
supports all the optional parameters oftfplot()
. Please see the help oftfplot()
for an exhaustive list.
- Output parameters:
Returns: The processed image data used in the plotting. Inputting this data directly to matshow()
or similar functions will create the plot. This is useful for custom post-processing of the image data.Return type: numpy.ndarray plotdgtreal(coef, a, M)
plots Gabor coefficients fromdgtreal()
. The parameters a and M must match those from the call todgtreal()
.
sgram¶
Module of spectrogram plotting
Ported from ltfat_2.1.0/gabor/sgram.m
-
ltfatpy.gabor.sgram.
sgram
(f, fs=None, tfr=1.0, wlen=None, nf=None, thr=None, fmax=None, xres=800, yres=600, norm='2', **kwargs)[source]¶ Spectrogram
- Input:
Parameters: - f (numpy.ndarray) – Analyzed signal
- fs (float) – Sampling rate in Hz of the analyzed signal
- tfr (float) – Set the ratio of frequency resolution to time resolution.
A value
tfr = 1.0
is the default. Settingtfr > 1.0
will give better frequency resolution at the expense of a worse time resolution. A value of0.0 < tfr < 1.0
will do the opposite. - wlen (int) – Window length. Specifies the length of the window
measured in samples. See help of
pgauss()
on the exact details of the window length (parameter width). - nf (bool) – If
True
, display negative frequencies, with the zero-frequency centered in the middle. For real signals, this will just mirror the upper half plane. This is standard for complex signals. - thr (float) – Keep only the largest fraction thr of the coefficients, and set the rest to zero.
- fmax (float) – Display fmax as the highest frequency. Default value
of
None
means to use the Nyquist frequency. - xres (int) – Approximate number of pixels along x-axis / time
- yres (int) – Approximate number of pixels along y-axis / frequency
- norm (string) – Window normalization. See parameter norm in the
help of
pgauss()
. - **kwargs –
sgram
supports all the optional parameters oftfplot()
. Please see the help oftfplot()
for an exhaustive list.
- Output parameters:
Returns: The image to be displayed as a matrix. Use this in conjunction with imwrite etc. These coefficients are only intended to be used by post-processing image tools. Numerical Gabor signal analysis and synthesis should always be done using the dgt()
,idgt()
,dgtreal()
andidgtreal()
functions.Return type: numpy.ndarray sgram(f)
plots a spectrogram of f using a Discrete Gabor Transform (dgt()
).Examples:
The
greasy()
signal is sampled using a sampling rate of 16 kHz. To display a spectrogram ofgreasy()
with a dynamic range of 90 dB, use:>>> from ltfatpy import sgram, greasy >>> from matplotlib.pyplot import show >>> _ = sgram(greasy()[0], 16000., dynrange=90) >>> show()
To create a spectrogram with a window length of 20 ms (which is typically used in speech analysis) use :
>>> from ltfatpy import sgram, greasy >>> from matplotlib.pyplot import show >>> fs = 16000. >>> _ = sgram(greasy()[0], fs, dynrange=90, wlen=round(20./1000.*fs)) >>> show()
tfplot¶
Module of time-frequency plotting
Ported from ltfat_2.1.0/gabor/tfplot.m
-
ltfatpy.gabor.tfplot.
tfplot
(coef, step, yr, fs=None, dynrange=None, normalization='db', tc=False, clim=None, plottype='image', colorbar=True, display=True, time='Time', frequency='Frequency', samples='samples', normalized='normalized')[source]¶ Plot coefficient matrix on the TF plane
- Input parameters:
Parameters: - coef (numpy.ndarray) – 2D coefficient array
- step (float) – Shift in samples between each column of coefficients
- yr (numpy.ndarray) – 2 elements vector containing the lowest and highest normalized frequency
- fs (float) – Sampling rate in Hz of the original signal
- dynrange (float) – Limit the dynamical range to dynrange by using a colormap in the interval [chigh-dynrange, chigh], where chigh is the highest value in the plot. The default value of None means to not limit the dynamical range. If both clim and dynrange are specified, then clim takes precedence.
- normalization (str) – String specifying the normalization of the plot, possible values are listed below
- tc (bool) – Time centering: if
True
, move the beginning of the signal to the middle of the plot. This is usefull for visualizing the window functions of the toolbox. - clim (tuple) – Use a colormap ranging from clim[0] to clim[1]. If both clim and dynrange are specified, then clim takes precedence.
- plottype (str) – String specifying the type of plot, possible values are listed below
- colorbar (bool) – If
True
, display the colorbar (this is the default) - display (bool) – If
True
, display the figure (this is the default). Usingdisplay=False
to avoid displaying the figure is usefull if you only want to obtain the output for further processing. - time (str) – Text customization: the word denoting time
- frequency (str) – Text customization: the word denoting frequency
- samples (str) – Text customization: the word denoting samples
- normalized (str) – Text customization: the word denoting normalized
- Output parameter:
Returns: The processed image data used in the plotting. Inputting this data directly to matshow()
or similar functions will create the plot. This is usefull for custom post-processing of the image data.Return type: numpy.ndarray tfplot(coef, step, yr)
will plot a rectangular coefficient array on the TF-plane.tfplot
is not meant to be called directly. Instead, it is called by other plotting routines to give a uniform display format.- Possible values for normalization:
'db'
Apply to the coefficients. This makes it possible to see very weak phenomena, but it might show too much noise. A logarithmic scale is more adapted to perception of sound. This is the default. 'dbsq'
Apply to the coefficients. Same as the 'db'
option, but assume that the input is already squared.'lin'
Show the coefficients on a linear scale. This will display the raw input without any modifications. Only works for real-valued input. 'linsq'
Show the square of the coefficients on a linear scale. 'linabs'
Show the absolute value of the coefficients on a linear scale. - Possible values for plottype:
'image'
Use imshow to display the plot. This is the default. 'contour'
Do a contour plot. 'surf'
Do a surface plot. 'pcolor'
Do a pcolor plot.
See also
sgram()
,plotdgt()
,plotdgtreal()
,plotwmdct()
,plotdwilt()