Joint source localization and separation method for acoustic sources

Abstract

A method is provided for acoustic source direction of arrival estimation and acoustic source separation, via spatial weighting of the dictionary based display of the steered response function calculated for a certain number of directions from spherical harmonic decomposition coefficients obtained from microphone array recordings of the sound field. The usage of spatial band limited functions of plane waves to represent more complex directional maps of the sound field constitutes the algorithm. These functions are calculated for pre-defined directions on an analysis surface (such as a sphere). The directions of arrival of sound sources are calculated with the same method in order to group source estimates to localize sound sources. Thereby, directions of arrival can be obtained from the recordings of the sound sources captured by means of a microphone array and following this, sound sources can be separated by using this direction information or predetermined source arrival directions.

Claims

1. A method run by a computer for an estimation of an arrival direction from one or more sound source mixtures, and a separation of sound sources, comprising the following processing steps; obtaining a spherical harmonic decomposition of one or more digital sound signal data from a plurality of microphones or sensors and/or a sound field as an input from an interface, in the case that the input is sound data from the plurality of microphones or sensors, carrying out a spherical harmonic decomposition of the sound data and providing a time-frequency representation of spherical harmonic decomposition coefficients, in the case that the input is the spherical harmonic decomposition coefficients, providing the time-frequency representation of the spherical harmonic decomposition coefficients, forming spatial filters, wherein the spatial filters have a predetermined selectivity at predetermined directions, forming beams from the spherical harmonic decomposition coefficients by using the spatial filters, obtaining a directional map of an amplitude of the sound field by steering the beams in a given number of directions, showing the directional map as a combination of a limited number of directional elements by using a redundant series of non-orthogonal template vectors and/or matrices, obtaining a usage frequency distribution of the non-orthogonal template vectors and/or matrices used in the time-frequency representation, calculating sound arrival directions from the usage frequency distribution, weighting the time-frequency representation depending on the sound arrival directions to obtain weighted representations, obtaining time-frequency transforms from the weighted representations, and determining and obtaining separated sound sources by carrying out inverse time frequency transforms.

2. The method according to claim 1, wherein values used for weighting are exemplified from a directional function having a single global maximum.

3. The method according to claim 2, wherein values used for weighting are adapted according to the sound arrival directions.

4. The method according to claim 2, wherein template series and/or matrices are formed of band limited functions.

5. The method according to claim 2, wherein template series and/or matrices are exemplified from direction localized functions.

6. The method according to claim 2, wherein template series and/or matrices are exemplified from real valued functions.

7. The method according to claim 1, wherein values used for weighting are adapted according to the sound arrival directions.

8. The method according to claim 1, wherein template series and/or matrices are formed of band limited functions.

9. The method according to claim 1, wherein template series and/or matrices are exemplified from direction localized functions.

10. The method according to claim 1, wherein template series and/or matrices are exemplified from real valued functions.

11. A method run by a computer for a separation of sound sources from a mixture of two or more sound sources, comprising the following processing steps; obtaining a spherical harmonic decomposition and sound arrival directions, of one or more digital sound signal data from a plurality of microphones or sensors and/or a sound field as an input from an interface, in the case that the input is sound data from the plurality of microphones or sensors, carrying out a spherical harmonic decomposition of the sound data and providing a time-frequency representation of spherical harmonic decomposition coefficients, in the case that the input is the spherical harmonic decomposition coefficients, providing the time-frequency representation of the spherical harmonic decomposition coefficients, forming spatial filters, wherein the spatial filters have a predetermined selectivity at predetermined directions, forming beams from the spherical harmonic decomposition coefficients by using the spatial filters, obtaining a directional map of an amplitude of the sound field by steering the beams in a given number of directions, showing the directional map as a combination of a limited number of directional elements by using a redundant series of non-orthogonal template vectors and/or matrices, weighting the time-frequency representation using a function to obtain weighted representations, wherein the function depends on a direction, obtaining time-frequency transforms from the weighted representations, and determining and obtaining separated sound sources by carrying out inverse time frequency transforms.

12. The method according to claim 11, wherein values used for weighting are exemplified from a directional function having a single global maximum.

13. The method according to claim 11, wherein values used for weighting are adapted according to the sound arrival directions.

14. The method according to claim 11, wherein template series and/or matrices are formed of band limited functions.

15. The method according to claim 11, wherein template series and/or matrices are exemplified from direction localized functions.

16. The method according to claim 11, wherein template series and/or matrices are exemplified from real valued functions.

17. A method run by a computer for an estimation of arrival directions of one or more sound sources, comprising the following processing steps; obtaining a spherical harmonic decomposition of one or more digital sound signal data from a plurality of microphones or sensors and/or a sound field as an input from an interface, in the case that the input is sound data from a plurality of microphones or sensors, carrying out a spherical harmonic decomposition of the sound data and providing a time-frequency representation of spherical harmonic decomposition coefficients, in the case that the input is the spherical harmonic decomposition coefficients, providing the time-frequency representation of the spherical harmonic decomposition coefficients, forming spatial filters, wherein the spatial filters have a predetermined selectivity at predetermined directions, forming beams from the spherical harmonic decomposition coefficients by using the spatial filters, obtaining a directional map of an amplitude of the sound field by steering the beams in a given number of directions showing the directional map as a combination of a limited number of directional elements by using a redundant series of non-orthogonal template vectors and/or matrices, obtaining a usage frequency distribution of the non-orthogonal template vectors and/or matrices used in the time-frequency representation, calculating sound arrival directions from the usage frequency distribution.

18. The method according to claim 17, wherein template series and/or matrices are formed of band limited functions.

19. The method according to claim 17, wherein template series and/or matrices are exemplified from direction localized functions.

20. The method according to claim 17, wherein template series and/or matrices are exemplified from real valued functions.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1 is a flow diagram of the localization and separation of sound sources.

(2) FIG. 2 is the flow diagram of the separation method.

(3) FIG. 3 is the flow diagram of the localization method.

(4) FIG. 4 shows the directional map obtained using steered response function that can be obtained from a single time-frequency bin.

(5) FIGS. 5A-5C show some dictionary elements that can be used in expressing the response function.

(6) FIG. 6 shows the neighborhood relations (related to the clustering method for different atoms) of the peaks in the histogram.

(7) FIG. 7 graphically shows the directional response obtained for different κ values of the Von Mises function and the directional response of maximum directivity (max DF) beamforming.

(8) The figures need not be scaled and details that are not critical for a clear understanding of the present invention may have been omitted. Apart from this, elements that are at least substantially identical or those that at least substantially have the same functions, have been shown with the same reference number.

DETAILED DESCRIPTION OF THE INVENTION

(9) In this detailed description, the preferred embodiments of the invention are described such that they do not have any limiting effect but have been provided to further describe the subject matter.

(10) The invention comprises two different algorithms for the localization and the separation of sound sources. These algorithms can be used together or independently from each other. The block diagram showing the flow of the disclosed invention is shown in FIG. 1.

(11) FIG. 2 shows the block diagram of the source separation method. The inputs are sound source positions and microphone array recordings and the outputs are the separated sound files. The details of the different steps of the algorithms are given below. A. Calculation of spherical harmonic decomposition coefficients: Harmonic series can be calculated using microphone array recordings and the positions of microphones that such arrays comprise. Harmonic series are used to define the sound field around the microphone array using spherically or cylindrically periodic functions. The disclosed method can also directly use the spherical harmonic decomposition of the sound field. In the case that such an input is present, this step does not need to be carried out. B. Time-frequency transform: Each of the spherical harmonic coefficient series that are to be processed, is expressed with a suitable invertible representation in time frequency domain. The procedures in further steps are carried out separately for each time-frequency bin. As the procedures in step A are linear, they can also be carried out in reverse order. C. Beamforming: The signals to be used in the next step are calculated for each time-frequency bin by means of steering a maximum directivity factor beam in a limited number of directions that are radially outward from the origin at which the spherical harmonic coefficients are obtained. This is achieved by weighting the spherical harmonic decomposition coefficients appropriately. The parameter that the algorithm uses is the number of directions at which the beam would be steered. D. Creation of the dictionary atoms at the determined directions: For a plane wave, the directional response of the beam with the maximum directivity can theoretically be described as a closed form function, as described below. In this step, the atoms to be used in the expression of the steered beamforming function are obtained by sampling this function on a sphere (or another analysis surface) at a finite number of directions. This process can not only be carried out offline in order to accelerate the method, but it can also be applied separately for each time-frequency bin at runtime based on the sound source directions obtained as a result of earlier analysis. E. Representation: This step involves the calculation of the representation of said beamforming results in an economical way according to certain criteria using the lowest number of atoms. The dictionary atoms mentioned above are used in this step. The result of this step is the calculation of complex or real valued coefficients for each of these atoms in the analyzed time-frequency bin by expressing the sound field as a linear sum of the previously calculated atoms in the specified directions. F. Directional weighting: The dictionary atoms determined in step D are spatially filtered using the predetermined sound source directions. For this process, the coefficient that is calculated for each atom whose direction is known, is multiplied with a directional gain that emphasizes the direction that is to be separated. Here, it is possible to use a weighting function defined in closed form in order to calculate this directional gain. It is also possible to carry out directional weighting adaptively. A directionally weighted beamform can be obtained using the weighted coefficients and corresponding atoms for each time-frequency bin. G. Reconstruction: Separated sound sources are reconstructed in the time domain, by inverting the new time-frequency representations that are obtained in the previous step.

(12) FIG. 3 shows the block diagram of the positioning method. The above mentioned A, B, C, D, E steps are common to the two algorithms and the below mentioned additional steps are used only for source direction estimation. H. Formation of a directional histogram based on selected atoms: The statistical distribution of atoms used to express the steered beamform at a certain time range is formed with a histogram or another method. If a histogram is used, the number of bins shall be selected to be the same with the number of atoms in the dictionary. I. Clustering: The peak points of the distribution obtained as a result of the previous step are calculated. Direction of arrival can be estimated by using the neighborhood relations between the atoms that these peaks correspond to.

(13) The definitions that were generally expressed above, have been used as a solution embodiment with the below mentioned preferred parameters. The spherical harmonic decomposition of the sound field is obtained from recordings made with a Rigid Spherical Microphone Array. Short time Fourier transform is used as the time-frequency transform. The Legendre impulse functions whose details are given below are sampled on the sphere to generate dictionary atoms. Orthogonal Matching Pursuit algorithm is used in the representation stage and maximum directivity factor beamforming is used for calculating steered beams. Von Mises function that is defined on the sphere is used for position dependent weighting. The distribution for direction of arrival estimation is obtained by using a histogram. In the preferred embodiment, the order of time-frequency transform and spherical harmonic decomposition has been swapped which leads to equivalent results due to the linearity of the concerned operations.

(14) Short-Time Fourier Transform: Each of the signals obtained from the microphone array is transformed into the time-frequency domain by means of a short time Fourier transform. Although any kind of window function and length can be used for this process, in the preferred embodiment a 2048 sample Hann window has been used with 50% overlap.

(15) The Calculation of Spherical Harmonic Decomposition: In this step the spherical harmonic decomposition for each time-frequency bin is calculated as follows:

(16) $p_{\mathrm{nm}}\left(k\right)=\underset{i=1}{\overset{M}{.Math.}}{\gamma }_i{p\left({\Omega }_i\right)\left[Y_n^m\left({\Omega }_i\right)\right]}^{*}$

(17) Here the M is the number of microphones, γ.sub.i is the related quadrature spherical weights, the k is the time-frequency bin index that has been obtained by using short time Fourier transform, Ω.sub.i=(θ.sub.i, ϕ.sub.i) is the position of the microphone on the spherical surface. Spherical harmonic function, Y.sub.n.sup.m is defined as follows:

(18) $Y_n^m\left(\Omega \right)=\sqrt{\frac{2n+1}{4\pi }\frac{\left(n-m\right)!}{\left(n+m\right)!}}{P}_n^m\left(\cos \theta \right)e^{\mathrm{im}\varphi }$

(19) Maximum directivity beamforming: This process is also known as the plane wave decomposition. It can be calculated as follows using spherical harmonic coefficients:

(20) $y\left(\Omega ,k\right)=\underset{n=0}{\overset{N}{.Math.}}\underset{m=-n}{\overset{n}{.Math.}}\frac{p_{\mathrm{nm}}\left(k\right)}{4\pi i^n{b}_n\left({\mathrm{kr}}_a\right)}{Y}_n^m\left(\Omega \right)$

(21) Wherein Ω=(θ, ϕ) is the steering direction of the maximum directivity factor beam, j.sub.n(.), h.sub.n.sup.(2)(.), j.sub.n′(.), and h.sub.n.sup.(2)′(.) are the spherical Bessel and Hankel functions, and the first-order derivatives thereof, r.sub.a is the radius of the spherical microphone, and frequency equalization function is given as:

(22) $b_n\left(\mathrm{kr}\right)=j_n\left(\mathrm{kr}\right)-\frac{j_n^{\prime }\left({\mathrm{kr}}_a\right)}{h_n^{{\left(2\right)}^{\prime }}\left({\mathrm{kr}}_a\right)}{h}_n^{\left(2\right)}\left(\mathrm{kr}\right)$

(23) Plane Wave Legendre Impulse Function Definitions at the Determined Directions: Maximum directivity factor beamform for a limited number of S plane wave is defined as given below:

(24) $y_S\left({\Omega }_q,k\right)=\underset{s=1}{\overset{S}{.Math.}}{a}_s\left(k\right)\phi \left({\Omega }_q.Math.{\Omega }_s\right)$

(25) Wherein

(26) $\phi \left(\Omega .Math.{\Omega }_s\right)=\frac{N+1}{4\pi }\left[\frac{P_{N+1}\left(\cos {\Theta }_s\right)-P_N\left(\cos {\Theta }_s\right)}{P_1\left(\cos {\Theta }_s\right)-P_0\left(\cos {\Theta }_s\right)}\right],$

(27) is the Legendre impulse with a maximum at Ω.sub.s=(θ.sub.s, ϕ.sub.s). This function is sampled at a finite number of points on the sphere to obtain the atoms in the dictionary used in Orthogonal Matching Pursuit algorithm in the following step.

(28) Orthogonal Matching Pursuit: Orthogonal matching pursuit is an iterative method used to express steered response function in a given time-frequency bin using a small number of dictionary atoms.

(29) As such, the steered response function at the given time-frequency bin can be expressed using a suitable selection of dictionary elements. The algorithm flow is as follows: 1. Maximum directivity factor beam is steered to calculate the steered response function at different directions covering the entire sphere for the analyzed time-frequency bin resulting in a directional map of the sound field for the given time-frequency bin. 2. The vector formed of these values is multiplied with the matrix comprising dictionary atoms and the atom corresponding to the highest value in the resulting vector is selected. 3. The approximation obtained using this atom is subtracted from the vector and a residual vector is formed. 4. The residual vector is multiplied with the matrix comprising dictionary atoms and the atom corresponding to the highest value in the resulting vector is selected. 5. The third and the fourth steps are repeated until the norm of the residual vector falls below a predetermined threshold value. 6. The coefficients of the approximation comprising a linear combination of atoms are obtained by using the Least Squares algorithm.

(30) For example the steered response function in FIG. 4, can be obtained by using only the 1st and 2nd atoms of the dictionary atoms given in FIGS. 5A-5C. The third atom is not used.

(31) Forming a Directional Histogram: The histogram calculated after finding the atoms that adequately express the steered response function by means of the orthogonal pursuit algorithm, shows how frequently these atoms are used in a given period of time.

(32) Histogram Clustering and Source Localization: Source localization is based on a clustering principle based on the neighborhood relations of the directions of local maxima points in the histogram. The neighborhood relations of the positions is side information, and the directions where the sources are located are calculated by averaging the directions that the clustered positions are facing. The outputs of this stage are the components and the directions of the sound sources in the environment. The neighborhood relations of the peaks in the histogram is shown in FIG. 6. Accordingly Group 1 is comprised of P7, P13; Group 2 is comprised of P6, P21 and P22.

(33) Directional Weighting: The source directions that have been calculated and the linear weights corresponding to these directions are used at this stage. In the preferred embodiment of the invention, the linear weights corresponding to each atom is weighted by using Von Mises Functions with a mean in the direction of the desired sound source evaluated at the center direction of that atom. The spatial filter obtained by means of weighting by the Von Mises function is shown in FIG. 7, for different density parameters (κ). The maximum directivity factor beam is also shown for comparison. The κ value determines the spatial selectivity of the Von Mises function. When this value is small, it causes the method to filter its input at a wider directional range and increasing this value results in a sharper beam with higher selectivity resulting in more accurate separation of sources. In this step, a complex value is obtained for each of the sound sources that are to be separated at each time-frequency bin.

(34) Inverse Short-Time Fourier Transform: The new time-frequency representations obtained for each of the each sound sources are transformed back into the time domain using the inverse short-time Fourier transform to obtain the separated source signals.