Free-Field TDOA-AOA Sound Source Localization Using Three Soundfield Microphones

Accurate source localization is an important problem in many research areas as well as practical applications in wireless communications and acoustic signal processing. This paper presents a passive three-dimensional sound source localization (SSL) method that employs a geometric conﬁguration of three soundﬁeld microphones. Two methods for estimating the angle of arrival (AOA) and time difference of arrival (TDOA) are proposed based on Ambisonics A and B format signals. The closed-form solution for sound source location estimation based on two TDOAs and three AOAs is derived. The proposed method is evaluated by simulations and physical experiments in our anechoic chamber. Simulations demonstrate that the estimation method can theoretically obtain Cramér-Rao lower bound for a small Gaussian noise present in AOA and TDOA observations. Investigation on the uncertainty of TDOA and AOA measurements depending on the length of measurement interval is also conducted. Experimental results in terms of RMSE indicate that the proposed solution can be used to accurately ﬁnd a 3D position of the sound source in free-ﬁeld environment. Performance evaluation regarding the number of estimation steps shows that higher accuracy can be achieved by longer observations of stationary sound source.


I. INTRODUCTION
The sound source localization (SSL) is an essential step in a wide range of audio/acoustic-based applications. Nowadaysconcerned research topics on SSL are ranging from detection of the speaker position in human-computer interaction [1] or smart video conferencing [2], robot movement in an unknown environment [3], [4], search and rescure [5] to advance military applications such as localization of a sniper [6] and medium-range aircraft localization [7]. Besides, SSL is usually used as a necessary preprocessing step before the enhancement of an acoustic signal from a particular location [8].
In general, the task of source localization can be viewed either as an active localization scenario, where transmitters actively emit signals to illuminate the target of interest while the target location is inferred by collecting reflections, wherein the scenario of passive localization, receiving sensors can only collect signals emitted from the source. In many The associate editor coordinating the review of this manuscript and approving it for publication was Manuel Rosa-Zurera. areas, only a passive source location is considered, where the signals usually do not carry the time information about their transmission. The sound localization as a type of passive localization refers to the problem of estimating the position from which a sound signal originates concerning the microphone array geometry. In this case, a localization system is unable to directly measure the time of arrival (TOA) between the source location and receiving sensor, but instead, only the difference between times when the different sensors receive the signals can be measured.
Various methods for SSL have been proposed where all methods can be grouped by their efforts to detect sound source either in 2D space [9], [10] or in 3D space [3], [11]. Fundamentally, there are two main approaches to finding a source with respect to recorded audio signals. Both approaches are mainly based on estimating the time difference of arrival (TDOA) obtained by using various configurations of microphone arrays, such as linear array [12], circular array [13], or distributed array [14] and different cross-correlation algorithm to estimate time lag between microphones. The first approach aims to maximize the VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ steered response power (SRP) of the output of a delay and sum beamformer [15]. This direct approach performs an exhaustive search in the whole SRP space to find a sound location, which is found to be computationally expensive.
In contrast, the indirect approach is based on using estimated TDOAs measurements, where the sound source is ascertained by addressing criteria such as the hypercone fitting problem [16]. Although results obtained by the indirect method are more prone to error as they are more sensitive to background noise, and reflections oppose to direct methods, the main advantage of indirect TDOA based approach is that it can be effectively used in distributed microphone network since it only requires for TDOA values to be transmitted and not raw sound signal data [17].

A. RELATED WORK
In many research papers, TDOA is a standard measurement used for passive localization [18], [19]. Source localization based solely on TDOA measurements demanded a specific number of measuring devices, i.e., a minimum number of three and four sensors are required to locate an unknown target respectively in 2D and 3D space. The principle of TDOA based approach is to estimate source position from the intersection of hyperbolic arcs and surfaces respectively for 2D and 3D cases. In the case of near-field applications where the source range to the sensor array baseline ratio is not large, the resulting intersection can be obtained by solving a set of nonlinear equations. On the contrary, in the far-field applications, resulting intersection produces a low location estimate since the hyperbolic arcs/surfaces become almost parallel to one another. Given the evidence that the accuracy of the source position estimate degrades when the source moves sufficiently far away from the sensors resulted in two different localization problems. Recently, unified near/fatfield TDOA based localization was proposed in [20]. The proposed approach consisted of two formulations for the unified localization problem in two different coordinate systems. The first formulation is of a nonlinear non-convex weighted least squares optimization based on the modified polar representation for the source position, and the other is the non-convex fractional programming formulation using the conventional Cartesian coordinates of the source position as the optimization variable. Besides the TDOA-based approach, a source location can also be calculated from the AOA measurements and its derivative. Obtaining AOA often involves a sensor that is equipped with an array of receivers; thus, it elevates a requirement of synchronization between different sensors since each produces angle by itself. In [21], the authors proposed a solution that can attain Cramér-Rao lower bound under mild conditions. Another 3D bearingonly localization is proposed by [22], where authors achieved a significant reduction in bias and root-mean-square error using a pseudo-linear estimator. In general, the passive source localization problem is not trivial since direct relationships between the position of a source and the measurements are complex, and procedures for solving equations for the TDOA and AOA methods are hard because of nonlinearity. One of the more recent research directions is in combining TDOA and AOA measurements, where algebraic manipulations allow transforming the relationships to the linear form and lowering energy consumption [23]. There are several advantages with the hybrid TDOA-AOA approach, such as improved localization performance [24], [25], reduced number of sensors required [26] and it can minimize the occurrence of ghost targets which is typical for localization approach with individual TDOA measurements [27]. Many studies have been carried out proposing a different solution to source localization [26], [28]- [31]. The most straightforward method to source localization is an exhaustive search in a feasible solution region, which is a time-consuming and inefficient solution for real-time application. In general, the maximum likelihood (ML) estimator is introduced to estimate the location since it is asymptotically efficient.
However, the aftermentioned approaches are computationally expensive, and it is hard to find a closed-form solution, or a closed-form solution does not exist at all. One of the solutions is to linearize equations with a recursive approach such as Taylor-series. Nevertheless, these numerical search techniques can converge to an optimal solution only if the ML function is convex. The numerical methods are prone to error since they depend on the right initial position guess, and thus it is difficult to guarantee its global convergence and calculating time. Their iterative nature does not make them very suitable for real-time applications. To improve the robustness and reduce the complexity, a closedform solution is required. A linear least-squares estimator with the closed-form solution, called pseudo linear estimator (PLE), was proposed in [22]. Although this approach is less computationally demanding, the estimated source position is biased because of the correlation between system and measurement. Recently, authors [32] proposed a new method for localization. The proposed method represents a simple algebraic solution that does not suffer from the local convergence problem. However, this method also has a larger bias since localization accuracy is affected by the deviation of the noise correlation matrix. Another fine localization solution, based on generalized trust-region subproblem technique is proposed by [33], were authors analyzed the necessary optimal conditions of squared range difference least square cost function. To reduce the bias of the estimator, a few different methods where proposed. In [34] BR-PLE method is proposed for reducing bias. Authors in [35] proposed a solution to find a rotation angle that could reduce the bias of the estimator since they demonstrated that estimator performance was sensitive to origin rotation. [36] proposed another solution for reducing estimator bias in the presence of sensors error. The proposed method introduces a quadratic constraint so that estimator expectation cost function can attain the minimum value at a true position and thus to achieve Cramer-Rao lower bound. By analyzing different approaches, one could conclude that for the most methods, an estimation bias arises from the least-squares techniques. One of the solutions is to construct a new pseudo linear system such as [27], where a new weight least-square method is proposed. Authors claimed that the novel structured total least squares method could reduce estimation bias when the target is outside the convex hull formed by measurement sensors. The total least squares estimator was also considered in [37], where it has demonstrated improved accuracy over least-square solutions. Another practical solution to position estimation was proposed by [38]. The proposed closed-form method was based on converting time-delay measurements to angular information, but it couldn't achieve the CRLB performance. Recently, an efficient closed-form solution has been proposed for passive source location using only two stations [26]. In their work, a new relationship between hybrid TDOA and AOA measurements and unknown source positions was constructed. The theoretical simulation shows that the proposed solution can achieve CRLB for Gaussian noise, where bias compared to variance can be ignored.
By analyzing the aforementioned different source localization approaches, some general remarks regarding localization performance can be made. Investigating different types of input data for source estimation, one can observe that using hybrid TDOA and AOA measurements reduces the number of sensors required, which is especially important for a reallife scenario where the basic assumption about line-of-sight (LOS) between unknown source position and sensor may not be respected. In many the non-line-of-sight (NLOS) environments, signals emitted from the source are often inaccessible by all measurement sensors, and one of the solutions is to use wireless sensor networks (WSN) in which received signals by each sensor are transmitted to the fusion center for source localization. Authors in [39] proposed a distributed NLOS cooperative localization algorithm. The proposed localization algorithm employs the multiplicative convex model based on the physical mechanism of the NLOS propagation to achieve robustness in changing environments. A TDOA-based cooperative localization approach for mixed LOS/NLOS conditions is proposed by [40]. For the location of the multiple stationary target nodes, the authors formulated a non-convex robust weighted least squares problem (RWLS). To efficiently solve RWLS, the semidefinite relaxation technique is used to transform RWLS into a convex mixed semidefinite and second-order cone programming problem. Authors in [41] presented the energy-based localization solution in WSNs using received signal strength (RSS) and received signal strength difference (RSSD). In the proposed solution, RSSD is based on transforming the nonlinear and non-convex objective functions into a convex optimization problem via relaxation and semidefinite programming. Another mixed semidefinite and second-order cone relaxation for source localization in 3D WSN was proposed by [42]. The proposed target node localization is based on hybrid RSS-AOA measurements in both noncooperative and cooperative WSNs, where the authors proposed new LS estimator to reduce the implementation costs. From conducted analysis authors concluded that proposed RSS-AOA based estimator is more suitable for large-scale cooperative WSNs compared with AOA-TDOA based estimators. In general, centralized algorithms suffer from high computation complexity, network traffic bottleneck, and as such, are not recommended in scenarios where each sensor node cannot get raw measurements directly. To address such problem authors in [43] proposed a completely decentralized localization approach based on augmented Lagrangian methods and alternating direction method of multipliers (ADMM). Discussing performance comparison in terms of computational efficiency, it is visible that the unknown source position can be directly calculated from different sets of geometric relationships. Although computational complexity for this kind of approach is low, it frequently does not perform sufficiently when measurement error exists. Since TDOA and AOA observations have measurement error, the location of a source is often estimated by propagating these errors trough the computation, where iterative nonlinear minimization is required for an optimal solution. The most straightforward approaches for handling measurement noise are iterative algorithms based on initial position estimate, obtained, for example, by the Gauss-Newton method, which has high computational requirements. The maximum likelihood estimator is asymptotically efficient, but it requires a good initial guess. To avoid the need for good initial position guess characteristic for iterative ML approach, different closed-form source estimation methods are proposed. A linear LS approach is an alternative approach, which can achieve CRLB but reports large estimation bias. In general, various closed-form solutions have been proposed, each designed to reduce biases, or to work with a different number of sensors. To asses the performance of localization estimation, most state-of-the-art studies perform benchmark in terms of CRLB. Analyzing literature, one can see that the weighted least square estimator is superior in comparison to LS solutions in terms of lower bias. Also, it can be concluded that even both WSL and ML approaches can achieve CRLB, the WSL approach is more computationally attractive and suitable for real-time application and does not really on the initial guess. From the performance comparison between different methods, one can conclude that while most state-of-the-art closed-form solutions achieve CRLB, each approach is designed to work under restricted measurement error.
In the context of room acoustics, it is difficult to establish measurement error as the measurement microphone captures not only the direct-path component of the source signal but also the multipath component caused by reflections. The multipath component, together with the background noise, can lead to distortion of the time delay information from received signals, and thus it can degrade the localization performance. To address key challenges of the realistic environment such as room reverberation, background noise, and sound interference, different methods to compute the TDOAs across various combinations of pairs of spatially separated microphones were proposed [44]- [47]. Recent work also suggested that deep learning can successfully be applied for modeling rooms acoustic [48]. For instance, in [49] authors employed deep model for SSL, where it was shown that deep learning-based system achieved higher accuracy under low SNR conditions in comparison with cross-correlation phase transform (GCC-PHAT) method. Authors in [17] proposed a novel learning approach for SSL based on TDOA estimation, where coordinates of a sound source were defined as functions of TDOA. In their work, pre-recorded sound measurements and their corresponding source locations were used to train the multilevel B-Splines based learning model. A new dataset for learning-based SSL was proposed by [50] which contained different acoustic events recorded in the anechoic chamber, where the anechoic chamber environment was used to verify the feasibility of the proposed baseline model. Authors in [51] addressed SSL for indoor environments with high reverberation and low signal-to-noise ratio. The authors proposed a novel sound source localization method using a probabilistic neural network for the classification of 3D space clusters. Another learning-based approach was demonstrated in [52], where authors evaluated the Structural Sparse Bayesian Learning model with signals recorded in an anechoic chamber with one reflective plate. Moreover, it can be seen that many researchers validated their SSL methods by using recorded audio signals simulating free field conditions or experimental analysis was carried out inside a semi-anechoic and anechoic chamber. The anechoic chamber provides a good simulation of the outdoor conditions due to the low level of reflections. This environment provides the possibility to precisely control the conditions and to measure the levels of sound events and noise, which is substantial for many experiments. Kotus et al. [53] performed multiple sound sources localization in the anechoic chamber, where different methods for obtaining the direction of arrival were tested. In the work [54], the authors proposed a modified cross-correlation algorithm to obtain a more reliable measurement of time difference of arrival in the reverberation environment. They performed a triangulation procedure using calculated TDOA values to obtain a sound source position in 3D space. Another three-dimensional method for SSL was proposed by Ding et al. [55]. The authors performed theoretical simulation accompanied by experimental results in the anechoic chamber. They proposed the use of a planar microphone array combined with a beamforming technique to obtain the location of the point sound source.

B. OUTLINE
This paper presents a passive sound source localization method using three soundfield microphone stations. The geometric configuration of three soundfield microphones can be employed to obtain two TDOA and three AOA measurements concerning unknown source position. A closed-form mathematical solution for SSL estimation is presented. Results are given in terms of RMSE, where it was shown that simulation of the sound source estimation algorithm can reach Cramér-Rao lower bound for small Gaussian noise presented in mea-surements. Theoretical simulation is supported by the experimental analysis conducted in our anechoic chamber [56]. In this work, TDOA and AOA measurements are directly obtained by exploiting the A and B format of the soundfield microphones, which can achieve small measurement deviations. This paper also demonstrates the value of using a soundfield microphone for the SSL task because AOA can be easily obtained due to the configuration of the soundfield microphone capsules.
The rest of the paper is organized as follows. Section II describes the mathematical formulation of the SSL algorithm. Section III presents our TDOA and AOA estimation methods using soundfield microphones. The simulation and experimental results in the anechoic chamber are given in Section IV, and Section V concludes the paper.

II. 3D SOURCE LOCATION ESTIMATION
This section explains a novel sound source localization method based on three soundfiled microphone stations. First, we establish a geometrical relationship between an unknown position of the sound source and known positions of sensors, over obtained AOA and TDOA observations. Second, we define a weighted least squares estimator, where the sum of the squared residuals is minimized with respect to the error measurement vector. Measurement error is modeled as a covariance matrix containing measurement uncertainties. In experimental phase those measurement error parameters are obtained through testing in the anechoic chamber [56].
Here we presented the theoretical formulation of the localization scenario where a geometric configuration of three ∈ R 3 is used to estimate the position of the single sound source s = [s x , s y , s z ] T ∈ R 3 in 3D space. Assuming that each station on position m i can determine bearing angles of received sound wave transmitted from the unknown source position s, geometric relationship between s and m i can be expressed trough observed AOAs by nonlinear equation: where θ i ∈ (−π, π) and φ i ∈ (0, π/2) form an AOA pair corresponding to azimuth and elevation angles in the righthand coordinate system. Let r i be the Euclidean distance between the source and microphone m i : and r 1 is the true distance from m 1 to the unknown source, then all TDOA observations with respect to m i are given as: where τ i1 is the range difference, t i1 is the actual TDOA measurement obtained as the propagation lag between m i and m 1 and v is the speed of signal propagation.
To estimate the unknown source position s using true TDOA and AOA values, we define vector µ = [τ 21 , τ 31 , κ T 1 , κ T 2 , κ T 3 ] T ∈ R 8 containing two TDOA observations with respect to m 1 and three AOA observations produced from positions m 1 , m 2 and m 3 . In practice, the AOA and TDOA are noisy observations where we assume that true values are influenced by additive Gaussian noise as: When expressing noisy TDOA and AOA in vector form results in:μ whereμ is actual noisy measurement vector µ and ε = [ε 21 , For localization scenario shown in Fig 1. we construct a unit norm vector and to satisfy the given geometric relation: By constructing the matrix G i in order to have columns that are orthonormal basis of the plane orthogonal to b i , it follows that G T i G i = I 2×2 and G T i b i = 0 2 where I 2×2 and 0 2 are identity matrix and zero vector. Let G i ∈ R 3×2 be: By performing left multiplication on (6) with G T i we get the following expression: Equation (9) can be rewritten with respect to m 1 as where r i is the geometric distance from m i to s. Substituting (3) as r i in (10) we have 1 and after using (6) as (10) becomes: A closed form solution was obtained in [26] for two stations and one source, while in this work, we expand the solution to work with three stations. Using (8) and (11) together in matrix form gives: (15) Equation (12) represents an ideal condition which doesn't hold in practice since there is a measurement error ε in h and G matrix. To analyze the influence of the error ε for a given measurement vectorμ in (5) on localization accuracy, ε will be introduced in (1) and (3) where the true geometric values of TDOA and AOA will be expressed in terms of their noisy observations (4). We can approximate (12) up to first order noise terms which gives: whereĥ andĜ are matrices containing noisy measurements instead of true TDOA and AOA values as their counterparts h and G. The formulation of T matrix is shown in (18) - (20), as shown at the bottom of the next page. From (5) it follows that T in (16) is zero-mean Gaussian with covariance matrix Q.
To calculate a weighted Least Square estimate of s from (16), the sum of squared residuals is minimized with the respect to error measurement vector ε.
where W = T QT T is weight matrix. The expression (17) is based on the implicit assumption that the measurement errors are uncorrelated with each other and that TDOA and AOA observations have corresponding σ RD and σ AOA uncertainties. This is insured by modeling Q as a matrix, where the diagonal elements are given in form of two TDOA and three AOA variances, setting the rest of off-diagonal entries to null. VOLUME 8, 2020

III. TDOA-AOA ESTIMATION METHOD
In this section, we explain the proposed approach for estimating TDOA and AOA values, computed from the real sound measurements collected by three soundfield microphone stations. Given that the solution (17) performs minimization of the sum of squared residuals with respect to measurement error, it was essential to determine observation uncertainties σ RD and σ AOA of the proposed TDOA-AOA estimation method in free field conditions. In this work, to ensure universal reference regarding methods uncertainties, all measurements are made in a small anechoic chamber [56] which facilitates free-field conditions where no reverberation of the sound source is present. Fig. 2 displays our measurement setup with three soundfield microphones and a loudspeaker representing a single sound source inside an anechoic chamber. Soundfiled microphone uses 4 subcardioid capsules mounted as close as possible to form a tetrahedron. Each soundfield microphone can be viewed as four symmetrical receivers positioned on the surface of the sphere where it can produce two distinct sets of audio signals called A-format and B-format. A-format consists of 4 signals coming from each microphone capsule arranged as shown in Fig. 3. The B-format signals comprise a truncated spherical harmonic decomposition of the sound field. They correspond to the sound pressure and the three components of the pressure gradient at a point in space. The transformation from A to B-format can be easily performed by knowing the measurement values of the individual capsules in A-format. A linear system of equations displayed below can be used for format conversion: where p w is sound pressure signal at the microphone position, p x is the sound velocity in the direction back and forth, p y is the sound velocity information in the direction left to right, and p z is the sound velocity information in the direction up and down. Additionally signal filtering can be used to compensate for inequalities between individual capsules. To validate estimation algorithm (17) using real measurements, this paper proposes two methods for calculating TDOA and AOA values from observed sound signal emitted from unknown source location s.

A. TDOA ESTIMATION
The system is composed of an unknown sound source and three soundfield microphones. If the sound source transmits a signal at time T = 0, the microphones will sense the 87754 VOLUME 8, 2020 signals at the unknown times T1, T2, and T3. Following (3), the difference between TOAs can be measured as follows: where the microphone at position m 1 is taken as a reference.
Note that for the case of three microphones T 32 value could be additionally calculated to introduce redundancy in the system. To calculate the time lag between two microphones, in this work we used the generalized cross-correlation phase transform (GCC-PHAT) algorithm as follows:  (17) assumes that the measurement vector is corrupted with small Gaussian noise defined by covariance matrix Q, pointing to the need of investigating the uncertainty of our TDOA estimation method. In order to analyze our TDOA based estimation performance given measurement duration t, we examined uncertainty reduction when the length of the sound signal used for estimation increases. From the results shown in Fig. 4 it is visible that by choosing measurement duration t > 300 ms will provide stable TDOA observations (RD < 5cm).

B. AOA ESTIMATION
In this work, we propose a method for calculating AOA estimation by exploiting B-format signals of the soundfield microphone. The AOA estimation technique is based on The procedure is performed by calculating the power of each velocity signal for a given integration interval t, where each velocity signals corespondent to X, Y, Z component of B-format. Next, normalization is applied by dividing each p component with R = p 2 x + p 2 y + p 2 z where p x , p y and p z represents the total power over t for each discrete axis. To produce a pair of AOA observations from the calculated directional vector d, similar  expression as (1) is used: where θ i and φ i are angles concerning microphone internal origin. It should be noted here that (24) is prone to the error if the axis of each soundfield microphone used for AOA estimation are not perfectly aligned with the axis of global geometry. Following similar logic as in previous subsection III-A, we investigated the dependency between the directivity error and integration time i.e. duration of the sound signal used for calculating AOA observations. From the Fig.5 it is visible that for t > 400 ms our methods starts to stabilize, providing estimation results with standard deviation of σ AOA < 0.5 • . Based on our findings in terms of dependency between σ AOA and t and by following the conclusion regarding σ RD and t from subsection III-A, we have chosen duration of the signal measurement to be t = 500ms in all our subsequent experiments.

IV. LOCALIZATION RESULTS AND PERFORMANCE EVALUATION
This section presents the performance study regarding the localization accuracy of the proposed SSL approach in the free field conditions. First, a simulation is made to demonstrate a margin of localization error regarding uncertainty in TDOA-AOA observations. To alleviate the dependency on the particular geometry, TDOA-AOA values are calculated by sampling the angles between three microphone stations from a uniform distribution. In total, three scenarios are conducted to analyze how measurement error affects the deviation of localization results compared to theoretical CRLB. Second, the experimental validation of the simulation results is made in the anechoic chamber containing one sound source and three soundfield microphones. Anechoic environment, however far from being realistic, gives the possibility to precisely control the environmental conditions and to accurately determine the measurement deviations of the proposed TDOA-AOA estimation method, which is substantial in this experiment. If the experiment was carried out in a reverberant room, the room acoustics would influence the estimation of TDOA-AOA values and thus obtained measurement uncertainties would not make a universal reference. Identifying measurement uncertainties σ RD and σ AOA , proposed solution (17) can be tested using estimated TDOA-AOA values calculated from real sound source.

A. SSL SIMULATION UNDER GAUSSIAN NOISE
This subsection presents performance analysis of the proposed estimation method where it will be shown that for the given localization scenario, the proposed method can reach CRLB for Gaussian noise present in TDOA and AOA measurements. For the known source position s, CRLB matrix is defined by [57]: where FIM is the Fisher information matrix given as (26) for zero-mean Gaussian with a covariance matrix Q ∈ R 8×8 . Taking partial derivative of µ with respect to s x , s y and s z yields: where: values d x,i = s x − m x,i and respectively d y,i and d z,i are the difference from source position s and m i for a given direction. The theoretical value of CRLB presents the best possible accuracy that estimator can achieve given small measurement error. To adequately analyze the performance of the proposed method, evaluation is made in terms of the difference between the theoretical value of square root CRLB and root mean square error (RMSE) of the estimation algorithm for a given measurement error variance σ 2 AOA , σ 2 RD . In this work, measurement error covariance matrix Q is modeled as: In total, three different simulations for a sound source localization scenario is made. Firstly, we inspect the performance of the proposed estimation method when the uncertainty of TDOA measurement increases. Secondly, we perform a similar test by increasing the standard deviation σ AOA . And for the third simulation, a difference of RMSE from theoretical root CRLB will be analyzed in a scenario where the range of the source s is increasing from the origin.
For an input simulation data, TDOA and AOA observations are modeled as (4), where τ i1 and κ i correspond to true geometrical values calculated from (1) and (2), while ε i1 and i are corresponding zero-mean Gaussian noise defined by error variance σ 2 AOA , σ 2 RD . We should also note that measurement noise between each TDOA and AOA observation is uncorrelated.
Each simulation results are reported in terms of RMSE which is defined as RMSE(s) = L l=1 ||s l − s|| 2 /L, wherë s l is the estimated source position after l-th iteration and s is true position of the source. L = 500 is the number of estimation runs for a given σ RD , σ AOA and s. Note that weight matrix W in (17) depends on actual source position via (18) which is not know in advance. We can define W as an identity matrix I ∈ R 8×8 to obtain an initial position estimatess 1 from (17) which in this case reduces (17) back to LS estimation. For the rest of the simulation, the matrix W is computed using previously estimated source positions l−1 together with the current AOA observation. Simulation procedure respects geometric configuration shown in Fig. 1 where three stations are placed 1 m apart forming vertices of an equilateral triangle 87756 VOLUME 8, 2020  with the center at the origin and lying on the XY plane. For each simulation, the calculated RMSE value is given by averaging results obtained for ten different rotations of the equilateral triangle around Z-axis.
Results of the first simulation are shown in Fig. 6 For a fixed parameters σ AOA = 0.5 • and s = [160, 40,40] T cm, reported RMSE closely follows theoretical CRLB given increase in σ RD . Fig. 7 illustrates the performance of the second simulation where σ AOA varies from [0.5 • − 2.5 • ] while parameter σ RD is fixed at 10 cm leaving source at position s = [160, 40,40] T cm. The third simulation shows the performance evaluation depending on the distance from the source for the fixed σ RD = 10 cm and σ AOA = 0.5 • . For this case the source moves away from the origin in accordance with the expression s = a [50,50,50] T where a ranges from 1 to 5. From Fig. 8 it is clear that proposed method provides near CRLB performance.

B. EXPERIMENTAL EVALUATION
In this section, we present experimental results on sound source localization in 3D space. Following theoretical simulation given in subsection IV-A, three TetraMic soundfield microphones were placed 1 m apart forming vertices of an equilateral triangle. Arta measurement software was used as the white noise generator, which was then converted to the analog audio signal by MOTU STAGE B-16 audio interface, amplified with Pioneer SA-8500 II sound amplifier and then played through Magnat 145 801 12 loudspeaker, which represented the source. The TetraMic is connected to audio interface via four XLR-M connectors coming from the PPA receiver. For each soundfield microphone, the four-channel signal is digitized with the B-16 audio interface and recorded with Reaper software in A-format. To obtain B-format from four-channel A-format in this work, we used an official software VVMic which also performs the corrections using the calibration file, provided with the specific microphone. In the experiment we used 48 different sound source locations. Due to space limitations of our anechoic chamber all measurements were conducted regarding s x > 0 by following measurement plan shown in Fig. 9, where after performing measurements for each 12 different sound source locations, an equilateral triangle holding TetraMic microphones was rotated around Z-axis for 90 • .
For each source location s i white noise signal of 18s duration was recorded. As concluded in the previous section, t = 500 ms interval was chosen for measurements' duration which resulted in 36 observations for each source position s i . VOLUME 8, 2020 Following the procedure explained in section III, for each segment we calculated two TDOAs from A-format and three AOAs from B-format. Measured observations were first used in (17) to get the initial position estimates 1 of the source which is needed to calculate the weight matrix W . After obtaining an initial position estimate, method (17) was used 35 more times where for each position update, source estimation from the previous step was used to calculate W which was used in conjunction with new TDOA and AOA observations to estimate new sound source position in 3D space. Since the computational step (17) is performed for each position update, it is interesting to analyze the efficiency of (17) in terms of computation time given its corresponding TDOA-AOA values. In this work, all processing procedures were carried out on the platform of Matlab 2019 on a 64-bit PC with the computational capability of Intel(R) Core(TM) i5 CPU @ 3.80 GHz and memory of 32.00 GB 400 MHz DDR4. In the average running time of our algorithm (17) for initial position estimation is 25.4ms and 3.9 ms for each additional position update.
To show the effectiveness of the proposed approach, Fig. 10 illustrates the estimation results for 4 different sound source locations, where reported RMSE results were 11.17, 15.56 18.90, and 10.47 cm respectively for each location. It is interesting to note that the estimate (shown in cyan) distant from the rest of the estimation group corresponds to the initial position guess, where the true source location is indicated by the coordinates label. Narrow estimation spread is clearly visible.
It was also interesting to analyze the dependency between reported RMSE and the number of iterations. From Fig. 11 it is visible that by using more iterations, (17) will provide better results. For longer measurement intervals, we can achieve even higher accuracy. The obtained result indicates that the proposed procedure can be used to locate the stationary sound source with excellent precision in the free-field. We plot the experimental cumulative distribution function (CDF) of positioning error RMSE in Fig.12 to show the localization performance of the proposed weighted least square estimator (17). CDF in Fig. 12    From the presented results, it is evident that the proposed measurement method provides accurate estimations of azimuth and elevation angles of the unknown sound source in the anechoic chamber. In particular, free field conditions play a major role in the acoustic measurements and sound perception experiments, as the results are influenced only by the direct-path component of the sound source and not by multipath component caused by room reflections. Given that the experimental evaluation in the SSL task was performed using real sound measurements, the above results indicate the potential for further research, where it can be expected that by employing the more robust method for estimating TDOA, such as method in [54], the results can translate to real-world cases. Although it is undeniable that the room acoustics would influence the SSL recognition results, conducting experimental evaluation in a controlled environment provides universal reference to localization performance. Since most state-of-the-art TDOA-based SSL either do not have experimental analysis (where results are given through simulation) or experimental evaluation is performed in non reproducible environmental conditions, a direct comparison is often inaccessible. Regardless of the results obtained under certain conditions, our comparison with state-of-the-art is made to demonstrate the performance of the proposed SSL approach. For instance, the state-of-the-art method in [17] reported a localization performance of 0.3m in terms of RMSE concerning the typical room experiment. Authors also compared themselves with commonly used SRP-HAT method, and they showed that their approach outperformed the SRP-HAT method, which reported a localization error of 0.76m. It is to expect that evaluation results obtained in the free-field will outperform results obtained in a real-life environment [17], and from comparison, it can be seen that the proposed TDOA-AOA approach, based on three soundfiled microphones, achieves 0.09 m lower localization error compared to Multilevel B-Splines-Based Learning Approach and 0.55 m compared to SRP-HAT method [17]. Given the fact that the SSL results presented are 30% and 72% lower in terms of RMSE indicate the effectiveness of the proposed approach in free-filed with great potential to operate under realistic conditions.

V. CONCLUSION
This paper proposes a sound source localization method in 3D space using a geometric configuration of three microphone stations. The closed-form solution for estimating a sound source location based on two TDOAs and three AOAs is presented. In this work, a soundfield microphone is used as a measurement station. By exploiting A-format of the soundfield microphone, a pair-wise TDOA estimation using a general cross-correlation is obtained as the time lag between the signals from two microphones. The method for obtaining AOA measurements is proposed, based on a calculating directional vector derived from B-format. The investigation of the impact of the signal measurement interval length on TDOA and AOA estimation performance is made. The proposed method was evaluated by simulations and physical experiments in our anechoic chamber. The simulation results of the sound source localization reported that (17) reaches the theoretical CRLB performance regarding a small Gaussian noise presented in measurements. Since (17) can be viewed as iterative algorithm, an investigation between accuracy and number of iteration steps was made, concluding that with a longer recording of the analyzed signals from unknown source position, we can achieve even higher accuracy. The obtained experimental results indicate that the proposed procedure can be used to locate the stationary sound source with good precision in free-field environment. In our future work, it is planned to extend the proposed method for multiple source localization and tracking, first in free-field anechoic chamber conditions, and later in different SNR and reverberation conditions. Special efforts will be put on problems tackling real-life scenarios, like detecting speech in noisy environments when distant recording is required, in order to increase speech recognition performance. He is currently involved as an Associate Professor and is leading the Laboratory of Electroacoustics. He has participated in several scientific and professional projects. His research interests include programming, programming languages, electroacoustics, geographic information systems, and simulations and virtual and augmented reality technologies. He is a member of the Association for Computing Machinery, National Geographic Society, and Croatian Acoustical Society. VOLUME 8, 2020