A sound target positioning method based on a uniform concentric circle microphone array

Abstract

The application discloses a sound target positioning method based on a uniform concentric circle microphone array, which can accurately extract a target sound signal in a noisy signal and calculate the direction of an outgoing wave. The method comprises the following steps: firstly, receiving a noisy signal based on the uniform concentric circle microphone array, pre-filtering and sampling the received signal; then, reading data and extracting the target sound signal; then, performing frequency spectrum analysis on the extracted signal to determine a peak frequency and performing re-filtering to reduce the bandwidth and improve positioning precision; finally, positioning the target sound by using an MVDR beam forming method. The method has the advantages of small data points, small calculation amount, reduced relevant hardware cost, high positioning precision and good stability.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of signal processing, in particular to a sound target positioning method based on a uniform concentric circle microphone array, which relates to beam forming, array signal processing and other theories. BACKGROUND

[0002] Sound is an important information for perceiving the world and an important way for interacting with the external environment. According to different sounds, the position information of the sound source can be judged, and specific meanings can be obtained. The traditional single microphone input is difficult to meet this application, therefore, the sound source positioning based on the microphone array has been rapidly developed, which has wide application in practical scenes such as road whistle detection positioning and mechanical equipment noise source positioning.

[0003] The commonly used sound source positioning algorithm can be divided into three categories: the first category is the method based on time difference of arrival. The algorithm firstly estimates the relative time delay between each microphone, and then determines the position of the sound source by using the estimated time delay and the geometric structure of the array. The second category is the method based on beam forming. The basic idea is to filter, weight and sum the sound signals received by the microphone to form a beam, and then search the possible position of the sound source to guide the beam, and the direction with the maximum output power of the beam is the direction of the target sound source. The third category is the method based on high resolution spectrum estimation. The method divides the array received data into two mutually orthogonal subspaces through mathematical decomposition: the signal subspace consistent with the array flow space of the signal source, and the noise subspace orthogonal to the signal subspace. By using the orthogonal characteristics of the two subspaces, a "needle-shaped" spatial spectrum peak is constructed.

[0004] In the existing sound source positioning algorithm, the positioning error of the positioning algorithm based on time difference of arrival is large under low signal-to-noise ratio, and the method based on high resolution spectrum estimation generally has the disadvantage of large amount of calculation, which is not suitable for real-time positioning. SUMMARY

[0005] The present application overcomes the above-mentioned shortcomings of the prior art, and provides a sound target positioning method based on a uniform concentric circle microphone array, which can accurately extract the target sound signal in the noisy signal and calculate the direction of the wave according to the extracted signal.

[0006] The sound target positioning method based on the uniform concentric circle microphone array of the present application comprises the following steps:

[0007] Step 1: define the coordinate system as follows, define the array plane as the XY plane, the Z axis direction perpendicular to the array center; the roll angle The roll angle is defined as the angle between the line connecting the sound source and the array center and the YZ plane, and the pitch angle θ is defined as the angle between the line connecting the sound source and the array center and the XZ plane.

[0008] The signals are received by a uniform concentric circular microphone array, assuming that the array has N microphones, and a narrowband plane wave with a center frequency f0is incoming from a direction Θ with a sound speed c in air. θ s and θ s are the yaw and pitch angles of the incoming signal, respectively, and satisfy -90°≤θ ≤90°. The signal received by the i-th microphone can be expressed as

[0009]

[0010] where s(k) represents the original signal, v i (k) represents additive noise, k represents the sampling point, and τ i (Θ s ) is the time delay of the signal incoming from Θ s to the i-th microphone relative to a selected reference point. Let the spatial coordinates of the i-th microphone be r i = (x i , y i , 0), i = 1, 2, …, N, then τ i (Θ s ) can be calculated according to the following formula

[0011]

[0012] Let the array received signal vector and noise vector be as follows (the superscript "T" represents transposition, and the subscripts "1, 2, …, N" represent the 1st, 2nd, …, Nth elements in the array, respectively)

[0013] X(k) = [x1(k) x2(k) … x N (k)] T (3)

[0014] V(k) = [v1(k) v2(k) … v N (k)] T (4)

[0015] Define the N × 1 dimensional array direction vector a(Θ

[0016]

[0017] Therefore, the signal received by the array can be expressed by the matrix as

[0018] X(k) = a(Θ s )s(k) + V(k) (6)

[0019] where the direction vector a(Θ s ) can be specifically expressed as

[0020]

[0021] The number of array elements N of the uniform concentric circular microphone array used above is 32, which consists of 4 uniformly distributed concentric circles, and the number of array elements on the 4 concentric circles from outside to inside is 12, 10, 6 and 4 respectively.

[0022] Step 2: Pre-filtering is performed on the received signal, which has two effects, one is to eliminate noise such as direct current component generated by electronic equipment, and the other is to filter out non-target sound. The pre-filtering uses a 4th order Butterworth band-pass filter, and since the typical frequency range of the target sound source to be positioned in this example is between 400-5000Hz, the passband range of the band-pass filter is set to 300-6000Hz.

[0023] Step 3: Read data and extract target sound signal. Extraction is mainly performed by calculating the short-time average amplitude and setting an appropriate amplitude threshold.

[0024] First, the signal with target sound signal is windowed and framed, assuming that the signal to be processed is x(k), which is divided into D frames with a frame shift length of T and a window function of w(k) with a length of L, then the dth frame data can be represented as

[0025] x d (k) = w(k)x(k + dT), 0≤k≤L-1, d = 1, 2, …, D (8)

[0026] Among them, the selected window function is the Hanning window, and its function expression is as follows

[0027]

[0028] The short-time average amplitude is used to determine whether each frame of signal is a target sound signal, and the specific steps are as follows: calculate the short-time average amplitude of each frame of signal, if the short-time average amplitude of a certain frame of signal is greater than the set threshold, it is determined as a target sound signal frame, otherwise it is an invalid frame or a noise frame.

[0029] The short-time average amplitude is defined as follows

[0030]

[0031] Among them, subscript "d" represents the dth frame, d = 1, 2, …, D, x d (k) represents the dth frame signal, and L represents the number of points in each frame of signal. Finally, each frame of target sound signal is recombined to extract the target signal from the original signal.

[0032] However, in actual situation, the received signal may be affected by short strong noise, so that the short-time average amplitude of a frame containing strong noise signal also meets the threshold condition and is considered as a target signal frame. To solve this problem, the continuity of the frame is judged, as follows. First, all signal frames (may contain strong noise) meeting the threshold condition are extracted from the signal according to the above rule and sorted. Then, the signal frames are grouped according to the index order, and the signal frames with continuous indexes are grouped respectively to obtain a plurality of groups of continuous signal frames. The group with the largest number of frames or the number of frames greater than a threshold is considered as the target signal, wherein the threshold can be set according to the actual target sound source time.

[0033] Step 4: Perform FFT analysis on the extracted target sound signal, determine the peak frequency, and filter according to the peak frequency to further remove noise.

[0034] Suppose the peak frequency determined after FFT analysis of the target signal is f max , and according to the passband of (f max -200, f max +400) Hz, the 4th order Butterworth band-pass filter is used for re-filtering to further remove irrelevant noise, reduce the signal bandwidth, and improve the positioning accuracy of the algorithm.

[0035] Step 5: Use the Minimum Variance Distortion Response (MVDR) beamforming method for positioning, calculate the output power of beamforming, and through traversing the roll angle and pitch angle, obtain the angle corresponding to the peak value of the output power spectrum, which is the estimated target sound direction.

[0036] The output of general beamforming can be represented as

[0037] Y(k) = W H X(k) (11)

[0038] where W is the weighting vector of the beamformer, and the superscript "H" represents the conjugate transpose. The power of the beamformer output is

[0039] P = E[|(Y(k))| 2 ] = W H E[X(k)X H (k)]W = W H R X W (12)

[0040] In the formula, E[·] represents the expectation operation, and R X is the covariance matrix of the array received signal and R X = E[X(k)X H(k)], for an array of N microphones, assuming the array receives K samples of the signal, the covariance matrix can be estimated by

[0041]

[0042] The principle of MVDR beamforming is to minimize the total output power of the array while constraining the amplitude response in the desired direction, which can be expressed as the following optimization problem

[0043]

[0044] The above optimization problem can be solved by using the Lagrange multiplier method. Let the objective function be

[0045] L(W) = W H R X W + l[W H a(Θ s )-1] (15)

[0046] where l is the Lagrange multiplier. Taking the partial derivative of the above equation with respect to W, we have

[0047]

[0048] Setting the above equation to zero, we can obtain the optimal weight vector of the MVDR beamformer for receiving the desired signal from direction Θ s as

[0049]

[0050] where μ is a proportionality constant. Note that the constraint W H a(Θ s ) = 1 can be equivalently written as a H (Θ s ) W = 1, and multiplying both sides of equation (16) by a H (Θ s ) simultaneously, we can obtain

[0051]

[0052] That is, the optimal weight vector of the MVDR beamformer is

[0053]

[0054] Substituting equation (19) into equation (12), we can obtain the power spectrum of the MVDR beamforming as

[0055]

[0056] Scan each angle in the observation space, and calculate the power spectrum value according to formula (20), wherein the angle value corresponding to the maximum value is the azimuth of the signal.

[0057] The estimated covariance matrix R in the actual situation X There is a certain error, and the diagonal loading method can be used to correct it and improve the stability of the subsequent algorithm. The diagonal loading method used in the application is as follows: first, the uncorrected covariance matrix R is estimated according to formula (13) X , and then the standard deviation σ of the diagonal elements of the covariance matrix and the mean value μ0 of the diagonal elements are calculated, and the correction coefficient α is

[0058] α = min <σ, μ0> (21)

[0059] Finally, the diagonal loaded covariance matrix R is obtained

[0060]

[0061] Where I is a unit matrix with the same dimension as R X .

[0062] The data length used by the MVDR beam forming method for positioning is 0.1s, and the target sound signal obtained in the foregoing steps is grouped according to the 0.1s length, and the positioning result is calculated for each group of data, and the final positioning calculation result is obtained by comparison and optimization method, thereby improving the positioning accuracy.

[0063] The application has the advantages that: the algorithm has less data points and less calculation amount, and can reduce the cost of related hardware; the positioning algorithm adopts the MVDR beam forming method and combines the diagonal loading technology, and has high positioning accuracy and good stability. BRIEF DESCRIPTION OF DRAWINGS

[0064] Figure 1 is a flowchart of the method of the application;

[0065] Figure 2 is a schematic diagram of the microphone array positioning experiment of the application;

[0066] Figure 3 is a schematic diagram of the uniform concentric circle microphone array of the application;

[0067] Fig. 4(a) and Fig. 4(b) are time domain waveform comparison diagrams before and after pre-filtering of a certain channel of the application, wherein Fig. 4(a) is a time domain waveform diagram before filtering, and Fig. 4(b) is a time domain waveform diagram after filtering;

[0068] Figure 5 is a flowchart of target sound recognition and extraction of the application;

[0069] Figure 6 is a flow chart of the MVDR beamforming positioning of the present application;

[0070] Fig. 7(a) and Fig. 7(b) are positioning experimental result graphs of the present application, wherein Fig. 7(a) is a positioning two-dimensional azimuth spectrum graph, and Fig. 7(b) is a positioning two-dimensional overhead view;

[0071] Figure 8 is a multiple positioning error result graph of the present application. DETAILED DESCRIPTION

[0072] The specific implementation of the positioning algorithm of the present application is further described below in combination with the drawings, and the present application includes but is not limited to the following embodiments.

[0073] Referring to Figure 1 , the whistle detection positioning method based on the uniform concentric circular microphone array described in the present embodiment includes three stages: data acquisition, target sound extraction, and MVDR beamforming positioning. It should be noted that the present algorithm aims to calculate the yaw angle and pitch angle of the target sound wave, as shown in Figure 2 is a schematic diagram of the experimental arrangement of the microphone array positioning, and the experiment is carried out in a full anechoic chamber. In the diagram, point O is the center position of the microphone array, which is arranged to face the sound source, and point S is the position where the sound source is placed. The sound signals emitted by the sound source are target sounds of different frequencies with background noise, which are recorded in advance, respectively represent the yaw angle and the pitch angle of the wave. In the present embodiment, the two angles of the sound source are calculated through multiple measurements , which are about 46.98° and 39.44°, respectively. The algorithm implementation process is described below by taking a target sound signal sample as an example.

[0074] The data acquisition stage is mainly completed by using the self-made uniform concentric circular microphone array. The schematic diagram of the uniform concentric circular microphone array is shown in Figure 3 , the diameter of which is 300 mm, the spacing between each circular ring is 75 mm, the number of array elements is 32, and all the array elements are in the same plane (XY plane). The sampling rate of the data acquisition is set to 65536 Hz; the pre-filtering uses a 4th order Butterworth band-pass filter, and the passband range is set to 300-6000 Hz according to the frequency range of the target sound. The pre-filtering is used to eliminate the direct current component noise and filter out other background noise. As shown in Fig. 4, the time domain waveform comparison diagram before and after pre-filtering of a certain channel is shown, wherein Fig. (a) is the time domain waveform diagram before filtering, and Fig. (b) is the time domain waveform diagram after filtering. As shown in the diagram, the pre-filtering can better eliminate the noise such as the direct current component.

[0075] The target sound extraction stage is realized by the following steps, and it is noted that this stage is processed in a certain channel and then mapped to the remaining channels in order to reduce the amount of calculation. The present embodiment is taken as an example of channel 7 to perform the following processing.

[0076] First, the signal with the target sound is framed and windowed, and the short-time average amplitude of each frame is calculated. The target signal frame is extracted by setting the amplitude threshold for subsequent positioning calculation. However, in the actual situation, the background noise is complex, and the influence of short-time strong noise may exist, so that the signal frame determined according to the short-time average amplitude contains strong noise. Therefore, in view of this problem, the frame continuity is judged.

[0077] The frame continuity is judged as follows: first, all signal frames (may contain strong noise) satisfying the threshold condition are extracted from the signal according to the above rule and sorted; then, the signal frames are grouped according to the index order, and the signal frames with continuous indexes are grouped respectively to obtain a plurality of groups of continuous signal frames; among the plurality of groups of continuous signal frames, the group with the maximum number of frames or the number of frames greater than the frame number threshold is considered as the target signal, otherwise it is discarded. The frame number threshold is set according to the frame length and the target signal duration.

[0078] Finally, the signal frames satisfying the above two conditions are recombined into the target signal and input to the subsequent algorithm for positioning calculation. The specific process of the target sound recognition and extraction can be referred to Figure 5 .

[0079] The MVDR beam forming positioning stage:

[0080] Step one: FFT analysis is performed on the target signal extracted in the foregoing to determine the peak frequency f max of the target signal, and the target signal is filtered again according to the frequency (f max -200, f max +400) Hz using a fourth-order Butterworth band-pass filter to filter out the remaining frequency components that interfere with positioning, reduce the bandwidth to approximate a narrow-band signal, and improve the positioning accuracy.

[0081] Step two: the minimum variance distortion response (Minimum Variance Distortion Response, MVDR) beam forming method is used for positioning, and the output power spectrum of beam forming is calculated. The angle corresponding to the peak value of the output power spectrum is obtained by traversing the roll angle and the pitch angle, that is, the estimated target sound arrival direction. The power spectrum can be calculated as follows

[0082]

[0083] where R X is the signal covariance matrix, and a(Θ) is the steering vector. The angle represents the direction of the incoming signal wave.

[0084] First, estimate the covariance matrix R of the horn signal X(k) according to the following formula. X , where K is the number of signal sampling points.

[0085]

[0086] The estimated covariance matrix R in actual practice X There will be some error, which can be corrected by diagonal loading to improve the stability of subsequent algorithms. The diagonal loading method used in this invention is as follows: First, the uncorrected covariance matrix is ​​estimated according to equation (2), denoted as R. X Then, the standard deviation σ of the diagonal elements and the mean μ0 of the diagonal elements of the covariance matrix are calculated, with the correction coefficient α being...

[0087] α=min<σ,μ0> (3)

[0088] The final covariance matrix after diagonal loading is:

[0089]

[0090] Where I is the identity matrix.

[0091] Next, calculate the covariance matrix after diagonal loading. inverse matrix Then, a guide vector a(Θ) is constructed based on the array geometry. Next, Θ is scanned in steps at a certain angle. The MVDR power spectrum value P is calculated according to equation (1). Finally, a peak search is performed on the P value, and the Θ corresponding to the peak point is the sound source angle calculated by the localization algorithm. The specific process of MVDR beamforming localization can be found in [reference needed]. Figure 6 As shown.

[0092] The sound source angle position set in this case experiment is: The MVDR beamforming positioning calculation uses data with a duration of 0.1s (approximately 6554 data points).

[0093] The experimental results are shown in Figure 7. The calculated angle value is... Figure 7(a) is a spatial orientation spectrum, and Figure 7(b) is a top view of the positioning. As can be seen from the figures, the positioning accuracy of the method of the present invention is high, with a positioning error of less than 3°.

[0094] To verify the stability of the positioning method of this invention, multiple positioning experiments were conducted by dividing the target acoustic data into segments with a duration of 0.1 seconds. Figure 8The angle error of multiple positioning experiments (the rest of the experimental samples are similar) is given, and as shown in the figure, the error of multiple positioning of the method is basically consistent, and the robustness is good.

[0095] The content described in the embodiments of the present specification is only a list of implementation forms of the inventive concept, and the protection scope of the present application should not be regarded as being limited to the specific forms stated in the embodiments, and the protection scope of the present application also extends to equivalent technical means that can be thought of by those skilled in the art according to the inventive concept.

Claims

1. A method for acoustic target localization based on a uniform concentric circular microphone array, characterized in that, Includes the following steps: Step 1: Define the coordinate system as follows: The array plane is defined as the XY plane, with the Z-axis perpendicular to the array center; the roll angle... The pitch angle θ is defined as the angle between the line connecting the sound source and the center of the array and the YZ plane, and the elevation angle θ is defined as the angle between the line connecting the sound source and the center of the array and the XZ plane. Using a uniform concentric circular microphone array to receive signals, assuming the array has N microphones, and a narrowband plane wave with a center frequency of f0 is present in space, with a direction angle of [missing information]. The speed of sound in air incident on the array is c; θ s Let be the azimuth and elevation angles of the incident signal, respectively, and satisfy . -90°≤θ s ≤90°; at this time, the signal received by the i-th microphone can be represented as: Where s(k) represents the original signal; v i (k) represents additive noise; k represents the sampling point; τ i (Θ S ) is from Θ S The time delay relative to a selected reference point when the original signal in the direction is incident on the i-th microphone; let the spatial coordinates of the i-th microphone be r. i =(x i ,y i If ,0),i=1,2,…N, then τ i (Θ S The following formula can be used for calculation: The received signal vector and noise vector of the array are represented by the following two equations, where the superscript "T" indicates transpose, and the subscripts "1,2,...,N" represent the 1st, 2nd...Nth array elements, respectively: X(k)=[x1(k) x2(k)…x N (k)] T (3) V(k)=[v1(k) v2(k)…v N (k)] T (4) Define the direction vector of an N×1 dimensional array: Therefore, the signal received by the array can be represented by a matrix as follows: X(k)=a(Θ S )s(k)+V(k) (6) Wherein, the direction vector a(Θ) S Specifically, it can be expressed as: The uniform concentric circle microphone array used above has N=32 array elements, consisting of 4 uniformly distributed concentric circles. The number of array elements on the 4 concentric circles from the outside to the inside are 12, 10, 6 and 4 respectively. Step 2: Pre-filter the received signal. The pre-filter has two functions: first, to eliminate noise generated by electronic devices, and second, to filter out non-target sounds. The pre-filter uses a 4th-order Butterworth bandpass filter. Since the typical frequency range of the target sound source to be located is between 400 and 5000 Hz, the passband range of the bandpass filter is set to 300 to 6000 Hz. Step 3: Read the data and extract the target acoustic signal; this is mainly done by calculating the short-time average amplitude and setting an appropriate amplitude threshold. First, the signal containing the target sound source is windowed and framed. Assuming the signal to be processed is x(k), it is divided into D frames with a frame shift length of T and a window function w(k) of length L. Then, the data of the d-th frame can be represented as follows: x d (k)=w(k)x(k+dT),0≤k≤L-1,d=1,2…,D (8) The selected window function is the Hanning window, and its function expression is as follows: The short-time average amplitude is used to determine whether each frame of signal is the target signal. The specific steps are as follows: calculate the short-time average amplitude of each frame of signal. If the short-time average amplitude of a frame of signal is greater than a set threshold, it is determined to be a target signal frame; otherwise, it is an invalid frame or a noise frame. The short-time average amplitude is defined as follows: Where the subscript "d" represents the d-th frame, d = 1, 2, ..., D, x d (k) represents the signal of the d-th frame, and L represents the number of signal points in each frame; finally, the target signal can be extracted from the original sound signal by recombining the target signals of each frame. However, in reality, the received signal may contain brief periods of strong noise, causing the short-time average amplitude of a frame containing strong noise to also meet the threshold condition and be considered a target signal frame. To address this issue, a frame continuity determination method is proposed, as follows: First, extract all signal frames that meet the threshold condition from the signal according to the above rules and sort them; then, group the signal frames according to their index order, grouping consecutively indexed signal frames into several groups of continuous signal frames; among these groups of continuous signal frames, the group with the largest number of frames or the number of frames greater than a certain threshold is selected as the target signal, where the threshold can be set based on the actual duration of the target sound source. Step 4: Perform FFT analysis on the extracted target acoustic signal to determine its peak frequency and filter it based on the peak frequency to further remove noise; Assume the peak frequency determined by FFT analysis of the target signal is f. max According to (f max -200,f max The +400 Hz passband is filtered again using a fourth-order Butterworth bandpass filter to further remove the influence of irrelevant noise, reduce the signal bandwidth, and improve the positioning accuracy of the algorithm. Step 5: Use the minimum variance distortionless response beamforming method for positioning, calculate the output power of the beamforming, and obtain the angle corresponding to the peak value of the output power spectrum by traversing the azimuth and elevation angles, which is the estimated target acoustic azimuth. The output of a typical beamforming beam can be expressed as Y(k) = W H X(k) (11) Where W is the weighting vector of the beamformer, and the superscript "H" indicates the conjugate transpose; the output power of the beamformer is P=E[|Y(k))| 2 ]=W H E[X(k)X H (k)]W=W H R X W (12) In the formula, E[·] represents the expectation operation, and R X Let R be the covariance matrix of the received signal from the array. X =E[X(k)X H [k] For an array of N microphones, assuming the number of sampling points for the received signal is K, the covariance matrix can be estimated by the following formula: The principle of MVDR beamforming is to minimize the total output power of the array while keeping the signal amplitude response constant in the desired direction. This can be represented by the following optimization problem: The above optimization problem can be solved using the Lagrange multiplier method; let the objective function be: L(W)=W H R X W+l[W H a(Θ S )-1] (15) In the formula, l is a Lagrange multiplier; taking the partial derivative of the above formula with respect to W, we get: Setting the above expression to zero, we can obtain the received direction Θ. S The optimal weight vector for the desired signal in the MVDR beamformer is: In the formula, μ is a proportionality constant; note the constraint W H a(Θ S ) = 1 can be equivalently written as a H (Θ S W = 1, multiply both sides of equation (16) by a on the left. H (Θ S From this, we can obtain: That is, the optimal weight vector for the MVDR beamformer is: Substituting equation (19) into equation (12), the power spectrum of MVDR beamforming can be obtained. Scan each angle in the observation space and calculate the power spectrum value according to equation (20), where the angle value corresponding to the maximum value is the azimuth of the signal; The estimated covariance matrix R in actual practice X There will be some error. A diagonal loading method is used to correct it and improve the stability of the subsequent algorithm. The diagonal loading method is as follows: First, the uncorrected covariance matrix is ​​estimated according to equation (13) and denoted as R. X Then, the standard deviation σ of the diagonal elements and the mean μ0 of the diagonal elements of the covariance matrix are calculated, and the correction coefficient α is: α=min<σ,μ0> (21) The final covariance matrix after diagonal loading is: Where I is related to R X Identity matrices of the same dimension; The MVDR beamforming method uses a data duration of 0.1s for positioning. The target acoustic signals obtained in the previous steps are grouped according to a duration of 0.1s, and the positioning result is calculated for each group of data. The final positioning calculation result is obtained by comparison and averaging optimization method, thereby improving the positioning accuracy.

Images