MeanFrequency - Maple Help

SignalProcessing

 MeanFrequency
 calculate the mean frequency of a 1-D signal

 Calling Sequence MeanFrequency( data, samplerate ) MeanFrequency( data, samplerate, frequencyrange) MeanFrequency( data, frequencies ) MeanFrequency( 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 meanfrequency.

Description

 • The MeanFrequency command takes a 1-D rtable or list of data and computes mean frequency based on the remaining parameters and options. Other statistics, like band power 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).
 • As the underlying implementation of the SignalProcessing package is a module, it is also possible to use the form SignalProcessing:-MeanFrequency to access the command from the package. For more information, see Module Members.
 • The MeanFrequency command is not thread safe.

Examples

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

Real-Valued Signal

Signal

 • To create a 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)$

Mean Frequency in Time Domain

 • Clearly, the signal is dominated by the frequencies 500, 1000, and 1500. To confirm this, we can use MeanFrequency in bands around these values:
 > $\mathrm{MeanFrequency}\left(X,r,400\mathrm{Unit}\left(\mathrm{Hz}\right)..600\mathrm{Unit}\left(\mathrm{Hz}\right)\right)$
 ${499.813729680681604}$ (5)
 > $\mathrm{MeanFrequency}\left(X,r,900\mathrm{Unit}\left(\mathrm{Hz}\right)..1100\mathrm{Unit}\left(\mathrm{Hz}\right)\right)$
 ${999.808430284599694}$ (6)
 > $\mathrm{MeanFrequency}\left(X,r,1.4\mathrm{Unit}\left(\mathrm{kHz}\right)..1.6\mathrm{Unit}\left(\mathrm{kHz}\right)\right)$
 ${1499.87447320330216}$ (7)

Mean Frequency in Frequency Domain

 • We can output the power spectral density and frequencies Vectors from MeanFrequency so that they are only calculated once:
 > $P,F≔\mathrm{MeanFrequency}\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}$ (8)
 • Computation in the frequency domain, of course, gives the same results as above in the time domain:
 > $\mathrm{MeanFrequency}\left(P,F,\mathrm{variety}=\mathrm{psd},400\mathrm{Unit}\left(\mathrm{Hz}\right)..600\mathrm{Unit}\left(\mathrm{Hz}\right)\right)$
 ${499.813729680681604}$ (9)
 > $\mathrm{MeanFrequency}\left(P,F,\mathrm{variety}=\mathrm{psd},900\mathrm{Unit}\left(\mathrm{Hz}\right)..1100\mathrm{Unit}\left(\mathrm{Hz}\right)\right)$
 ${999.808430284599694}$ (10)
 > $\mathrm{MeanFrequency}\left(P,F,\mathrm{variety}=\mathrm{psd},1.4\mathrm{Unit}\left(\mathrm{kHz}\right)..1.6\mathrm{Unit}\left(\mathrm{kHz}\right)\right)$
 ${1499.87447320330216}$ (11)

Tempering Endpoints

 • Consider the following signal:
 > $n≔{2}^{13}:$
 > $a≔0:$
 > $b≔1:$
 > $\mathrm{fs}≔\mathrm{evalhf}\left(\frac{n-1}{b-a}\right):$
 > $T≔\mathrm{Vector}\left(1..n,i↦\mathrm{evalhf}\left(\frac{i-1}{\mathrm{fs}}\right),'\mathrm{datatype}'='\mathrm{float}\left[8\right]'\right):$
 > $\mathrm{\alpha }≔0:$
 > $\mathrm{\beta }≔2\mathrm{\pi }\cdot 1000:$
 > $g≔\mathrm{unapply}\left(\mathrm{cos}\left(\mathrm{\alpha }t\right)+\mathrm{sin}\left(\mathrm{\beta }t\right),t\right)$
 ${g}{≔}{t}{↦}{1}{+}{\mathrm{sin}}{}\left({2000}{\cdot }{\mathrm{\pi }}{\cdot }{t}\right)$ (12)
 > $X≔\mathrm{map}\left['\mathrm{evalhf}'\right]\left(g,T\right):$
 > $\mathrm{Periodogram}\left(X,'\mathrm{samplerate}'=\mathrm{fs},'\mathrm{powerscale}'="dB/Hz"\right)$
 • There are two peaks in the frequencies, 0 and 1000. Interestingly, the peak for 0 is so sharp that only the first element of the FFT captures it:
 > $\mathrm{\mu }\left[1\right],K\left[1\right]≔\mathrm{MeanFrequency}\left(X,\mathrm{fs},0...,'\mathrm{variety}'='\mathrm{signal}','\mathrm{output}'=\left['\mathrm{meanfrequency}','\mathrm{indices}'\right]\right)$
 ${{\mathrm{\mu }}}_{{1}}{,}{{K}}_{{1}}{≔}{199.975709831543128}{,}\left[{1}{,}{4097}\right]$ (13)
 > $\mathrm{\mu }\left[2\right],K\left[2\right]≔\mathrm{MeanFrequency}\left(X,\mathrm{fs},1.0..,'\mathrm{variety}'='\mathrm{signal}','\mathrm{output}'=\left['\mathrm{meanfrequency}','\mathrm{indices}'\right]\right)$
 ${{\mathrm{\mu }}}_{{2}}{,}{{K}}_{{2}}{≔}{999.976201585510694}{,}\left[{2}{,}{4097}\right]$ (14)
 • Moreover, tempering of the endpoints also has a noticeable effect:
 > $\mathrm{\mu }\left[3\right]≔\mathrm{MeanFrequency}\left(X,\mathrm{fs},'\mathrm{variety}'='\mathrm{signal}','\mathrm{temperendpoints}'='\mathrm{true}'\right)$
 ${{\mathrm{\mu }}}_{{3}}{≔}{333.298268326632410}$ (15)

Compatibility

 • The SignalProcessing[MeanFrequency] command was introduced in Maple 2021.