BandPower - Maple Help

SignalProcessing

 BandPower
 calculate the band power of a 1-D signal

 Calling Sequence BandPower( data ) BandPower( data, samplerate ) BandPower( data, samplerate, frequencyrange ) BandPower( data, frequencies ) BandPower( data, frequencies, frequencyrange )

Parameters

 data - 1-D rtable or list of data, representing either a signal or power spectral density (PSD). samplerate - (optional) Positive numeric value for the sampling rate. frequencies - (optional) 1-D list or rtable of numeric frequency components. frequencyrange - (optional) Range of numeric frequency values used to filter (i.e. band-limit) the data. If omitted, all frequencies are considered.

Options

 • variety: (optional) Specifies the type of data passed. The options are signal and psd (for a power spectral density). The default is signal.
 • frequencyunit: (optional) Specifies the underlying unit of frequency used when computing statistics. Either of the forms algebraic and Unit(algebraic) are accepted. The default is Unit(Hz).
 • fftnormalization: (optional) One of none, symmetric, or full, indicates the normalization to be applied when using the Fast Fourier Transform (FFT). The default is symmetric.
 • temperendpoints: (optional) Either true or false, specifies whether the power spectrum is to be tempered at the endpoints. The default is false.
 • window: (optional) Either a list, name, or string, specifies the windowing command to be applied to the signal when variety=signal. The default is "none" (for no windowing to be applied).
 • windownormalization (optional) Either true or false, indicates if the windowing function is to be normalized. The default is true.
 • method: (optional) The method of numerical integration to be used to find band power. The options are leftendpoint, rightendpoint, simpson, and trapezoid. The default is trapezoid, which is equivalent to finding the area under the linear interpolation through all the points.
 • output : (optional) The type of output. The supported options are:
 – bandpower: Returns a value of type float[8] for the band power.
 – entropy: Returns a value of type float[8] for the spectral entropy.
 – meanfrequency: Returns a value of type float[8] for the mean frequency.
 – frequencies: Returns a Vector, of float[8] datatype and length the same as signal, containing the frequencies.
 – indices: Returns a list of the form $\left[\mathrm{\alpha },\mathrm{\beta }\right]$, where $\mathrm{\alpha }$ and $\mathrm{\beta }$ are, respectively, the lower and upper indices of the frequencies Vector that correspond to the frequencies in frequencyrange.
 – psd: Returns a vector of float[8] datatype containing the power spectral density.
 – record: Returns a record with the previous options.
 – list of any of the above options: Returns an expression sequence with the corresponding outputs, in the same order.
 The default output is bandpower.

Description

 • The BandPower command takes a 1-D rtable or list of data and computes band power based on the remaining parameters and options. Other statistics, like mean frequency and spectral entropy, can be computed and returned, along with Vectors for the power spectral density and frequencies.
 • If data is an rtable of type AudioTools:-Audio, the sample rate in hertz ($\mathrm{Hz}$) is inferred from the attributes. Should samplerate also be passed, it will be overridden.
 • Values for/in samplerate, frequencyrange, and frequencies can include units for frequency, and those without an explicit unit are assumed to be in terms of frequencyunit. An error will be thrown if a unit which cannot be converted to hertz is passed. Caution: Passing long vectors with elements having units requires time-consuming conversions. It is recommended that the frequencies Vector consist only of the numeric values, with units indicated by frequencyunit.
 • If window is passed as a list, the first element provides the name of the windowing command, and any remaining terms are passed as options to the command.
 • The value of window, when not passed as a list, should be the name or string, with or without the Window suffix, that corresponds to the windowing command. For example, to use a Hamming window, you can pass window=Hamming or window="HammingWindow". In both cases, the command SignalProcessing[HammingWindow] will be used internally. Similarly, you can pass window=["Exponential",0.5] or window=[ExponentialWindow,0.5] to use SignalProcessing[ExponentialWindow] with parameter value 0.5.
 • To apply a window to a Vector $V$ of length $n$, the window is first applied to another Vector $W$ of size $n$ and filled with ones, and then $V$ is multiplied element-wise by $W$. When windownormalization=true, $W$ is first normalized with respect to its Root Mean Square (RMS).
 • Two or more points are always required in data to compute band power, mean frequency, and spectral entropy in the frequency domain.
 • When variety=psd, the elements of data must be non-negative, the elements of frequencies must be strictly increasing, and the left end of frequencyrange must be less than the right end.
 • A value for samplerate is required and used when variety=signal.
 • The frequencies rtable is required and used only when variety=psd. In this case, the rtable must be of the same size as data.
 • Input rtables (data and frequencies) are converted to Vectors of either float[8] or complex[8] datatype, and an error will be thrown if this is not possible. For this reason, it is most efficient for input rtables to already be Vectors having the appropriate datatype. Specifically:
 • frequencies should be float[8].
 • When variety=signal, data should be float[8] (for real-valued data) or complex[8] (for complex-valued data).
 • When variety=psd, data should be float[8].
 • The band power of signal $X$ with length $n$ is computed as the area beneath the graph of the power spectral density $P$ of $X$ versus the frequencies $F$. To this end, first calculate the FFT of $X$. If $X$ is real-valued, we take $Y$ to be the first half of the FFT (due to symmetry). The end result will be scaled appropriately, giving the same band power as if the truncation had not been performed. When $X$ has complex values, on the other hand, a full FFT is used to define $Y$. We take the half-size to be the following:

$m=⌊\frac{n}{2}⌋+1$

 Second, compute the power spectral density $P$ from $Y$:

$\mathrm{P__i}=\frac{c{\left|\mathrm{Y__i}\right|}^{2}}{r},i=1..k$

 where $k$ is the size of $Y$ (either $m$ or $n$), $r$ is the sample rate, and $c=2$ when $X$ is real-valued and $c=1$ when $X$ is complex-valued. To create the Vector $F$ of frequencies, we consider the two cases. When $X$ is real-valued:

$\mathrm{F__i}=\frac{\left(i-1\right)r}{n},i=1..m$

 When $X$ is complex-valued:

$\mathrm{F__i}=\frac{\left(i-1\right)r}{n},i=1..n$

 Caution: Maple uses the convention of the frequencies ranging from $0$ to $T$, where $0 and $\left[0,T\right]$ is the interval over which the signal is sampled. The other common convention is to use frequencies in the range $\left[-\frac{T}{2},\frac{T}{2}\right]$. Usually, $T$ is the period of a periodic function.
 • The mean frequency of signal $X$ is computed as the average of the frequencies $F$ weighted by the power spectral density $P$:

$\mathrm{MeanFrequency}\left(X\right)=\frac{{\sum }_{i=1}^{k}\mathrm{F__i}\mathrm{P__i}}{{\sum }_{i=1}^{k}\mathrm{P__i}}$

 where $k$ is the size of $F$ (either $m$ or $n$). When all elements of $P$ are zero, then the mean frequency is HFloat(undefined).
 • The spectral entropy of signal $X$ is computed using the formula from Shannon Information Theory applied to the power spectral density $P$ treated as a probability distribution:

$\mathrm{SpectralEntropy}\left(X\right)=-\left({\sum }_{i=1}^{k}\mathrm{Q__i}{\mathrm{log}}_{2}\left(\mathrm{Q__i}\right)\right)$

 where

$Q=\frac{P}{{\sum }_{i=1}^{k}\mathrm{P__i}}$

 and $k$ is the size of $P$ (either $m$ or $n$). When all elements of $P$ are zero, then the spectral entropy is HFloat(undefined). Note that another convention for spectral entropy is to divide the above by ${\mathrm{log}}_{2}\left(k\right)$, which quantifies the maximum spectral entropy of white noise.
 • When temperendpoints=true, the endpoints of the power spectrum, determined by the method described above, are halved. Note that the signal size must be three or more for tempering.
 • The frequencyrange option is used to restrict the computations for statistics to a specific frequency band. It can be in the form a..b (to restrict frequencies to those between a and b), a.. (to restrict frequencies to those that are no less than a), ..b (to restrict frequencies to those that are no more than b), or .. (to allow all frequencies). The endpoints of frequencyrange are used to determine indices $\mathrm{\alpha }$ and $\mathrm{\beta }$, with the band power estimated as the the numeric integral of ${P}_{\mathrm{\alpha }..\mathrm{\beta }}$ over ${F}_{\mathrm{\alpha }..\mathrm{\beta }}$, the mean frequency estimated as the mean of ${F}_{\mathrm{\alpha }..\mathrm{\beta }}$ with weights ${P}_{\mathrm{\alpha }..\mathrm{\beta }}$, and the spectral entropy computed as the entropy of ${Q}_{\mathrm{\alpha }..\mathrm{\beta }}$ where $Q$ is $P$ normalized as a probability distribution. In the simplest case, where $a$ and $b$ are both provided with $\mathrm{min}\left(F\right)\le a$, $a, and $b\le \mathrm{max}\left(F\right)$, we take $\mathrm{\alpha }$ to be the largest index $i$ such that ${F}_{i}\le a$, and $\mathrm{\beta }$ to be the smallest index $j$ such that $b\le {F}_{j}$.
 • Band powers for a disjoint union of intervals typically add up to close to the total. For example, BandPower(X,0.01,..50)+BandPower(X,0.01,50..) should be close to (but not necessarily equal, due to discretization and numerical errors) BandPower(X).
 • The BandPower command is not thread safe.
 • As the underlying implementation of the SignalProcessing package is a module, it is also possible to use the form SignalProcessing:-BandPower to access the command from the package. For more information, see Module Members.

Examples

 > $\mathrm{with}\left(\mathrm{SignalProcessing}\right):$

Real-Valued Signal

Signal

 • For an illustrative example involving a real-valued signal, first define the sampling rate and a Vector of times:
 > $r≔5000$
 ${r}{≔}{5000}$ (1)
 > $T≔\mathrm{Vector}\left(\left[\mathrm{seq}\right]\left(0..1,\frac{1}{r}\right),\mathrm{datatype}=\mathrm{float}\left[8\right]\right)$
 ${T}{≔}\begin{array}{c}\left[\begin{array}{c}{0.}\\ {0.000200000000000000}\\ {0.000400000000000000}\\ {0.000600000000000000}\\ {0.000800000000000000}\\ {0.00100000000000000}\\ {0.00120000000000000}\\ {0.00140000000000000}\\ {0.00160000000000000}\\ {0.00180000000000000}\\ {⋮}\end{array}\right]\\ \hfill {\text{5001 element Vector[column]}}\end{array}$ (2)
 • Second, sample the signal at these times:
 > $g≔t↦\mathrm{cos}\left(1000\cdot \mathrm{\pi }\cdot t\right)+3\cdot \mathrm{cos}\left(2000\cdot \mathrm{\pi }\cdot t\right)+2\cdot \mathrm{cos}\left(3000\cdot \mathrm{\pi }\cdot t\right)$
 ${g}{≔}{t}{↦}{\mathrm{cos}}{}\left({1000}{\cdot }{\mathrm{\pi }}{\cdot }{t}\right){+}{3}{\cdot }{\mathrm{cos}}{}\left({2000}{\cdot }{\mathrm{\pi }}{\cdot }{t}\right){+}{2}{\cdot }{\mathrm{cos}}{}\left({3000}{\cdot }{\mathrm{\pi }}{\cdot }{t}\right)$ (3)
 > $X≔\mathrm{map}\left[\mathrm{evalhf}\right]\left(g,T\right)$
 ${X}{≔}\begin{array}{c}\left[\begin{array}{c}{6.}\\ {1.11803398874990}\\ {-3.73606797749979}\\ {-1.11803398874990}\\ {0.736067977499790}\\ {0.}\\ {0.736067977499789}\\ {-1.11803398874989}\\ {-3.73606797749979}\\ {1.11803398874989}\\ {⋮}\end{array}\right]\\ \hfill {\text{5001 element Vector[column]}}\end{array}$ (4)

Periodogram

 > $\mathrm{Periodogram}\left(X,\mathrm{samplerate}=r,\mathrm{powerscale}="dB/Hz"\right)$

Band Power

 • Intuitively, we expect the overall band power to be dominated by the frequencies 500, 1000, and 1500, in the proportion 1:9:4. The frequencies can be gleaned from the standard expression $A\mathrm{cos}\left(2\mathrm{\pi }ft\right)$, and the proportions follow from the squares of the respective amplitudes.
 • To confirm this, first find the overall band power:
 > $p\left[0\right]≔\mathrm{BandPower}\left(X,r,\mathrm{variety}=\mathrm{signal}\right)$
 ${{p}}_{{0}}{≔}{7.00579884023191113}$ (5)
 • Second, we output the power spectral density and frequencies Vectors from BandPower so that they are only calculated once:
 > $P,F≔\mathrm{BandPower}\left(X,r,\mathrm{variety}=\mathrm{signal},\mathrm{output}=\left[\mathrm{psd},\mathrm{frequencies}\right]\right)$
 ${P}{,}{F}{≔}\begin{array}{c}\left[\begin{array}{c}{2.87942411523574}{×}{{10}}^{{-6}}\\ {2.87943139176523}{×}{{10}}^{{-6}}\\ {2.87945322153172}{×}{{10}}^{{-6}}\\ {2.87948960569330}{×}{{10}}^{{-6}}\\ {2.87954054530436}{×}{{10}}^{{-6}}\\ {2.87960604232246}{×}{{10}}^{{-6}}\\ {2.87968609931651}{×}{{10}}^{{-6}}\\ {2.87978071917066}{×}{{10}}^{{-6}}\\ {2.87988990586231}{×}{{10}}^{{-6}}\\ {2.88001366228741}{×}{{10}}^{{-6}}\\ {⋮}\end{array}\right]\\ \hfill {\text{2501 element Vector[column]}}\end{array}{,}\begin{array}{c}\left[\begin{array}{c}{0.}\\ {0.999800039992002}\\ {1.99960007998400}\\ {2.99940011997600}\\ {3.99920015996801}\\ {4.99900019996001}\\ {5.99880023995201}\\ {6.99860027994401}\\ {7.99840031993601}\\ {8.99820035992802}\\ {⋮}\end{array}\right]\\ \hfill {\text{2501 element Vector[column]}}\end{array}$ (6)
 • Of course, we can pass the original signal using the variety=signal option to compute band powers and mean frequencies for various frequency ranges, but each call would compute internally the power spectral density and frequencies Vectors.
 • Now, find the band powers in a narrow band around each of the frequencies:
 > $p\left[1\right]≔\mathrm{BandPower}\left(P,F,400..600,\mathrm{variety}=\mathrm{psd}\right)$
 ${{p}}_{{1}}{≔}{0.502979149085503674}$ (7)
 > $p\left[2\right]≔\mathrm{BandPower}\left(P,F,900..1100,\mathrm{variety}=\mathrm{psd}\right)$
 ${{p}}_{{2}}{≔}{4.50079569444223004}$ (8)
 > $p\left[3\right]≔\mathrm{BandPower}\left(P,F,1400..1600,\mathrm{variety}=\mathrm{psd}\right)$
 ${{p}}_{{3}}{≔}{1.99657149902250763}$ (9)
 • The total of these individual band powers is indeed very close to the overall band power:
 > $\mathrm{abs}\left(\frac{p\left[1\right]+p\left[2\right]+p\left[3\right]}{p\left[0\right]}-1\right)$
 ${0.0007782836}$ (10)
 • Units are also accepted:
 > $\mathrm{BandPower}\left(P,F,..1.250\mathrm{Unit}\left(\mathrm{kHz}\right),\mathrm{variety}=\mathrm{psd}\right)$
 ${5.00764465432853001}$ (11)
 • BandPower can also return the indices for the frequencies Vector corresponding to the provided frequency range:
 > $\mathrm{BandPower}\left(P,F,1.25\mathrm{Unit}\left('\mathrm{kHz}'\right)..,\mathrm{variety}=\mathrm{psd},\mathrm{output}=\mathrm{indices}\right)$
 $\left[{1251}{,}{2501}\right]$ (12)

Complex-Valued Signal

 • Consider now the following complex-valued signal:
 > $r≔5000$
 ${r}{≔}{5000}$ (13)
 > $T≔\mathrm{Vector}\left(\left[\mathrm{seq}\right]\left(0..1,\frac{1}{r}\right),\mathrm{datatype}=\mathrm{float}\left[8\right]\right)$
 ${T}{≔}\begin{array}{c}\left[\begin{array}{c}{0.}\\ {0.000200000000000000}\\ {0.000400000000000000}\\ {0.000600000000000000}\\ {0.000800000000000000}\\ {0.00100000000000000}\\ {0.00120000000000000}\\ {0.00140000000000000}\\ {0.00160000000000000}\\ {0.00180000000000000}\\ {⋮}\end{array}\right]\\ \hfill {\text{5001 element Vector[column]}}\end{array}$ (14)
 > $g≔t↦5\cdot \mathrm{cos}\left(1000\cdot \mathrm{\pi }\cdot t\right)-3\cdot I\cdot \mathrm{sin}\left(2000\cdot \mathrm{\pi }\cdot t\right)$
 ${g}{≔}{t}{↦}{5}{\cdot }{\mathrm{cos}}{}\left({1000}{\cdot }{\mathrm{\pi }}{\cdot }{t}\right){-}{3}{\cdot }{I}{\cdot }{\mathrm{sin}}{}\left({2000}{\cdot }{\mathrm{\pi }}{\cdot }{t}\right)$ (15)
 > $X≔\mathrm{Vector}\left(\mathrm{map}\left(g,T\right),\mathrm{datatype}=\mathrm{complex}\left[8\right]\right)$
 ${X}{≔}\begin{array}{c}\left[\begin{array}{c}{5.}{-}{0.}{}{I}\\ {4.04508497163362}{-}{2.85316954903757}{}{I}\\ {1.54508497109448}{-}{1.76335575608094}{}{I}\\ {-1.54508497304513}{+}{1.76335575807213}{}{I}\\ {-4.04508497283919}{+}{2.85316954827701}{}{I}\\ {-5.}{-}{2.46124380671923}{×}{{10}}^{{-9}}{}{I}\\ {-4.04508497042806}{-}{2.85316954979814}{}{I}\\ {-1.54508496914383}{-}{1.76335575408976}{}{I}\\ {1.54508497499578}{+}{1.76335576006332}{}{I}\\ {4.04508497404476}{+}{2.85316954751644}{}{I}\\ {⋮}\end{array}\right]\\ \hfill {\text{5001 element Vector[column]}}\end{array}$ (16)
 • First, determine the power spectral density and frequencies Vectors, so that they are only calculated once:
 > $P,F≔\mathrm{BandPower}\left(X,r,\mathrm{variety}=\mathrm{signal},\mathrm{output}=\left[\mathrm{psd},\mathrm{frequencies}\right]\right)$
 ${P}{,}{F}{≔}\begin{array}{c}\left[\begin{array}{c}{9.99801302228631}{×}{{10}}^{{-7}}\\ {9.98772081349862}{×}{{10}}^{{-7}}\\ {9.97759106040794}{×}{{10}}^{{-7}}\\ {9.96762346817629}{×}{{10}}^{{-7}}\\ {9.95781775321601}{×}{{10}}^{{-7}}\\ {9.94817364346192}{×}{{10}}^{{-7}}\\ {9.93869087798189}{×}{{10}}^{{-7}}\\ {9.92936920722300}{×}{{10}}^{{-7}}\\ {9.92020839121018}{×}{{10}}^{{-7}}\\ {9.91120820794676}{×}{{10}}^{{-7}}\\ {⋮}\end{array}\right]\\ \hfill {\text{5001 element Vector[column]}}\end{array}{,}\begin{array}{c}\left[\begin{array}{c}{0.}\\ {0.999800039992002}\\ {1.99960007998400}\\ {2.99940011997600}\\ {3.99920015996801}\\ {4.99900019996001}\\ {5.99880023995201}\\ {6.99860027994401}\\ {7.99840031993601}\\ {8.99820035992802}\\ {⋮}\end{array}\right]\\ \hfill {\text{5001 element Vector[column]}}\end{array}$ (17)
 • The overall power can be found in either the time domain or the frequency domain:
 > $\mathrm{BandPower}\left(X,r,\mathrm{variety}=\mathrm{signal}\right)$
 ${17.0015986811105506}$ (18)
 > $\mathrm{BandPower}\left(P,F,\mathrm{variety}=\mathrm{psd}\right)$
 ${17.0015986811105506}$ (19)
 • We can consult the plot of the power spectral density versus the frequencies to see where the power is concentrated:
 > $\mathrm{dataplot}\left(F,P,\mathrm{view}=\left[0..\mathrm{max}\left(F\right),0..\mathrm{max}\left(P\right)\right],\mathrm{labels}=\left[\mathrm{Frequency},\mathrm{Power}\right],\mathrm{title}="Power Spectral Density",\mathrm{font}=\left[\mathrm{Verdana},15\right]\right)$
 > $\mathrm{BandPower}\left(P,F,400..600,\mathrm{variety}=\mathrm{psd}\right)$
 ${6.24851796795285619}$ (20)
 > $\mathrm{BandPower}\left(P,F,900..1100,\mathrm{variety}=\mathrm{psd}\right)$
 ${2.24883444827089463}$ (21)

Audio Signal

 • The signal passed to BandPower can also come directly from an audio file:
 > $\mathrm{with}\left(\mathrm{AudioTools}\right):$
 > $\mathrm{file}≔\mathrm{cat}\left(\mathrm{kernelopts}\left(\mathrm{datadir}\right),\mathrm{kernelopts}\left(\mathrm{dirsep}\right),"audio",\mathrm{kernelopts}\left(\mathrm{dirsep}\right),"ViolinThreePosVibrato.wav"\right)$
 ${\mathrm{file}}{≔}{"/maple/cbat/active/169462/data/audio/ViolinThreePosVibrato.wav"}$ (22)
 > $\mathrm{Violin}≔\mathrm{ToMono}\left(\mathrm{Read}\left(\mathrm{file},\mathrm{samples}=1000..60000\right)\right)$
 ${\mathrm{Violin}}{≔}\left[\begin{array}{cc}{"Sample Rate"}& {44100}\\ {"File Format"}& {\mathrm{PCM}}\\ {"File Bit Depth"}& {16}\\ {"Channels"}& {1}\\ {"Samples/Channel"}& {59001}\\ {"Duration"}& {1.33789}{}{s}\end{array}\right]$ (23)
 > $\mathrm{Periodogram}\left(\mathrm{Violin}\right)$
 • Here, let's find the band power for frequencies between 0 and 5000, and the band power for frequencies above 5000:
 > $p\left[1\right]≔\mathrm{BandPower}\left(\mathrm{Violin},..5000\right)$
 ${{p}}_{{1}}{≔}{0.000360588871864467189}$ (24)
 > $p\left[2\right]≔\mathrm{BandPower}\left(\mathrm{Violin},5000..\right)$
 ${{p}}_{{2}}{≔}{0.0000138401583488221083}$ (25)
 • Units can also be supplied for frequencies:
 > $\mathrm{BandPower}\left(\mathrm{Violin},..5\mathrm{Unit}\left(\mathrm{kHz}\right)\right)$
 ${0.000360588871864467189}$ (26)
 • As percentages:
 > $p\left[0\right]≔\mathrm{BandPower}\left(\mathrm{Violin}\right)$
 ${{p}}_{{0}}{≔}{0.000374426425028005923}$ (27)
 > $q\left[1\right]≔\frac{p\left[1\right]}{p\left[0\right]}\cdot 100$
 ${{q}}_{{1}}{≔}{96.30433319}$ (28)
 > $q\left[2\right]≔\frac{p\left[2\right]}{p\left[0\right]}\cdot 100$
 ${{q}}_{{2}}{≔}{3.696362603}$ (29)

Compatibility

 • The SignalProcessing[BandPower] command was introduced in Maple 2021.
 • For more information on Maple 2021 changes, see Updates in Maple 2021.