MISCELLANEOUS TECHNICAL ARTICLES BY Dr A R COLLINS

JavaScript Spectrum Analyser

Spectrum analyzer basics

A Spectrum Analyzer is designed to capture and digitize time varying signals such as the output from a microphone. It then performs a Fast Fourier Transform (FFT) on the captured data and displays the resulting frequency spectrum.

Here is a simulated Spectrum Analyzer, it accepts a mathematical expression in JavaScript, describing the input signal waveform as a function of time. This signal optionally has an anti-aliasing filter applied and may also be windowed with Hanning or Kaiser-Bessel window. The FFT is then calculated. The output spectral data are complex values representing the RMS (root mean square) magnitude and average phase of the frequency components that fall within the bandwidth of each frequency bin. These may be plotted as Real and Imaginary components, or Magnitude-Phase or Power. The Power display plots the magnitude squared value. The logarithmic plots of magnitude or power convert their respective values to dB re 1V rms, both return the same since magnitude in dB is 20*log of the value, and power in dB is 10*log of value.

The Anti-aliasing filter was designed using the Filter Design page, more detail is presented below. The time domain windows and the FFT are written in JavaScript, all the computations are done within the browser. The data are plotted using the Cango canvas graphics library.

Simulated input signal

The input waveform function may be any suitable JavaScript expression of the single variable \(t\), representing time. The time variable must be written as 't' and the time data will be generated by using the JavaScript eval function for 512 equally spaced values of t from 0 to 511*dt, where dt = 1/Fs. For brevity, wrapper functions for the Math library methods are provided, Math.PI may be entered as PI, Math.random() may be entered as random() and so on.

The sample frequency may also be changed as required, it has been arbitrarily limited to values between 1mHz and 1MHz.

The following JavaScript expressions may be copied and pasted to the input as examples:

  square pulse: ((t>0.03)&&(t<0.04))?1:0
  triangle:     2*(abs(64*t%1-0.5)-0.25)
  square wave:  (sin(2*PI*80*t)>0)?1:-1
  amplitude modulation: (1+0.3*sin(2*PI*60*t))*sin(2*PI*400*t)
  sine wave + noise:  random()+random()-1+0.3*sin(2*PI*270*t)

With the expression entered, check the anti-aliasing 'on' or 'off' and select the type of windowing to be applied. Windowing is recommended for continuous signals, and not for transients whose duration is within the 512 sample periods used in for generating the spectrum. Press the 'Generate Time Date' button to see the time domain data. The resulting spectrum of the input can then be viewed in various formats, Magnitude-Phase, Real-Imaginary or Power Spectrum by pressing the corresponding button. For the Magnitude-Phase and Power Spectrum plots, logarithmic or linear scaling are available by clicking the corresponding checkbox.

A measuring marker is available for the frequency domain plots. Press 'Mkr' and a small square will appear on the spectrum waveform, rotating the knob will move the marker left or right and the corresponding waveform frequency and amplitude will be displayed at the top of the plot.

Time Domain
Frequency Domain

Notes on JavaScript spectrum analyzer

Aliasing

A real time varying signal has a spectrum that is symmetric around 0Hz. If the spectrum has a component line at 300Hz then there will be a corresponding frequency line at -300Hz. Digitizing a real, time domain signal at a sampling frequency \(F_s\) generates a time domain waveform that is a series of delta functions (spikes) at intervals of \(dT = 1 /F_s\) the envelope of these spikes reproduces the original waveform. The spectrum of this digital series is a series of copies of the original spectrum. Spectrum analysers display only the positive frequencies of the copy centred on 0Hz.

If the sampling frequency is too low, then the negative frequency side of each spectrum will overlap the positive frequency side of the adjacent copy of spectrum. This overlap results in spurious frequency components appearing in the waveform spectrum. This effect is called aliasing.

For the frequency spectrum of a waveform to be accurately described by the FFT, the input signal should not have any frequency components higher than half the sample frequency. Any spectral components of the input waveform with frequency higher than \(F_s / 2 \), known as the Nyquist rate, appear to fold back into the lower frequency part of the spectrum. They are just the negative frequency aliases of the signal's frequency components that lie above the Nyquist rate from the next highest frequency copy of the waveform's spectrum.

Anti-aliasing Filters

Spectrum Analysers avoid displaying aliased frequencies by using a low-pass, 'Anti-aliasing' filter, which is designed to remove all the frequency components of the input signal that have frequencies above the Nyquist rate before the signal is sampled. Ideally the anti-aliasing filter would pass all frequencies up to \(F_s / 2 \) and then its response would drop to zero after that. Filters can't be made infinitely steep, so the cutoff must start at a lower frequency and have a steep enough roll off to so that signal above the Nyquist rate is sufficiently attenuated so that when aliased back it is below the dynamic range of the display at the cut-off frequency of the filter.

Anti-aliasing filters typically start about 80% of the Nyquist frequency and all frequencies above this value are considered corrupt. Thus only the first 400 bins of a 512 bin frequency spectrum are usually displayed. In the early days of digital signal analysis Analogue-to-Digital Converters (ADCs) were much slower than today and they set the upper limit of sampling rate. Over sampling a signal to reduce the complexity and tight tolerance on the multi-pole anti-aliasing filters was out of the question. Modern ADCs are so fast that the waveform may be over sampled at frequencies many times the displayed frequency and simple 2 or 3 pole anti-aliasing filters are all that is required to produce a digital signal with the low frequency part of the spectrum free from aliases. A second digital filter can then do the hard anti-aliasing filtering required prior to sub-sampling at the desired value of \(F_s\). The advantage of using the digital filter is that it has a perfectly reproducible response within the mathematical precision used, unlike the components of multi-pole analogue filter. As an example, using 128 times over-sampling, a simple, low tolerance two pole filter which attenuates aliased signals by 12dB per octave, will give 96dB rejection at the Nyquist rate.

JavaScript anti-aliasing filter

Anti-aliasing filters must be applied before the input signal is digitized, since it is the digitizing that causes the aliasing. This JavaScript Spectrum Analyser has its input created as discrete samples so there is no opportunity to use an analogue anti-aliasing filter. The anti-aliasing filter used here over-samples by a factor of 16, but instead of an analogue filter, it relies of the nature of the synthesized signals effectively providing a 1 pole (6dB per octave) low pass filter. The idea behind this is that all repetitive waveforms, have harmonics that fall in amplitude by at least \(1/n\), i.e. the 10th harmonic will have amplitude 1/10th of the fundamental. Since signals to be studied in this JavaScript Spectrum Analyser should have their fundamental frequency in the displayed spectrum, i.e. below 0.4*\(F_s\) (80% of the Nyquist rate), thus the amplitude of the first harmonic that aliases below 0.8 of the Nyquist rate will be (16*\(F_s\))/(0.8*\(F_s / 2\)), the 40th harmonic. This has amplitude 20*log(1/40) = -32dB with respect to the fundamental. Signals with lower frequency fundamentals will have greater rejection. This JavaScript spectrum analyser has a dynamic range limited to 40dB very little aliasing should be noticed if the anti-aliasing filter is used.

The first stage of the anti-aliasing filter delivers a digital representation of the waveform which has a spectrum with aliased components typically suppressed by 40dB for frequencies up to half the requested sample rate. The second phase of the anti-aliasing filter is to re-sample the waveform at the requested sample rate. This will cause the spectral components above the Nyquist rate to be aliased into the displayed spectrum, so a second digital anti-aliasing filter is used prior to the sub-sampling. This filter is a multi-pole (order 161) low pass filter with cut off frequency at 0.8*\(F_s / 2\) and very steep fall off. The filters transition band extends to the frequency 1.2*\(F_s / 2\), which will be aliased to 0.8*\(F_s / 2\). Signal level at this point should be below 40dB. This filter's coefficients were calculated using the Kaiser-Bessel windowing technique and the filter is applied to the over sampled digital waveform with the output simultaneous sub-sampled by a factor of 16 to the generate the requested sample rate.

To demonstrate, examine the spectrum of a 115Hz square wave with and without the anti-aliasing filter.

Set the sample frequency to be 2000Hz, this gives a Nyquist rate of 1000Hz, select anti-aliasing filter 'off', check the 'K-B' window then copy and paste the following input function which generates a 115Hz square wave.

(sin(2*PI*115*t)>0)?2:-2

Press 'Generate Time Data' and then press 'Power' to display the power spectrum, then check the 'log' display mode. As expected the frequency components of a square wave consist of every odd harmonic, 155, 345 (=3*115), 575 (=5*115). The harmonics above the Nyquist rate appear folded back into the displayed spectrum, they are actually the negative frequency harmonics from the copy of the spectrum centred at the sample frequency 2000Hz. The 11th negative frequency harmonic is 11*-115=-1265 this appears at 2000-1265=735, the 13th harmonic appears at 505, the 15th at 275, the 17th at 45Hz. It gets worse, the 19th harmonic from the spectrum copy centred on -2000Hz appears at 185Hz (-2000+19*115) and so on. Figure 1a shows the messy result.

Now check the anti-aliasing filter 'On' and press 'Generate Time Data', then press 'Power' and check the 'Log' checkbox. The anti-aliasing filter removes all the aliases and leaves only the harmonics that lie in the displayed frequency range, 115, 345 and 575Hz. This anti-aliased spectrum is shown in Fig 1b.

Aliasing

Figure 1a. Square wave spectrum demonstrating aliasing.

Aliasing

Figure 1b. Same waveform as Fig 1a but with anti-alias filter on.

Windowing

The Discrete Fourier Transform (DFT) is an approximation to the Fourier transform. It takes in a time record of fixed length T, and the output spectrum is the Fourier transform of this data block replayed end on end. Therefore, if the original signal was longer than T, then the repetition of a fixed length section will produce discontinuities at the abrupt transition from the signal approaching the back end to the restart of the signal at the front. If we were listening to the repetitive replay of the effective input signal, these discontinuities would sound like 'clicks' every T seconds.

Rectangular Window

The only part of the input signal processed is the part that fits in the 'window' starting at 0 of length T seconds. The time record has effectively been multiplied by the rectangular window function w(t) given by:

  w(t) = 0,  t<0 and t>T;
  w(t) = 1   0<t<T;

The Fourier transform of signals with sharp discontinuities contain broad frequency components and generate many side lobes of prominent spectral lines. This spread of energy reduces the peak value of spectral lines. If the waveform happens to join seamlessly at the ends then the power will all land within the appropriate frequency bin without this artificial splatter. The analyser measurement of peak spectral level of a narrow band signal will thus vary as a function of its frequency, depending on how smoothly the front and back ends of the windowed signal join. This reduction in peak spectral level is called 'scalloping loss'.

Scalloping can be demonstrated by measuring the peak spectral level of a sine wave signal with small changes to the frequency. Try an 80Hz sine wave input sampled at \(F_s\) = 2048Hz. This \(F_s\) value will give output frequency bins spaced RBW of 4Hz:

Input Function: 20*sin(2*PI*80*t)
Sample Frequency: 2048
Window: Rect

The logarithmic power display shows the signal has peak power 23.0dBV. Edit the '80' Hz frequency in the input function to be 82Hz, press 'Generate Time Data' and show the log power spectrum. The signal now falls half way between two frequency bins, 80Hz and 84Hz, the situation with most spectral splatter and largest peak power loss. The peak has gone down to 19.1dB, a loss of 3.9dB. Try other frequencies, say 83Hz, you should observe several dB variation in peak signal and varying amount of spectral splatter. The scalloping loss for the rectangular window varies from 0 to -3.9dB.

To minimize this variation caused by cutting off the time record sharply at the beginning and end, the input time data can be multiplied by a smoothly varying window function that will taper the signal down to zero at the beginning and end of the time record.

Hanning window

The Hanning window is simply one period of a sinusoid having a length equal to the time record, T, but its lifted to start and stop at 0 and its scaled to have a central maximum of 2.

The effect of the Hanning versus the Rectangular window may be seen by repeating the measurements of the sine waves but having the 'Hann' windowing button checked. The power at 80Hz is still at 23dB but the peak level at 82Hz should be about 21.6dB, a smaller peak loss than the rectangular window and greatly reduced spectral splatter. The maximum scalloping loss for the Hanning window is reduced to 1.4dB.

These improvements in measurement quality by using the Hanning window applies to continuous type inputs only. For transient signals whose full extent fits within the length of the input time record, there will be no discontinuities at the beginning and end so the rectangular windowing will have no effect on the output spectrum and the FFT gives an accurate Fourier transform. If a Hanning window is applied to this type of signal there will be little effect on the spectrum and if the compensation for the expected windowing loss is applied the answer would be less accurate.

Kaiser-Bessel Window

The Kaiser-Bessel window is another window shape designed to perform the same function as the simpler Hanning window but it has some advantages. The window generation code takes a parameter 'alpha', a number with a useful range of 1.5 to 4.0. This parameter can be used to set the maximum height of the side lobe level the window will generate. This is a very useful characteristic particularly in limiting the effects of aliasing. The alpha value here is preset to 2.5. The maximum scalloping loss for K-B window with alpha 2.5 is 1.2dB

Resolving spectral features

The display of spectral features on a Spectrum Analyzer is greatly enhanced by the use of windowing. As a demonstration try the following input signal, comprising to sine waves one at 162Hz and the other at 192Hz. Copy and paste the waveform expression into the Function input box.

10*sin(2*PI*162*t)+0.4*sin(2*PI*192*t)

Again use 2048Hz as the sampling frequency. Generate the time data with the Rect window selected. Observe the log power display, spectral splatter masks the lower level 192Hz signal. Next select the Hann window and press 'Generate Time Data'. Now when the log power is displayed, the peaks are easily resolved.

Fourier transform

The Fourier transform converts information from the time domain into the frequency domain or vice versa. No information is gained or lost in transforming from one domain to the other. The idea of Fourier analysis is to present the information in such a way that it is easy to interpret and facilitate solutions to problems.

The Fourier transform of a time signal \(a(t)\) defines the complex spectrum \(A(f)\) and is given by:

$$A(f) = \int^{+\infty}_{-\infty} a(t) e^{\displaystyle{-i2\pi ft}}dt$$

Discrete Fourier Transform

The approximation to the Fourier transform used in modern digital signal processing is called the Discrete Fourier Transform (DFT). The DFT operates on a series of data points obtained by sampling the input signal discrete intervals. The input series is of a fixed length, say \(N\) points. The DFT will generate a series of \(N\) frequency domain points which are a uniformly spaced samples of the spectrum of the input signal if it were repeated infinitely end to end.

A general expression for the forward DFT is given by:

$$A_n = \frac{1}{N}\sum^{N-1}_{k=0}a_k e^{-\displaystyle\frac{i2\pi kn}{N}} \qquad\text{for}\ n=0\ ...\ N\!-\!1$$

A typical application may be to determine the spectrum of a time domain signal. If a DFT takes an input of length \(N\), then input series \(a_n\) is generated by sampling the input signal \(N\) times at time intervals of \(dt\). The total length of the input record will be \(T=N*dt\). If the input signal is defined as \(f(t)\) then the input series will be:

$$a_0 = f(0),\ \ a_1 = f(dt),\ \ a_2 = f(2dt), \; ... \; a_N = f(N dt)$$

The sampling frequency, \( F_s \) is then defined as:

$$F_s = \frac{1}{dt}$$

The output series \( A_n \) represents the frequency spectrum of the input. The frequency spacing of the discrete points in the transformed data is termed the frequency resolution, or resolution bandwidth \(dF\), given by:

$$dF = \frac{1}{T}$$

where \(T\) is the length of the time domain record.

$$dF = \frac{F_s}{N}$$

The frequency data represent both positive and negative frequency components of the input signal. The frequencies generated range from \( -F_s/2+dF \) Hz through \(0\) Hz (DC) to \(F_s/2\) Hz.

The Fourier transform and the DFT may be equally well applied to waveforms domains other than time, for example if the input series samples distance, then the output data represents a spatial frequency response. A scenario often used in image processing and antenna beam forming.

Fast Fourier Transform

The Fast Fourier Transform (FFT) is just an algorithm which computes the DFT with a greatly reduced number of arithmetic operations compared to a direct computation.

The FFT algorithm used in this page is written in JavaScript, it is a radix2, in place, complex FFT. The transform size may be any power of 2, i.e. 8, 16, 32, ... 2048 etc. A length of 512 is used here. The DC value plus the 256 positive frequency points are plotted in the frequency spectrum displays. The full 512 time domain data points are plotted in the time domain display.