Fast deconvolution acoustic source localization method based on array sparsity and broadband comprehensive processing

Abstract

This invention discloses a fast deconvolution sound source localization method based on array sparsity and broadband integrated processing, comprising the following steps: selecting microphones to form a sparse array, transforming the regular focusing plane into an irregular grid coordinate; calculating the cross-spectrum matrix of the sparse array based on the coordinates, and obtaining the function beamforming output through eigenvalue decomposition and power exponentiation; constructing a frequency domain deconvolution model based on its convolution relationship with the sound source distribution and the power-law spread function, and solving the sound source intensity distribution matrix using fast Fourier transform and iterative optimization; estimating the noise floor and transforming it into a signal-to-noise ratio matrix, traversing frequency points to obtain the maximum signal-to-noise ratio and fusing it into a broadband imaging matrix; introducing a time smoothing factor to weight and fuse the inter-frame imaging matrix to complete the localization; this method reduces the load through array sparsity and fast deconvolution, maintains accuracy through function beamforming and Gaussian filtering, achieves multi-frequency point fusion adaptive imaging, and suppresses flicker through inter-frame smoothing, adapting to complex scenes and realizing low-power, high-precision, and robust sound source localization.

Description

Technical Field

[0001] This invention relates to the field of acoustic imaging technology in the field of sound source localization, specifically to a fast deconvolution sound source localization method based on array sparsity and broadband integrated processing. Background Technology

[0002] Acoustic imaging technology has advanced to the point where fast deconvolution algorithms based on microphone arrays offer significant performance advantages over traditional beamforming and deconvolution algorithms due to their balance of high resolution and computational efficiency. However, these algorithms generally work in narrowband environments and require prior knowledge of the sound source frequency and bandwidth. When faced with real-world sound sources with unknown spectrum and bandwidth, and low signal-to-noise ratio, traditional algorithms can lead to incomplete imaging, location drift, and weak source suppression if the source frequency and bandwidth cannot be accurately determined. Although function beamforming and deconvolution techniques have significantly improved narrowband imaging performance, a comprehensive imaging mechanism that can adaptively fuse effective information across the entire frequency band without requiring prior spectral knowledge is still lacking in broadband applications. Furthermore, the high computational cost of large arrays also limits their engineering deployment. Summary of the Invention

[0003] To address the shortcomings of existing technologies, this invention provides a fast deconvolution sound source localization method based on array sparsity and broadband integrated processing, aiming to solve the problems in the background technology.

[0004] To achieve the above objectives, the present invention provides the following technical solution: a fast deconvolution sound source localization method based on array sparsity and broadband integrated processing, comprising the following steps:

[0005] Step S1: Divide the microphone array into several clusters using k-means clustering, and select the microphones closest to the centroid in each cluster to form a sparse array;

[0006] Step S2: Convert the regular focusing plane into an irregular focusing plane, and calculate the coordinates of the irregular grid by sampling the sine values ​​of the angles at equal intervals;

[0007] Step S3: Calculate the cross-spectral matrix of the sparse array based on the irregular grid coordinates and perform eigenvalue decomposition. Introduce the power exponent to obtain the function beamforming output.

[0008] Step S4: Based on the convolution relationship between the function beamforming output, the sound source distribution, and the power-law spread function, construct a frequency domain deconvolution model;

[0009] Step S5: Use Fast Fourier Transform to convert the convolution in the model into wavenumber domain multiplication, and solve for the sound source intensity distribution matrix through zero-padding and iterative optimization;

[0010] Step S6: Estimate the noise floor of the sound source intensity distribution matrix and set a threshold, then convert it into a sound source distribution signal-to-noise ratio matrix;

[0011] Step S7: For the signal-to-noise ratio matrix of the sound source distribution, traverse all frequency points and take the maximum signal-to-noise ratio of each grid point to fuse them into a broadband integrated imaging matrix;

[0012] Step S8: Introduce a time smoothing factor to weighted fuse the imaging matrices of the current frame and the previous frame to obtain the final acoustic imaging result and complete the sound source localization.

[0013] Furthermore, the specific process of step S1 is as follows:

[0014] Step S1.1: First, determine the location information of all microphones to form a location matrix containing the total number of microphones, where the coordinates of each microphone correspond to a row of data in the matrix;

[0015] Step S1.2: The k-means clustering algorithm is used to divide the entire microphone array into several clusters. The clustering process aims to minimize the sum of squares within each cluster. The cluster assignment of the microphones is adjusted and the centroid of the clusters is updated iteratively until the clustering results are stable.

[0016] Step S1.3: Calculate the distance from all microphones in each cluster to the cluster centroid, select the microphone closest to the centroid as the representative of the cluster, and finally form a sparse array from the representative microphones of all clusters.

[0017] Furthermore, the specific process of step S2 is as follows:

[0018] Step S2.1: Determine the coordinates of the grid points in the regular focusing plane, and the vertical distance between the focusing grid plane and the microphone array plane;

[0019] Step S2.2: By calculating the distance from the grid point in the focusing plane to the reference point, and the trigonometric function relationship between this distance and the horizontal and vertical coordinates, the corresponding angle sine value is obtained;

[0020] Step S2.3: Within the set range of angle sine values, perform discretization sampling at equal intervals to obtain a series of angle sine values ​​for sampling points;

[0021] Step S2.4: Based on the sine values ​​of the angles obtained from these samples, and combined with the vertical distance between the focusing plane and the microphone array plane, calculate the specific coordinates of each grid point in the irregular focusing plane, and complete the transformation from a regular focusing plane to an irregular focusing plane.

[0022] Furthermore, the specific process of step S3 is as follows:

[0023] Step S3.1: Perform Fourier transform on the sound signal recorded by each microphone in the sparse array to obtain the pressure amplitude data corresponding to different frequencies of each microphone, and then construct a column matrix containing the number of microphones after sparse based on the pressure amplitude data;

[0024] Step S3.2: Multiply the column matrix by its conjugate transpose to obtain the cross-spectral matrix of the sparse array signal;

[0025] Step S3.3: Perform eigenvalue decomposition on the cross-spectral matrix to obtain a unitary matrix composed of the eigenvectors of the cross-spectral matrix and a diagonal matrix composed of the eigenvalues;

[0026] Step S3.4: Introduce a positive integer power exponent and combine it with the constructed guide vector to calculate the function beamforming output; the guide vector is calculated from the microphone array position, sound velocity, and frequency through an exponential function.

[0027] Furthermore, the specific process of step S4 is as follows: Based on the function beamforming output obtained in step S3, analyze its relationship with the sound source intensity distribution and the power-law spread function: the function beamforming output is represented as the spatial convolution of the sound source intensity distribution and the power-law spread function.

[0028] Based on the fact that the output of function beamforming is represented as the spatial convolution of the sound source intensity distribution and the power-law spread function, a deconvolution model in the frequency domain is constructed.

[0029] Furthermore, the specific process of step S5 is as follows:

[0030] Step S5.1: Perform zero-filling processing on the function beamforming output, the sound source intensity distribution matrix, and the power-law spread function;

[0031] Step S5.2: Calculate the Fast Fourier Transform of the power-up spread function after zero-filling, and the sum of the absolute values ​​of the results;

[0032] Step S5.3: Initialize the zero-filled sound source intensity distribution matrix as a zero matrix, and then perform iterative calculations using iterative formulas: First, calculate the function beamforming estimation result for the current iteration based on the sound source intensity distribution matrix of the current iteration, the Fourier transform result of the power-law spread function, and the Gaussian regularized filter function; then, update the sound source intensity distribution matrix by combining the zero-filled function beamforming output.

[0033] Step S5.4: After the iteration is completed, remove the zero-filling part to obtain the final estimated sound source intensity distribution matrix.

[0034] Furthermore, the specific process of step S6 is as follows:

[0035] Step S6.1: Convert the final estimated sound source intensity distribution matrix obtained in step S5 into column vector form, where each column vector contains the sound source intensity data of all grid points at the corresponding frequency;

[0036] Step S6.2: Sort the sound source intensity column vectors corresponding to each frequency in ascending order; select the minimum value of a preset ratio according to the length of the sorted vectors, calculate the average value of these minimum values, and use it as the noise floor estimate at the corresponding frequency.

[0037] Step S6.3: Set a threshold based on the noise floor estimate and process the sound source intensity matrix at each frequency: if the sound source intensity at a certain grid point is less than or equal to the noise floor estimate, then set the value of that point to 0; if it is greater than the noise floor estimate, then divide the sound source intensity at that point by the noise floor estimate, and finally obtain the sound source distribution signal-to-noise ratio matrix at each frequency.

[0038] Furthermore, the specific process of step S7 is as follows: For the source distribution signal-to-noise ratio matrix at all frequencies obtained in step S6, for each spatial grid point, traverse the signal-to-noise ratio values ​​corresponding to the grid point at all frequencies, and select the largest signal-to-noise ratio value as the final imaging intensity value of the grid point; integrate the maximum signal-to-noise ratio values ​​of all grid points to form a broadband integrated imaging matrix.

[0039] Furthermore, the specific process of step S8 is as follows: a time smoothing factor is introduced to perform inter-frame weighted fusion processing on the broadband integrated imaging matrix obtained in step S7; during fusion, the value of each grid point in the broadband integrated imaging matrix of the current frame is multiplied by the time smoothing factor, and the value of the corresponding grid point in the imaging matrix of the previous frame is multiplied by 1 minus the time smoothing factor. The two results are added together to obtain the final acoustic imaging result of the current frame, thereby completing the sound source localization.

[0040] A non-volatile computer storage medium storing computer-executable instructions that execute a fast deconvolution sound source localization method based on array sparsity and broadband integrated processing.

[0041] Compared with existing technologies, the present invention has the following advantages:

[0042] (1) This method significantly reduces the computational load by sparsifying the k-means clustering array while preserving the core characteristics of the array, and improves computational efficiency by combining it with a fast deconvolution algorithm, thus solving the problem of high deployment costs for large array projects; it increases the dynamic range by using function beamforming and uses a frequency domain deconvolution model with Gaussian filtering to suppress noise, thus ensuring positioning accuracy; it uses a broadband integration mechanism that fuses the extreme values ​​of signal-to-noise ratio at multiple frequency points, which does not require prior frequency band information and adaptively synthesizes clear images, thus avoiding weak source flooding and positioning drift; it then uses inter-frame weighted smoothing to suppress cloud image flicker and simultaneously optimizes the weak source identification capability, thus achieving low-power, high-precision, and robust sound source positioning, which is suitable for complex real-world scenarios with unknown spectrum and low signal-to-noise ratio.

[0043] (2) This invention introduces function beamforming and combines it with power exponent to improve the dynamic range of imaging and make it easier to distinguish the differences between sound sources of different intensities; it constructs a frequency domain deconvolution model, converts the convolution into a wavenumber domain product through fast Fourier transform, and combines zero-filling and iterative optimization with Gaussian regularization filtering to suppress high wavenumber noise, accurately solves the sound source intensity distribution matrix, and ensures positioning accuracy; the broadband integration mechanism uses multi-frequency signal-to-noise ratio extreme value fusion to select the frequency components with the most prominent energy and the best signal-to-noise ratio for each spatial grid point, avoiding the imaging bias and weak source submersion problems caused by unknown spectrum and bandwidth in traditional algorithms, and achieving sound source imaging with clear focus and accurate intensity.

[0044] (3) This invention does not require prior determination of the sound source frequency and bandwidth. By independently executing the fast deconvolution algorithm across the entire frequency band, combined with noise basis estimation and signal-to-noise ratio calculation, it adaptively fuses effective information across the entire frequency band, which can stably process actual sound sources with unknown spectrum and low signal-to-noise ratio, and solve the localization drift problem of traditional algorithms. The introduction of a time smoothing factor for inter-frame weighted fusion effectively suppresses the flicker of the imaging cloud map, improves the visual stability of the video output, and enhances the application adaptability in dynamic sound source scenarios.

[0045] (4) The present invention estimates the noise floor and sets a threshold for the sound source intensity distribution matrix, converts the sound source intensity matrix into a sound source distribution signal-to-noise ratio matrix, highlights the difference between weak sources and noise, enhances the identification of weak sources, and ensures that weak sources can be effectively identified in scenarios where multiple sound sources coexist and there are differences in strength. Attached Figure Description

[0046] Figure 1 This is a flowchart of the method of the present invention;

[0047] Figure 2 This is an acoustic imaging image of the present invention in an environment with a signal-to-noise ratio of -30dB;

[0048] Figure 3 This is an acoustic imaging image of the present invention in an environment with a signal-to-noise ratio of -25dB;

[0049] Figure 4 This is an acoustic imaging image of the present invention in an environment with a signal-to-noise ratio of -20dB;

[0050] Figure 5 This is an acoustic imaging image of the present invention in an environment with a signal-to-noise ratio of -15dB;

[0051] Figure 6 This is an acoustic imaging image of the present invention under an environment with a signal-to-noise ratio of -10dB. Detailed Implementation

[0052] like Figure 1 As shown, the present invention provides a technical solution: a fast deconvolution sound source localization method based on array sparsity and broadband integrated processing, comprising the following steps:

[0053] Step S1: Divide the microphone array into several clusters using k-means clustering, and select the microphones closest to the centroid in each cluster to form a sparse array (while reducing subsequent computation and ensuring array spatial coverage).

[0054] Step S2: Convert the regular focusing plane into an irregular focusing plane, and calculate the coordinates of the irregular grid by sampling the sine values ​​of the angles at equal intervals;

[0055] Step S3: Calculate the cross-spectral matrix of the sparse array based on the irregular grid coordinates and perform eigenvalue decomposition. Introduce the power exponent to obtain the function beamforming output (to improve the imaging dynamic range).

[0056] Step S4: Based on the convolution relationship between the function beamforming output, the sound source distribution, and the power-law spread function, construct a frequency domain deconvolution model;

[0057] Step S5: Use Fast Fourier Transform to convert the convolution in the model into wavenumber domain multiplication, and solve for the sound source intensity distribution matrix through zero-padding and iterative optimization;

[0058] Step S6: Estimate the noise floor of the sound source intensity distribution matrix and set a threshold, then convert it into a sound source distribution signal-to-noise ratio matrix (to enhance weak source identification).

[0059] Step S7: For the signal-to-noise ratio matrix of the sound source distribution, traverse all frequency points and take the maximum signal-to-noise ratio of each grid point to fuse them into a broadband integrated imaging matrix;

[0060] Step S8: Introduce a time smoothing factor to weighted fuse the imaging matrices of the current frame and the previous frame to obtain the final acoustic imaging result and complete the sound source localization.

[0061] The specific process of step S1 is as follows:

[0062] Step S1.1: First, determine the location information of all microphones, forming a location matrix containing the total number of microphones, where the coordinates of each microphone correspond to a row of data in the matrix; microphone location matrix. Represented as:

[0063] ;

[0064] In the formula, Indicates the first The x-coordinate of each microphone in a plane coordinate system , Indicates the total number of microphones; Indicates the first The ordinate of each microphone in a plane coordinate system; Represents a real matrix.

[0065] Step S1.2: Use the k-means clustering algorithm to divide the entire microphone array into... The clustering process aims to minimize the intra-cluster sum of squares (the sum of the squares of the distances from all microphones in a cluster to the centroid of that cluster, where the centroid is the average of the coordinates of all microphones within the cluster). The clustering is iteratively adjusted by reassigning microphones to the cluster containing the nearest centroid and updating the centroids of the clusters until the clustering results stabilize.

[0066] Specifically, cluster allocation and centroid update are iteratively performed by minimizing the intra-cluster sum of squares using the objective function:

[0067] The objective function is expressed as:

[0068] ;

[0069] In the formula, Indicates the first Clusters, ; Represents the microphone's position vector; Indicates the number of clusters to be divided; Indicates the first The centroid vector of each cluster;

[0070] Iterative cluster allocation is represented as:

[0071] ;

[0072] In the formula, Indicates the first The first iteration A cluster; Indicates the first The position vectors of each microphone; Indicates the first During the nth iteration The centroid vector of each cluster; Indicates the first During the nth iteration The centroid vector of each cluster;

[0073] The centroid update is represented as:

[0074] .

[0075] Step S1.3: Calculate the distance from all microphones within each cluster to the centroid of that cluster. Choose the microphone closest to the center of mass. As representatives of this cluster, the final array is composed of representative microphones from all clusters. , can be represented as:

[0076] ;

[0077] ;

[0078] ;

[0079] ;

[0080] In the formula, This represents an index array used to record the index of the representative microphone for each cluster in the original array; This indicates that all columns are selected.

[0081] The specific process of step S2 is as follows:

[0082] Step S2.1: Determine the coordinates of the grid points in the regular focusing plane, and the vertical distance between the focusing grid plane and the microphone array plane.

[0083] Step S2.2: By calculating the distance from the grid point in the focusing plane to the reference point, and the trigonometric function relationship between this distance and the horizontal and vertical coordinates, the corresponding angle sine value is obtained.

[0084] Step S2.3: Within the set range of angle sine values ​​(including both horizontal and vertical directions), perform discretization sampling at equal intervals to obtain a series of angle sine values ​​for sampling points.

[0085] Step S2.4: Based on the sine values ​​of the angles obtained from these samples, and combined with the vertical distance between the focusing plane and the microphone array plane, calculate the specific coordinates of each grid point in the irregular focusing plane, and complete the transformation from a regular focusing plane to an irregular focusing plane.

[0086] Wherein, the coordinates of the grid points in the focal plane are The vertical distance between the focusing plane and the microphone array plane is ,pass , ,exist and Perform equidistant divergence sampling within the range:

[0087] ;

[0088] In the formula, This represents the distance from a grid point in the focal plane to the reference point; This represents the vertical component, which is combined with practical calculations. Sure; The sine value representing the angle in the horizontal direction; The sine value representing the angle in the vertical direction; , These represent the minimum and maximum ranges of the sine value of the horizontal angle, respectively. , These represent the minimum and maximum ranges of the sine value of the angle in the vertical direction, respectively. Indicates the horizontal direction. The sine value of the angle at each sampling point; Indicates the vertical direction of the first The sine value of the angle at each sampling point; Indicates the number of sampling points in the horizontal direction; Indicates the number of sampling points in the vertical direction;

[0089] The final result is the irregular grid coordinates:

[0090] ;

[0091] In the formula, Indicates the first Line number Column of irregular grid points; Indicates the first Line number The vertical coordinates of the irregular grid points.

[0092] The specific process of step S3 is as follows:

[0093] Step S3.1: Perform Fourier transform on the sound signal recorded by each microphone in the sparse array to obtain the pressure amplitude data corresponding to different frequencies of each microphone, and then construct a column matrix containing the number of sparsed microphones based on the pressure amplitude data.

[0094] Step S3.2: Multiply the column matrix by its conjugate transpose to obtain the cross-spectral matrix of the sparse array signal.

[0095] Step S3.3: Perform eigenvalue decomposition on the cross-spectral matrix to obtain a unitary matrix composed of eigenvectors of the cross-spectral matrix and a diagonal matrix composed of eigenvalues.

[0096] Step S3.4: Introduce a positive integer power exponent and combine it with the constructed guide vector (the guide vector is calculated by a specific exponential function based on parameters such as microphone array position, sound velocity, and frequency) to calculate the output result of the function beamforming. This result can reflect the power of the sound source at different spatial locations.

[0097] Wherein, the cross-spectral matrix is ​​defined. for:

[0098] ;

[0099] ;

[0100] In the formula, This indicates that after sparsification, each microphone has a frequency... A column matrix composed of pressure amplitude data at a given location; Indicates the first A sparsified microphone at frequency Pressure amplitude at the location, , This indicates the number of microphones after sparsification; express The conjugate transpose of .

[0101] Among them, eigenvalue decomposition is performed on the cross-spectral matrix: , This represents a unitary matrix composed of the eigenvectors of the cross-spectral matrix. express The conjugate transpose of .

[0102] Among them, positive integer exponents are introduced. Calculation function beamforming output:

[0103] ;

[0104] ;

[0105] ;

[0106] ;

[0107] In the formula, Represents the set of positive integers; Indicates that the function beam is formed at the location The output at that location (reflecting the power of the sound source at that location). Indicates the position of the focusing plane in beamforming; Indicates the location where the beam is formed. The weight vector at the location (used to weight the microphone signal to focus on the position) ); express The conjugate transpose of; Indicates position The steering vector, determined by parameters such as the microphone array position, sound velocity, and frequency, describes the propagation characteristics from the sound source to the microphone array. express The conjugate transpose of; Indicates position The guide vector, Indicates the location of the sound source; express The conjugate transpose of; Indicates the location of the sound source Intensity distribution at the location; Indicates the speed of sound; Indicates the first The position vector of each microphone; Indicates wave number; It represents the imaginary unit.

[0108] Among them, location Power estimate at the location Represented as:

[0109] ;

[0110] In the formula, A column matrix representing the pressure amplitude of the microphone array (composed of the pressure amplitude of each microphone); express The conjugate transpose of .

[0111] Among them, in position The cross-spectral matrix at point is represented as:

[0112] .

[0113] Among them, location Power estimate at the location Represented as:

[0114] .

[0115] The specific process of step S4 is as follows: Based on the function beamforming output obtained in step S3, analyze its relationship with the sound source intensity distribution and the power-law spread function. It is found that the function beamforming output can be expressed as the spatial convolution of the sound source intensity distribution and the power-law spread function (the power-law spread function is the result of the power operation of the point spread function, and it is only related to the position of the observation point relative to the sound source, that is, it has spatial shift invariance, and the point spread function at the center of the focusing plane is used as the spatial shift-invariant point spread function). Based on this convolution relationship, a deconvolution model in the frequency domain is constructed (this deconvolution model relates the function beamforming output, the sound source intensity distribution, and the power-law spread function).

[0116] Among them, the relationship between the function beamforming output obtained in step S3 and the sound source intensity distribution and the power-law spread function is analyzed and expressed as follows:

[0117] ;

[0118] ;

[0119] In the formula, This represents the vector indicating the center position of the focusing plane. Represents the power-increasing point spread function;

[0120] Constructing a frequency domain deconvolution model:

[0121] ;

[0122] In the formula, The matrix representation of the sound source intensity distribution; express The corresponding matrix representation; This represents the Fourier transform operator; This represents the inverse Fourier transform operator; This represents the Hadamard product operator in the frequency domain.

[0123] The specific process of step S5 is as follows:

[0124] Step S5.1: Perform zero-padding on the function beamforming output, the sound source intensity distribution matrix, and the power-law spread function to ensure that the dimensions of the zero-padding matrix meet the computational requirements (at least twice the number of rows and columns of the original sound source intensity distribution matrix); the zero-padding function beamforming output is represented as follows: The sound source intensity distribution matrix after zero filling is expressed as follows: The zero-filled power-law spread function is expressed as: .

[0125] Step S5.2: Calculate the Fast Fourier Transform of the power-up spread function after zero-filling, and the sum of the absolute values ​​of the results.

[0126] Step S5.3: Initialize the zero-filled sound source intensity distribution matrix as a zero matrix, and then perform iterative calculations using the iterative formula: First, based on the sound source intensity distribution matrix of the current iteration, the Fourier transform result of the power-law spread function, and the Gaussian regularized filter function (used to suppress high wavenumber noise, whose parameters are related to the cutoff wavenumber, which is determined by parameters such as the grid spacing distance of the focusing plane), calculate the function beamforming estimation result of the current iteration; then, combine the zero-filled function beamforming output to update the sound source intensity distribution matrix (ensuring that the updated matrix value is non-negative).

[0127] Step S5.4: After the iteration is completed, remove the zero-filling part to obtain the final estimated sound source intensity distribution matrix.

[0128] Among them, the calculation is:

[0129] ;

[0130] ;

[0131] In the formula, This represents the wavenumber domain result of the zero-filled power-law point spread function after Fourier transform; This represents the sum of the absolute values ​​of the spread functions at the power-up point in the wavenumber domain;

[0132] Let the zero-filling sound source intensity matrix be in the 0th iteration (initial state). Zero matrix Through the iterative formula:

[0133] ;

[0134] ;

[0135] In the formula, Indicates the first Function beamforming estimation results at the next iteration; , They represent the first sequence The sound source intensity distribution matrix after zero-filling in the next iteration; This represents the Gaussian regularization filter function; , These represent the wavenumber components (two-dimensional components in the frequency domain) in the horizontal and vertical directions, respectively.

[0136] Solve for the sound source intensity distribution matrix, and After removing zero-padding, the final estimated sound source intensity distribution matrix is ​​obtained.

[0137] Among them, the Gaussian regularization filter function Represented as:

[0138] ;

[0139] In the formula, Indicates the cutoff wavenumber. This represents parameters related to filtering characteristics. This indicates the distance between grid points in the focused plane.

[0140] The specific process of step S6 is as follows:

[0141] Step S6.1: Convert the final estimated sound source intensity distribution matrix obtained in step S5 into column vector form, where each column vector contains the sound source intensity data of all grid points at the corresponding frequency.

[0142] Step S6.2: Sort the sound source intensity column vectors corresponding to each frequency in ascending order; select the minimum value of a certain proportion (the proportion value is between 0 and 1) according to the length of the sorted vectors, calculate the average value of these minimum values, and use it as the noise floor estimate at that frequency.

[0143] Step S6.3: Set a threshold based on the noise floor estimate and process the sound source intensity matrix at each frequency: if the sound source intensity at a certain grid point is less than or equal to the noise floor estimate, then set the value of that point to 0; if it is greater than the noise floor estimate, then divide the sound source intensity at that point by the noise floor estimate, and finally obtain the sound source distribution signal-to-noise ratio matrix at each frequency.

[0144] In step S5, the noise floor threshold is estimated using the final estimated sound source intensity distribution matrix. First, the final estimated sound source intensity distribution matrix... Vectorization is as follows:

[0145] ;

[0146] In the formula, This represents the final estimated sound source intensity distribution matrix. The vectorized column vector contains the sound source intensity data of all grid points at the corresponding frequency; Indicates transpose;

[0147] right Sort in ascending order to get:

[0148] ;

[0149] ;

[0150] In the formula, Indicates to The column vector obtained after sorting in ascending order; Indicates the sorted order of the first... The sound source intensity value of each element; Indicates the total number of grid points;

[0151] Pick The noise basis estimate at that frequency point is obtained by averaging the minimum values.

[0152] ;

[0153] ;

[0154] In the formula, Represents frequency The noise basis estimate is as follows; This indicates the proportion of selecting the minimum value; Indicates rounding down; Indicates the sorted order of the first... The sound source intensity value of each element;

[0155] Then, based on the noise floor threshold, the sound source distribution matrix is ​​converted into a sound source distribution signal-to-noise ratio matrix:

[0156] ;

[0157] In the formula, Represents frequency Next Line number Elements of the signal-to-noise ratio matrix representing the sound source distribution at column grid points; Represents frequency Below, the first in the focused planar mesh Line number The sound source intensity value corresponding to the grid point.

[0158] The specific process of step S7 is as follows: For the source distribution signal-to-noise ratio (SNR) matrix obtained in step S6 at all frequencies, for each spatial grid point, iterate through the SNR values ​​corresponding to that grid point at all frequencies, and select the largest SNR value as the final imaging intensity value of that grid point; integrate the maximum SNR values ​​of all grid points to form a broadband integrated imaging matrix, represented as:

[0159] ;

[0160] In the formula, Indicates the first In the frame broadband synthetic imaging matrix, the first Line number The final image intensity value of the column grid points; Indicates the first frequency points The signal-to-noise ratio matrix of the lower sound source distribution Line number Signal-to-noise ratio values ​​for grid points; This indicates the total number of frequency points participating in broadband convergence.

[0161] The specific process of step S8 is as follows: Introduce a time smoothing factor with a value between 0 and 1. The broadband composite imaging matrix obtained in step S7 is subjected to inter-frame weighted fusion processing. During fusion, the value of each grid point in the current frame's broadband composite imaging matrix is ​​multiplied by a time smoothing factor, and the value of the corresponding grid point in the previous frame's imaging matrix is ​​multiplied by 1 minus the time smoothing factor. The two results are added together to obtain the final acoustic imaging result of the current frame, thereby completing the sound source localization, as shown below:

[0162] ;

[0163] In the formula, Indicates the first The first frame in the final acoustic imaging result Line number The imaging intensity values ​​of the grid points; Indicates the first In the frame broadband synthetic imaging matrix, the first Line number The final imaging intensity value of the grid points.

[0164] A non-volatile computer storage medium storing computer-executable instructions that execute a fast deconvolution sound source localization method based on array sparsity and broadband integrated processing.

[0165] Figures 2-6 Acoustic imaging images were obtained in a low signal-to-noise ratio environment of -30dB to -10dB using the method of the present invention. The test environment was as follows: a 19KHz sound source was placed at the center of two noise sources, with the sound source position at the origin and a vertical distance of 1m from the array plane.

[0166] from Figure 2 As can be seen, even in an environment with extremely low signal-to-noise ratio, the core location of the sound source can still be clearly presented, with no obvious imaging blurring or positional shift. Background noise interference is effectively suppressed, and the distinction between weak sources and noise is significant, proving the anti-interference capability and positioning effectiveness of this invention under extremely poor noise conditions.

[0167] from Figure 3As can be seen, the sound source imaging focusing is further improved, the main lobe width is narrow, the side lobe suppression effect is excellent, there is no situation where the side lobe obscures the weak source or interferes with the localization, the sound source position accurately corresponds to the actual origin position, the imaging stability is strong, and the high resolution advantage of the method is demonstrated in medium and low signal-to-noise ratio scenarios.

[0168] from Figure 4 As can be seen, the imaging clarity and contrast are improved, the sound source outline is clearer, there is no imaging distortion (such as stretching or shifting), and the background noise floor is further reduced, verifying the actual effect of this method in suppressing spatial displacement variation of the point spread function through irregular focusing plane and noise suppression by Gaussian regularization filtering.

[0169] from Figure 5 As can be seen, the sound source localization accuracy remains stable with no localization drift, the imaging cloud map shows no flickering, and the details of weak sources (if there are multiple sources or weak source components) are clearly distinguishable, indicating the effective role of the broadband integration strategy (multi-frequency signal-to-noise ratio extreme value fusion) and the inter-frame smoothing mechanism in this signal-to-noise ratio range.

[0170] from Figure 6 As can be seen, the imaging effect is optimal, with high accuracy and strong focusing of sound source localization, significant main-side lobe ratio, and no noise redundancy information interference, which fully meets the engineering-level localization requirements. This further proves the robustness and practicality of the method in the entire low signal-to-noise ratio range (-30dB to -10dB), and stable imaging can be achieved without prior sound source frequency and bandwidth information.

[0171] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A fast deconvolution sound source positioning method based on array sparsity and broadband comprehensive processing, characterized in that, Comprising the following steps: Step S1: divide the microphone array into several clusters by k-means clustering, and select the microphone closest to the cluster center in each cluster to form a sparse array; Step S2: convert the regular focusing plane to an irregular focusing plane, and calculate the irregular grid coordinates by angle sine value equidistant sampling; Step S3: calculate the cross-spectral matrix of the sparse array based on the irregular grid coordinates and perform eigenvalue decomposition, and introduce a power index to obtain the function beamforming output; Step S4: according to the convolution relationship between the function beamforming output and the sound source distribution and the up-power point spread function, a frequency domain deconvolution model is constructed; Step S5: the convolution in the model is converted to a wavenumber domain product by using fast Fourier transform, and the sound source intensity distribution matrix is solved by zero padding and iterative optimization; Step S6: estimate the noise basis of the sound source intensity distribution matrix and set a threshold to convert it to a sound source distribution signal-to-noise ratio matrix; Step S7: for the sound source distribution signal-to-noise ratio matrix, traverse all frequency points, and take the maximum signal-to-noise ratio of each grid point to fuse into a wideband comprehensive imaging matrix; Step S8: introduce a time smoothing factor to weight and fuse the imaging matrix of the current frame and the last frame to obtain the final acoustic imaging result and complete the sound source positioning; The specific process of step S2 is as follows: Step S2.1: determine the coordinates of the grid points in the regular focusing plane, and the vertical distance between the focusing grid plane and the microphone array plane; Step S2.2: by calculating the distance from the grid point in the focusing plane to the reference point, and the trigonometric function relationship between the distance and the horizontal and vertical coordinates, the corresponding angle sine value is obtained; Step S2.3: within the set range of angle sine value, discrete sampling is carried out in an equidistant manner to obtain a series of angle sine values of sampling points; Step S2.4: according to the angle sine values obtained by sampling, combined with the vertical distance between the focusing plane and the microphone array plane, the specific coordinates of each grid point in the irregular focusing plane are calculated, and the conversion from the regular focusing plane to the irregular focusing plane is completed; The specific process of step S3 is as follows: Step S3.1: Fourier transform is performed on the sound signals recorded by each microphone in the sparse array to obtain the pressure amplitude data corresponding to different frequencies for each microphone, and then a column matrix containing the number of sparse microphones is constructed based on the pressure amplitude data; Step S3.2: the cross-spectral matrix of the sparse array signal is obtained by multiplying the column matrix with its conjugate transpose matrix; Step S3.3: the eigenvalue decomposition is performed on the cross-spectral matrix to obtain a unitary matrix composed of the eigenvectors of the cross-spectral matrix and a diagonal matrix composed of the eigenvalues; Step S3.4: a positive integer power index is introduced, and the function beamforming output is calculated combined with the constructed steering vector; the steering vector is calculated by the microphone array position, sound speed, frequency, and exponential function; The specific process of step S4 is as follows: according to the function beamforming output obtained in step S3, analyze the relationship between it and the sound source intensity distribution and the up-power point spread function: the function beamforming output is represented as the spatial convolution of the sound source intensity distribution and the up-power point spread function; Based on the function beamforming output represented as the spatial convolution of the sound source intensity distribution and the elevated point spread function, a deconvolution model in the frequency domain is constructed.

2. The fast deconvolution acoustic source localization method based on array sparsity and broadband synthesis processing according to claim 1, characterized in that: The specific process of step S1 is as follows: Step S1.1: First, the position information of all microphones is determined to form a position matrix containing the total number of microphones, wherein the coordinates of each microphone correspond to a row of data in the matrix; Step S1.2: The k-means clustering algorithm is used to divide the entire microphone array into several clusters, and the clustering process aims to minimize the within-cluster sum of squares; the cluster assignment of the microphones is adjusted and the cluster centroids are updated through continuous iteration until the clustering result is stable; Step S1.3: The distance of all microphones in each cluster to the cluster centroid is calculated, and the microphone closest to the centroid is selected as the representative of the cluster, and finally the representative microphones of all clusters form the sparse array.

3. The fast deconvolution acoustic source localization method based on array sparsity and broadband synthesis processing according to claim 2, characterized in that: The specific process of step S5 is as follows: Step S5.1: The function beamforming output, the sound source intensity distribution matrix, and the elevated point spread function are zero-padded; Step S5.2: The fast Fourier transform result of the zero-padded elevated point spread function is calculated, as well as the sum of the absolute values of the result; Step S5.3: The zero-padded sound source intensity distribution matrix is initialized to a zero matrix, and then iterative calculation is performed through the iterative formula: first, the function beamforming estimation result of the current iteration is calculated according to the sound source intensity distribution matrix of the current iteration, the Fourier transform result of the elevated point spread function, and the Gaussian regularization filter function; then the sound source intensity distribution matrix is updated combined with the zero-padded function beamforming output; Step S5.4: After the iteration is completed, the zero-padded part is removed to obtain the final estimated sound source intensity distribution matrix.

4. The fast deconvolution acoustic source localization method based on array sparsity and broadband synthesis processing according to claim 3, characterized in that: The specific process of step S6 is as follows: Step S6.1: The final estimated sound source intensity distribution matrix obtained in step S5 is converted into a column vector form, and each column vector contains the sound source intensity data of all grid points at the corresponding frequency; Step S6.2: The sound source intensity column vector corresponding to each frequency is sorted in ascending order; according to the length of the sorted vector, the minimum values of the preset proportion are selected, and the average value of these minimum values is calculated as the noise floor estimate value at the corresponding frequency; Step S6.3: According to the noise floor estimate value, a threshold is set to process the sound source intensity matrix at each frequency: if the sound source intensity of a grid point is less than or equal to the noise floor estimate value, the value of the point is set to 0; if it is greater than the noise floor estimate value, the sound source intensity of the point is divided by the noise floor estimate value, and finally the sound source distribution signal-to-noise ratio matrix at each frequency is obtained.

5. The fast deconvolution acoustic source localization method based on array sparsity and broadband synthesis processing according to claim 4, characterized in that: The specific process of step S7 is as follows: For the sound source distribution signal-to-noise ratio matrices at all frequencies obtained in step S6, for each spatial grid point, the signal-to-noise ratio values corresponding to the grid point at all frequencies are traversed, and the maximum signal-to-noise ratio value is selected as the final imaging intensity value of the grid point; all the maximum signal-to-noise ratio values of the grid points are integrated to form a broadband comprehensive imaging matrix.

6. The fast deconvolution acoustic source localization method based on array sparsity and broadband synthesis processing according to claim 5, characterized in that: The specific process of step S8 is: introducing a time smoothing factor, and performing inter-frame weighted fusion processing on the wideband comprehensive imaging matrix obtained in step S7; during the fusion, multiplying the value of each grid point in the current frame wideband comprehensive imaging matrix by the time smoothing factor, multiplying the value of the corresponding grid point in the previous frame imaging matrix by 1 minus the time smoothing factor, adding the two results to obtain the final acoustic imaging result of the current frame, so as to complete the sound source positioning.

7. A non-transitory computer storage medium storing computer-executable instructions, the computer-executable instructions comprising instructions for: The computer executable instructions perform the fast deconvolution sound source positioning method based on array sparsity and wideband comprehensive processing according to any one of claims 1-6.

Images