Electronic Journal «Technical Acoustics» http://www.ejta.org
2006, 7
Patrick J. Vitarius1, Don A. Gregory1*, John T. Wiley2, Valentin Korman3
lThe University of Alabama in Huntsville, Department of Physics, Huntsville, Alabama 35899
2Madison Research Corporation, 401 Wynn Drive, Huntsville, Alabama 35805
3Marshall Space Flight Center, Marshall Space Flight Center, Alabama 35812
Sampling rate error in acoustic measurements
Received 01.03.2006, published 04.04.2006
Two aspects of data collection and analysis are considered in the paper: the transfer function of the detector and the sampling method used on the data the detector reports. Following a brief look at transfer function theory, a simple model is constructed which shows the effect of sampling time dependent functions (acoustic or otherwise) at different rates. The average value of a time dependent parameter (pressure for example) is calculated to illustrate the analysis method. Four different type functions were chosen to represent the parameter: sinusoidal, pseudo-sinusoidal, asymmetric triangular, and random. The results illustrate the important role played by sampling rate when analyzing time dependent data.
INTRODUCTION
Data from any measurement is affected and altered by everything that stands between the actual event and the eventual recorded data; all of these components may be jointly called ‘the detector’. For example, suppose a pressure reading is obtained using a transducer with a long response time [1]. If such a sensor is used to measure a rapid, transient pressure change, then the shape of the data curve will differ significantly from that of the true pressure spike. The sensor will still produce a measurement because a pressure change will have occurred, but exactly what that measurement means will be debatable. It is not the peak pressure but a progressive time average that will be seen and often the peak will be orders of magnitude larger than the average. While the area under the two curves may be nearly the same, the fact that the actual peak pressure was much higher than the recorded peak pressure may lead to a system being under-engineered or exposed to pressures too high for its components. The response time of the sensor must be much faster than the event in order to accurately capture the data. Thus the temporal qualities of the detector, as well as the detection process itself, are both of considerable importance.
1. THEORY
Transfer functions describe what happens to a signal as it is transmitted through a system. The term is used extensively in acoustic signal processing and optics applications [2]. If the system is linear, as most optical imaging systems and electrical systems are when operating under normal conditions, the total output from the system can be written in terms of the effect
Corresponding author, e-mail: [email protected]
the system had on individual input signals that made up the total input. Let the input and output signals be given by gi and go, and let the system operator be designated as S :
The input function can be rewritten by applying the sifting property of Dirac delta functions [3]:
The operator S only acts on functions of t, and the integral is over all values of t, so the order of these two can be reversed:
which is independent of the input function g(t). This function h(t, t) is called the impulse response function. With h thus defined, Equation (3) can be written
For systems or components whose characteristics do not change implicitly with time, h(t, t) is in fact only a function of the difference between t and t, i.e.,
This is immediately recognized to be the convolution of gj(t) and h(t). From the convolution theorem, the Fourier transform of a convolution is equal to the product of the Fourier transforms of its constituents. Thus
where Go(a>) and Gt(a) are the Fourier transforms of the output and input functions and H(a>) is actually also a Fourier transform: it is the Fourier transform of the impulse response function. This is the response of the system to a delta function input [5]. The concern in this investigation is in developing a technique for calculating the transfer function based on accurately sampling the raw data from the sensor. If the input could be true white noise, which ideally contains equal amounts of all frequencies, then the transfer function could easily be determined. A truly flat spectrum input signal is impossible to create, and in this
go (t) = S{gi (t)}.
(1)
(2)
(3)
and the operator acting on a delta function input can be defined as [4]
h(t, t) = S{5(t -t)},
(4)
go (t) = j gi (T)h(t,T)dT
(5)
h(t, t) = h(t - t) ,
(6)
so
go (t) = j gi (T)h(t -T)dT
(7)
Go (©) = H (ffl)G (©)
(8)
or
(9)
investigation it was found that a sum of equal amplitude sinusoidal input functions served as a reasonable substitute. All of the output functions, one for each input, were then summed and Fourier transformed. The Fourier transform of a sum is the sum of the Fourier transforms, so that the equation for H(a>) is still valid. Examples are given below for fabricated inputs just to illustrate the technique.
2. SIMULATIONS
Consider an input function formed by summing sine functions of many different frequencies (Figure 1a). In frequency (Fourier transform) space, the power spectrum (magnitude squared) of this function is a sum of delta functions (a comb function), as shown in Figure 1b. Now consider a representative transfer function consisting of a tapered filter with an element of added noise (Figure 2a). This may be represented mathematically as the sum of a Gaussian function centered at zero and a low-amplitude random signal. The output signal spectrum will be the product of the Fourier transform of the input signal and the transfer function. The output signal itself is obtained by Fourier transforming the output signal spectrum. The subtle but normally important differences between the definitions of the Fourier transform and the power spectrum are not applicable in this case because of the functions chosen. The magnitude squared of the Fourier transform of the product of the signal in Figure 1a and the transfer function of Figure 2(a) is still a series of delta functions as shown in Figure 2(b). Comparing the input function of Figure 1(a) to the output function of Figure 2(c), the effect of the transfer function is apparent. The general shapes are the same but the high frequency ripple in the input has been lost. This is not noise elements that have been removed from the signal, but rather high frequency components of the signal itself that have been lost. The introduction of noise also results in a loss in symmetry in the minor peaks and in a shift of the locations of the major peaks from Figure 1(a).
(a)
0 100 200 300 frequency 500
(b)
Figure 1. (a) A sum of sine functions of many different equally spaced frequencies representing a temporally dependent input.
(b) The Fourier transform of the function given in (a)
6--------.-------.--------.--------.-------
5 - -
4 - -
3 - -
2 - -
1 - -
0 -------L-------L--------L-1 I I il .--lL
0 100 200 300 frequency 500
(b)
Figure 2.
(a) An example transfer function of the system, versus frequency.
(b) Magnitude squared of the Fourier transform of the output function, obtained by multiplying the transfer function by the Fourier transform of the input function (Equation 8).
(c) The output signal versus time, obtained by Fourier transforming the function in (b) above
(c)
3. DATA ANALYSIS
Digitally analyzing any analog signal first requires that signal to be sampled. Accurately calculating the transfer function depends on sampling the time dependent input and output signals. The classic sampling theorem of Nyquist is routinely quoted and is in fact often sufficient for sinusoidal periodic data streams [6]. If the signal is periodic but not sinusoidal, the rules are still strictly applicable, but the application is not always useful [7]. The Nyquist sampling theorem states that at least two points per period are necessary to reproduce a sinusoidal function. Fourier analysis is based on the notion that any function can be written in terms of sine and cosine basis functions.
It seems reasonable then to Fourier transform the time dependent signal and look for the highest non-negligible frequency present and sample at a frequency of at least twice that value. The problem comes from the fact that, at these high frequencies, there is practically no signal, so it is pointless to sample the signal at a high frequency when the amount of information at that frequency is very small. A series of functions, not necessarily representing realistic data, is chosen here to illustrate what can happen when data is not sampled at the correct rate.
(a)
0.8 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.6
If a physical parameter, such as pressure, is being monitored by an analog device, such as a pressure transducer, and the output voltage is digitized, then it is very important to know at what rate the digitization takes place. The true time dependence is not always displayed but rather an average is often computed and displayed. In fact, all sensors have characteristic time constants that must be accounted for in data interpretation. The averaging is done digitally or electronically with the use of a low pass filter. This average value is then taken as the value of the pressure, when in reality it probably should not be used, at least not for important decision making purposes.
The functions chosen here for illustrating the analysis are simple and are constructed to intentionally have a zero average value. The functions are then sampled at various frequencies and the average of the sampled values computed. The difference in this value from zero is an indication of the error that would be inherent in the average. The errors are expressed in terms of the unit amplitude of the periodic signal, or the unit standard deviation of the random signal. There are some surprises in the results and some things are as expected. The high frequency sampling produces the lowest error as expected, but there are other frequencies that can produce low errors as well. A random signal was also chosen for analysis. The signal was constructed so as to have zero mean, just as the periodic signals mentioned previously.
To illustrate the importance of sampling, four representative functions were chosen. These functions would be functions of time if they represented real data from, for example, a pressure transducer. The first function is a single-frequency sinusoidal function. The second function looks sinusoidal but it is not. It is actually constructed from parabolas. The third function is a normal triangular periodic function except that it is slightly skewed in one direction. The final function is Gaussian pseudo-noise drawn from a normal distribution using Marsaglia's ziggurat algorithm [8]. The periodic functions each have a fundamental frequency of 0.2 Hertz and unit amplitude; the normal distribution has a standard deviation of one unit. A sketch of these four functions is given in Figure 3 below.
Sinusoidal Signal Pseudo-Sinusoidal Signal
0 5 10 0 5 10
Time (sec) Time (sec)
Figure 3. Representative functions used in the model: sinusoidal, pseudo-sinusoidal, skewed triangular, and random
4. SIMULATION RESULTS
The results of the sampling investigation for the three representative functions are shown in Figure 4 below. All four functions in Figure 3 are designed to have an average value of zero over sufficiently large time intervals. The plots in the figure represent the error, expressed in terms of the unit amplitude (or unit standard deviation for the random function), in the calculated average of the function versus the number of samples taken of the function. In figure 4(a), as expected, the error is generally large for sampling less than one per cycle. However there are regions between 0 and 1 sample per cycle where the error is quite small. The error is very large at 1 sample per cycle. This actually is expected because the first sample was chosen to occur at the first maximum of the function, and for 1 sample per cycle, each successive sample chosen would also be a maximum. The average calculated using these samples should have the largest possible error. This first sample rule (phasing) was used for all sampling cases. Where the sampling actually begins can significantly affect the result, particularly in the low sampling frequency range.
It was expected that a sampling rate of 2 per cycle would produce a low error but this was not the case. The Nyquist criteria really only applies to truly sinusoidal functions, or to the highest frequency contained in nonsinusoidal functions. The error for the pseudo-sinusoidal function (figure 4b) looks similar to that of the skewed triangular function in the 0 to 1 sample per cycle range. This function was similar enough to a true sinusoid however, that when the function was sampled at the Nyquist frequency (2 samples per cycle), the error was very small. The 2 samples per period criteria actually corresponds to sampling at one of the lowest frequencies in the spectrum of the triangular function (figure 4c), not at the highest.
The random function sampling was different from the previous two in that it was not periodic or sinusoidal and thus could not be analyzed in terms of number of samples per cycle. It was examined however just for general completeness. The total number of random data points used in constructing the random function was 100,000 and the sampling rate was done in terms of samples per 5000 of these original data points. The general behavior is the same as the earlier functions. As the number of sampled points per 5000 increases, the error decreases as expected. There are zeroes and near zeroes in the result but they occur randomly, even in the low sampling rate region.
Figure 4. Error (in the same units as amplitude) versus number of samples per cycle for the sinusoid (a), pseudo-sinusoid (b), asymmetric triangle (c), and random function (d).
The mean error for the last plot is versus the number of samples taken per 5000 data points since the function is not periodic. The error for the Gaussian noise function exceeds unity at
some very low sample rates
Figure 5 explores the importance of phase in the offset of the initial term. In figure 5(a), it can be seen that the maximum positive errors for the sinusoidal function occur for zero initial phase offset, and the minimum negative errors occur at a phase offset of n. Figure 5(b) explores the structure of one peak in greater detail. Figure 5(c) examines the phase dependence for the asymmetric triangle function. The peak that occurs at two samples per cycle is particularly interesting, as it represents the overlap of one sample per cycle on the fundamental frequency and half a sample per cycle on the first harmonic.
a)
b)
c)
Figure 5. Contour plots show the dependence of the mean value error upon the phasing term for the sinusoidal function (a) and a close-up of one peak (b), and the asymmetric triangle function (c). Appreciable errors at larger sampling rates are encountered in the asymmetric triangle function, which has non-negligible higher harmonics
CONCLUSIONS
The simple models created in this investigation have drawn attention to some of the data analysis procedures utilized in acoustic systems testing. Transfer functions, even simple linear ones, must be determined carefully if data collected from an experiment is to be believed. Determining transfer functions accurately depends on careful sampling of the raw laboratory generated data using the actual sensor to be used in the field.
It is obvious that the sampling rate is not a trivial decision, and that it must be determined based on the requirements of the test being done. For a true sinusoidal function the Nyquist criterion is sufficient to reproduce the function with zero error as expected. It was found however that the condition of periodicity alone is not sufficient for Nyquist sampling to reproduce the function exactly. The sinusoidal function, and two periodic, but not sinusoidal, functions were examined and errors in calculating the average value of the function were found to be as large as 8 percent of one amplitude, even if the sampling was 3 per cycle. The random function showed serious sensitivity to sampling as well. An error of about 20 percent of one standard deviation is seen for sampling rates of about 5 per 5000 points.
REFERENCES
1. Halliday D., Resnick R., Walker J. Fundamentals of Physics, John Wiley and Sons, 5th ed., 1997.
2. Ono T., Yamaguchi K., Suzuki H., Kawaura J., Ohashi, M. A New Error Index for the Determination of the Number of Averages in Transfer Function Estimation. Appl. Acoustics, vol. 28, issue 1, 1989, p. 21.
3. Poularikas A., Seely S. Properties of Delta Functions. Signals and Systems, 2nd ed., Krieger Publishing, Boston, 1991, pp. 690-694.
4. Goodman J. Fourier Optics. McGraw-Hill Publishing, New York, 2nd ed., 1996, pp. 2122.
5. Hecht, E. Optics. Addison-Wesley Publishing, 4th ed., 2002.
6. Randall R., Tech B. Application of B&K Equipment to Frequency Analysis. Bruel & Kjaer, Denmark, 2nd ed., 1977.
7. Yamaguchi S., Kato Y. A Statistical Study for Determining the Minimum Sample Size Leq Estimation of Periodic Nonstationary Random Noise. Appl. Acoustics, vol. 32, issue 1, 1991, p 35.
8. Moler Cleve. Ziggurat Algorithm Generates Normally Distributed Random Numbers. Matlab News and Notes, The Mathworks, Spring, 2001.