Abstract

Based on the fast Fourier wave superposition spectrum method, a new equivalent source method (ESM) with a sparse sampling technique is proposed. First, the equivalent source intensities are expanded on a rectangular virtual surface using a bidirectional Fourier series, resulting in a semi-analytic and half-numerical acoustic pressure expression. The Fourier coefficients result in good sparsity for continuous acoustic pressures from structural vibration sources, and the proposed sparse sampling method can further reduce correlation in the measurement matrix. Better results can be obtained by solving the l1 norm optimization problem. Finally, the method was verified using several examples. The proposed method offers two main advantages compared with the traditional compressive equivalent source method: (1) the unknown source intensity vector is expanded into a bidirectional Fourier series, thereby transforming an unknown source intensity vector into a sparse Fourier coefficient vector, which has better sparsity; (2) the proposed method constructs a random sampling matrix, which is expanded into a sparse sampling matrix by random distribution, thereby improving the reconstruction accuracy of planar near-field acoustic field compared with the traditional random position sampling method reducing correlation in the transfer matrix.

1. Introduction

Near-field acoustic holography (NAH) is a powerful technique for identifying and localizing acoustic sources and visualizing acoustic fields [1]. With the development of NAH technology, several new algorithms have been introduced including the spatial Fourier transform (SFT) [2], boundary element method (BEM) [3], and equivalent source method (ESM) [4]. Among them, the ESM has been widely studied and applied, mainly owing to its simple principle and strong adaptability [5]. According to the principle of the ESM, the acoustic field radiated by a vibrating body of arbitrary shape can be approximated as the superposition of the acoustic field generated by a series of virtual equivalent sources distributed inside the structure. The acoustic field is reconstructed using the weight coefficients of each equivalent source based on measured acoustic pressure on a holographic surface [6]. Since noise in measured data cannot be avoided, the least-squares method based on the l2 norm is usually used with the traditional ESM to obtain a stable solution. In practical applications of NAH, due to the measurement conditions and cost limitations, the acoustic reconstruction problem usually involves solving a set of underdetermined equations. The least-squares method based on the l2 norm is used to solve the undetermined system of equations, which makes it difficult to obtain an ideal reconstruction result.

In recent years, compressed sensing (CS) has been widely used in signal and image processing [7] and has become a hot topic of research [8, 9]. The novel compressive sampling technique uses the sparsity of the signal to solve the underdetermined system of equations. The sparse solution is easier to obtain because the reconstruction algorithm is based on the l0 or l1 norm instead of the traditional least-squares method based on the l2 norm. Thus, good reconstruction accuracy is guaranteed with fewer measurement points. In 2012, Chardon et al. [10] introduced CS into NAH for the first time and experimentally demonstrated higher accuracy compared with the traditional method when reducing the measurement points. Fernandez and Xenaki [11] proposed the compressive equivalent source method (C-ESM) by combining CS with ESM and analyzed the influence of column correlation in the transfer matrix on the reconstruction result [12]. On the basis of the C-ESM, Hu et al. [13] studied the reconstruction of sparse sound field through the sparse basis obtained by singular value decomposition of the transfer matrix; a fast sparse acoustic field reconstruction method was proposed to combine the Bayesian compressed sensing with the sparse basis function in the following study [14]. In addition to the C-ESM, there are many methods combining compressed sensing with equivalent source method, for example, compressed fused model equivalent source method (CFMESM) [15], compressed velocity-mode equivalent source method (CVMESM) [16], and fused total generalized variation (FTGV) [17]. The particle velocity maps will have sharp peaks in CVMESM [15]. The mode scaling factor was set in CFMESM, and the need for the parameter is a weakness of the method [16]. FTGV combines sparsity in the second derivatives with sparsity in the amplitudes, but FTGV must be solved using the CVX toolbox [17]. In previous studies, sampling points (microphones) are often randomly arranged to satisfy noncorrelation between the observation matrix and the sparse basis. Random sampling can reduce correlation in the measurement matrix through random positioning. Moreover, since previous methods are based on the assumed sparsity of the equivalent source intensities, it is difficult to obtain accurate prior information of the acoustic source when the type of acoustic source is unknown, and the sparsity requirement is not guaranteed. Using corresponding sparse sampling to the transform matrix, we can not only randomly select the sampling position to reduce correlation in the measurement matrix but also use the random coefficient matrix to further reduce correlation.

In this study, a planar acoustic field reconstruction method based on the fast wave superposition spectrum and the sparse sampling matrix is proposed. The bidirectional Fourier series expansion of the equivalent source intensity on a rectangular virtual surface is used to establish a semi-analytical and semi-numerical Fourier series expansion form, and the sparsity of the Fourier coefficient vector is analyzed. To reduce the correlation of the measurement matrix further, a sampling matrix with random spatial positioning and random summation coefficients is introduced. Finally, simulations were performed to compare the proposed method and the traditional C-ESM for planar acoustic field reconstruction.

2. Planar Acoustic Field Reconstruction Based on Fast Fourier Wave Superposition Spectrum

2.1. Fourier Series Expansion of Acoustic Pressure

From the wave superposition method (WSM) [4], the acoustic pressure at any field point outside the acoustic source can be superposed by a series of virtual equivalent sources in the virtual surface inside the acoustic source. Thus, the acoustic pressure can be expressed using the following integral:where is the position vector of the equivalent source, is the strength of the equivalent source at , and is Green’s function. As shown in Figure 1, the position coordinates of field point and equivalent source point in the rectangular coordinate system are and , respectively. To reconstruct the planar acoustic field, the virtual surface is arranged in a rectangular domain of size .

Expanding the equivalent source intensity in (1) into a Fourier series along the x-axis and y-axis of the selected virtual surface , is given, and it becomes the following:where the Fourier series coefficient is as follows: .

Substituting equation (2) into equation (1), the pressure field becomes

In practical calculations, the summation range of the series in equation (3) needs to be truncated; that is,

2.2. Fast Fourier Wave Superposition Spectrum Method

For convenience, the integral term in equation (4) can be expressed as follows:

Dividing the integration range and in equation (5) into equal divisions of M and N, respectively, it becomes

Using the trapezoidal formula to numerically solve equation (5), it becomes

The two-dimensional discrete Fourier transform (DFT) and inverse DFT are defined as [18] follows:where the value of samples is as follows: .

It can be seen from equation (8) that equation (7) is the standard two-dimensional DFT. When M and N are positive integer powers of 2, can be quickly calculated using the discrete fast Fourier transform.

It is worth mentioning that when the DFT is used to calculate , the range of m and n is and , respectively. However, the summation ranges in equation (4) are and , respectively, and only need to be adjusted by the periodicity of with respect to and .

The fast Fourier wave superposition spectrum of acoustic pressure at any point can be obtained using equations (4) and (7):

The number of summation truncation terms requires and .

For acoustic field reconstruction in NAH, measured holographic acoustic pressure data are substituted into equation (9) and written in matrix form:where is the column vector of acoustic pressure on the holographic surface at sampling point and is the column vector composed of coefficients of the Fourier series of equivalent source intensity on the virtual surface . is the transfer matrix of between the acoustic pressure on the holographic sampling surface and the Fourier coefficient of the equivalent source intensity on the virtual surface . Then, the least-squares method based on the l2 norm is used to solve equation (10):where is the regularization parameter and is the l2 norm.

After solving the Fourier coefficient column vector of equation (10), the acoustic pressure of any reconstructed surface in the acoustic field can be obtained as follows:where is the transfer matrix between the acoustic pressure on any reconstructed surface and the Fourier coefficient of the equivalent source strength on the virtual surface .

From the above derivation, the equivalent source intensity on the virtual surface can be expanded into a Fourier series and the acoustic pressure expression can be obtained using the fast Fourier wave superposition spectrum. Since this formula is semi-analytical and semi-numerical, the solution accuracy will be higher than that of the traditional ESM [19]. To improve the calculation results, the traditional least-squares method based on the l2 norm must satisfy and equation (10) must be a set of overdetermined equations.

In practical applications of NAH, the number of holographic measuring points is often insufficient to meet the above requirements due to measurement conditions and cost constraints. Therefore, equation (10) is generally a set of underdetermined equations, making it difficult to obtain satisfactory results using the least-squares method based on the l2 norm. However, if the sparsity of the Fourier coefficient vector of equivalent source intensity is known, sparse sampling can be carried out on the holographic surface and the underdetermined equations can be solved using the l0 or l1 norm methods of CS, which can obtain a better solution result.

According to the theory of structural dynamics, the dynamic response (vibration velocity) of a structure under an accidental load excitation is mainly composed of a superposition of low-order modes. Then, the Rayleigh integral shows that the radiated acoustic field is also mainly concentrated in the low-order modes. In other words, the contribution of structural vibration to the acoustic field is mainly concentrated in the low-order modes, whereas the high-order modes are similar to the evanescent modes in the SFT. The higher the order, the faster the attenuation; that is, the higher-order becomes smaller and smaller. Therefore, provided that the number of summation truncation terms and of equation (9) is large enough, the acoustic field represented by this formula will contain both low-order modes with strong radiation ability and evanescent modes. Thus, the Fourier coefficient vector must be a sparse vector with certain sparsity. As the number of summation truncation terms and increases, the sparsity is strengthened. Therefore, the proposed fast Fourier wave superposition spectrum method can directly solve the acoustic field reconstruction using sparse sampling combined with the l0 or l1 norm of CS.

To intuitively explain the sparsity of the Fourier coefficient vector , a simple analysis of the sparsity of the Fourier coefficient vector of this kind of continuous vibration acoustic source is presented. The simulation calculations were performed with the vibration of a simply supported plate as the acoustic source and compared with the results of the C-ESM [12]. The dimensions of the simply supported plate were 0.5 m × 0.5 m × 0.03 m (length × width × height). A harmonic excitation force of 1 N was applied at the center of the plate, and the excitation point was located at the center of the simply supported plate. The sampled acoustic pressure on the holographic surface was calculated by the Rayleigh integral method to obtain the theoretical acoustic pressure at the measuring point [20]. The Gaussian white noise with a signal-to-noise ratio (SNR) of 20 dB was added to simulate the actual acoustic pressure. The holographic surface was 0.1 m above the simply supported plate, the reconstruction surface was 0.05 m above the simply supported plate, and the virtual equivalent source surface was arranged on the simply supported plate. Sizes of the holographic surface, reconstruction surface, and virtual equivalent source surface are consistent with those of the simply supported plate. In the C-ESM, 81 measuring points (a 9 × 9 uniform array) were set on the holographic surface and the number of equivalent sources was 441. According to the relationship between the number of summation truncation terms and the sound wavelength, the number of summation and truncation terms in the x and y directions was and , respectively. According to equation (9), the number of integral segments in the x and y directions of the virtual surface was set as . The equation was solved by basis pursuit de-noising (BPDN) based on the l1 norm in the SPGL1 toolbox in MATLAB [21].

The amplitudes obtained using the C-ESM and the method in this study with a sampling rate of 1000 Hz are presented in Figure 2. The virtual source intensity obtained by the C-ESM has large amplitudes at most of the source intensity sequence points and the sparsity is poor, as shown in Figure 2(a). The Fourier coefficient vector of this study has only a small number of large amplitudes and the rest are close to zero, as shown in Figure 2(b), which indicates that the Fourier coefficient vector has stronger sparsity for continuous structural vibration sources such as the simply supported plate. The larger the number of truncated terms and , the stronger the sparsity. Therefore, compared with the C-ESM, the vector has better sparsity and can be better solved using the l1 norm.

3. Principle of Compressed Sensing and Sparse Sampling Matrix

3.1. Principle of Compressed Sensing

Assuming an N-dimensional signal can be linearly represented by a set of basis vectors and defining a matrix composed of basis vectors as , then the signal can be expressed as follows:where is the decomposition coefficient. If there are only K nonzero values in the decomposition coefficient vector and , the sparsity of the coefficient vector is K and is the sparse basis matrix of signal . The signal is sampled M times based on CS through the measurement matrix , which is not related to the sparse basis matrix, to obtain the measurement value:where is a matrix of . Since the decomposition coefficient vector has some sparsity, it can be optimally solved using the l0 norm of CS:

However, the solution of equation (15) belongs to the NP-hard problem, which means that the nondeterministic polynomial (NP) cannot be solved by an exact algorithm, and an effective approximation algorithm for such problems must be sought. To obtain the correct solution, it is necessary to exhaustively enumerate combinations, which will consume a lot of calculation time. Usually, the l1 norm is used instead of the l0 norm to obtain

3.2. Construction of Sparse Sampling Matrix

In a previous study, the C-ESM reduced correlation in the measurement matrix (transfer matrix) by randomly arranging the sensor positions [12]. In this study, a random sampling matrix is constructed to reduce the correlation. There is a random matrix of : , where is the number of actual sampling points (microphones). Extending this random matrix to a sampling matrix of and randomly distributing column vectors in to satisfies , that is . Zero column vectors of matrix S are the positions of no sampling points. Multiplying both sides of the sampling matrix in equation (10), it becomes

Let and ;where is the column vector of and is the matrix of . Equation (18) is a set of underdetermined equations that can be solved using the l1 norm minimization method, that is,

In contrast to the random arrangement of sampling points used in the traditional methods [12], this study proposes a new sparse sampling method. A random matrix is constructed and expanded into the sampling matrix. However, the method of selecting the random sampling matrix is important. Commonly used random sampling matrices are as follows:(1)Random matrix: construct a random matrix of , with values of the array elements evenly distributed between 0 and 1. Measurement matrix P can be obtained by normalizing the column vector.(2)Gaussian matrix [22]: construct a matrix of such that each element of the matrix independently satisfies a Gaussian distribution with mean 0 and variance (strongly random, not related to most orthogonal bases).(3)Bernoulli matrix [23]: construct a matrix of such that every element in the matrix obeys a Bernoulli distribution, independently. Compared with the Gaussian matrix, the elements of the Bernoulli matrix are , which are easier to implement and store in practical applications.(4)Circulant matrix [24]: first generate an dimensional random vector, circulate it times to construct the remaining row vector, and finally normalize the column vector to obtain the measurement matrix .(5)Part Fourier matrix [25]: first generate an orthogonal matrix of , then randomly select row vectors of the orthogonal matrix, and finally normalize the column vectors of to obtain the measurement matrix .

In this study, the above five random sampling matrices are used, and sparse sampling matrices with random distributions are constructed. Compared with the traditional randomly selected sampling positions for reducing correlation in the measurement matrix, the randomly distributed coefficient matrix ensures that the sparse sampling matrix also has the randomness of the summation coefficients of the acoustic pressure at the measurement points. Therefore, the sparse sampling matrix contains not only the randomness of the position but also the randomness of the summation coefficients of the acoustic pressure at measurement points, which further reduces correlation in the sparse sampling matrix.

4. Simulation Analysis

4.1. Simply Supported Plate Acoustic Source under Central Excitation

Since the C-ESM has been compared with FTGV, CVMESM, and CFMESM in reference [26], o reduce the amount of computation, the proposed method was only compared with the C-ESM method. To verify the correctness of the proposed method and the accuracy of the acoustic field reconstruction, a simply supported plate was used as a vibration source for simulation calculations. The size of the simply supported steel plate was 0.5 m × 0.5 m, and the thickness was 0.003 m. The plate was driven by a harmonic excitation force of 1 N, and the excitation point was located at the center of the plate. A theoretical radiation acoustic field was calculated using the Rayleigh first integral [20]. The holographic surface was located 0.05 m above the plate, and its dimensions were consistent with the acoustic source surface. The reconstruction surface was located 0.05 m from the holographic surface, and the equivalent source surface was arranged on the acoustic source surface.

Dimensions of the reconstruction surface and equivalent source surface were consistent with the holographic surface, and the center of all three surfaces was on the z-axis with the center of the plate. In simulations, the actual measured acoustic pressure on the holographic surface was obtained by adding 30 dB of the Gaussian white noise to the theoretically calculated acoustic pressure. For the C-ESM, the number of virtual equivalent source points was 625 and the number of holographic measurement points was 81. The sampling points were randomly selected from 625 measurement points on the 25 × 25 holographic grid. The distribution of sampling point locations is shown in Figure 3. When the calculations were performed using the method of this study, the microphone array on the measurement surface was consistent with [12]. According to the relationship between the number of summation truncation terms and the sound wavelength, the number along the x and y directions of summation truncation terms was and , respectively. According to equation (9), the sampling number along the x and y directions of the virtual surface was .

The reconstruction error is defined as follows:where is the reconstructed acoustic pressure and is the analytic expression of acoustic pressure.

Figure 4 shows the reconstructed acoustic pressure and theoretical acoustic pressure obtained using the C-ESM or proposed method with five different sparse sampling matrices at frequencies of 500 Hz, 1500 Hz, and 2500 Hz. When the frequency is 500 Hz, all reconstructed pressures match the analytical acoustic pressure, except for the proposed method, which combines part of the Fourier matrix. As the frequency increases, the C-ESM can no longer match the analytical acoustic pressure well, whereas the five different sparse sampling matrices are sufficient for obtaining better results with the proposed method. It can also be seen that the proposed method results better reconstruction effects when the Gaussian matrix, Bernoulli matrix, and circular measurement matrix were used as the random sampling measurement matrix and was in better agreement with the analytic expression of acoustic pressure.

Figure 4 is only the reconstruction results at three specific frequencies: f = 500 Hz, 1500 Hz, and 2500 Hz. To analyze the reconstruction results at f = 100∼3000 Hz, Figure 5 shows the reconstruction error curves of the C-ESM and the proposed method with five different random sampling measurement matrices in the frequency band from 100 to 3000 Hz. The reconstruction error of the proposed method is lower than that of the C-ESM. Meanwhile, it is difficult to ensure stable reconstruction results in the computational frequency band using either the partial Fourier coefficient matrix or random matrix as the random sampling measurement matrix. The proposed method obtains good results with the Gaussian matrix, Bernoulli matrix, and circular measurement matrix, and the error curves were almost coincident.

4.2. Simple Supported Plate Acoustic Source under Eccentric Excitation

To compare the acoustic field reconstruction of the proposed method and the C-ESM for a simply supported plate acoustic source under eccentric excitation, the parameters of the simply supported plate in Section 3.1 were used. The excitation point was located on the top right of the plate (0.375, 0.375), as shown in Figure 6. The holographic surface was set at 0.06 m above the plate, and its size was consistent with the acoustic source surface. The reconstruction surface was 0.04 m above the plate, and the equivalent source surface was 0.01 m above the acoustic source surface. The sizes of the reconstruction surface and the equivalent source surface were consistent with those of the holographic surface.

In simulations, the same measured acoustic pressure on the holographic surface was used as in the previous example of Section 3.1, which was obtained by adding 30 dB of the Gaussian white noise to the theoretically calculated acoustic pressure. When the C-ESM is used for calculation, the number of virtual equivalent source points was 484 and the number of holographic measuring points was 64. The number of sampling points was randomly selected from the 484 measuring points on the 22 × 22 holographic grid (intervals of 0.024 m). The method used for selecting the random sampling matrix in this study was the same as that used in the C-ESM.

Because of the particularity of central excitation, it is more universal to choose eccentric excitation as random excitation. Figure 7 shows the theoretical acoustic pressures of the reconstructed surface at frequencies of 500 Hz, 1500 Hz, and 2000 Hz and acoustic pressure contours reconstructed using 64 sampling points using either the C-ESM or the proposed method with five different sparse sampling matrices. When f = 500 Hz, the reconstructed acoustic pressures of all methods are in good agreement with the analytical acoustic pressure and the error is less than 10%. As the frequency increases, the C-ESM no longer matches the analytical acoustic pressure well, but the proposed method has high accuracy with five different sparse sampling matrices. Using the Gaussian matrix, Bernoulli matrix, and circular measurement matrix as the random sampling measurement matrix, the proposed method can identify richer details of the acoustic field information.

Figure 8 shows the reconstruction error curves of the C-ESM and the proposed method with five different sparse sampling measurement matrices in the frequency band of 100∼2000 Hz. The reconstruction error of the proposed method is lower than that of the C-ESM. Meanwhile, the reconstruction errors using the partial Fourier coefficient matrix or random matrix are higher than the other three methods; the proposed method obtains good results with the Gaussian matrix, Bernoulli matrix, and circular measurement matrix.

According to Figures 5 and 8, the errors of the C-ESM and the proposed method with five different sparse sampling measurement matrices are less than 20% in the 100∼1400 Hz frequency range, and they are all acceptable error percentages. Due to different excitations, the calculation accuracy of five different sparse sampling measurement matrices is also different, the reconstruction results of the Gaussian matrix, Bernoulli matrix, and circular measurement matrix are better than that of the C-ESM, and the reconstruction results of the partial Fourier coefficient matrix or random matrix are similar to that of the C-ESM. However, the error of the C-ESM under central excitation exceeds 20%, the errors of the Gaussian matrix and circular measurement matrix are about 10∼16% in the 1400∼3000 Hz frequency range, and the errors of the Gaussian matrix and circular measurement matrix are acceptable error percentages. Because of the construction characteristic of the Gaussian matrix and circular measurement matrix, satisfactory results could be obtained even in high frequencies. The proposed method with the Gaussian matrix and circular measurement matrix can apply in further studies. In conclusion, the reconstruction results of the Gaussian matrix and circular measurement matrix are the best in five different sparse sampling measurement matrices, and it is further studied in other applications of near-field acoustic holography.

5. Conclusions

(1)Based on the ESM, the source intensity can be expanded into a bidirectional Fourier series on the rectangular virtual surface, and a semi-analytical and semi-numerical expression of acoustic pressure can be obtained. The proposed method has higher accuracy compared with the traditional numerical discrete source methods. Furthermore, since this method converts the source intensity vector into the sparse Fourier coefficient vector, the vector to be solved also has stronger sparsity and can be more effectively solved using the l1 norm.(2)Correlation in the measurement matrix is reduced by constructing a random sampling matrix, and the sparse sampling matrix is introduced by expanding it into the sampling matrix. The sparse sampling matrix not only includes the randomness of the sampling position, but also includes the randomness of the sum coefficient of the acoustic pressure at the measuring point, which further reduces correlation in the coefficient matrix.(3)The acoustic field reconstruction results of the traditional C-ESM and the proposed method for a simply supported plate source were compared. The simulation results show that the reconstruction error of the proposed method is lower than that of the C-ESM, and this method has higher accuracy with the Gaussian matrix, Bernoulli matrix, and circular measurement matrix. In particular, the errors of the Gaussian matrix and circular measurement matrix under central excitation are about less than 16% in the 1400∼3000 Hz frequency range, which are at least 4% lower than the C-ESM.(4)In this study, center excitation and arbitrary eccentric excitation are used for the excitation of a simply supported plate. The sound field reconstruction of a simply supported plate under other arbitrary excitations can be regarded as the superposition of the excitations in this study. Therefore, this method is of great significance for planar sound sources, but for the sound field reconstruction of nonplanar sound sources, especially the sound field reconstruction of rotating structures, this method can be further studied (Table 1).

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

The research was supported by the National Natural Science Foundation of China under grant no. 51775121 and Doctor Start-up Fund of Guangxi University of Science and Technology (No. 03210126). The supports are gratefully acknowledged.