The received Doppler signal of a stationary sensor, as emitted by a transiting acoustic source, is used to estimate source motion parameters, including speed, closest distance, rest frequency, and closest point of approach (CPA) time. First, the instantaneous frequency, amplitude, and CPA time are accurately estimated by the polynomial chirplet transform of the Doppler signal. Thereafter, the three other source motion parameters are obtained with a simplified amplitude-weighted nonlinear least squares method. The proposed scheme is successfully applied to the analysis of the characteristics of a moving truck.
1. Introduction
The frequency of a pure tone emitted by a transiting source, when received by a stationary observer, changes with time because of the acoustical Doppler effect. This variation in instantaneous frequency (IF) can be used for source motion parameter estimation, which is a significant process in acoustics and has attracted much research attention.1–5 The basic principles involved are measuring the temporal variation in the nonlinear IF of the Doppler signal and then minimizing the sum of the squared deviations of the IF estimates from their predicted values. Therefore, central to the success of this frequency-based motion parameter estimation scheme is the need for an accurate IF estimate. Various IF estimation methods for the Doppler signal have been used. A common approach is to locate the peak of the time–frequency distribution (TFD), such as short-time Fourier transform (STFT), Wigner–Ville distribution (WVD), and polynomial WVD.6–9 However, these methods have some drawbacks, such as their low energy concentration on the IF law of the Doppler signal, appearance of inner artifacts, and intensive computation.
In this study, we consider the polynomial chirplet transform (PCT), which is a parameterized time–frequency analysis method developed by extension of the conventional chirplet transform.10 By replacing the chirplet kernel that has a linear IF law in the chirplet transform with a new kernel that has a polynomial nonlinear IF law, the PCT can generate a TFD that satisfies the energy concentration for the Doppler signal without the inner artifacts. The accurate estimation of the IF and its corresponding amplitude can be effectively achieved through a signal-dependent iterative procedure, the first iteration of which is the STFT. Meanwhile, the closest point of approach (CPA) time is directly estimated through a search of the minimum of the instantaneous frequency rate (IFR), which is obtained with the polynomial kernel coefficients of the PCT. As the polynomial kernel coefficients are estimated by the iterative procedure based on the STFT, the PCT is more robust than the polynomial phase estimation method that obtains the phase through Hilbert transform.6,11 It is also different from the multilinear function class methods, which are need to design a complicated signal kernel and are suitable for the polynomial phase signal of low order (usually less than 5).12–14 With the obtained IF, amplitude, and CPA time, the other motion parameters, such as the source speed, closest distance, and rest frequency, can be estimated with a simplified amplitude-weighted nonlinear least squares method. This proposed scheme is applied to the analysis of the characteristics of a moving truck. These features include the four motion parameters, the cylinder firing rate, engine firing rate, number of cylinders, and tachometer reading.
2. Method
2.1 Doppler shift model
The geometry of the relative motion between a stationary sensor and a moving source is shown in Fig. 1. The source is traveling in a straight line at a constant speed , so at time , the source is at the CPA to the stationary sensor at a separation distance . The source emits an acoustic tone of constant frequency (i.e., the rest frequency), and the speed of sound in air is . Because of Doppler effect, the IF of the signal received by the sensor at time is given by2
where , , , and . The IF of the Doppler signal is a function with respect to the four parameters , which completely specify the source trajectory. Therefore, these motion parameters, or equivalently , can be estimated by measurement of the temporal variation in the IF and then the minimization of the sum of the squared deviations of the IF estimates from their predicted values over a period of time.
Geometry of a stationary sensor and a source traveling in a straight line at a constant speed.
Geometry of a stationary sensor and a source traveling in a straight line at a constant speed.
2.2 IF estimate based on PCT
As mentioned earlier, the IF law of the Doppler signal is a nonlinear function of time. The conventional chirplet transform, which is a parameterized time–frequency analysis method particularly designed for the analysis of chirp-like signals with the linear IF law, cannot track the evolution versus time of the IF of the Doppler signal. To improve the efficacy of the chirplet transform in analyzing signals with a nonlinear IF law, a modified version known as PCT is defined as15
where is the analytical signal; stands for a nonnegative, symmetric, and normalized real window, usually taken as the Gaussian function ; and and denote the nonlinear frequency rotating operator and the frequency shift operator, whose respective expressions are
Here, are the polynomial kernel parameters.
Weierstrass approximation theorem16 states that any continuous function on a closed and bounded interval can be uniformly approximated by a polynomial of a finite-order . Therefore, the Doppler signal can be expressed in the form of a polynomial phase signal, i.e., , where is the slowly varying amplitude, and denote the polynomial phase coefficients. Substituting the polynomial form of the Doppler signal into Eq. (2), we can obtain
When the kernel parameter vector exactly matches the polynomial phase coefficients , i.e., of the Doppler signal, Eq. (4) is reduced to
where is the polynomial form of the IF of the Doppler signal, is the Fourier transform of , and denotes the convolution operator with respect to . Obviously, the frequency resolution of the PCT at all specific moments is equal to the frequency bandwidth of the Gaussian window . Therefore, the PCT can produce a high-quality TFD with an energy concentration on the IF law of the Doppler signal by giving a proper polynomial kernel parameter vector , which can be determined by a signal-dependent iterative procedure with the termination condition
where is a predetermined threshold, and denotes the estimated IF in each iteration. The iterative procedure for the IF estimates is described as follows:
Initialization: Set as a specific value. Set the polynomial order . Set , and set the window size.
Estimation: Calculate the PCT with , and locate the peak in the obtained TFD of the PCT. Then, approximate the TFD peak data with a polynomial function of order with the least square method and denote the coefficients as .
Recursion: Calculate the termination criterion ; if , update with and repeat step (2); otherwise, go to step (4).
Quit: Take the IF trajectory as the obtained -order polynomial function of time.
The predetermination for the order has been discussed in Ref. 14. Following the IF estimates, the peak values on the estimated IF law of the TFD are the estimated amplitude . Because the initialized kernel parameter vector , the first iteration of the PCT is the STFT. The PCT can be viewed as a generalized STFT. Therefore, taking the Chirp z-transform17 instead of the Fourier transform in the PCT of the narrow-band Doppler signal can speed up the entire procedure.
The IFR is the derivative of the IF, so it can be also expressed in the following polynomial form:
The IFR has a minimum at the CPA time, i.e., when , . Therefore, searching the minimum of the IFR estimation, which is obtained as , can facilitate the direct estimation of the CPA time.
The performance of the PCT-based IF and CPA time estimation method is assessed by a test on the numerically generated Doppler signal, whose IF has the following parameters: , , , , and . The sampling frequency is 600 Hz and the duration of the Doppler signal is 30 s. The signal-to-noise ratio (SNR) is −5 dB, and it is defined as the ratio of average power of the Doppler signal and white Gaussian noise. When the PCT-based method is used to estimate the Doppler shift, the size of the Gaussian window is predetermined as 4500, the order of polynomial kernel K = 12, and the threshold . For the purpose of comparison, Figs. 2(a)–2(c) depict the TFDs of the Doppler signal, which are generated by the STFT (first iteration of the PCT), PCT, and WVD. All three TFDs could render the inherent time-frequency pattern of the Doppler signal. However, the TFD generated by the PCT has the best concentration, the TFD by the STFT has a spreading concentration, and the TFD by the WVD has inner artifacts that disperse the concentration. Figure 2(d) shows the estimated IF curves of the three methods, as well as the real Doppler shift (denoted by red line). Corresponding to the respective TFDs, the IF estimated by the PCT (denoted by green line) overlaps the real Doppler shift, the IF by the STFT (denoted by blue line) cannot follow the nonlinear variation of the Doppler shift specially when the source approaches to the CPA, and the IF by the WVD (denoted by black line) is biased and sensitive to the noise. Figure 2(e) clearly shows that the IFR estimated by the PCT (denoted by green dot line) is almost the same as the real IFR (denoted by red line), and the estimated CPA time is 14.96 s that is approximately equal to the actual value 15 s. In conclusion, the PCT has been demonstrated to be effective for the analysis of Doppler signal, even when the SNR is low.
(Color online) The analysis of the numerical Doppler signal by using the PCT, STFT, and WVD: (a) TFD generated by the STFT, (b) TFD generated by the PCT, (c) TFD generated by the WVD, (d) the IF estimation curves and real IF, and (e) the IFR curve and CPA time estimated by the PCT.
(Color online) The analysis of the numerical Doppler signal by using the PCT, STFT, and WVD: (a) TFD generated by the STFT, (b) TFD generated by the PCT, (c) TFD generated by the WVD, (d) the IF estimation curves and real IF, and (e) the IFR curve and CPA time estimated by the PCT.
2.3 Motion parameter estimation
After the IF estimates are obtained with the amplitude, the motion parameters , or equivalently , can be estimated with a nonlinear least squares method, i.e., by minimization of the sum of the squared deviations of the IF estimates from their predicted values. The CPA time has been accurately estimated, so the amplitude-weighted NLS estimates of are specifically given by
where is the IF estimate at time , and is the number of IF estimates. The weighted factor is used to highlight the accurate IF estimates in the duration when the transiting source approaches the CPA, and the SNR is relatively high. According to the linear form of the IF in Eq. (1), for a fixed , Eq. (8) is solved as
where , , and . Substituting Eqs. (9) and (10) into Eq. (8) reduces the three-dimensional problem to a one-dimensional minimization problem
or equivalently, the parameter is solved by minimizing the following cost function:
Substitution of the estimated value into Eqs. (9) and (10) obtains the corresponding parameters and . Then, the other motion parameters can be calculated with the formulas , , . This approach uses an accurately estimated CPA time in advance, so it not only simplifies the solution procedure of the nonlinear least squares method, but it also improves the robustness and precision of the motion parameter estimation.
3. Application
The proposed scheme is applied to the analysis of a radiated acoustic noise from a truck with six cylinders. The truck travels along a straight road with a constant speed of 20 km/h. It has a four-stroke diesel engine and a tachometer reading of 2350 revolutions per minute (rpm). The output of the 35 m far microphone is sampled at the rate of 12 000 Hz, and 30 s of data are recorded during transit of the truck.18 The speed of sound propagation in the air is a constant 347 m/s during data recording.
Figure 3 shows the results of the acoustic data processing. The spectrogram of the recorded signal in Fig. 3(a) shows there are a harmonic series with a fundamental frequency of about 20 Hz, which corresponds to cylinder firing rate. Among the harmonic series, the strongest line corresponding to the engine firing rate is observed at a frequency of about 120 Hz. Therefore, the number of cylinders is given as six through dividing engine firing rate by cylinder firing rate. This number is confirmed by visual inspection because every sixth harmonic of the cylinder firing rate is a strong line. Because of Doppler effect, the IF of the engine firing rate spectral line is observed to change with time. The parameters of PCT-based method are set as the same as the previous numerical simulation. Before the criterion reaches the predetermined threshold, three iterations of PCT have been implemented. Figure 3(b) shows the generated TFDs for the first (the top subgraph) and third iterations (the bottom subgraph). The concentration of the time–frequency representation of the third iteration is much better than that of the first, which is essentially produced by the STFT. Correspondingly, Fig. 3(c) shows that the IF estimated by the PCT (denoted by black circle) is more accurate than the IF by the STFT (denoted by blue line). The normalized amplitude estimates are presented in Fig. 3(d). Figure 3(e) depicts the IFR curve according to the estimated coefficients of the polynomial kernel, and the estimated CPA time . With the use of the estimated IF and amplitude, the source speed, closest distance, and rest frequency (i.e., engine firing rate) are estimated as , , and , respectively. Then, the estimated cylinder firing rate is known as 19.8 Hz, and the estimated tachometer reading is equal to 2 × cylinder firing rate × 60 = 2373 rpm.
(Color online) Data analysis of a moving truck with the proposed scheme: (a) spectrogram of the recording signal, (b) TFDs of the first and last iterations of the PCT-based IF estimate method, (c) IF curves obtained by the STFT and PCT, (d) the estimated amplitude, and (e) the CPA time estimated by the IFR.
(Color online) Data analysis of a moving truck with the proposed scheme: (a) spectrogram of the recording signal, (b) TFDs of the first and last iterations of the PCT-based IF estimate method, (c) IF curves obtained by the STFT and PCT, (d) the estimated amplitude, and (e) the CPA time estimated by the IFR.
Table 1 shows the compared results of the PCT and two existing schemes, which use the STFT and WVD, respectively, to estimate the IF of the Doppler signal and then solve the nonlinear least squares problem by Gauss–Newton algorithm. All the three methods provide acceptable estimates for the characteristic parameters of the moving truck. However, the estimation results of the PCT are more accurate than those of the STFT and WVD.
Compared results of the source parameter estimation obtained by the proposed and existing schemes.
. | Actual value . | PCT . | STFT . | WVD . |
---|---|---|---|---|
(km/h) | 20 | 21.4 | 22.1 | 21.9 |
(m) | 35 | 35.3 | 38.2 | 37.7 |
TRa (rpm) | 2350 | 2373 | 2375 | 2375 |
. | Actual value . | PCT . | STFT . | WVD . |
---|---|---|---|---|
(km/h) | 20 | 21.4 | 22.1 | 21.9 |
(m) | 35 | 35.3 | 38.2 | 37.7 |
TRa (rpm) | 2350 | 2373 | 2375 | 2375 |
TR denotes the tachometer reading.
4. Conclusion
The acoustical Doppler shift is used to estimate the parameters of a moving source. Central to the success of this frequency-based passive acoustic source motion parameter estimation scheme is the need for an accurate IF estimate. We consider the PCT method, which produces a high-quality TFD, with the energy concentrated on the IF law of the Doppler signal. Together with the accurate IF estimate, the CPA time is estimated through minimization of the IFR, which is also obtained by the PCT. With the known CPA time, estimating the other parameters with a simplified amplitude-weighted nonlinear least squares method is easy. The proposed scheme is validated through an analysis of the characteristics of a moving truck.
Acknowledgments
This work was supported by the National Natural Science Foundation of China under Grant No.11274253.