Multi-element microphone array sound source localization method based on multistage signal preprocessing and subspace spectrum optimization

Abstract

The invention relates to a multi-element microphone array sound source localization method based on multistage signal preprocessing and subspace spectrum optimization, and belongs to the technical field of acoustic detection. Aiming at the problems of poor noise immunity, weak multi-sound-source resolution capability and low calculation efficiency of the existing sound source positioning technology, a triple signal preprocessing and subspace collaborative optimization scheme is provided; firstly, incoherent noise is suppressed through phase coherent filtering, a signal is reconstructed through principal component analysis, and phase deviation is calibrated through fundamental frequency; then constructing a guiding matrix and decomposing a noise subspace, and extracting a coarse positioning result; and finally, high-precision angle optimization is realized based on a chaos initialization differential evolution algorithm, and the efficiency is improved by combining a dynamic search range and an early stop mechanism. According to the method, the anti-interference capability in a low signal-to-noise ratio environment is remarkably enhanced, the problems of missing detection and false detection during dense distribution of multiple sound sources are effectively solved, meanwhile, the positioning precision and the real-time performance are considered, and the method is suitable for acoustic fault detection of complex scenes such as power transmission line inspection.

Description

Technical Field

[0001] The present invention belongs to the technical field of acoustic detection and relates to a multi-microphone array sound source localization method based on multi-stage signal preprocessing and subspace spectrum optimization. Background Art

[0002] The safe and stable operation of power systems is highly dependent on the efficient detection of hidden faults in transmission lines. Traditional manual inspection methods suffer from low efficiency, high risk, and a high rate of missed detections. The limitations of manual inspection are particularly pronounced in complex terrain or harsh environments. Acoustic source localization technology, which captures the specific acoustic signals generated by faults such as partial discharge, can accurately pinpoint the location of faults and is becoming a research hotspot in the industry.

[0003] However, the performance of sound source localization systems is highly dependent on the design of the sensor array and the robustness of the signal processing algorithm. The array geometry directly impacts localization accuracy and noise immunity. Existing research has shown that appropriate array design (such as circular arrays) can improve robustness in environments with complex background noise and multipath interference.

[0004] Among the current mainstream sound source localization methods, algorithms based on time delay difference (TDOA) are widely used due to their simplicity and strong real-time performance. However, they require solving a system of nonlinear equations, which significantly increases computational complexity in multi-sensor scenarios. They are also sensitive to noise, with performance significantly degrading in low signal-to-noise ratio environments. The MUSIC algorithm, based on high-resolution spectral estimation, can improve localization accuracy, but it relies on the orthogonality of the noise subspace and is susceptible to interference from false peaks in low signal-to-noise ratio environments. Missed or false detections are a prominent problem when multiple sound sources are densely distributed.

[0005] In recent years, the introduction of intelligent optimization algorithms and neural networks has provided new directions for sound source localization. For example, particle swarm optimization (PSO) has been used to optimize localization results, but it is prone to falling into local optima, and its convergence speed is significantly affected by the initial population. Methods based on convolutional neural networks (CNNs) offer excellent noise immunity, but require large amounts of training data, struggle to adapt to dynamic environments, and lack real-time performance. To address these limitations, a variety of improvements have been proposed. Research has shown that differential evolution (DE) algorithms, while introduced to the optimization field due to their global search capabilities, are prone to premature convergence and search stagnation in complex problems. To address this, researchers have proposed strategies such as dynamic parameter adjustment and chaotic initialization. For example, through subpopulation division and parameter adaptive design, these strategies have significantly improved convergence speed and accuracy. Furthermore, signal preprocessing techniques (such as phase coherence filtering and principal component analysis) have been shown to effectively suppress noise and enhance the robustness of subsequent localization algorithms.

[0006] Existing methods still have the following bottlenecks: insufficient signal preprocessing and limited noise suppression capabilities of traditional methods, resulting in positioning failure in low signal-to-noise ratio scenarios; low algorithm robustness, and the MUSIC algorithm is prone to missed detections or false detections when multiple sound sources are densely distributed; it is difficult to balance computational efficiency and accuracy, and intelligent optimization algorithms either have excessively high computational complexity or premature convergence, making it difficult to strike a balance between real-time performance and accuracy.

[0007] To address the above challenges, the present invention adopts a two-dimensional multi-element circular microphone array to adapt to the characteristics of wide bandwidth and high signal-to-noise ratio in transmission line detection scenarios, and proposes a sound source localization method that integrates improved signal preprocessing and subspace spectrum optimization. Through triple preprocessing of phase coherent filtering, principal component analysis (PCA) reconstruction and phase calibration, the signal noise resistance is significantly improved; based on signal subspace decomposition, coarse positioning results are extracted, and high-precision angle optimization is achieved through an adaptive optimization algorithm, forming a collaborative optimization mechanism from coarse positioning to high-precision optimization. Compared with traditional schemes, this method effectively avoids local optimal problems through adaptive parameter strategies and chaotic initialization design, while maintaining high robustness in complex noise environments, providing an efficient and reliable technical solution for transmission line fault detection. Summary of the Invention

[0008] In view of this, the purpose of the present invention is to provide a multi-microphone array sound source localization method based on multi-level signal preprocessing and subspace spectrum optimization, which has high accuracy in locating high-frequency narrowband sound sources such as leakage current, and through triple preprocessing of phase coherence filtering, principal component analysis (PCA) reconstruction and phase calibration, it enhances the anti-interference ability in low signal-to-noise ratio environment, and uses the subspace spectrum optimization algorithm to improve the local optimal problem in the traditional subspace spectrum method, so as to adapt to the detection needs of complex noise environment and multiple sound source scenes in transmission line inspection.

[0009] In order to achieve the above object, the present invention provides the following technical solutions:

[0010] A multi-microphone array sound source localization method based on multi-stage signal preprocessing and subspace spectrum optimization is characterized by comprising the following steps:

[0011] S1: Establish a model for microphone array to collect sound source signals and generate time-domain discrete signals;

[0012] S2: Multi-stage signal preprocessing including phase coherent filtering, PCA reconstruction and phase calibration;

[0013] S3: Construct a steering matrix based on the preprocessed signal, extract the coarse positioning results through signal subspace decomposition, and achieve high-precision angle optimization through an adaptive optimization algorithm, forming a collaborative optimization mechanism from coarse positioning to high-precision optimization.

[0014] Furthermore, in S1, a simulation model for a multi-element microphone array to receive high-frequency narrowband signals is constructed, specifically including the following steps:

[0015] S11: Simulate the complex signal generated by the sound source to the microphone array: Assume N source The pitch angle θ of the sound source i and azimuth φ i , combined with N mic The position information of the microphone array composed of microphones can be used to obtain the phase delay τ of the microphone to each sound source direction. m,i And construct the steering vector a m,i :

[0016]

[0017]

[0018] Where c is the speed of sound, x m and y m is the Cartesian coordinate of the mth microphone, λ is the wavelength, a m,j is the steering vector of the i-th sound source at the m-th microphone, and all a m,i Constructing a steering matrix

[0019] First, define the signal-to-noise ratio (SNR):

[0020]

[0021] Using the sound pressure level SPL and the reference sound pressure P ref =20μpa, the amplitude amp is obtained by inversely solving the relationship and generating the complex Gaussian signal Q s :

[0022] amp = 2 × 10 -5 10 SPL / 20

[0023] Composition Q s The row vector q s(i) for:

[0024]

[0025] N is a simplified expression of Gaussian distribution, amp i is the amplitude of the i-th sound source.

[0026] S12: According to the signal power P true The noise power is calculated by adding the signal-to-noise ratio (SNR) and generating the noise matrix through the complex Gaussian distribution:

[0027]

[0028] P true,i and amp i is the power and amplitude of the i-th sound source, SNR is the signal-to-noise ratio (in dB), is the noise matrix, and L is the number of sampling points of the microphone array.

[0029] Combined with the steering matrix A, the model of the array signal Y is obtained:

[0030] Y=AQ s +N

[0031] in is the sound source signal matrix.

[0032] Furthermore, in S2, the signal collected by the microphone is preprocessed by phase coherence weighting, PCA signal reconstruction, and fundamental frequency phase alignment calibration method, which specifically includes the following steps:

[0033] S21: Original signal acquisition and preprocessing initialization:

[0034] Multi-channel signal acquisition. Using N array channels installed on the UAV platform, the acoustic signal of duration L is collected in the form of complex sampling to form a two-dimensional complex signal matrix:

[0035]

[0036] Where Y is the acquired two-dimensional complex signal matrix; y n (t) is the complex sound signal collected by the nth channel at time t; N is the number of channels; L is the sampling length of each channel (number of time points); C is the complex space;

[0037] Data normalization and outlier cleaning: Perform amplitude normalization on the acquired signal and remove outliers (such as NaN, Inf, etc.) to provide clean input for subsequent filtering steps.

[0038] S22: Channel weighted filtering based on phase coherence:

[0039] Extract the phase information of the complex signal. Extract the phase information of each channel in the signal at each moment and construct the phase matrix;

[0040] Φ∈R N×L ,Φ n,t =∠Y n,t

[0041] Where Φ is a two-dimensional phase matrix that stores the phase information of all channels at each moment; Φ n,t is the phase of the nth channel at the tth moment; Y n,t is an element in the complex signal matrix; ∠· is the operator for taking the complex phase;

[0042] Calculate the phase offset variance. For the phase differences of all channels, perform variance estimation based on a sliding window to measure signal cooperativity; first perform de-averaging on the phase values ​​of each column (at each moment) to eliminate the overall offset:

[0043]

[0044] Where, ΔΦ n,t is the phase difference after removing the mean value; Φ i,t is the phase value of the i-th channel at the t-th moment; is the average phase of all channels at that moment; Φ n,t is the original phase;

[0045] Then, at each time point, the variance of the phase shift ΔΦ is estimated using a sliding window to obtain the phase consistency measure

[0046] Construct phase coherence weighting coefficients. Using the phase variance obtained above, generate a corresponding weighting factor matrix W∈R at each time point N×L , used to adaptively suppress low-coherence channels:

[0047]

[0048] Where W n,t is the phase weighting factor of channel n at time t; α is the empirical adjustment factor (controlling the influence of variance on weight); ε is a very small positive number used to prevent zero division and numerical instability;

[0049] Perform weighted filtering on the original signal. Apply a dynamic weighting matrix to the original complex signal matrix and perform element-wise multiplication on the values ​​at each channel-time point to enhance coherence and suppress incoherent components.

[0050] S23: Signal reconstruction based on principal component analysis:

[0051] Construct the Hankel matrix for each channel. Before signal reconstruction, the filtered time series signal of each channel is first converted into a Hankel matrix to model the local dynamic pattern of the signal in the time domain;

[0052] Perform singular value decomposition (SVD) on the Hankel matrix. Decomposition yields principal components and singular value vectors:

[0053] H=U·S·V T

[0054] Where H is the Hankel matrix; U∈R m×m is the left singular vector matrix; V∈R n×n is the right singular vector matrix; S∈R m×nis a diagonal matrix whose diagonal elements s1≥s2≥…≥s r ≥0 is a singular value;

[0055] Determine the number of principal components to retain, reconstruct the signal, and restore it to a time series form. Reconstruct the matrix using principal components and restore it to a time series form by diagonal averaging.

[0056] S24: Channel calibration based on fundamental frequency phase alignment:

[0057] Perform a fast Fourier transform (FFT) on the reconstructed signal. Perform a one-dimensional fast Fourier transform on the reconstructed time domain signal of all channels and convert it to the frequency domain to extract the phase information of specific frequency components:

[0058] Y i (f) = F[y i (t)],i=1,2,...,N

[0059] Where y i (t) is the reconstructed signal of the i-th channel; F[·] is the Fourier transform operator; Y i (f) is the frequency domain response of the i-th channel at frequency f;

[0060] Select the first channel as the phase reference. Extract the fundamental frequency phase of the first channel as the reference. All subsequent channels will use this as the alignment target to eliminate relative phase errors.

[0061] Calibrate the phase deviation of each channel and uniformly apply the fundamental frequency alignment phase correction to the entire spectrum to ensure the alignment of the main frequency terms without affecting the signal frequency domain energy distribution;

[0062] Perform an inverse transform to obtain the corrected time domain signal. Perform an inverse fast Fourier transform (IFFT) on each channel to restore the calibrated frequency domain signal to a real-valued time domain waveform.

[0063]

[0064] Where, is the real-valued time domain signal output after calibration; Y i calib (f) is the complex signal after the frequency domain phase is aligned; F -1 [·] is the inverse Fourier transform operator; Re{·} is the real part operation (because it is still a complex number after the transformation).

[0065] At the same time, the results are amplitude normalized to prevent output signal overflow or imbalance.

[0066] Furthermore, in S3, a steering matrix is ​​constructed based on the preprocessed signal, the coarse positioning result is extracted through signal subspace decomposition, and high-precision angle optimization is achieved through an adaptive optimization algorithm, forming a collaborative optimization mechanism from coarse positioning to high-precision optimization.

[0067] Coarse positioning specifically includes the following steps:

[0068] S31: Calculate the covariance matrix based on the noisy time domain signal received by the microphone array.

[0069] R=(Y·Y′) / L

[0070] Where Y is the microphone array signal obtained after multi-stage signal preprocessing in step S2, L is the number of sampling points of the microphone array; Y′ is the conjugate transpose of Y; and R is the covariance matrix of the signal.

[0071] S32: Covariance matrix R = EΛE H Perform eigendecomposition, where Λ is the diagonal eigenvalue matrix, E=[E s ,E n ], E s 、E n They are signal subspace and noise subspace respectively. Use descending order of eigenvalues ​​and select the NN s The eigenvectors form the noise subspace matrix E n .

[0072] S33: Based on the steering matrix a(θ,φ) and the noise subspace E n , calculate the spatial spectrum power P(θ,φ):

[0073]

[0074] Where θ is the elevation angle, i.e. the angle between the received signal vector and the array plane; φ is the azimuth angle, i.e. the horizontal direction angle.

[0075] S34: Perform robust peak detection and screening on the spatial spectrum power P(θ,φ) to obtain a preliminary positioning result set (θ,φ).

[0076] Extract P by comparing neighborhood extreme values candidate Candidate peak power set:

[0077] P candidate ={(θ i ,φ i ,P i )|P i >P(θ±Δθ,φ±Δφ)}

[0078] Then adopt a dual screening mechanism combining adaptive threshold and spatial distance constraint:

[0079] Relative power threshold:

[0080] P(θ i ,φ i )>α·max(P),α∈[0.2,0.5]

[0081] Where, P(θ i ,φ i ) is the spatial spectrum power in the candidate peak direction; max(P) is the global maximum value of the spatial spectrum power; α is the adaptive threshold coefficient used to filter low-power false peaks, and the empirical value range is [0.2, 0.5].

[0082] Angle separation constraint:

[0083] Δθ ij =arccos(sinθ i sinθ j +cosθ i cosθ j cos(φ i -φ j ))>β,β∈[3°,10°]

[0084] Where Δθ ij are two candidate peak directions (θ i ,φ i ) and (θ j ,φ j ) between them, the angle is calculated using the vector dot product in the spherical coordinate system to ensure that the sound source directions corresponding to adjacent peaks do not overlap, and β is the minimum resolution interval;

[0085] High-precision optimization specifically includes the following steps

[0086] S35: Dynamic search range calculation. Based on the preliminary positioning result set (θ, φ), calculate the upper and lower limits [lb, ub] of the search range of the current sound source:

[0087] δ=min(base_delta,0.5×min(angular_distance))

[0088] [θ min ,φ min ]=[θ c -δ,φ c -δ]

[0089] [θ max ,φ max ]=[θ c +δ,φ c +δ]

[0090] lb=[θmin ,φ min ]

[0091] ub=[θ max ,φ max ]

[0092] Where base_delta is the basic search radius, angular_distance is the minimum spherical angular distance between the current sound source and other sound sources, min is the minimum value operation, and the minimum value of base_delta and 0.5×angular_distance is taken as the dynamically adjusted search radius δ to ensure that the search areas do not overlap in multi-sound source scenarios; θ c 、φ c are the current sound source coordinates respectively.

[0093] S36: Chaos initialization population. Use Tent chaos map to generate the initial population, set the population size to NP, and define the chaotic sequence The initial value x1 is a random number, and the subsequent sequence recursion generation satisfies:

[0094]

[0095] Map the chaotic sequence to the search space and obtain the coordinates of the initial population pop:

[0096] pop=lb+x i ×(ub-lb)

[0097] S37: Dynamic parameter adjustment. Update the mutation factor F and crossover probability CR in each iteration:

[0098]

[0099] Where gen is the current iteration number; maxGen is the total iteration number; F0 is the basic variation factor; CR0 is the basic crossover probability; F adjust , CR adjust In order to dynamically adjust the coefficients and control the rate of change of the mutation factor F and the crossover probability CR during iteration, F decreases with iteration to shift to local development, and CR increases to enhance convergence.

[0100] S38: Mutation and boundary processing. For each individual, the DE / currenttobest / 1 strategy is used to generate the mutant individual:

[0101] mutant=x i +F×(x best -x i )+F×(x rl -x i )

[0102] Where x i is the current individual coordinate; x best is the current optimal individual coordinate; x r1 are the coordinates of a randomly selected individual.

[0103] Reflection method is used to correct individuals who cross the boundary:

[0104]

[0105] Where x(d) is the θ or φ coordinate of the out-of-bounds individual; lb(d) and ub(d) are the upper and lower limits of the search range θ or φ.

[0106] S39: Crossover and selection. Perform binomial crossover to generate trial individuals:

[0107]

[0108] Where trial(d) is the θ or φ coordinate of the trial individual. By comparing the random index rand(1) with the crossover probability CR, the parameters of the mutant individual or the original individual are selected for mixing.

[0109] Elite retention strategy: If the fitness of the experimental individual P(θ,φ) is better (that is, the power matching is higher), the experimental individual replaces the original individual. The individual fitness calculation formula is as follows:

[0110]

[0111] In the formula, the steering vector F(θ,φ) represents the phase delay characteristics of the sound source direction (θ,φ), F H is the conjugate transpose of the steering vector F(θ,φ); y(k) is the signal vector received by the microphone array at time k; L is the number of sampling points of the microphone array.

[0112] 310: Early stopping mechanism. If there is no improvement in the optimal fitness for 10 consecutive generations, the iteration is terminated early.

[0113] Furthermore, in S2, the phase coherent weighted filtering is characterized by:

[0114] The phase variance is calculated using a sliding window with a window length of T. w Meet T w =0.1L, where L is the total number of signal sampling points;

[0115] The weighting factor is calculated as Where α = 0.5, ε = 10 -6 .

[0116] Furthermore, in S3, the adaptive optimization algorithm is an improved differential evolution algorithm, which is characterized by:

[0117] The population is initialized using the Tent chaotic map, and the chaotic sequence satisfies:

[0118]

[0119] Among them, x i is the i-th value of the chaotic sequence, and the initial value x1 is a random number in the interval [0,1];

[0120] The dynamic adjustment formula of the mutation factor F and the crossover probability CR is:

[0121]

[0122] Among them, F0=0.8, CR0=0.6, gen is the current number of iterations, and maxGen is the total number of iterations.

[0123] Furthermore, in S3, the dynamic search range calculation is characterized by:

[0124] The basic search radius is Δ base =10°;

[0125] The final search radius δ is calculated as follows:

[0126] δ=min(Δ base ,0.5·min(Δθ ij ))

[0127] where Δθ ij is the current sound source direction (θc,φ c ) and other sound source directions (θ j ,φ j ), the minimum spherical angle distance between the two is calculated as follows:

[0128] Δθ ij =arccos(sinθ c sinθ j +cosθ c cosθ j cos(φ c -φ j )).

[0129] Furthermore, the method further includes post-processing S4:

[0130] Perform spatial clustering on the optimized sound source angles (θ, φ) to eliminate isolated points;

[0131] The positioning result is corrected by energy-weighted average, and the formula is:

[0132]

[0133] Among them, P i is the spatial spectrum power of the i-th candidate sound source.

[0134] The beneficial effects of the present invention are:

[0135] (1) The triple preprocessing process of phase coherence filtering, principal component analysis reconstruction and fundamental frequency phase calibration can deeply suppress environmental noise and inter-channel phase deviation, enhance signal purity and consistency, and make sound source feature extraction more reliable in low signal-to-noise ratio scenarios.

[0136] (2) The coarse positioning results based on signal subspace decomposition guide the adaptive optimization algorithm to perform high-precision angle search, forming a collaborative mechanism of global coarse screening and local fine adjustment. This not only avoids the pseudo-peak interference caused by the destruction of the orthogonality of the noise subspace in the traditional algorithm, but also overcomes the problem of premature convergence of the intelligent optimization algorithm, thus achieving stable and high-precision positioning in complex sound field environments.

[0137] (3) Through dynamic search range calculation and chaos initialization strategy, the optimization space of different sound sources is adaptively adjusted to ensure that the search areas do not overlap when multiple sound sources are densely distributed. Combined with the dual screening mechanism of spatial spectrum power, the risk of missed detection and false detection is significantly reduced.

[0138] (4) The signal preprocessing stage reconstructs and effectively compresses the redundant data size through principal component analysis. The optimization algorithm adopts dynamic parameter adjustment and early stopping mechanism to significantly reduce the number of iterations while ensuring positioning accuracy, meeting the stringent real-time requirements of scenarios such as drone inspections.

[0139] (5) The phase calibration step can eliminate the phase inconsistency problem caused by microphone array hardware errors and installation offsets. The guidance matrix construction method is compatible with a variety of array configurations and is suitable for sound source positioning needs in multiple scenarios such as power transmission lines and industrial equipment, providing universal technical support for the detection of hidden defects in equipment.

[0140] Other advantages, objects, and features of the present invention will be described in part in the following description and, in part, will be apparent to those skilled in the art upon examination of the following description or may be learned from practice of the present invention. The objects and other advantages of the present invention may be realized and obtained through the following description. BRIEF DESCRIPTION OF THE DRAWINGS

[0141] In order to make the purpose, technical solutions and advantages of the present invention more clear, the present invention will be described in detail below with reference to the accompanying drawings, in which:

[0142] Figure 1 This is a flow chart of the sound source localization method of the present invention;

[0143] Figure 2A microphone array configuration;

[0144] Figure 3 This is the block diagram of the multi-stage signal preprocessing algorithm;

[0145] Figure 4 Block diagram of the subspace spectrum optimization algorithm;

[0146] Figure 5 The MUSIC algorithm prediction graph for the number of simulated sound sources is 15;

[0147] Figure 6 The prediction graph of the subspace spectrum optimization algorithm for the number of simulated sound sources is 15;

[0148] Figure 7 Comparison of positioning errors between the subspace spectrum optimization algorithm and the MUSIC algorithm under different SNRs. DETAILED DESCRIPTION

[0149] The following describes the embodiments of the present invention by means of specific examples, and those skilled in the art can easily understand other advantages and effects of the present invention from the contents disclosed in this specification. The present invention can also be implemented or applied through other different specific embodiments, and the details in this specification can also be modified or changed in various ways based on different viewpoints and applications without departing from the spirit of the present invention. It should be noted that the illustrations provided in the following embodiments are only schematic illustrations of the basic concept of the present invention, and the following embodiments and features in the embodiments can be combined with each other without conflict.

[0150] Among them, the accompanying drawings are only for illustrative purposes and represent only schematic diagrams rather than actual pictures, and should not be understood as limiting the present invention. In order to better illustrate the embodiments of the present invention, some parts of the accompanying drawings may be omitted, enlarged or reduced, and do not represent the dimensions of actual products. For those skilled in the art, it is understandable that some well-known structures and their descriptions may be omitted in the accompanying drawings.

[0151] The same or similar numbers in the drawings of the embodiments of the present invention correspond to the same or similar parts; in the description of the present invention, it should be understood that if there are terms such as "upper", "lower", "left", "right", "front", "back", etc. indicating directions or positional relationships, they are based on the directions or positional relationships shown in the drawings. They are only for the convenience of describing the present invention and simplifying the description, and do not indicate or imply that the device or element referred to must have a specific direction, be constructed and operate in a specific direction. Therefore, the terms describing the positional relationship in the drawings are only used for illustrative purposes and cannot be understood as limiting the present invention. For ordinary technicians in this field, the specific meanings of the above terms can be understood according to specific circumstances.

[0152] See also Figures 1 to 7 , Figure 1It is a flow chart of the method of the present invention, as shown in Figure 1 As shown, the present invention provides a multi-microphone array sound source localization method based on multi-level signal preprocessing and subspace spectrum optimization, comprising the following steps:

[0153] S1: Construct a simulation model for a multi-microphone array receiving high-frequency narrowband signals, which includes the following steps:

[0154] S11: Simulate the complex signal generated by the sound source to the microphone array: Assume N source The pitch angle θ of the sound source i and azimuth φ i , combined with N mic The position information of the microphone array composed of microphones can be used to obtain the phase delay τ of the microphone to each sound source direction. m,i And construct the steering vector a m,i :

[0155]

[0156] Where c is the speed of sound, x m and y m is the Cartesian coordinate of the mth microphone, λ is the wavelength, a m,j is the steering vector of the i-th sound source at the m-th microphone, and all a m,i Constructing a steering matrix

[0157] First, define the signal-to-noise ratio (SNR):

[0158]

[0159] Using the sound pressure level SPL and the reference sound pressure P ref =20μpa, the amplitude amp is obtained by inversely solving the relationship and generating the complex Gaussian signal Q s :

[0160] amp = 2 × 10 -5 10 SPL / 20

[0161] Composition Q s The row vector q s(i) for:

[0162]

[0163] N is a simplified expression of Gaussian distribution, amp i is the amplitude of the i-th sound source.

[0164] S12: According to the signal power P true The noise power is calculated by adding the signal-to-noise ratio (SNR) and generating the noise matrix through the complex Gaussian distribution:

[0165]

[0166] P true,i and amp i is the power and amplitude of the i-th sound source, SNR is the signal-to-noise ratio (in dB), is the noise matrix, and L is the number of sampling points of the microphone array.

[0167] Combined with the steering matrix A, the model of the array signal Y is obtained:

[0168] Y=AQ s +N

[0169] in is the sound source signal matrix.

[0170] S2: Preprocess the signal collected by the microphone through phase coherence weighting, PCA signal reconstruction, and fundamental frequency phase alignment calibration methods, specifically including the following steps:

[0171] S21: Original signal acquisition and preprocessing initialization:

[0172] Multi-channel signal acquisition. Using N array channels installed on the UAV platform, the acoustic signal of duration L is collected in the form of complex sampling to form a two-dimensional complex signal matrix:

[0173]

[0174] Where Y is the acquired two-dimensional complex signal matrix; y n (t) is the complex sound signal collected by the nth channel at time t; N is the number of channels; L is the sampling length of each channel (number of time points); C is the complex space;

[0175] Data normalization and outlier cleaning: Perform amplitude normalization on the acquired signal and remove outliers (such as NaN, Inf, etc.) to provide clean input for subsequent filtering steps.

[0176] S22: Channel weighted filtering based on phase coherence:

[0177] Extract the phase information of the complex signal. Extract the phase information of each channel in the signal at each moment and construct the phase matrix;

[0178] Φ∈R N×L ,Φ n,t =∠Y n,t

[0179] Where Φ is a two-dimensional phase matrix that stores the phase information of all channels at each moment; Φ n,t is the phase of the nth channel at the tth moment; Y n,tis an element in the complex signal matrix; ∠· is the operator for taking the complex phase;

[0180] Calculate the phase offset variance. For the phase differences of all channels, perform variance estimation based on a sliding window to measure signal cooperativity; first perform de-averaging on the phase values ​​of each column (at each moment) to eliminate the overall offset:

[0181]

[0182] Where, ΔΦ n,t is the phase difference after removing the mean value; Φ i,t is the phase value of the i-th channel at the t-th moment; is the average phase of all channels at that moment; Φ n,t is the original phase;

[0183] Then, at each time point, the variance of the phase shift ΔΦ is estimated using a sliding window to obtain the phase consistency measure

[0184] Construct phase coherence weighting coefficients. Using the phase variance obtained above, generate a corresponding weighting factor matrix W∈R at each time point N×L , used to adaptively suppress low-coherence channels:

[0185]

[0186] Where W n,t is the phase weighting factor of channel n at time t; α is the empirical adjustment factor (controlling the influence of variance on weight); ε is a very small positive number used to prevent zero division and numerical instability;

[0187] Perform weighted filtering on the original signal. Apply a dynamic weighting matrix to the original complex signal matrix and perform element-wise multiplication on the values ​​at each channel-time point to enhance coherence and suppress incoherent components.

[0188] S23: Signal reconstruction based on principal component analysis:

[0189] Construct the Hankel matrix for each channel. Before signal reconstruction, the filtered time series signal of each channel is first converted into a Hankel matrix to model the local dynamic pattern of the signal in the time domain;

[0190] Perform singular value decomposition (SVD) on the Hankel matrix. Decomposition yields principal components and singular value vectors:

[0191] H=U·S·V T

[0192] Where H is the Hankel matrix; U∈R m×mis the left singular vector matrix; V∈R n×n is the right singular vector matrix; S∈R m×n is a diagonal matrix whose diagonal elements s1≥s2≥…≥s r ≥0 is a singular value;

[0193] Determine the number of principal components to retain, reconstruct the signal, and restore it to a time series form. Reconstruct the matrix using principal components and restore it to a time series form by diagonal averaging.

[0194] S24: Channel calibration based on fundamental frequency phase alignment:

[0195] Perform a fast Fourier transform (FFT) on the reconstructed signal. Perform a one-dimensional fast Fourier transform on the reconstructed time domain signal of all channels and convert it to the frequency domain to extract the phase information of specific frequency components:

[0196] Y i (f) = F[y i (t)],i=1,2,...,N

[0197] Where y i (t) is the reconstructed signal of the i-th channel; F[·] is the Fourier transform operator; Y i (f) is the frequency domain response of the i-th channel at frequency f;

[0198] Select the first channel as the phase reference. Extract the fundamental frequency phase of the first channel as the reference. All subsequent channels will use this as the alignment target to eliminate relative phase errors.

[0199] Calibrate the phase deviation of each channel and uniformly apply the fundamental frequency alignment phase correction to the entire spectrum to ensure the alignment of the main frequency terms without affecting the signal frequency domain energy distribution;

[0200] Perform an inverse transform to obtain the corrected time domain signal. Perform an inverse fast Fourier transform (IFFT) on each channel to restore the calibrated frequency domain signal to a real-valued time domain waveform.

[0201]

[0202] Where, is the real-valued time domain signal output after calibration; Y i calib (f) is the complex signal after the frequency domain phase is aligned; F -1 [·] is the inverse Fourier transform operator; Re{·} is the real part operation (because it is still a complex number after the transformation).

[0203] At the same time, the results are amplitude normalized to prevent output signal overflow or imbalance.

[0204] S3: Construct a steering matrix based on the preprocessed signal, extract the coarse positioning results through signal subspace decomposition, and achieve high-precision angle optimization through an adaptive optimization algorithm, forming a collaborative optimization mechanism from coarse positioning to high-precision optimization.

[0205] Coarse positioning specifically includes the following steps:

[0206] S31: Calculate the covariance matrix based on the noisy time domain signal received by the microphone array.

[0207] R=(Y·Y′) / L

[0208] Where Y is the microphone array signal obtained after multi-stage signal preprocessing in step S2, L is the number of sampling points of the microphone array; Y′ is the conjugate transpose of Y; and R is the covariance matrix of the signal.

[0209] S32: Covariance matrix R = EΛE H Perform eigendecomposition, where Λ is the diagonal eigenvalue matrix, E=[E s ,E n ], E s 、E n They are signal subspace and noise subspace respectively. Use descending order of eigenvalues ​​and select the NN s The eigenvectors form the noise subspace matrix E n .

[0210] S33: Based on the steering matrix a(θ,φ) and the noise subspace E n , calculate the spatial spectrum power P(θ,φ):

[0211]

[0212] Where θ is the elevation angle, i.e. the angle between the received signal vector and the array plane; φ is the azimuth angle, i.e. the horizontal direction angle.

[0213] S34: Perform robust peak detection and screening on the spatial spectrum power P(θ,φ) to obtain a preliminary positioning result set (θ,φ).

[0214] Extract P by comparing neighborhood extreme values candidate Candidate peak power set:

[0215] P candidate ={(θ i ,φ i ,P i )|P i >P(θ±Δθ,φ±Δφ)}

[0216] Then adopt a dual screening mechanism combining adaptive threshold and spatial distance constraint:

[0217] Relative power threshold:

[0218] P(θ i ,φ i )>α·max(P),α∈[0.2,0.5]

[0219] Where, P(θ i ,φ i ) is the spatial spectrum power in the candidate peak direction; max(P) is the global maximum value of the spatial spectrum power; α is the adaptive threshold coefficient used to filter low-power false peaks, and the empirical value range is [0.2, 0.5].

[0220] Angle separation constraint:

[0221] Δθ ij =arccos(sinθ i sinθ j +cosθ i cosθ j cos(φ i -φ j ))>β,β=5°

[0222] Where Δθ ij are two candidate peak directions (θ i ,φ i ) and (θ j ,φ j ) and use the vector dot product in the spherical coordinate system to calculate the included angle to ensure that the sound source directions corresponding to adjacent peaks do not overlap. β = 5° is the minimum resolution interval and is an empirical value.

[0223] High-precision optimization specifically includes the following steps

[0224] S35: Dynamic search range calculation. Based on the preliminary positioning result set (θ, φ), calculate the upper and lower limits [lb, ub] of the search range of the current sound source:

[0225] δ=min(base_delta,0.5×min(angular_distance))

[0226] [θ min ,φ min ]=[θ c -δ,φ c -δ]

[0227] [θ max ,φ max ]=[θ c +δ,φ c +δ]

[0228] lb=[θ min ,φ min ]

[0229] ub=[θ max ,φ max ]

[0230] Where base_delta is the basic search radius, angular_distance is the minimum spherical angular distance between the current sound source and other sound sources, min is the minimum value operation, and the minimum value of base_delta and 0.5×angular_distance is taken as the dynamically adjusted search radius δ to ensure that the search areas do not overlap in multi-sound source scenarios; θ c 、φ c are the current sound source coordinates respectively.

[0231] S36: Chaos initialization population. Use Tent chaos map to generate the initial population, set the population size to NP, and define the chaotic sequence The initial value x1 is a random number, and the subsequent sequence recursion generation satisfies:

[0232]

[0233] Map the chaotic sequence to the search space and obtain the coordinates of the initial population pop:

[0234] pop=lb+x i ×(ub-lb)

[0235] S37: Dynamic parameter adjustment. Update the mutation factor F and crossover probability CR in each iteration:

[0236]

[0237] Where gen is the current iteration number; maxGen is the total iteration number; F0 is the basic variation factor; CR0 is the basic crossover probability; F adjust , CR adjust In order to dynamically adjust the coefficients and control the rate of change of the mutation factor F and the crossover probability CR during iteration, F decreases with iteration to shift to local development, and CR increases to enhance convergence.

[0238] S38: Mutation and boundary processing. For each individual, the DE / currenttobest / 1 strategy is used to generate the mutant individual:

[0239] mutant=x i +F×(x best -x i )+F×(x rl -x i )

[0240] Where x i is the current individual coordinate; x best is the current optimal individual coordinate; x r1 are the coordinates of a randomly selected individual.

[0241] Reflection method is used to correct individuals who cross the boundary:

[0242]

[0243] Where x(d) is the θ or φ coordinate of the out-of-bounds individual; lb(d) and ub(d) are the upper and lower limits of the search range θ or φ.

[0244] S39: Crossover and selection. Perform binomial crossover to generate trial individuals:

[0245]

[0246] Where trial(d) is the θ or φ coordinate of the trial individual. By comparing the random index rand(1) with the crossover probability CR, the parameters of the mutant individual or the original individual are selected for mixing.

[0247] Elite retention strategy: If the fitness of the experimental individual P(θ,φ) is better (that is, the power matching is higher), the experimental individual replaces the original individual. The individual fitness calculation formula is as follows:

[0248]

[0249] In the formula, the steering vector F(θ,φ) represents the phase delay characteristics of the sound source direction (θ,φ), F H is the conjugate transpose of the steering vector F(θ,φ); y(k) is the signal vector received by the microphone array at time k; L is the number of sampling points of the microphone array.

[0250] 310: Early stopping mechanism. If there is no improvement in the optimal fitness for 10 consecutive generations, the iteration is terminated early.

[0251] Experimental verification:

[0252] To verify the effectiveness of the present invention, 15 simulated sound source samples with mixed distributions of near-field and far-field were generated based on the simulation environment. The spatial coordinates of the sound sources covered the spatial range of the spherical coordinate system (θ: [0,90], φ: [0,360]), and the signal-to-noise ratio (SNR) range was set to -25dB to 0dB to simulate complex noise scenarios. After all samples collected signals through the designed multi-microphone array, multi-level signal preprocessing (phase coherent filtering → PCA reconstruction → phase calibration) was performed in sequence. Coarse positioning results were extracted based on signal subspace decomposition, and high-precision angle optimization was achieved through an adaptive optimization algorithm. The experiment compared the performance differences between the traditional MUSIC algorithm and the algorithm of the present invention under the same sound source distribution.

[0253] Figure 5 The MUSIC algorithm prediction diagram for 15 simulated sound sources is shown below. Figure 6 The prediction graph of the subspace spectrum optimization algorithm for 15 simulated sound sources shows that the subspace spectrum optimization algorithm can still accurately predict the sound source location in a multi-sound source environment. In comparison, the traditional MUSIC algorithm has obvious missed detections. Figure 7 The figure below compares the positioning error trends of the subspace spectrum optimization algorithm and the traditional MUSIC algorithm when the signal-to-noise ratio (SNR) is in the range of -25dB to 0dB. It is obvious that the positioning error of the subspace spectrum optimization algorithm is smaller than that of the MUSIC algorithm under the same SNR.

[0254] The above results show that through the collaborative design of signal preprocessing and subspace spectrum optimization algorithm, the present invention has high precision, strong noise resistance and real-time performance in complex sound source distribution and noise interference scenarios, providing reliable technical support for transmission line fault detection.

[0255] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not limiting. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present invention can be modified or replaced by equivalents without departing from the purpose and scope of the technical solutions, which should all be included in the scope of the claims of the present invention.

Claims

1. A multi-microphone array sound source localization method based on multi-stage signal preprocessing and subspace spectrum optimization, characterized by: The following steps are involved: S1: Establish a model for microphone array to collect sound source signals and generate time-domain discrete signals; S2: Multi-stage signal preprocessing including phase coherence filtering, principal component analysis (PCA) reconstruction and phase calibration; S3: Construct a steering matrix based on the preprocessed signal, extract the coarse positioning results through signal subspace decomposition, and achieve high-precision angle optimization through an adaptive optimization algorithm, forming a collaborative optimization mechanism from coarse positioning to high-precision optimization.

2. The method for sound source localization using a multi-microphone array based on multi-stage signal preprocessing and subspace spectrum optimization according to claim 1, wherein: In S1, a simulation model of a multi-element microphone array receiving a high-frequency narrowband signal is constructed, which specifically includes the following steps: S11: Simulate the complex signal generated by the sound source to the microphone array: Let N source The pitch angle θ of the sound source i and azimuth φ i , combined with N mic The position information of the microphone array composed of microphones is used to obtain the phase delay τ of the microphone to each sound source direction. m,i And construct the steering vector a m,i : Where c is the speed of sound, x m and y m is the Cartesian coordinate of the mth microphone, λ is the wavelength, a m,j is the steering vector of the i-th sound source at the m-th microphone, and all a m,i Constructing a steering matrix First, define the signal-to-noise ratio (SNR): Using the sound pressure level SPL and the reference sound pressure P ref =20μpa, the amplitude amp is obtained by inversely solving the relationship and generating the complex Gaussian signal Q s : amp=2×10 -5 ·10 SPL / 20 Composition Q s The row vector q s(i) for: N is a simplified expression of Gaussian distribution, amp i is the amplitude of the i-th sound source; S12: According to the signal power P true The noise power is calculated by adding the signal-to-noise ratio (SNR) and generating the noise matrix through the complex Gaussian distribution: P true,i and amp i is the power and amplitude of the i-th sound source, SNR is the signal-to-noise ratio, in dB, is the noise matrix, L is the number of sampling points of the microphone array; Combined with the steering matrix A, the model of the array signal Y is obtained: Y=AQ s +N in is the sound source signal matrix.

3. The multi-microphone array sound source localization method based on multi-stage signal preprocessing and subspace spectrum optimization according to claim 1, characterized in that: In the above S2, the signal collected by the microphone is preprocessed by phase coherence weighting, PCA signal reconstruction, and fundamental frequency phase alignment calibration method, which specifically includes the following steps: S21: Original signal acquisition and preprocessing initialization: Multi-channel signal acquisition: Using N array channels installed on the UAV platform, acoustic signals of duration L are collected in the form of complex sampling to form a two-dimensional complex signal matrix: Where Y is the acquired two-dimensional complex signal matrix; y n (t) is the complex sound signal collected by the nth channel at time t; N is the number of channels; L is the sampling length of each channel; C is the complex space; Data normalization and anomaly cleaning: normalize the amplitude of the collected signal and remove anomalies to provide clean input for subsequent filtering steps; S22: Channel weighted filtering based on phase coherence: Extract the phase information of the complex signal; extract the phase information of each channel in the signal at each moment and construct the phase matrix; Φ∈R N×L ,F n,t =∠Y n,t Where Φ is a two-dimensional phase matrix that stores the phase information of all channels at each moment; Φ n,t is the phase of the nth channel at the tth moment; Y n,t is an element in the complex signal matrix; ∠· is the operator for taking the complex phase; Calculate the phase offset variance; estimate the variance of the phase differences across all channels based on a sliding window to measure signal cooperativity; first remove the mean of each column, i.e., the phase value at each moment, to eliminate the overall offset: Where, ΔΦ n,t is the phase difference after removing the mean value; Φ i,t is the phase value of the i-th channel at the t-th moment; is the average phase of all channels at that moment; Φ n,t is the original phase; Then, at each time point, the phase offset ΔΦ is estimated using a sliding window to obtain the phase consistency measure and construct the phase coherence weighting coefficient; using the phase variance obtained above, a corresponding weighting factor matrix W∈R is generated at each time point. N×L , used to adaptively suppress low-coherence channels: Where W n,t is the phase weighting factor of channel n at time t; α is the empirical adjustment factor, which indicates the influence of the control variance on the weight; ε is a very small positive number used to prevent zero division and numerical instability; α∈[0.1,1.0],ε=10 -6 ; Perform weighted filtering on the original signal; apply a dynamic weighting matrix to the original complex signal matrix, and perform element-by-element multiplication on the values ​​at each channel-time point to achieve coherence enhancement and incoherent component suppression; S23: Signal reconstruction based on principal component analysis: Construct the Hankel matrix of each channel; before signal reconstruction, first convert the filtered time series signal of each channel into a Hankel matrix to model the local dynamic pattern of the signal in the time domain; Perform singular value decomposition (SVD) on the Hankel matrix; decomposition yields principal components and singular value vectors: H=U·S·V T Where H is the Hankel matrix; U∈R m×m is the left singular vector matrix; V∈R n×n is the right singular vector matrix; S∈R m×n is a diagonal matrix whose diagonal elements s1≥s2≥…≥s r ≥0 is a singular value; Determine the number of principal components to be retained, reconstruct the signal and restore it to a time series form; reconstruct the matrix through the principal components and restore it to a time series through diagonal averaging; S24: Channel calibration based on fundamental frequency phase alignment: Perform a fast Fourier transform (FFT) on the reconstructed signal; perform a one-dimensional fast Fourier transform on the reconstructed time domain signal of all channels and convert it to the frequency domain to extract the phase information of a specific frequency component: Y i (f)=F[y i (t)],i=1,2,...,N Where y i (t) is the reconstructed signal of the i-th channel; F[·] is the Fourier transform operator; Y i (f) is the frequency domain response of the i-th channel at frequency f; Select the first channel as the phase reference; extract its fundamental frequency phase using the first channel as the reference, and all subsequent channels will use this as the alignment target to eliminate relative phase errors; Calibrate the phase deviation of each channel and uniformly apply the fundamental frequency alignment phase correction to the entire spectrum to ensure the alignment of the main frequency terms without affecting the signal frequency domain energy distribution; Perform inverse transform to obtain the corrected time domain signal; perform inverse fast Fourier transform IFFT on each channel to restore the calibrated frequency domain signal to a real-valued time domain waveform Where, is the real-valued time domain signal output after calibration; Y i calib (f) is the complex signal after the frequency domain phase is aligned; F -1 [·] is the inverse Fourier transform operator; Re{·} is the real part operation; At the same time, the results are amplitude normalized to prevent output signal overflow or imbalance.

4. The method for sound source localization using a multi-microphone array based on multi-stage signal preprocessing and subspace spectrum optimization according to claim 1, wherein: In S3, a steering matrix is ​​constructed based on the preprocessed signal, a coarse positioning result is extracted by signal subspace decomposition, and high-precision angle optimization is achieved through an adaptive optimization algorithm, forming a collaborative optimization mechanism from coarse positioning to high-precision optimization; Coarse positioning specifically includes the following steps: S31: Calculate the covariance matrix based on the noisy time domain signal received by the microphone array. R=(Y·Y′) / L Where Y is the microphone array signal obtained after multi-stage signal preprocessing in step S2, L is the number of sampling points of the microphone array; Y′ is the conjugate transpose of Y; R is the covariance matrix of the signal; S32: Covariance matrix R = EΛE H Perform eigendecomposition, where Λ is the diagonal eigenvalue matrix, E=[E s ,E n ], E s 、E n They are signal subspace and noise subspace respectively; use descending order of eigenvalues ​​and select NN s The eigenvectors form the noise subspace matrix E n ; S33: Based on the steering matrix a(θ,φ) and the noise subspace E n , calculate the spatial spectrum power P(θ,φ): Where θ is the elevation angle, i.e. the angle between the received signal vector and the array plane; φ is the azimuth angle, i.e. the horizontal direction angle; S34: Perform robust peak detection and screening on the spatial spectrum power P(θ,φ) to obtain a preliminary positioning result set (θ,φ); Extract P by comparing neighborhood extreme values candidate Candidate peak power set: P candidate ={(θ i ,f i ,P i )∣P i >P(θ±Δθ,φ±Δφ)} Then adopt a dual screening mechanism combining adaptive threshold and spatial distance constraint: Relative power threshold: P(θ i ,f i )>α·max(P),α∈[0.2,0.5] Where, P(θ i ,φ i ) is the spatial spectrum power in the candidate peak direction; max(P) is the global maximum value of the spatial spectrum power; α is the adaptive threshold coefficient used to filter low-power false peaks, and the empirical value range is [0.2, 0.5]; Angle separation constraint: Dth ij =arccos(sinθ) i sinth j +cosθ i cosθ j cos(φ i -f j ))>β,β∈[3°,10°] Where Δθ ij are two candidate peak directions (θ i ,φ i ) and (θ j ,φ j ) between them, the angle is calculated using the vector dot product in the spherical coordinate system to ensure that the sound source directions corresponding to adjacent peaks do not overlap, and β is the minimum resolution interval; High-precision optimization specifically includes the following steps S35: Dynamic search range calculation; based on the preliminary positioning result set (θ, φ), calculate the upper and lower limits [lb, ub] of the search range of the current sound source: δ=min(base_delta,0.5×min(angular_distance)) [i min ,f min ]=[θ c -d,f c -d] [i max ,f max ]=[θ c +d,f c +d] lb=[θ min ,φ min ] ub=[θ max ,f max ] Where base_delta is the basic search radius, angular_distance is the minimum spherical angular distance between the current sound source and other sound sources, min is the minimum value operation, and the minimum value of base_delta and 0.5×angular_distance is taken as the dynamically adjusted search radius δ to ensure that the search areas do not overlap in multi-sound source scenarios; θ c 、φ c are the current sound source coordinates respectively; S36: Chaos initialization population; use Tent chaos map to generate the initial population, set the population size to NP, and define the chaotic sequence The initial value x1 is a random number, and the subsequent sequence recursion generation satisfies: Map the chaotic sequence to the search space and obtain the coordinates of the initial population pop: pop=lb+x i ×(ub-lb) S37: Dynamic parameter adjustment; update the mutation factor F and crossover probability CR in each iteration: Where gen is the current iteration number; maxGen is the total iteration number; F0 is the basic variation factor; CR0 is the basic crossover probability; F adjust , CR adjust To dynamically adjust the coefficients, the rate of change of the mutation factor F and the crossover probability CR during iteration is controlled. F decreases with iteration to shift to local development, and CR increases to enhance convergence. S38: Mutation and boundary processing; for each individual, the DE / currenttobest / 1 strategy is used to generate mutant individuals: mutant=x i +F×(x best -x i )+F×(x rl -x i ) Where x i is the current individual coordinate; x best is the current optimal individual coordinate; x r1 are the coordinates of randomly selected individuals; Reflection method is used to correct individuals who cross the boundary: Where x(d) is the θ or φ coordinate of the out-of-bounds individual; lb(d) and ub(d) are the upper and lower limits of the search range θ or φ; S39: Crossover and selection; perform binomial crossover to generate trial individuals: In the formula, trial(d) is the θ or φ coordinate of the test individual. By comparing the random index rand(1) with the crossover probability CR, the parameters of the mutant individual or the original individual are selected for mixing. Elite retention strategy: If the fitness of the experimental individual P(θ,φ) is better (that is, the power matching is higher), the experimental individual replaces the original individual. The individual fitness calculation formula is as follows: In the formula, the steering vector F(θ,φ) represents the phase delay characteristics of the sound source direction (θ,φ), F H is the conjugate transpose of the steering vector F(θ,φ); y(k) is the signal vector received by the microphone array at time k; L is the number of sampling points of the microphone array; 310: Early stopping mechanism; if there is no improvement in the optimal fitness for 10 consecutive generations, the iteration is terminated early.

5. The method for sound source localization using a multi-microphone array based on multi-stage signal preprocessing and subspace spectrum optimization according to claim 1, wherein: In S2, phase coherent weighted filtering is characterized by: The phase variance is calculated using a sliding window with a window length of T. w Meet T w =0.1L, where L is the total number of signal sampling points; The weighting factor is calculated as Where α = 0.5, ε = 10 -6 .

6. The method for sound source localization using a multi-microphone array based on multi-stage signal preprocessing and subspace spectrum optimization according to claim 1, characterized in that: In S3, the adaptive optimization algorithm is an improved differential evolution algorithm, which is characterized by: The population is initialized using the Tent chaotic map, and the chaotic sequence satisfies: Among them, x i is the i-th value of the chaotic sequence, and the initial value x1 is a random number in the interval [0,1]; The dynamic adjustment formula of the mutation factor F and the crossover probability CR is: Among them, F0=0.8, CR0=0.6, gen is the current number of iterations, and maxGen is the total number of iterations.

7. The method for sound source localization using a multi-microphone array based on multi-stage signal preprocessing and subspace spectrum optimization according to claim 1, wherein: In S3, dynamic search range calculation is characterized by: The basic search radius is Δ base =10°; The final search radius δ is calculated as follows: where Δθ ij is the current sound source direction (θc,φ c ) and other sound source directions (θ j ,φ j ), the minimum spherical angle distance between the two is calculated as follows: Dth ij =arccos(sinθ) c sinth j +cosθ c cosθ j cos(φ c -f j ))。 8. The method for sound source localization using a multi-microphone array based on multi-stage signal preprocessing and subspace spectrum optimization according to claim 1, characterized in that: The method further comprises post-processing S4: Perform spatial clustering on the optimized sound source angles (θ, φ) to eliminate isolated points; The positioning result is corrected by energy-weighted average, and the formula is: Among them, P i is the spatial spectrum power of the i-th candidate sound source.

Images