Front page|Spectrum - Spectral Analysis in Python (0.5.2)

5.1. Fourier Methods

5.1.1. Power Spectrum Density based on Fourier Spectrum

This module provides Periodograms (classics, daniell, bartlett)

Periodogram(data[, sampling, window, NFFT, ...]) The Periodogram class provides an interface to periodogram PSDs
DaniellPeriodogram(data, P[, NFFT, detrend, ...]) Return Daniell’s periodogram.
speriodogram(x[, NFFT, detrend, sampling, ...]) Simple periodogram, but matrices accepted.
WelchPeriodogram(data[, NFFT, sampling]) Simple periodogram wrapper of numpy.psd function.
speriodogram(x[, NFFT, detrend, sampling, ...]) Simple periodogram, but matrices accepted.

Code author: Thomas Cokelaer 2011

References:See [Marple]

Usage

You can compute a periodogram using speriodogram():

from spectrum import *
from pylab import plot
p = speriodogram(marple_data)
plot(p)

However, the output is not always easy to manipulate or plot, therefore it is advised to use the class Periodogram instead:

from spectrum import *
p = Periodogram(marple_data)
p()
p.plot()

This class will take care of the plotting and internal state of the computation. For instance, if you can change the output easily:

p.plot(sides='twosided')
class pdaniell(data, P, sampling=1.0, window='hann', NFFT=None, scale_by_freq=True, detrend=None)

Bases: spectrum.psd.FourierSpectrum

The pdaniell class provides an interface to DaniellPeriodogram

from spectrum import *
data = data_cosine(N=4096, sampling=4096)
p = pdaniell(data, 8, NFFT=4096)
p()
p.plot()

pdaniell Constructor

Parameters:
  • data (array) – input data (list or numpy.array)
  • P (int) – number of neighbours to average over.
  • sampling (float) – sampling frequency of the input data.
  • window (str) – a tapering window. See Window.
  • NFFT (int) – total length of the final data sets (padded with zero if needed; default is 4096)
  • scale_by_freq (bool) –
  • detrend (str) –
speriodogram(x, NFFT=None, detrend=True, sampling=1.0, scale_by_freq=True, window='hamming', axis=0)

Simple periodogram, but matrices accepted.

Parameters:
  • x – an array or matrix of data samples.
  • NFFT – length of the data before FFT is computed (zero padding)
  • detrend (bool) – detrend the data before co,puteing the FFT
  • sampling (float) – sampling frequency of the input data.
  • scale_by_freq
  • window (str) –
Returns:

2-sided PSD if complex data, 1-sided if real.

if a matrix is provided (using numpy.matrix), then a periodogram is computed for each row. The returned matrix has the same shape as the input matrix.

The mean of the input data is also removed from the data before computing the psd.

from pylab import *
from spectrum import *
data = data_cosine(N=1024, A=0.1, sampling=1024, freq=200)
semilogy(speriodogram(data, detrend=False, sampling=1024), marker='o')
grid(True)

[hires.png, pdf]

../_images/b8d1c9dd3f.png
from spectrum import *
from pylab import *
# create N data sets and make the frequency dependent on the time
N = 100
m = numpy.concatenate([data_cosine(N=1024, A=0.1, sampling=1024, freq=x) for x in range(1, N)]);
m.resize(N, 1024)
res = speriodogram(m)
figure(1)
semilogy(res)
figure(2)
imshow(res.transpose(), aspect='auto')

[hires.png, pdf]

../_images/51e6805f41_00.png

[hires.png, pdf]

../_images/51e6805f41_01.png

Todo

a proper spectrogram class/function that takes care of normalisation

class Periodogram(data, sampling=1.0, window='hann', NFFT=None, scale_by_freq=False, detrend=None)

Bases: spectrum.psd.FourierSpectrum

The Periodogram class provides an interface to periodogram PSDs

from spectrum import *
data = data_cosine(N=1024, A=0.1, sampling=1024, freq=200)
p = Periodogram(data, sampling=1024)
p()
p.plot(marker='o')

[hires.png, pdf]

../_images/cb48a2be44.png

Periodogram Constructor

Parameters:
  • data (array) – input data (list or numpy.array)
  • sampling (float) – sampling frequency of the input data.
  • window (str) – a tapering window. See Window.
  • NFFT (int) – total length of the final data sets (padded with zero if needed; default is 4096)
  • scale_by_freq (bool) –
  • detrend (str) –
WelchPeriodogram(data, NFFT=None, sampling=1.0, **kargs)

Simple periodogram wrapper of numpy.psd function.

Parameters:
  • A – the input data
  • NFFT (int) – total length of the final data sets (padded with zero if needed; default is 4096)
  • window (str) –
Technical documentation:
 

When we calculate the periodogram of a set of data we get an estimation of the spectral density. In fact as we use a Fourier transform and a truncated segments the spectrum is the convolution of the data with a rectangular window which Fourier transform is

W(s)= \frac{1}{N^2} \left[ \frac{\sin(\pi s)}{\sin(\pi s/N)} \right]^2

Thus oscillations and sidelobes appears around the main frequency. One aim of t he tapering is to reduced this effects. We multiply data by a window whose sidelobes are much smaller than the main lobe. Classical window is hanning window. But other windows are available. However we must take into account this energy and divide the spectrum by energy of taper used. Thus periodogram becomes :

D_k \equiv \sum_{j=0}^{N-1}c_jw_j \; e^{2\pi ijk/N}  \qquad k=0,...,N-1

P(0)=P(f_0)=\frac{1}{2\pi W_{ss}}\arrowvert{D_0}\arrowvert^2

P(f_k)=\frac{1}{2\pi W_{ss}} \left[\arrowvert{D_k}\arrowvert^2+\arrowvert{D_{N-k}}\arrowvert^2\right]        \qquad k=0,1,...,     \left( \frac{1}{2}-1 \right)

P(f_c)=P(f_{N/2})= \frac{1}{2\pi W_{ss}} \arrowvert{D_{N/2}}\arrowvert^2

with

{W_{ss}} \equiv N\sum_{j=0}^{N-1}w_j^2

from spectrum import *
psd = WelchPeriodogram(marple_data, 256)

[hires.png, pdf]

../_images/6dd3b03194.png
DaniellPeriodogram(data, P, NFFT=None, detrend='mean', sampling=1.0, scale_by_freq=True, window='hamming')

Return Daniell’s periodogram.

To reduce fast fluctuations of the spectrum one idea proposed by daniell is to average each value with points in its neighboorhood. It’s like a low filter.

\hat{P}_D[f_i]= \frac{1}{2P+1} \sum_{n=i-P}^{i+P} \tilde{P}_{xx}[f_n]

where P is the number of points to average.

Daniell’s periodogram is the convolution of the spectrum with a low filter:

\hat{P}_D(f)=   \hat{P}_{xx}(f)*H(f)

Example:

>>> DaniellPeriodogram(data,8)

if N/P is not integer, the final values of the original PSD are not used.

using DaniellPeriodogram(data, 0) should give the original PSD.

Correlogram PSD estimates

This module provides Correlograms methods

CORRELOGRAMPSD(X[, Y, lag, window, norm, ...]) PSD estimate using correlogram method.
pcorrelogram(data[, sampling, lag, window, ...]) The Correlogram class provides an interface to CORRELOGRAMPSD().

Code author: Thomas Cokelaer 2011

References:See [Marple]
CORRELOGRAMPSD(X, Y=None, lag=-1, window='hamming', norm='unbiased', NFFT=4096, window_params={}, correlation_method='xcorr')

PSD estimate using correlogram method.

Parameters:
  • X (array) – complex or real data samples X(1) to X(N)
  • Y (array) – complex data samples Y(1) to Y(N). If provided, computes the cross PSD, otherwise the PSD is returned
  • lag (int) – highest lag index to compute. Must be less than N
  • window_name (str) – see window for list of valid names
  • norm (str) – one of the valid normalisation of xcorr() (biased, unbiased, coeff, None)
  • NFFT (int) – total length of the final data sets (padded with zero if needed; default is 4096)
  • correlation_method (str) – either xcorr or CORRELATION. CORRELATION should be removed in the future.
Returns:

  • Array of real (cross) power spectral density estimate values. This is a two sided array with negative values following the positive ones whatever is the input data (real or complex).

Description:

The exact power spectral density is the Fourier transform of the autocorrelation sequence:

P_{xx}(f) = T \sum_{m=-\infty}^{\infty} r_{xx}[m] exp^{-j2\pi fmT}

The correlogram method of PSD estimation substitutes a finite sequence of autocorrelation estimates \hat{r}_{xx} in place of r_{xx}. This estimation can be computed with xcorr() or CORRELATION() by chosing a proprer lag L. The estimated PSD is then

\hat{P}_{xx}(f) = T \sum_{m=-L}^{L} \hat{r}_{xx}[m] exp^{-j2\pi fmT}

The lag index must be less than the number of data samples N. Ideally, it should be around L/10 [Marple] so as to avoid greater statistical variance associated with higher lags.

To reduce the leakage of the implicit rectangular window and therefore to reduce the bias in the estimate, a tapering window is normally used and lead to the so-called Blackman and Tukey correlogram:

\hat{P}_{BT}(f) = T \sum_{m=-L}^{L} w[m] \hat{r}_{xx}[m] exp^{-j2\pi fmT}

The correlogram for the cross power spectral estimate is

\hat{P}_{xx}(f) = T \sum_{m=-L}^{L} \hat{r}_{xx}[m] exp^{-j2\pi fmT}

which is computed if Y is not provide. In such case, r_{yx} = r_{xy} so we compute the correlation only once.

from spectrum import *
from pylab import *

psd = CORRELOGRAMPSD(marple_data, marple_data, lag=15)
f = linspace(-0.5, 0.5, len(psd))
psd = cshift(psd, len(psd)/2)
plot(f, 10*log10(psd/max(psd)))
axis([-0.5,0.5,-50,0])
grid(True)

[hires.png, pdf]

../_images/22de69f7c3.png

See also

create_window(), CORRELATION(), xcorr(), pcorrelogram.

class pcorrelogram(data, sampling=1.0, lag=-1, window='hamming', NFFT=None, scale_by_freq=True, detrend=None)

Bases: spectrum.psd.FourierSpectrum

The Correlogram class provides an interface to CORRELOGRAMPSD().

It returns an object that inherits from FourierSpectrum and therefore ease the manipulation of PSDs.

from spectrum import *
p = pcorrelogram(data_cosine(N=1024), lag=15)
p()
p.plot()
p.plot(sides='twosided')

[hires.png, pdf]

../_images/c25b524297.png

Correlogram Constructor

Parameters:
  • data (array) – input data (list or numpy.array)
  • sampling (float) – sampling frequency of the input data.
  • lag (int) –
  • window (str) – a tapering window. See Window.
  • NFFT (int) – total length of the final data sets (padded with zero if needed; default is 4096)
  • scale_by_freq (bool) –
  • detrend (str) –

5.1.2. Tapering Windows

This module contains tapering windows utilities.

Window(N[, name, norm]) Window tapering object
window_visu([N, name]) A Window visualisation tool
create_window(N[, name]) Returns the N-point window given a valid name
window_hann(N) Hann window (or Hanning). (wrapping of numpy.bartlett)
window_hamming(N) Hamming window
window_bartlett(N) Bartlett window (wrapping of numpy.bartlett) also known as Fejer
window_bartlett_hann(N) Bartlett-Hann window
window_blackman(N[, alpha]) Blackman window
window_blackman_harris(N) Blackman Harris window
window_blackman_nuttall(N) Blackman Nuttall window
window_bohman(N) Bohman tapering window
window_chebwin(N[, attenuation]) Cheb window
window_cosine(N) Cosine tapering window also known as sine window.
window_flattop(N[, mode, precision]) Flat-top tapering window
window_gaussian(N[, alpha]) Gaussian window
window_hamming(N) Hamming window
window_hann(N) Hann window (or Hanning). (wrapping of numpy.bartlett)
window_kaiser(N[, beta, method]) Kaiser window
window_lanczos(N) Lanczos window also known as sinc window.
window_nuttall(N) Nuttall tapering window
window_parzen(N) Parsen tapering window (also known as de la Valle-Poussin)
window_tukey(N[, r]) Tukey tapering window (or cosine-tapered window)

Code author: Thomas Cokelaer 2011

References:See [Nuttall], [Marple], [Harris]
class Window(N, name=None, norm=True, **kargs)

Window tapering object

This class provides utilities to manipulate tapering windows. Plotting functions allows to visualise the time and frequency response. It is also possible to retrieve relevant quantities such as the equivalent noise band width.

The following examples illustrates the usage. First, we create the window by providing a name and a size:

from spectrum import *
w = Window(64, 'hamming')

The window has been computed and the data is stored in:

w.data

This object contains plotting methods so that you can see the time or frequency response.

from spectrum.window import Window
w = Window(64, 'hamming')
w.plot_frequencies()

[hires.png, pdf]

../_images/f7d03e28f1.png

Some windows may accept optional arguments. For instance, window_blackman() accepts an optional argument called \alpha as well as Window. Indeed, we use the factory create_window(), which manage all the optional arguments. So you can write:

w = Window(64, 'blackman', alpha=1)

See also

create_window().

Constructor:

Create a tapering window object

Parameters:
  • N – the window length
  • name – the type of window, e.g., ‘Hann’
  • norm – normalise the window in frequency domain (for plotting)
  • kargs – any of create_window() valid optional arguments.

Attributes:

  • data: time series data
  • frequencies: getter to the frequency series
  • response: getter to the PSD
  • enbw: getter to the Equivalent noise band width.
N

Getter for the window length

compute_response(**kargs)

Compute the window data frequency response

Parameters:
  • norm – True by default. normalised the frequency data.
  • NFFT (int) – total length of the final data sets( 2048 by default. if less than data length, then NFFT is set to the data length*2).

The response is stored in response.

Note

Units are dB (20 log10) since we plot the frequency response)

data

Getter for the window values (in time)

enbw

getter for the equivalent noise band width. See enbw() function

frequencies

Getter for the frequency array

info()

Print object information such as length and name

mean_square

returns :math:` rac{w^2}{N}`

name

Getter for the window name

norm

Getter of the normalisation flag (True by default)

plot_frequencies(mindB=None, maxdB=None, norm=True)

Plot the window in the frequency domain

Parameters:
  • mindB – change the default lower y bound
  • maxdB – change the default upper lower bound
  • norm (bool) – if True, normalise the frequency response.
from spectrum.window import Window
w = Window(64, name='hamming')
w.plot_frequencies()

[hires.png, pdf]

../_images/e06573be1a.png
plot_time_freq(mindB=-100, maxdB=None, norm=True)

Plotting method to plot both time and frequency domain results.

See plot_frequencies() for the optional arguments.

from spectrum.window import Window
w = Window(64, name='hamming')
w.plot_time_freq()

[hires.png, pdf]

../_images/e15669a00d.png
plot_window()

Plot the window in the time domain

from spectrum.window import Window
w = Window(64, name='hamming')
w.plot_window()

[hires.png, pdf]

../_images/3b4c6c193b.png
response

Getter for the frequency response. See compute_response()

create_window(N, name=None, **kargs)

Returns the N-point window given a valid name

Parameters:

The following windows have been simply wrapped from existing librairies like NumPy:

The following windows have been implemented from scratch:

Todo

on request taylor, potter, Bessel, expo, rife-vincent, Kaiser-Bessel derived (KBD)

from pylab import plot, hold, legend
from spectrum import create_window

data = create_window(51, 'hamming')
plot(data, label='hamming')
hold(True)
data = create_window(51, 'kaiser')
plot(data, label='kaiser')
legend()

[hires.png, pdf]

../_images/a9a671fe04.png
from pylab import *
from spectrum import create_window

A = fft(create_window(51, 'hamming'), 2048) / 25.5
mag = abs(fftshift(A))
freq = linspace(-0.5,0.5,len(A))
response = 20*log10(mag)
mindB = -60
response = clip(response,mindB,100)
plot(freq, response)

[hires.png, pdf]

../_images/05a10f8316.png

See also

window_visu(), Window(), spectrum.dpss

enbw(data)

Computes the equivalent noise bandwidth

ENBW = N \frac{\sum_{n=1}^{N} w_n^2}{\left(\sum_{n=1}^{N} w_n \right)^2}

>>> from spectrum import *
>>> w = create_window(64, 'rectangular')

>>> enbw(w)
1.0

The following table contains the ENBW values for some of the implemented windows in this module (with N=16384). They have been double checked against litterature (Source: [Harris], [Marple]).

If not present, it means that it has not been checked.

name ENBW litterature
rectangular
triangle 1.3334 1.33
Hann 1.5001 1.5
Hamming 1.3629 1.36
blackman 1.7268 1.73
kaiser 1.7  
blackmanharris,4 2.004
riesz 1.2000 1.2
riemann 1.32 1.3
parzen 1.917 1.92
tukey 0.25 1.102 1.1
bohman 1.7858 1.79
poisson 2 1.3130 1.3
hanningpoisson 0.5 1.609 1.61
cauchy 1.489 1.48
lanczos 1.3  
window_bartlett(N)

Bartlett window (wrapping of numpy.bartlett) also known as Fejer

Parameters:N (int) – window length

The Bartlett window is defined as

w(n) = \frac{2}{N-1} \left(
\frac{N-1}{2} - \left|n - \frac{N-1}{2}\right|
\right)

from spectrum import window_visu
window_visu(64, 'bartlett')

[hires.png, pdf]

../_images/ce3175695d.png

See also

numpy.bartlett, create_window(), Window.

window_bartlett_hann(N)

Bartlett-Hann window

Parameters:N – window length

w(n) = a_0 + a_1 \left| \frac{n}{N-1} -\frac{1}{2}\right| - a_2 \cos \left( \frac{2\pi n}{N-1} \right)

with a_0 = 0.62, a_1 = 0.48 and a_2=0.38

from spectrum import window_visu
window_visu(64, 'bartlett_hann')

[hires.png, pdf]

../_images/6c7ab7040c.png
window_blackman(N, alpha=0.16)

Blackman window

Parameters:N – window length

a_0 - a_1 \cos(\frac{2\pi n}{N-1}) +a_2 \cos(\frac{4\pi n }{N-1})

with

a_0 = (1-\alpha)/2, a_1=0.5, a_2=\alpha/2 \rm{\;and\; \alpha}=0.16

When \alpha=0.16, this is the unqualified Blackman window with a_0=0.48 and a_2=0.08.

from spectrum import window_visu
window_visu(64, 'blackman')

[hires.png, pdf]

../_images/b5ad039628.png

Note

Although Numpy implements a blackman window for \alpha=0.16, this implementation is valid for any \alpha.

See also

numpy.blackman, create_window(), Window

window_blackman_harris(N)

Blackman Harris window

Parameters:N – window length

w(n) = a_0 - a_1 \cos\left(\frac{2\pi n}{N-1}\right)+ a_2 \cos\left(\frac{4\pi n}{N-1}\right)- a_3 \cos\left(\frac{6\pi n}{N-1}\right)

coeff value
a_0 0.35875
a_1 0.48829
a_2 0.14128
a_3 0.01168
from spectrum import window_visu
window_visu(64, 'blackman_harris', mindB=-80)

[hires.png, pdf]

../_images/81699181db.png
window_blackman_nuttall(N)

Blackman Nuttall window

returns a minimum, 4-term Blackman-Harris window. The window is minimum in the sense that its maximum sidelobes are minimized. The coefficients for this window differ from the Blackman-Harris window coefficients and produce slightly lower sidelobes.

Parameters:N – window length

w(n) = a_0 - a_1 \cos\left(\frac{2\pi n}{N-1}\right)+ a_2 \cos\left(\frac{4\pi n}{N-1}\right)- a_3 \cos\left(\frac{6\pi n}{N-1}\right)

with a_0 = 0.3635819, a_1 = 0.4891775, a_2=0.1365995 and 0_3=.0106411

from spectrum import window_visu
window_visu(64, 'blackman_nuttall', mindB=-80)

[hires.png, pdf]

../_images/2676a9b85b.png
window_bohman(N)

Bohman tapering window

Parameters:N – window length

w(n) = (1-|x|) \cos (\pi |x|) + \frac{1}{\pi} \sin(\pi |x|)

where x is a length N vector of linearly spaced values between -1 and 1.

from spectrum import window_visu
window_visu(64, 'bohman')

[hires.png, pdf]

../_images/a69334c105.png
window_cauchy(N, alpha=3)

Cauchy tapering window

Parameters:
  • N (int) – window length
  • alpha (float) – parameter of the poisson window

w(n) = \frac{1}{1+\left(\frac{\alpha*n}{N/2}\right)**2}

from spectrum import window_visu
window_visu(64, 'cauchy', alpha=3)
window_visu(64, 'cauchy', alpha=4)
window_visu(64, 'cauchy', alpha=5)

[hires.png, pdf]

../_images/b37e4a5474.png
window_chebwin(N, attenuation=50)

Cheb window

Parameters:N – window length
from spectrum import window_visu
window_visu(64, 'chebwin', attenuation=50)

[hires.png, pdf]

../_images/57d2493897.png

See also

scipy.signal.chebwin, create_window(), Window

window_cosine(N)

Cosine tapering window also known as sine window.

Parameters:N – window length

w(n) = \cos\left(\frac{\pi n}{N-1} - \frac{\pi}{2}\right) = \sin \left(\frac{\pi n}{N-1}\right)

from spectrum import window_visu
window_visu(64, 'cosine')

[hires.png, pdf]

../_images/a7cdb67406.png
window_flattop(N, mode='symmetric', precision=None)

Flat-top tapering window

Returns symmetric or periodic flat top window.

Parameters:
  • N – window length
  • mode – way the data are normalised. If mode is symmetric, then divide n by N-1. IF mode is periodic, divide by N, to be consistent with octave code.

When using windows for filter design, the symmetric mode should be used (default). When using windows for spectral analysis, the periodic mode should be used. The mathematical form of the flat-top window in the symmetric case is:

w(n) = a_0 
- a_1 \cos\left(\frac{2\pi n}{N-1}\right)
+ a_2 \cos\left(\frac{4\pi n}{N-1}\right)
- a_3 \cos\left(\frac{6\pi n}{N-1}\right)
+ a_4 \cos\left(\frac{8\pi n}{N-1}\right)

coeff value
a0 0.21557895
a1 0.41663158
a2 0.277263158
a3 0.083578947
a4 0.006947368
from spectrum import window_visu
window_visu(64, 'bohman')

[hires.png, pdf]

../_images/a69334c105.png
window_gaussian(N, alpha=2.5)

Gaussian window

Parameters:N – window length

\exp^{-0.5 \left( \sigma\frac{n}{N/2} \right)^2}

with \frac{N-1}{2}\leq n \leq \frac{N-1}{2}.

Note

N-1 is used to be in agreement with octave convention. The ENBW of 1.4 is also in agreement with [Harris]

from spectrum import window_visu
window_visu(64, 'gaussian', alpha=2.5)

[hires.png, pdf]

../_images/6d4896ce38.png

See also

scipy.signal.gaussian, create_window()

window_hamming(N)

Hamming window

Parameters:N – window length

The Hamming window is defined as

0.54 -0.46 \cos\left(\frac{2\pi n}{N-1}\right)
\qquad 0 \leq n \leq M-1

from spectrum import window_visu
window_visu(64, 'hamming')

[hires.png, pdf]

../_images/76c5bfde42.png

See also

numpy.hamming, create_window(), Window.

window_hann(N)

Hann window (or Hanning). (wrapping of numpy.bartlett)

Parameters:N (int) – window length

The Hanning window is also known as the Cosine Bell. Usually, it is called Hann window, to avoid confusion with the Hamming window.

w(n) =  0.5\left(1- \cos\left(\frac{2\pi n}{N-1}\right)\right)
\qquad 0 \leq n \leq M-1

from spectrum import window_visu
window_visu(64, 'hanning')

[hires.png, pdf]

../_images/a3d9b26b28.png

See also

numpy.hanning, create_window(), Window.

window_kaiser(N, beta=8.6, method='numpy')

Kaiser window

Parameters:
  • N – window length
  • beta – kaiser parameter (default is 8.6)

To obtain a Kaiser window that designs an FIR filter with sidelobe attenuation of \alpha dB, use the following \beta where \beta = \pi \alpha.

w_n = \frac{I_0\left(\pi\alpha\sqrt{1-\left(\frac{2n}{M}-1\right)^2}\right)} {I_0(\pi \alpha)}

where

  • I_0 is the zeroth order Modified Bessel function of the first kind.
  • \alpha is a real number that determines the shape of the window. It determines the trade-off between main-lobe width and side lobe level.
  • the length of the sequence is N=M+1.

The Kaiser window can approximate many other windows by varying the \beta parameter

beta Window shape
0 Rectangular
5 Similar to a Hamming
6 Similar to a Hanning
8.6 Similar to a Blackman
from pylab import plot, legend, hold, xlim
from spectrum import window_kaiser
N = 64
for beta in [1,2,4,8,16]:
    plot(window_kaiser(N, beta), label='beta='+str(beta))
    hold(True)
xlim(0,N)
legend()

[hires.png, pdf]

../_images/829e0ab17c.png
from spectrum import window_visu
window_visu(64, 'kaiser', beta=8.)

[hires.png, pdf]

../_images/7c9284a40a.png

See also

numpy.kaiser, spectrum.window.create_window()

window_lanczos(N)

Lanczos window also known as sinc window.

Parameters:N – window length

w(n) = sinc \left(  \frac{2n}{N-1} - 1 \right)

from spectrum import window_visu
window_visu(64, 'lanczos')

[hires.png, pdf]

../_images/bcf42c1af4.png
window_nuttall(N)

Nuttall tapering window

Parameters:N – window length

w(n) = a_0 - a_1 \cos\left(\frac{2\pi n}{N-1}\right)+ a_2 \cos\left(\frac{4\pi n}{N-1}\right)- a_3 \cos\left(\frac{6\pi n}{N-1}\right)

with a_0 = 0.355768, a_1 = 0.487396, a_2=0.144232 and a_3=0.012604

from spectrum import window_visu
window_visu(64, 'nuttall', mindB=-80)

[hires.png, pdf]

../_images/8f1ac5e438.png
window_parzen(N)

Parsen tapering window (also known as de la Valle-Poussin)

Parameters:N – window length

Parzen windows are piecewise cubic approximations of Gaussian windows. Parzen window sidelobes fall off as 1/\omega^4.

if 0\leq|x|\leq (N-1)/4:

w(n) = 1-6 \left( \frac{|n|}{N/2} \right)^2 +6 \left( \frac{|n|}{N/2}\right)^3

if (N-1)/4\leq|x|\leq (N-1)/2

w(n) = 2 \left(1- \frac{|n|}{N/2}\right)^3

from spectrum import window_visu
window_visu(64, 'parzen')

[hires.png, pdf]

../_images/3b163de63c.png
window_poisson(N, alpha=2)

Poisson tapering window

Parameters:N (int) – window length

w(n) = \exp^{-\alpha \frac{|n|}{N/2} }

with -N/2 \leq n \leq N/2.

from spectrum import window_visu
window_visu(64, 'poisson')
window_visu(64, 'poisson', alpha=3)
window_visu(64, 'poisson', alpha=4)

[hires.png, pdf]

../_images/73e558f5f1.png
window_poisson_hanning(N, alpha=2)

Hann-Poisson tapering window

This window is constructed as the product of the Hanning and Poisson windows. The parameter alpha is the Poisson parameter.

Parameters:
  • N (int) – window length
  • alpha (float) – parameter of the poisson window
from spectrum import window_visu
window_visu(64, 'poisson_hanning', alpha=0.5)
window_visu(64, 'poisson_hanning', alpha=1)
window_visu(64, 'poisson_hanning')

[hires.png, pdf]

../_images/483a2e9663.png
window_rectangle(N)

Kaiser window

Parameters:N – window length
from spectrum import window_visu
window_visu(64, 'rectangle')

[hires.png, pdf]

../_images/4d5a00169f.png
window_riemann(N)

Riemann tapering window

Parameters:N (int) – window length

w(n) = 1 - \left| \frac{n}{N/2}  \right|^2

with -N/2 \leq n \leq N/2.

from spectrum import window_visu
window_visu(64, 'riesz')

[hires.png, pdf]

../_images/d62f4b0c27.png
window_riesz(N)

Riesz tapering window

Parameters:N – window length

w(n) = 1 - \left| \frac{n}{N/2}  \right|^2

with -N/2 \leq n \leq N/2.

from spectrum import window_visu
window_visu(64, 'riesz')

[hires.png, pdf]

../_images/d62f4b0c27.png
window_tukey(N, r=0.5)

Tukey tapering window (or cosine-tapered window)

Parameters:
  • N – window length
  • r – defines the ratio between the constant section and the cosine section. It has to be between 0 and 1.

The function returns a Hanning window for r=0 and a full box for r=1.

from spectrum import window_visu
window_visu(64, 'tukey')
window_visu(64, 'tukey', r=1)

[hires.png, pdf]

../_images/f5901ef5c2.png

0.5 (1+cos(2pi/r (x-r/2))) for 0<=x<r/2

0.5 (1+cos(2pi/r (x-1+r/2))) for x>r/2

window_visu(N=51, name='hamming', **kargs)

A Window visualisation tool

Parameters:
  • N – length of the window
  • name – name of the window
  • NFFT – padding used by the FFT
  • mindB – the minimum frequency power in dB
  • maxdB – the maximum frequency power in dB
  • kargs – optional arguments passed to create_window()

This function plot the window shape and its equivalent in the Fourier domain.

from spectrum import window_visu
window_visu(64, 'kaiser', beta=8.)

[hires.png, pdf]

../_images/7c9284a40a.png