FunctionOne-sided power-density spectrum
Syntax C/C++#include <VFstd.h>
float VF_spectrum( fVector Spc, ui specsiz, fVector X, ui xsiz, fVector Win );
float VFb_spectrum( fVector Spc, ui specsiz, fVector X, ui xsiz, fVector Win, fVector Buf );
C++ VecObj#include <OptiVec.h>
T vector<T>::spectrum( const vector<T>& X, const vector<T>& Win );
T vector<T>::b_spectrum( const vector<T>& X, const vector<T>& Win, vector<T>& Buf );
Pascal/Delphiuses VFstd;
function VF_spectrum( Spc:fVector; specsiz:UIntSize; X:fVector; xsiz:UIntSize; Win:fVector ): Single;
function VFb_spectrum( Spc:fVector; specsiz:UIntSize; X:fVector; xsiz:UIntSize; Win, Buf:fVector ): Single;
CUDA function C/C++#include <cudaVFstd.h>
int cudaVF_spectrum( float *h_psdfNyq, fVector d_Spc, ui specsiz, fVector d_X, ui xsiz, fVector d_Win );
int cusdVF_spectrum( float *d_psdfNyq, fVector d_Spc, ui specsiz, fVector d_X, ui xsiz, fVector d_Win );
float VFcu_spectrum( fVector h_Spc, ui specsiz, fVector h_X, ui xsiz, fVector h_Win );
CUDA function Pascal/Delphiuses VFstd;
function cudaVF_spectrum( var h_psdfNyq:Single; d_Spc:fVector; specsiz:UIntSize; d_X:fVector; xsiz:UIntSize; d_Win:fVector ): IntBool;
function cusdVF_spectrum( d_psdfNyq:PSingle; d_Spc:fVector; specsiz:UIntSize; d_X:fVector; xsiz:UIntSize; d_Win:fVector ): IntBool;
function VFcu_spectrum( h_Spc:fVector; specsiz:UIntSize; h_X:fVector; xsiz:UIntSize; h_Win:fVector );
DescriptionThe data set X is analyzed for its power spectral density (PSD), i.e. the mean square amplitude. The result is stored in Spc. xsiz must be at least 2*specsiz, and specsiz has to be an integer power of 2. Internally, X is divided into (xsiz/specsiz)−1 segments and the average over the spectra of the individual segments is calculated. Each segment of length 2*specsiz yields the PSD for specsiz+1 frequencies (see VF_FFT). In order to keep specsiz an integer power of 2, there are only specsiz points stored in Spc and the last one, the PSD at the Nyquist frequency fNyquist = 0.5 / sampling_interval, is given as the return value of the function. It may either be neglected (by calling the function like a void function) or stored as the last element in Spc by calling the function as
Spc[specsiz] = VF_spectrum( Spc, specsiz, X, xsiz, Win );
in this case, Spc must have a length of specsiz+1.

Win is a window that is applied to the data segments. The size of the Win vector must be 2*specsiz. Within the VectorLib library, three functions are available that give suitable Windows: VF_Welch,   VF_Parzen, and VF_Hann. A square window (i.e. no windowing at all) is achieved by setting all elements of Win to 1.0 using VF_equ1. Use of the square window is not recommended here, though.

You may wish to test the quality of the calculated spectrum by applying Parseval's theorem (provided you called VF_spectrum as in the above example and stored the PSD for the Nyquist frequency):
1.0/xsize * VF_ssq( X, xsize ) must be about equal to VF_sum( Spc, specsiz+1 ).
If the deviation between both results is large, the sampling interval in X probably is too large.

Internally, VF_spectrum allocates and frees additional workspace memory. For repeated calls, this would be inefficient. In such a case, it is recommended to use VFb_spectrum instead. The size of Buf must be ≥ 2*xsiz + 2*specsiz. Additionally, Buf must be 128-bit (P8) or 256-bit (P9) aligned. This usually means you can only take a vector allocated by the VF_vector family as Buf.

For an example for VF_spectrum, see the demo program VDEMO.

Error handlingIf size is not a power of 2, VF_FFT (on which VF_spectrum is based) complains "Size must be an integer power of 2" and the program is aborted.
Return valuePSD at the Nyquist frequency
See alsoVF_FFT,   VF_xspectrum,   VF_convolve,   VF_autocorr,   VF_xcorr,   VF_filter

VectorLib Table of Contents  OptiVec home