HIGH DYNAMIC RANGE FOURIER ANALYSIS
This tutorial covers the Fourier spectral analysis capabilities of AutoSignal for that instance where you want to characterize very low power components within wide sense stationary signals.
Generating a Test Signal
Select the Generate Signal option in the Edit menu or Main toolbar.
For this tutorial, we will create a data set that contains four equal power reference sinusoids across the Nyquist range and three low power sinusoids placed between them. The low power components will be -40dB, -50dB, and -60dB below the reference sinusoids. White noise will then be added so that a noise floor exists at approximately -75dB. The object of this test is the accurate frequency and power estimation of these three low power components. A secondary focus will be the recovery of phase information for these components.
Click Read and select the file tutor1.sig from the Signals subdirectory.
The following signal expression is imported:
The X (time) values vary from 0 to 0.1 with a 0.0001 sample increment. The Nyquist frequency is thus 5000 (half the 10000 sampling frequency). The four reference sinusoids span the Nyquist range. The first low power sinusoid (-40dB) has 1% of the amplitude and 0.01% of the power of one of the reference sinusoids. The next low power sinusoid (-50dB) has 0.001% of the power of the references. The last of the test sinusoids (-60dB) has only 0.1% of the amplitude and 0.0001% of the power of the reference sinusoids. The 0.15 Gaussian Noise% creates a white noise floor at about -75dB. It should thus be possible to resolve even the -60dB component given a well designed spectral analysis.
AutoSignal compiles test signal expressions for fast evaluation. Certain rules thus apply for constructing signal expressions. Briefly, expressions that evaluate to constants can be assigned any variable name that is not a built-in function or constant. Expressions containing the X variable must be assigned the function names F1 through F25. The final expression must be assigned to Y.
Click OK to process the current signal.
An AutoSignal graph is presented containing the generated data.
An AutoSignal graph is presented containing the generated data.
Click OK to accept the generated data. Click Yes when asked to update the main data table with the revised data.
Note that the OK button is a green checkmark and the Cancel button is a red X.
The 1001 point generated data set is now displayed in the main status window. In AutoSignal, the main window is informational only. The various analyses and processing take place in separate procedures that are tailored to each specific option.
The Non-Windowed Fourier Spectrum
For high dynamic range Fourier spectral analysis, data tapering windows are almost always used. In order to see why, we will first explore a basic or non-windowed FFT.
Select the Fourier Spectrum option in the Spectral menu or Spectral toolbar. Be sure the default Best Exact N algorithm is selected and Nmin is 1001. Select the dB Norm plot. Be sure the component count (sig) is set to 4.
If necessary, click the Display Maxima button until frequency labels are shown.
The Best Exact N is the default FFT algorithm. This algorithm uses the fastest available FFT algorithm that exactly transforms the data. In this instance, the Prime Factor FFT algorithm is used since it is the fastest AutoSignal FFT procedure that exactly processes the data size of 1001. The Fourier spectrum that is produced readily detects the four reference peaks. The -40dB peak is observed. There is only a hint of the -50dB peak and no evidence of the -60 dB peak. Even though AutoSignal uses a cubic spline bin interpolation to refine the FFT peak frequencies, the reference frequencies of 1005, 2005, 3005, and 4005 are not exactly recovered.
Note that the four references consist of peaks that have substantially different shapes. The spectral peaks become sharper as the frequencies increase. This is not a product of frequency, but rather of spectral leakage to adjacent bins. The more closely a harmonic is centered within a spectral frequency bin, the sharper will be the FFT frequency response.
Click the List Data button.
Observe that the first peak's actual frequency of 1005 must be accommodated by FFT bins with frequencies of 999 and 1009. The fourth peak's actual frequency of 4005 can be mostly accommodated by a 4006 frequency FFT bin. Only when a frequency falls in the exact center of a bin can a harmonic be entirely represented by a single FFT channel.
Close the FFT Data Summary window and change Nmin to 8192 using the drop down list.
This increases the size of the FFT using zero padding. There are now approximately eight times as many frequencies within the spectrum. The bin-interpolated frequencies of the references are now accurately determined, but the three low power components are lost in the baseline. The spectral leakage increases since zero padding introduces a discontinuity in the data stream.
Enter 2000 in the Nmin field.
This produces an FFT spaced at exactly 5Hz intervals. The frequencies of all seven of the sinusoids are now exactly in the center of an FFT bin:
The spectral peaks are now impulse functions, but the spectral leakage introduced by zero padding is sufficient to mask the three lower power components. It is not possible to resolve the low power components without using a data tapering window.
Close the Fourier Spectrum window using the red X Cancel button. The current settings are not retained. You are returned to the main menu.
To preserve the current settings of a spectral option, exit the procedure with the green check OK button.
To discard the current settings, use the red X Cancel button to exit the procedure.
Windowed Fourier Spectra
A data tapering window reduces the magnitude of the data values to zero at each end of the data stream. While this increases the influence of data in the center of the data stream and diminishes the influence of data near the bounds, this does eliminate the discontinuity at the bounds of the data. If the data were to be viewed in a wraparound format, a continuity would be observed where the two ends meet.
The influence of adding a data tapering window is not neutral in the Fourier domain. The data window will have a frequency response that widens any spectral peak, thus reducing its resolution. Low power harmonics are also introduced. The benefit is a dramatic reduction in spectral leakage, enabling any measure of zero padding to be used.
In AutoSignal, a power normalization is made whenever a data window is used. The power of a windowed Fourier spectrum is identical to the power of its non-windowed counterpart. On the other hand, even with this conservation of power, the height of the spectral peaks will evidence some attenuation as a consequence of the windowing. When a data tapering window is used, it is still possible to judge relative power using a normalized dB scale, but it is no longer possible to directly read absolute amplitudes or powers using the spectral peak heights.
Select the Fourier Spectrum with Data Window option from the Spectral menu or toolbar. Select the cs4 BHarris min window, set Nmin to 8192, change the plot to dB Norm, and set the signal count in the sig field to 7.
Click the Display Maxima button until both frequencies and spectral magnitudes are displayed as labels.
This Blackman-Harris window is a popular one for high dynamic range Fourier analysis. The four reference peaks are accurately mapped for both frequency and power. The three low power peaks should also be well characterized for both frequency and power. The results you see may be slightly different due to the random noise component.
Data Tapering Windows
There is a significant art to selecting a data window. In AutoSignal, there are twenty fixed property windows, and 9 whose properties vary as a consequence of an adjustable parameter.
Select the Explore Data Tapering Windows option.
This option displays the frequency response of up to three data tapering windows. This is the best way to learn about the various data windows that AutoSignal offers.
The Mainlobe width is directly related to the frequency resolution. A window with a mainlobe width of 3 will have only 1/3 the frequency resolution of an unwindowed Fourier spectrum. The Sidelobe value is the extent in dB the largest sidelobe peak is below the primary peak. The higher the sidelobe value, the better the likelihood of resolving a low power peak in close frequency proximity to a higher power component. The Rolloff value is reported in dB/octave. The higher the rolloff value, the greater the chance for resolving low power components far from the primary spectral peak.
The primary tradeoff is resolution against dynamic range, and this is determined by the main lobe width. The secondary tradeoff is resolution of adjacent peaks versus widely-spaced peaks. This is determined by the sidelobe and rolloff characteristics.
In the graph, the white line represents the frequency response from a rectangular window (no tapering window is applied). The frequency resolution is good but the dynamic range is very limited. The blue curve is the frequency response for the csx maxRolloff window with an adjustable width of 4. The other extreme is represented in the green curve where a Chebyshev window with an adjustable sidelobe of 103 dB is set.
The Chebyshev window produces the sharpest spectral peak for a given sidelobe level. It is an excellent tapering window for resolving closely spaced components of different power levels. The csx maxRolloff window produces the fastest decay in the frequency domain making it possible to resolve very low power peaks located some distance from the primary spectral peaks. The Kaiser-Bessel window is a very good adjustable window that offers an effective compromise between the two extremes.
After exploring the different AutoSignal tapering windows, close this Explore Data Tapering Windows option using the green check OK button. These settings will be preserved. You are back in the Fourier Spectrum with Data Window option.
Click the Display Maxima button until no labels are displayed.
AutoSignal offers peak-based critical limits for ascertaining whether or not a signal's peak is statistically significant. The reference is normally white (normally-distributed) noise, although a red noise background can also be approximated with an order one AR (autoregressive) model. Red noise is present when the background power decreases with increasing frequency.
A 95% critical limit is that level where in only 1 of 20 separate random noise signals would the largest peak present achieve this height strictly due to random chance. A 99.9% critical limit is similarly that level where in only 1 of 1000 separate random noise sets would the largest peak attain this height.
AutoSignal's critical limits are based on extensive Monte Carlo trials using the actual algorithms and adjustments in the program. In this instance, separate critical limit models are used for each of AutoSignal's data windows. For adjustable windows, 3D surface models are used. Since this approach incorporates the peak height attenuation that occurs when data windows are used, the critical limits are accurate for all n, all data windows, and all data window adjustable parameters.
Click the Show Significance Levels button at the far right in the graph's toolbar.
Five different levels are plotted. By default, the 50% levels will be white, the 90% levels will be cyan, the 95% levels will be green, the 99% levels will be yellow, and the 99.9% levels will be red. In this data set, the four primary reference peaks are shown to be significant in excess of a 99.9% confidence. The three low power peaks are not statistically significant relative to white noise.
Extracting The Three Target Components By Eigendecomposition
Is there another way to extract quantitative information associated with the three low power peaks without using a windowed FFT? Is there some way to determine if the three low power peaks would be statistically significant if the four primary power peaks were mathematically removed? Indeed there is. In fact, it is possible to isolate any one of the three low power spectral components for independent analysis.
AutoSignal offers four local options in most of its procedures. These are options which alter the data table only for the duration of the existing procedure (the main data table is not modified). One of these local procedures is an Eigendecomposition option.
Eigendecomposition offers the means to filter signal components based on signal strength. Although the four reference peaks each have identical power, there is a wide gulf in power between these four peaks and the three low power components. Eigendecomposition cannot isolate one of the four reference peaks, but it is able to effectively isolate one or more of the low power peaks.
Click the Apply Eigendecomposition and Reconstruction to Local Copy of Data button. Be sure the CovM FB (forward-backward covariance matrix) algorithm is selected and that the Normalize box is checked. Right click the Order field and select an order of 60. Click the FFT display option.
The upper graph contains the normalized eigenvalues from the decomposition. The lower graph displays a windowed FFT of the reconstructed data using a Kaiser-Bessel window whose width is automatically determined. Initially, all of the eigenmodes are selected which results in an exact reconstruction of the input data.
Using the left mouse button, box the first eight eigenmodes in the upper graph.
The first eight eigenmodes will be shown circled and the following FFT spectra are displayed:
The three low power sinusoids have been completely removed from the reconstructed data. A significant measure of the Gaussian noise has also been filtered out.
Again using the left mouse button, box eigenmodes 9 and 10.
This recovers the -40dB component.
Box eigenmodes 11 and 12.
This isolates and reconstructs the -50dB component.
Box eigenmodes 13 and 14.
Eigendecomposition is even able to recover the -60dB component. For the purpose of this tutorial, we will reconstruct all three low power components.
Box eigenmodes 9 through 14.
This recovers the three low power components. The four reference sinusoids are completely removed.
Close the Eigenfiltering window using the OK button to update the Fourier spectrum with the filtered data. Set the signal count in the sig field to 3.
The Eigenfiltered Fourier Spectrum
The three low power peaks are now recovered with good frequency estimation accuracy and a respectable power estimation. The values you find may vary appreciably as a consequence of the random noise component, especially with the lowest power peak. Further, the first two spectral peaks are shown to be highly significant statistically relative to white noise having the same variance as this reconstructed signal. The -60dB peak cannot be considered statistically valid in this particular reconstruction.
Change the Window to None and Nmin to 1001 in the drop down box.
Absent the four major reference components, and their contributions to spectral leakage, the three lower power spectral peaks that comprise the filtered data are now recovered using a non-windowed FFT.
Parametric Sinusoidal Modeling
We will conclude this tutorial by addressing how accurate the amplitude, frequency, and phase values are from the filtered FFT. We will then perform linear and non-linear parametric sinusoidal fits to see if these estimates can be refined.
Change the Window back to cs4 BHarris min and again set Nmin to 8192.
Click the Numeric Summary button. Be sure the first four items are checked in the Format menu. Inspect the Frequency Analysis section of the report.
If the original data (before eigenfiltering) were free of noise, the reference and low power components would have the following characteristics:
Note that the total power in the three low power components is only 0.0028% of the total signal power. The Frequency Analysis thus reflects not only the estimation accuracy of the FFT, but also the accuracy of the eigendecomposition in extracting the components associated with this tiny fraction of the overall power.
All of the information in the Frequency Analysis is drawn directly from the FFT. In this instance, accurate frequencies are realized for the first two spectral peaks. The frequency for the lowest power peak is not estimated as accurately. The phase values are quite good for all three peaks. The amplitudes are lower, however, as the power of each component is only about half of that known to be in the signal. The results you see may vary significantly due to the random noise component.
Inspect the Sine Component Fit section:
This section presents a linear sinusoidal fit that uses the estimated frequencies from the FFT in a least-squares procedure that determines the amplitudes and phases. The amplitudes and phases thus derive from this separate time-domain curve-fitting procedure. In this case, the two components are now very accurately characterized. The least-power component is not estimated as accurately. Note also the very high r² correlation coefficient value. An r² of 1.0 is a perfect fit. The fit statistics are important in confirming that a harmonic model is valid in representing signal data.
Close the Numeric Summary window and click the Non-Linear Optimization button. Accept the default settings and click the OK button. At the conclusion of the iterative fit, click the Review Fit button.
Due to the high frequency spectral content, the time domain sinusoidal fit plots are not particularly helpful.
Click the Numeric Summary button in the Non-Linear Optimization Review. Inspect the fit results.
|r² Coef Det
||DF Adj r²
||Fit Std Err
When the frequencies are also allowed to vary in the non-linear sinusoidal fit, the fit is further refined and even the least power component is respectably estimated.
Close the Non-Linear Optimization Numeric Summary window. Click the OK button to close the Non-Linear Optimization window. Click OK to close the windowed FFT procedure.