file: d:\jon\wpfiles\wrkshop 6-1-01 Talk for Single dish workshop, Arecibo, June, 2001 ed. 10-29

Spectrometry and Autocorrelation

Jon Hagen, NAIC

Abstract

The spectral power density, at a frequency f0, of an electrical signal is simply the average power per unit bandwidth when the signal is passed through a lossless bandpass filter centered at f0, in the limit that the bandwidth goes to zero. In radio astronomy, however, it is common to use a digital processing method whereby uniformly-spaced samples of the signal are used to estimate an average of the signal's autocorrelation function. This function is then Fourier transformed to produce an estimate of the spectral power density. In the limit that the averaging time goes to infinity, this is an equivalent definition of the spectral power density. The equivalence is reviewed below.

1. Operational Definitions of Spectral Density

Just as the light entering an optical telescope contains a spectrum of optical frequencies (colors), the electromagnetic waves intercepted by a radio telescope also span a band of radio frequencies. Analyzing, i.e. displaying, the spectrum allows us to see spectral lines and study other spectral structures, as well as measure their Doppler shifts. Light is commonly analyzed using narrow transmission filters. In radio work we do not have to put filters over the aperture of the radio telescope. The antenna converts the incident E field (volts/meter vs. time) into a voltage at the antenna terminals (volts vs. time), and we can install an electrical bandpass filter anywhere along the line connecting the antenna to the receiver, or at any point in the receiver before the detector. If we use a tunable narrow bandpass filter and a square-law detector or a bank of filters, each with its own square-law detector, we can measure the incoming power vs. frequency by averaging the output power of each detector. (Remember that power is proportional to the square of voltage). If a lossless bandpass filter centered at f has a bandwidth of B Hertz and its power output is P watts, the input spectral density at f0 must be S(f) = P/B watts/Hz. This is an operational "filter" definition of power spectral density. The basic bandpass filter/squarer/averager is shown below in Figure 1.

Figure 1. Basic filter channel used to estimate spectral power density at one frequency.

2. Spectral Power Density Defined in Terms of the Signal's Time History

Since the signal from the antenna is simply a time-varying voltage, the spectral power density function must also be definable directly in terms of observed sequences of voltage vs. time. Such a definition is needed in order to implement digital processors for spectral analysis. Let us jump ahead and first present the definition indirectly, in the form of a recipe for doing spectral analysis via the signal's autocorrelation function. This function, R( ), is defined as the time average of the "lagged products" V(t)V (t- ), i.e.

R( ) = < V(t) V(t- )> 1)

where the brackets < > denote averaging over t. If we have sequences of receiver output voltage values, we can calculate an estimate this function, R( ). We can also build hardware autocorrelators that operate in real time to simultaneously calculate estimates of R( ) for many successive values of .

The spectral power density is the Fourier Transform of the autocorrelation function:

2)

In autocorrelation spectrometry, a hardware autocorrelator does the heavy computing work. After a suitable integration time, it dumps the lagged product averages into a general purpose computer which, with a single Fourier transform, converts this estimated R( ) into an estimate of the spectral density function. The frequency resolution is inversely proportional to max , the maximum delay for which R is estimated.

Note that the inverse of the transform in Equation 2 is given by

3)

To see that Equation 2 (which is often used as the definition of spectral density) agrees with our operational understanding of spectral density, the following relation can be derived (see Ref. 1, for example) by substituting Equation 1) into Equation 2):

4) where 5)

In engineering terms, Equations 4 and 5 can be seen as a frequency conversion (the multiplication by ej2 ft in Equation 5) which shifts the spectrum so that the frequency in question becomes dc, followed by a "box-car" lowpass filter (the integration from -T to T), which passes only a narrow frequency range near dc. Successive evaluations of the integral in Equation 5 produce a sequence of random numbers (random and with zero mean, no matter how large T becomes), equivalent to a sequence of output voltages from a narrow bandpass filter. These are squared and averaged in Equation 4 to estimate the power output of the filter. Finally, the T to infinity limit makes the analyzing filter infinitesimally narrow around the frequency f. These operations are shown below in Figure 2. Note that Equations 4 and 5 provide a recipe for FFT spectral analysis, except that a finite T would be used; its value determines the spectral resolution. Likewise, when the circuit of Figure 2 is used for the channels of an analog filter bank, the lowpass filters have finite width.

Figure 2. Circuit equivalent to the operations specified by Equations 4 and 5.

3. Equivalent circuit using one bandpass filter rather than two lowpass filters

While analog filters and filter banks are sometimes implemented at base band, as shown in Figure 2, it is more common to use bandpass filters. Consider the circuit of Figure 3, which is equivalent to the circuit of Figure 2:

Figure 3. An alternate circuit, equivalent to the circuit of Figure 2

The circuitry in the dotted box is equivalent to a bandpass filter.

To see the equivalence, note that the baseband outputs (signals than extend down to dc) of the two low pass filters are converted back up to the original center frequency, f, to produce the signals I cos(2 ft) and Q sin(2 ft), where I(t) and Q(t) denote the slowly varying outputs of the lowpass filters. These two signals are added to produce I cos(2 ft)+ Q sin(2 ft). The square of this voltage contains the term I2 cos2(2 ft) + Q2 sin2(2 ft), whose average, ½(< I2 > + <Q2 >), is the same as the output of the circuit of Figure 2. The other terms in the square contain double frequency components that average to zero. The components inside the dotted box in Figure 3 are equivalent to the simple bandpass filter of Figure 1, as can been seen by considering the overall response at any single frequency.

4. Properties of the Autocorrelation Function

Given that the Fourier transform of the autocorrelation function is, indeed, the power spectral density, consider the following properties: The "zero-lag" value of the autocorrelation function, R(0), is just the average square of the voltage, i.e. the power. No lag will have an absolute value greater than the value of the zero lag. The autocorrelation function, R( ), is real and even (as is the spectral density function) and contains all the information available in the signal. Higher moments of the received signal (which is a Gaussian process, as explained below), contain no additional information, because they can be expressed in terms of the second moment function, R( ). For example,

<V(t)V(t- ) V(t- ')V(t- ''> =

<V(t)V(t- )>< V(t- ')V(t- ''> + <V(t)V(t- ')>< V(t- )V(t- ''> + <V(t)V(t- '')>< V(t- )V(t- ')>. 5)

5. Gaussian Nature of the Signal

In radio astronomy, the incoming signal will always be a random process whose amplitude has a zero mean (i.e. no dc) and whose amplitude probability distribution is Gaussian. That is, if we take a series of independent voltage samples (samples taken with adequate spacing between them) and plot a histogram, normalized to have unit area, the result will be

6)

where 2 is the average square of the voltage. This always-Gaussian amplitude distribution results from the fact that the signal voltage is always the superposition of voltages from a great number of individual physical radiators. The Central Limit Theorem states that if a random process is a sum of random processes with any given amplitude distribution or distributions, then the sum process will have a Gaussian amplitude distribution. You probably have a computer program that can generate random numbers, xi, uniformly distributed between -1 and 1. You can easily show that the average value of x2 is 1/3. Sum 100 of these numbers to create a new random variable, y. Perform this sum a few thousand times and plot a histogram of the values of y. You will see that it has a Gaussian shape. The point is that you could have predicted this Gaussian by calculating the expected value of y2 (the variance of y). To calculate <y2 >, write <y2 > as a double sum over i and j. The only terms that have non-zero expectations are the terms where i=j. There are 100 such terms, so <y2 > = 100<x2 > = 100 x 1/3 = 33.3. Use this value, 2 = 33.3, to plot the Gaussian distribution function of Equation 6 to see how closely it agrees with the histogram.

Since our signals always have a Gaussian amplitude distribution, it is not absolutely necessary that we average the square of the voltage if we need to estimate the power. We might average the absolute value of the voltage, since its expectation (using Equation 6) is given by (2/ )1/2 . Or we might just measure the fraction of the time that the absolute value of the voltage is greater than some threshold, , since the expectation of this fraction is easily calculated (again, using Equation 6) as a function of . In this case, we are able to make an accurate estimate of the power, even though the samples of the signal have been quantized very coarsely to only three levels. As you might expect, it turns out that the autocorrelation function can also be accurately estimated when the signal has been sampled by a digitizer of very coarse precision. The extreme case is two-level (one-bit) sampling, where only the polarity of the voltage is measured. In this case

8)

This equation, containing both the trigonometric sin function and the algebraic sign function was first derived by J.H. VanVleck (2), who investigated the spectrum of clipped noise during WWII. Jamming transmitters could be much more efficient if they were allowed to clip (limit) the peaks of the transmitted noise. VanVleck showed that communication and radar receivers could be jammed as effectively with clipped Gaussian noise as with undistorted Gaussian noise.

The advantage of using coarse quantization is hardware simplicity. The disadvantage is the addition of quantization noise, requiring longer integrations to produce estimates of R with statistical fluctuations no greater than if precise many-bit quantization had been used. Another minor disadvantage, in the case of one-bit sampling, is that all amplitude information is lost, so the power, R(0), must be estimated separately.

6. Hardware Autocorrelators

Evaluating a lagged product requires only a multiplier and averager (adder/accumulator). If the signal has one-bit quantization, the one-bit x one-bit multiplier can be a simple exclusive-or logic gate and the adder/accumulator can be just a digital counter. A multi-channel autocorrelator, as shown in Figure 4, consists of many multiplier/averager "lag units" plus a digital delay line so that each lag unit receives the present ('prompt") signal and a delayed signal. The delay line provides successively delayed versions of the prompt signal.

Figure 4. Autocorrelator Block Diagram

Three-level quantization(the digitizer output values can be assigned the values -1, 0 and +1) requires only a little more complexity and introduces less quantization noise than two-level quantization. Autocorrelator architecture is very simple, especially when coarse quantization is used. It can be expanded in several ways so that chips containing N lag units can be arranged to implement correlators with any number of lags and with overall effective clock speeds that are multiples of the chip's maximum clock speed. The FFT lacks this flexibility. Nevertheless, at some high number of spectral channels, the intrinsic higher efficiency of the FFT algorithm (NlogN operations vs. N2 operations, where N is the number of spectral bins) should make FFT spectrometers more economical than autocorrelation spectrometers.

References:

1. Thomas, John B. "An Introduction to Statistical Communication Theory", John Wiley & Sons, Inc. 1969

2. Van Vleck, J.H. and D. Middleton, "The spectrum of clipped noise," Proc. IEEE, 54, 2-19, 1966.

Appendix - Autocorrelation function and spectrum of a pure sine wave

The case of an infinitely narrow spectral line, a sine wave, requires some care. Let the input voltage be a sine wave at frequency f0: V(t)= sin(2 f0t). From Equation 1, it is readily shown that R( ) = ½ cos(2 f0 ). The power, R(0), is therefore ½, which was to be expected, since

< sin2(2 f0t)> = ½. If we calculate the spectral density by substituting this autocorrelation function into Equation 2, we get infinity; the line is infinitesimally narrow so the spectral density at that point goes to infinity. To approach this infinity in an orderly way, we would give our sine wave some spectral width, via modulation controlled by an adjustable parameter. Then we would calculate the autocorrelation function and its transform, the spectral density. Finally we would adjust the parameter to narrow the line and find that the spectrum approaches a sum of two delta functions. This result is easy to guess: since the total power is ½, the spectrum must be given by S(f) = (-f0)/4 + (f0)/4. If we substitute this spectrum in Equation 2, we do indeed get R = ½ cos(2 f0 ).