# Research on Ambiguity Resolution Algorithm by Quaternion Based on Acoustic Vector Sensor

Authors: Guibao Wang; Xinkuan Wang; Lanmei Wang; Xiangyu Wang
Journal: Journal of Sensors
 (2020)

---

## Abstract

The increase in element spacing can increase the aperture of the array and improve its resolution performance. However, phase ambiguity will occur when the array element interval is larger than the minimum half wavelength of the incident signal. The three acoustic velocity components of the acoustic vector are ingeniously constructed into a new kind of quaternions because of the special structure of the acoustic vector sensor array, and the rough estimation of the direction of arrival (DOA) is obtained using the rotation relationship between the subarray steering vectors corresponding to quaternion data. The rough estimate is used to resolve the phase ambiguity of the spatial phase difference between the array elements, and the high-precision DOA estimation of the signal can be obtained. Simulation results show that the method is effective.


---

## Body

## 1. Introduction

The use of electromagnetic wave for communication or target location has become mainstream. However, the acoustic wave has more advantages than the electromagnetic wave in some special circumstances. For example, sound waves can travel around obstacles, such as mountains, rivers, and forests, and the sound source detection system can effectively detect the information of the target. The loss of acoustic wave propagation is small in water. Thus, the wave can spread for a long distance. Sonar is a key equipment in underwater communication. The acoustic vector sensor is the core part of new sonars. A single acoustic vector sensor not only can output the acoustic pressure information of the scalar sensor but also can measure the acoustic velocity information of each point directly and synchronously [1, 2]. The acoustic vector sensor consists of an acoustic pressure component sensor and a speed sensor in three-speed directions. The vibration speed sensor may be an acceleration, a displacement, or a sound pressure gradient sensor [3, 4]. In 2010, Wu et al. obtained the near-field source model of the acoustic vector sensor by using the wave equations of the acoustic and particle velocity fields [5]. In 2012, Song et al. proposed to place the components of the three speed directions at the coordinate origin and an acoustic pressure component on one of the coordinate axes; the former is used to estimate the angle of the sound source, and the latter is used to estimate the distance information of the sound source [6]. In the array direction finding, the resolution of the angle with large array spacing will improve correspondingly. However, when the array spacing increases to more than half the wavelength of the incident signal, the pattern of the array will show grating lobes, which will lead to the ambiguity problem of parameter estimation [7, 8]. Many ways can be used to solve the ambiguity problem. Wang et al. proposed a virtual baseline method of resolving ambiguity based on the virtual ring [9]. Wong and Zoltowski proposed a uniform plane array method for solving ambiguity [10]. A nonuniform array method for solving ambiguity problems was proposed in [11], and the problem of ambiguity resolution by virtual array expansion was studied in [12]. In fact, the general mode of the ambiguity-solving method is to obtain two different arrays through array arrangement, virtual reconstruction, or data processing. One nonsparse array without ambiguity can enable rough parameter estimation, and the other sparse array with phase ambiguity can realise accurate parameter estimation. The accurate and unambiguous estimations of the direction of arrival (DOA) are obtained by using the rough estimation of the nonsparse array to solve the phase ambiguity of the sparse array.In recent years, array signal processing technology based on quaternion [13] has attracted increasing attention. In [14], the singular value decomposition of the quaternion matrix was studied, and a polarisation wave separation algorithm based on the real signal quaternion model was proposed. In [15, 16], the quaternion and biquaternion music algorithms for two- and three-component vector sensor arrays were proposed, respectively. In [17], the DOA estimation for the six-component electromagnetic vector sensor array was obtained using a quad-quaternion music operation. In [18], a quaternion ESPRIT algorithm based on the translation-invariant cross dipole array was proposed. These studies show that the parameter estimation algorithm based on the multivariate number has improved model error tolerance by using the constraint relationship between the different imaginary parts of the multivariate number and the strong constraint condition between the orthogonal vectors of multivariate number. In this work, a quaternion–ESPRIT algorithm for solving the direction-finding ambiguity based on the acoustic vector is proposed. This algorithm is different from the previous ambiguity-solving method. The DOA is roughly estimated using the rotation invariant relationship between the quaternion subarrays and is used to solve the phase ambiguity of the spatial steering vector. The phase ambiguity number vector is attained using the relationship between the estimated value of the phase difference between sparse elements with phase ambiguity and the rough value without ambiguity but with low accuracy. The high-accuracy estimated value without ambiguity is obtained in the end. This method has a small amount of computation and can enable high-precision DOA estimation without ambiguity.


## 2. Acoustic Vector Array Model

A uniform circular array is composed ofM elements; R is the radius of the circular array, as shown in Figure 1. Phase ambiguity occurs when the interval between the adjacent elements is greater than the half wavelength, that is,Δd>λmin/2, where λmin is the minimum wavelength of the incident signal. The array elements are the acoustic vector sensors; each of these elements is composed of three orthogonal velocity sensors along the x-axis, y-axis, and z-axis, which induce the velocity component of the incoming wave signal in the x-axis, y-axis, and z-axis directions, respectively.Figure 1
Uniform circular array geometry.Assuming thatK narrow-band, zero-mean, and stationary signals are incident from the far-field onto the uniform circular array shown in Figure 1, the incident direction of the kth signal source is θk,ϕk, k=1,⋯,K, and φn is the angular coordinates of the nth element. The velocity components along the x-, y- and z-axes of the kth unit energy signal in the Cartesian coordinate system can be expressed as:
(1)uk=vxkvykvzk=sinθkcosϕksinθksinϕkcosθk.


## 3. Quaternion Method

### 3.1. Mathematical Model of Quaternion Method

Thex-axis velocity component vxk, y-axis velocity component vyk, and z-axis velocity component vzk on the same array element are expressed as quaternion models as Bk:
(2)Bk=vzk+ivxk+jvyk=cosθk+isinθkcosϕk+jsinθksinϕk=cosθk+isinθkejϕk=cosθk1+itanθkejϕk,where i and j are the imaginary parts of a quaternion.The spatial steering vector of the array consisting of the phase difference between theM array elements and the reference origin on the circumference is as follows:
(3)qθk,ϕk=ej2πR/λksinθkcosϕk−φ1,⋯,ej2πR/λksinθkcosϕk−φn⋯,ej2πR/λksinθkcosϕk−φM.In the formula (3) above, φn=2πn−1/M,n∈1,⋯,M represents the angular coordinates of the nth element on the circumference, and λk is the wavelength of the kth incident signal.The received data of the array att time can be expressed as:
(4)X1t=∑k=1KBk⊗qθk,ϕkSkt+N1t=C1St+N1t,whereSkt is the kth incident signal, ⊗ represents the Kronecker product, C1 is the steering vector matrix at time t, St is the K incident signals, and N1t is the additive white Gaussian noise at time t.Then, the received data of the array att+ΔT time can be represented as:
(5)X2t=X1t+ΔT=∑k=1KBk⊗qθk,ϕkej2πfkΔTSkt+N2t=C2St+N2t=C1ΘSt+N2t,where ΔT is the interval between two snapshots, fk is the frequency of the kth signal, C2 is the steering vector matrices at time t+ΔT, N2t is the additive white Gaussian noise at time t+ΔT, and
(6)Θ=ej2πf1ΔT⋱ej2πfKΔT.So the whole received data can be written as:(7)Zt=X1tX2t=C1C2St+N1tN2t=CSt+Nt,where C=C1C2 is the whole steering vector matrix, and Nt=N1tN2t is the whole additive white Gaussian noise.


### 3.2. Quaternion Ambiguity Resolution Algorithm

The quaternion correlation matrix of the received data isR(8)R=EZZH=CESSHCH+ENNH=CRsCH+σ2I,where E· denotes the mathematical expectation, ·H represents the conjugate transpose, Rs=ESSH, ENNH=σ2I, σ2 indicates the white noise power, and I is the identity matrix.The signal subspaceVs and noise subspace Vn can be obtained by quaternion eigendecomposition operations of the data correlation matrix R, V1 and V2 are the first M rows and the last M rows of Vs, respectively. According to the subspace theory, it can be got as follows:
(9)V1=C^1TV2=C^2T=C^1Θ^T,(10)V1#V2T−1=T−1Θ^,where T is the nonsingular matrix; the matrices C^1, C^2, and Θ^ are the respective estimates of C1, C2, and Θ; ·# denotes the pseudoinverse operator; and ·−1 represents the inverse operation. Let ψ=V1#V2=V1HV1−1V1HV2, then formula (10) can be expressed as:
(11)ψT−1=T−1Θ^,where the eigenvalue of ψ constitutes the estimation Θ^ of Θ, which contains the frequency information of the source, and T−1 is the corresponding eigenvector. The matrix factorisation of the quaternion is different from that of the long vector, that is, the eigendecomposition operation of the quaternion matrices cannot be performed directly. After the quaternion matrix is obtained, it is usually used to extract one real part and the three imaginary parts of the quaternion matrix and reconstruct the complex term to form the complex adjoint matrix. The properties of quaternion matrices are obtained by studying the properties of the plural adjoint matrix. Therefore, the corresponding information should be extracted from the quaternion matrix ψ. Accordingly, the estimate of the steering vector C1 is obtained as follows:
(12)C^1=V1T^−1,where T^−1 is the estimation of T−1.


#### 3.2.1. Rough Estimate of DOA

The eigenvalueΘ^ of ψ can be expressed as follows:
(13)Θ^=ej2πf^1ΔT⋱ej2πf^KΔT,wheref^1,⋯,f^K are the frequencies of the corresponding signals.The frequency of thekth signal source can be estimated as follows:
(14)f^k=argΘ^kk2πΔT,where arg· denotes the operation of the extracting phase and Θ^kk is the kth row and kth column element of the matrix Θ^.According to formulas (2) and (12), the steering vector C^1 can be written as another form:
(15)C^1=C^1x+iC^1y.Among them,C^1x and C^1y satisfy the following relationships:
(16)C^1x=C^1yΦ^1,Φ^1=C^1y#C^1x=C^H1yC^1y−1C^H1yC^1x,where
(17)Φ^1=tanθ^1ejϕ∧1⋱tanθ^Kejϕ∧K.According to formula (17), the rough estimate of DOA without ambiguity is obtained.
(18)θ^k=tan−1∣Φ^1kk∣ϕ^k=argΦ^1kk,where Φ^1kk is the kth row and kth column element of the matrix Φ^1.


#### 3.2.2. Ambiguity Resolution Algorithm

According to the unambiguous rough estimate ofθ^k and ϕ^k, the rough phase difference estimation Ω^1 can be obtained:
(19)Ω^1=ΔW⋅Γ^1kΓ^2k=W2−W1⋅Γ^1kΓ^2k,where ΔW=W2−W1, W1=W1:M−1,:, W2=W2:M,:, Γ^1k=sinθ^ksinϕ^k, Γ^2k=sinθ^kcosϕ^k, W1:M−1,: represents the submatrix extracted from the first row to M row of W, and W2:M,: represents the submatrix extracted from the second row to M row of W, with
(20)W=2πR/λksinφ1cosφ1⋮⋮sinφMcosφM.Let(21)g^1k=C^11:M−1,k,g^2k=C^12:M,k,where C^11:M−1,k and C^12:M,k denote vectors composed of elements from row 1 to row M and from row 2 to row M in the column k of the matrix C^1, respectively.From formula (21), the estimate of spatial steering vectors Δq^k can be written as:
(22)Δq^k=g^2k./g^1k,where ./ denotes the dot division which is the division corresponding elements of two vectors.By taking the phase on both sides of equation (22), the estimated value containing the ambiguous phase difference can be expressed as:
(23)Ω^2=argΔq^k.The ambiguous vectorPomp can be solved by the following formula:
(24)Pomp=argminPΩ^2+2πP−Ω^1The expressionargmin· represents the vector composed of the smallest integer satisfying the condition and P is the possible vector of ambiguity numbers.By substituting the obtained ambiguous vector into formula (25), the accurate unambiguous phase difference Ω^3 can be obtained as follows:
(25)Ω^3=Ω^2+2πPomp=ΔW⋅Γ^1k′Γ^2k′.From (25), the accurate DOA estimate without ambiguity is estimated as follows:
(26)θ^k′=arcsinΓ^1k′2+Γ^2k′2,ϕ^k′=arctanΓ^1k′Γ^2k′,Γ^2k′≥0,ϕ^k′=π+arctanΓ^1k′Γ^2k′,Γ^2k′<0.


## 3.1. Mathematical Model of Quaternion Method

Thex-axis velocity component vxk, y-axis velocity component vyk, and z-axis velocity component vzk on the same array element are expressed as quaternion models as Bk:
(2)Bk=vzk+ivxk+jvyk=cosθk+isinθkcosϕk+jsinθksinϕk=cosθk+isinθkejϕk=cosθk1+itanθkejϕk,where i and j are the imaginary parts of a quaternion.The spatial steering vector of the array consisting of the phase difference between theM array elements and the reference origin on the circumference is as follows:
(3)qθk,ϕk=ej2πR/λksinθkcosϕk−φ1,⋯,ej2πR/λksinθkcosϕk−φn⋯,ej2πR/λksinθkcosϕk−φM.In the formula (3) above, φn=2πn−1/M,n∈1,⋯,M represents the angular coordinates of the nth element on the circumference, and λk is the wavelength of the kth incident signal.The received data of the array att time can be expressed as:
(4)X1t=∑k=1KBk⊗qθk,ϕkSkt+N1t=C1St+N1t,whereSkt is the kth incident signal, ⊗ represents the Kronecker product, C1 is the steering vector matrix at time t, St is the K incident signals, and N1t is the additive white Gaussian noise at time t.Then, the received data of the array att+ΔT time can be represented as:
(5)X2t=X1t+ΔT=∑k=1KBk⊗qθk,ϕkej2πfkΔTSkt+N2t=C2St+N2t=C1ΘSt+N2t,where ΔT is the interval between two snapshots, fk is the frequency of the kth signal, C2 is the steering vector matrices at time t+ΔT, N2t is the additive white Gaussian noise at time t+ΔT, and
(6)Θ=ej2πf1ΔT⋱ej2πfKΔT.So the whole received data can be written as:(7)Zt=X1tX2t=C1C2St+N1tN2t=CSt+Nt,where C=C1C2 is the whole steering vector matrix, and Nt=N1tN2t is the whole additive white Gaussian noise.


## 3.2. Quaternion Ambiguity Resolution Algorithm

The quaternion correlation matrix of the received data isR(8)R=EZZH=CESSHCH+ENNH=CRsCH+σ2I,where E· denotes the mathematical expectation, ·H represents the conjugate transpose, Rs=ESSH, ENNH=σ2I, σ2 indicates the white noise power, and I is the identity matrix.The signal subspaceVs and noise subspace Vn can be obtained by quaternion eigendecomposition operations of the data correlation matrix R, V1 and V2 are the first M rows and the last M rows of Vs, respectively. According to the subspace theory, it can be got as follows:
(9)V1=C^1TV2=C^2T=C^1Θ^T,(10)V1#V2T−1=T−1Θ^,where T is the nonsingular matrix; the matrices C^1, C^2, and Θ^ are the respective estimates of C1, C2, and Θ; ·# denotes the pseudoinverse operator; and ·−1 represents the inverse operation. Let ψ=V1#V2=V1HV1−1V1HV2, then formula (10) can be expressed as:
(11)ψT−1=T−1Θ^,where the eigenvalue of ψ constitutes the estimation Θ^ of Θ, which contains the frequency information of the source, and T−1 is the corresponding eigenvector. The matrix factorisation of the quaternion is different from that of the long vector, that is, the eigendecomposition operation of the quaternion matrices cannot be performed directly. After the quaternion matrix is obtained, it is usually used to extract one real part and the three imaginary parts of the quaternion matrix and reconstruct the complex term to form the complex adjoint matrix. The properties of quaternion matrices are obtained by studying the properties of the plural adjoint matrix. Therefore, the corresponding information should be extracted from the quaternion matrix ψ. Accordingly, the estimate of the steering vector C1 is obtained as follows:
(12)C^1=V1T^−1,where T^−1 is the estimation of T−1.


### 3.2.1. Rough Estimate of DOA

The eigenvalueΘ^ of ψ can be expressed as follows:
(13)Θ^=ej2πf^1ΔT⋱ej2πf^KΔT,wheref^1,⋯,f^K are the frequencies of the corresponding signals.The frequency of thekth signal source can be estimated as follows:
(14)f^k=argΘ^kk2πΔT,where arg· denotes the operation of the extracting phase and Θ^kk is the kth row and kth column element of the matrix Θ^.According to formulas (2) and (12), the steering vector C^1 can be written as another form:
(15)C^1=C^1x+iC^1y.Among them,C^1x and C^1y satisfy the following relationships:
(16)C^1x=C^1yΦ^1,Φ^1=C^1y#C^1x=C^H1yC^1y−1C^H1yC^1x,where
(17)Φ^1=tanθ^1ejϕ∧1⋱tanθ^Kejϕ∧K.According to formula (17), the rough estimate of DOA without ambiguity is obtained.
(18)θ^k=tan−1∣Φ^1kk∣ϕ^k=argΦ^1kk,where Φ^1kk is the kth row and kth column element of the matrix Φ^1.


### 3.2.2. Ambiguity Resolution Algorithm

According to the unambiguous rough estimate ofθ^k and ϕ^k, the rough phase difference estimation Ω^1 can be obtained:
(19)Ω^1=ΔW⋅Γ^1kΓ^2k=W2−W1⋅Γ^1kΓ^2k,where ΔW=W2−W1, W1=W1:M−1,:, W2=W2:M,:, Γ^1k=sinθ^ksinϕ^k, Γ^2k=sinθ^kcosϕ^k, W1:M−1,: represents the submatrix extracted from the first row to M row of W, and W2:M,: represents the submatrix extracted from the second row to M row of W, with
(20)W=2πR/λksinφ1cosφ1⋮⋮sinφMcosφM.Let(21)g^1k=C^11:M−1,k,g^2k=C^12:M,k,where C^11:M−1,k and C^12:M,k denote vectors composed of elements from row 1 to row M and from row 2 to row M in the column k of the matrix C^1, respectively.From formula (21), the estimate of spatial steering vectors Δq^k can be written as:
(22)Δq^k=g^2k./g^1k,where ./ denotes the dot division which is the division corresponding elements of two vectors.By taking the phase on both sides of equation (22), the estimated value containing the ambiguous phase difference can be expressed as:
(23)Ω^2=argΔq^k.The ambiguous vectorPomp can be solved by the following formula:
(24)Pomp=argminPΩ^2+2πP−Ω^1The expressionargmin· represents the vector composed of the smallest integer satisfying the condition and P is the possible vector of ambiguity numbers.By substituting the obtained ambiguous vector into formula (25), the accurate unambiguous phase difference Ω^3 can be obtained as follows:
(25)Ω^3=Ω^2+2πPomp=ΔW⋅Γ^1k′Γ^2k′.From (25), the accurate DOA estimate without ambiguity is estimated as follows:
(26)θ^k′=arcsinΓ^1k′2+Γ^2k′2,ϕ^k′=arctanΓ^1k′Γ^2k′,Γ^2k′≥0,ϕ^k′=π+arctanΓ^1k′Γ^2k′,Γ^2k′<0.


## 3.2.1. Rough Estimate of DOA

The eigenvalueΘ^ of ψ can be expressed as follows:
(13)Θ^=ej2πf^1ΔT⋱ej2πf^KΔT,wheref^1,⋯,f^K are the frequencies of the corresponding signals.The frequency of thekth signal source can be estimated as follows:
(14)f^k=argΘ^kk2πΔT,where arg· denotes the operation of the extracting phase and Θ^kk is the kth row and kth column element of the matrix Θ^.According to formulas (2) and (12), the steering vector C^1 can be written as another form:
(15)C^1=C^1x+iC^1y.Among them,C^1x and C^1y satisfy the following relationships:
(16)C^1x=C^1yΦ^1,Φ^1=C^1y#C^1x=C^H1yC^1y−1C^H1yC^1x,where
(17)Φ^1=tanθ^1ejϕ∧1⋱tanθ^Kejϕ∧K.According to formula (17), the rough estimate of DOA without ambiguity is obtained.
(18)θ^k=tan−1∣Φ^1kk∣ϕ^k=argΦ^1kk,where Φ^1kk is the kth row and kth column element of the matrix Φ^1.


## 3.2.2. Ambiguity Resolution Algorithm

According to the unambiguous rough estimate ofθ^k and ϕ^k, the rough phase difference estimation Ω^1 can be obtained:
(19)Ω^1=ΔW⋅Γ^1kΓ^2k=W2−W1⋅Γ^1kΓ^2k,where ΔW=W2−W1, W1=W1:M−1,:, W2=W2:M,:, Γ^1k=sinθ^ksinϕ^k, Γ^2k=sinθ^kcosϕ^k, W1:M−1,: represents the submatrix extracted from the first row to M row of W, and W2:M,: represents the submatrix extracted from the second row to M row of W, with
(20)W=2πR/λksinφ1cosφ1⋮⋮sinφMcosφM.Let(21)g^1k=C^11:M−1,k,g^2k=C^12:M,k,where C^11:M−1,k and C^12:M,k denote vectors composed of elements from row 1 to row M and from row 2 to row M in the column k of the matrix C^1, respectively.From formula (21), the estimate of spatial steering vectors Δq^k can be written as:
(22)Δq^k=g^2k./g^1k,where ./ denotes the dot division which is the division corresponding elements of two vectors.By taking the phase on both sides of equation (22), the estimated value containing the ambiguous phase difference can be expressed as:
(23)Ω^2=argΔq^k.The ambiguous vectorPomp can be solved by the following formula:
(24)Pomp=argminPΩ^2+2πP−Ω^1The expressionargmin· represents the vector composed of the smallest integer satisfying the condition and P is the possible vector of ambiguity numbers.By substituting the obtained ambiguous vector into formula (25), the accurate unambiguous phase difference Ω^3 can be obtained as follows:
(25)Ω^3=Ω^2+2πPomp=ΔW⋅Γ^1k′Γ^2k′.From (25), the accurate DOA estimate without ambiguity is estimated as follows:
(26)θ^k′=arcsinΓ^1k′2+Γ^2k′2,ϕ^k′=arctanΓ^1k′Γ^2k′,Γ^2k′≥0,ϕ^k′=π+arctanΓ^1k′Γ^2k′,Γ^2k′<0.


## 4. Long Vector Method

### 4.1. Mathematical Model of the Long Vector Method

The received data of the acoustic velocity subarray in thex-axis, y-axis, and z-axis directions can be expressed as:
(27)Xxt=AxSt+Nxt,Xyt=AySt+Nyt,Xzt=AzSt+Nzt,where St is the incident signal; Nxt, Nyt, and Nzt are, respectively, the additive white Gaussian noise along the x-axis, y-axis, and z-axis directions; and Ax, Ay, and Az are the subarray steering vector along the x-axis, y-axis, and z-axis directions.Ax, Ay, and Az can be further expressed as follows:
(28)Ax=axθ1,ϕ1,⋯,axθk,ϕk,⋯,axθK,ϕK,Ay=ayθ1,ϕ1,⋯,ayθk,ϕk,⋯,ayθK,ϕK,Az=azθ1,ϕ1,⋯,azθk,ϕk,⋯,azθK,ϕK,where axθk,ϕk=vxkqθk,ϕk, ayθk,ϕk=vykqθk,ϕk, and azθk,ϕk=vzkqθk,ϕk, with qθk,ϕk the spatial steering vector of the array; its expression is the same as formula (3).The relationship between the three subarray steering vectors is:(29)Ax=AzΦx,Ay=AzΦy,where
(30)Φx=tanθ1cosϕ1⋱tanθkcosϕk⋱tanθKcosϕKΦy=tanθ1sinϕ1⋱tanθksinϕk⋱tanθKsinϕK.The received data of the acoustic velocity subarray along thex-axis, y-axis, and z-axis directions are arranged in sequence to form a long vector; then, the received data of the array at t and t+ΔT time can be expressed as:
(31)X¯1t=∑k=1Kuk⊗qθk,ϕkSkt+N1t=C¯1St+N1t,X¯2t=X¯1t+ΔT=∑k=1Kuk⊗qθk,ϕkej2πfkΔTSkt+N2t=C¯1ΘSt+N2t=C¯2St+N2t,where C¯1=Ax;Ay;Az is an array steering vector at time t, ⊗ denotes the Kronecker product, Skt is the kth element of St, N1t is the receiving noise of the whole array at time t, N2t is the receiving noise of the whole array at time t+ΔT, Nxt is the noise received by the x-axis vibration velocity subarray, Nyt and Nzt are the noises received by the y-axis and z-axis vibration velocity subarrays, respectively, and the expression of uk is the same as formula (1). C¯2=C¯1Θ is an array steering vector at time t+ΔT.In summary, the whole received data can be written as:(32)Z¯t=X¯1tX¯2t=C¯1C¯2St+N1tN2t=C¯St+N¯t,whereC¯=C¯1C¯2, N¯t=N1tN2t.


### 4.2. Long Vector Ambiguity Resolution Algorithm

The plural received data correlation matrix is:(33)R¯=EZ¯Z¯H=C¯ESSHC¯H+EN¯N¯H=C¯RsC¯H+σ2I.The signal and noise subspaceUs and Un can be obtained by the complex domain eigendecomposition of the data correlation matrix R¯. The signal subspace is divided into front and back halves, and the estimation of signal frequency fk and steering vector C¯^1 can be obtained by the ESPRIT method. The process is similar to Section 3.2.where(34)C¯^1=A^xA^yA^z.


#### 4.2.1. Rough Estimate of DOA

According to formulas (29) and (34), the following estimates can be obtained:
(35)Φ^x=A^Z#A^x,Φ^y=A^Z#A^y.The rough estimate value of DOA without ambiguity is given by:(36)θ~k=tan−1Φ^x2k,k+Φ^y2k,k,ϕ~k=tan−1Φ^yk,kΦ^xk,kΦ^yk,k≥0π+tan−1Φ^yk,kΦ^xk,kΦ^yk,k<0,where Φ^xk,k and Φ^yk,k are the kth row and kth column element of Φ^x. The ambiguity resolution algorithm is the same as in Section 3.2.2.


## 4.1. Mathematical Model of the Long Vector Method

The received data of the acoustic velocity subarray in thex-axis, y-axis, and z-axis directions can be expressed as:
(27)Xxt=AxSt+Nxt,Xyt=AySt+Nyt,Xzt=AzSt+Nzt,where St is the incident signal; Nxt, Nyt, and Nzt are, respectively, the additive white Gaussian noise along the x-axis, y-axis, and z-axis directions; and Ax, Ay, and Az are the subarray steering vector along the x-axis, y-axis, and z-axis directions.Ax, Ay, and Az can be further expressed as follows:
(28)Ax=axθ1,ϕ1,⋯,axθk,ϕk,⋯,axθK,ϕK,Ay=ayθ1,ϕ1,⋯,ayθk,ϕk,⋯,ayθK,ϕK,Az=azθ1,ϕ1,⋯,azθk,ϕk,⋯,azθK,ϕK,where axθk,ϕk=vxkqθk,ϕk, ayθk,ϕk=vykqθk,ϕk, and azθk,ϕk=vzkqθk,ϕk, with qθk,ϕk the spatial steering vector of the array; its expression is the same as formula (3).The relationship between the three subarray steering vectors is:(29)Ax=AzΦx,Ay=AzΦy,where
(30)Φx=tanθ1cosϕ1⋱tanθkcosϕk⋱tanθKcosϕKΦy=tanθ1sinϕ1⋱tanθksinϕk⋱tanθKsinϕK.The received data of the acoustic velocity subarray along thex-axis, y-axis, and z-axis directions are arranged in sequence to form a long vector; then, the received data of the array at t and t+ΔT time can be expressed as:
(31)X¯1t=∑k=1Kuk⊗qθk,ϕkSkt+N1t=C¯1St+N1t,X¯2t=X¯1t+ΔT=∑k=1Kuk⊗qθk,ϕkej2πfkΔTSkt+N2t=C¯1ΘSt+N2t=C¯2St+N2t,where C¯1=Ax;Ay;Az is an array steering vector at time t, ⊗ denotes the Kronecker product, Skt is the kth element of St, N1t is the receiving noise of the whole array at time t, N2t is the receiving noise of the whole array at time t+ΔT, Nxt is the noise received by the x-axis vibration velocity subarray, Nyt and Nzt are the noises received by the y-axis and z-axis vibration velocity subarrays, respectively, and the expression of uk is the same as formula (1). C¯2=C¯1Θ is an array steering vector at time t+ΔT.In summary, the whole received data can be written as:(32)Z¯t=X¯1tX¯2t=C¯1C¯2St+N1tN2t=C¯St+N¯t,whereC¯=C¯1C¯2, N¯t=N1tN2t.


## 4.2. Long Vector Ambiguity Resolution Algorithm

The plural received data correlation matrix is:(33)R¯=EZ¯Z¯H=C¯ESSHC¯H+EN¯N¯H=C¯RsC¯H+σ2I.The signal and noise subspaceUs and Un can be obtained by the complex domain eigendecomposition of the data correlation matrix R¯. The signal subspace is divided into front and back halves, and the estimation of signal frequency fk and steering vector C¯^1 can be obtained by the ESPRIT method. The process is similar to Section 3.2.where(34)C¯^1=A^xA^yA^z.


### 4.2.1. Rough Estimate of DOA

According to formulas (29) and (34), the following estimates can be obtained:
(35)Φ^x=A^Z#A^x,Φ^y=A^Z#A^y.The rough estimate value of DOA without ambiguity is given by:(36)θ~k=tan−1Φ^x2k,k+Φ^y2k,k,ϕ~k=tan−1Φ^yk,kΦ^xk,kΦ^yk,k≥0π+tan−1Φ^yk,kΦ^xk,kΦ^yk,k<0,where Φ^xk,k and Φ^yk,k are the kth row and kth column element of Φ^x. The ambiguity resolution algorithm is the same as in Section 3.2.2.


## 4.2.1. Rough Estimate of DOA

According to formulas (29) and (34), the following estimates can be obtained:
(35)Φ^x=A^Z#A^x,Φ^y=A^Z#A^y.The rough estimate value of DOA without ambiguity is given by:(36)θ~k=tan−1Φ^x2k,k+Φ^y2k,k,ϕ~k=tan−1Φ^yk,kΦ^xk,kΦ^yk,k≥0π+tan−1Φ^yk,kΦ^xk,kΦ^yk,k<0,where Φ^xk,k and Φ^yk,k are the kth row and kth column element of Φ^x. The ambiguity resolution algorithm is the same as in Section 3.2.2.


## 5. Performance Analysis of the Algorithm

### 5.1. Scope of Application Analysis

The quaternion ESPRIT parameter estimation method is convenient for estimating the roughly estimated angle of arrival of the signal. It ingeniously combines three velocity sensor subarrays into a quaternion subarray by using the rotation relationship between the quaternion subarrays. Comparing formulas (18) and (36) reveals that the quaternion method has evident advantages in calculating the rough estimation of the angle of arrival. The rotation invariant relationship between direct and delayed sampling data can be used to estimate the parameters for narrow-band far-field signals with different frequencies; in this scenario, no requirement is imposed for the array configuration. The translation invariant relationship between subarrays should be used for parameter estimation of the same-frequency signal; this approach is only suitable for uniform array and results in array aperture loss. For uniform L-array, the eigendecomposition should be performed using the translation-invariant relationship in the x- and y-axis directions to obtain the spatial steering vector in these directions.


### 5.2. Computational Complexity Analysis

The following is a comparison of the computational cost of the algorithms based on the two models. Given that the amount of algorithmic computation depends on hardware and software, we only focus on the computational complexity of the covariance matrix of the algorithm. Without loss of generality, an array ofM acoustic velocity vector sensors is considered. In data processing, N snapshots are used to estimate the data covariance matrix. The data covariance matrices based on the quaternion and plural model are expressed as follows:
(37)RS=1N∑n=1NZsnZsnH=1N∑n=1MRsn,RL=1N∑n=1NZLnZLnH=1N∑n=1NRLn.Each of the array elements is composed of three orthogonal velocity sensors along thex-, y-, and z-axis directions. Quaternion is Zsn∈HM, Rsn contains M2quaternion elements, storage matrix Rsn requires 4M2 machine memory units, and the multiplication of two quaternions is equivalent to 16 real multiplications and 12 real additions. The operation of the matrix Rsn requires 16M2 real multiplications and 12M2 real additions. Thus, 16M2N real multiplications, 4M2 real divisions, and 12M2N+N‐14M2=16M2N‐4M2 real additions are required for the operation of the matrix RS. The plural is ZLn∈C3M, RLn contains 9M2 complex units, storage matrix RLn requires 9M2 machine memory units, and the multiplication of two plurals is equivalent to 4 real multiplications and 2 real additions. The operation of the matrix RLn requires 36M2 real multiplications and 18M2 real additions. Thus, the operation of the matrix RL requires 36M2N real multiplications, 18M2real divisions, and 18M2N+N‐19M2=27M2N‐9M2 real additions. Under the same array conditions, the storage units needed for the calculation of the quaternion covariance matrix and the real number addition, multiplication, and division operations are less than those for the long-vector method.


## 5.1. Scope of Application Analysis

The quaternion ESPRIT parameter estimation method is convenient for estimating the roughly estimated angle of arrival of the signal. It ingeniously combines three velocity sensor subarrays into a quaternion subarray by using the rotation relationship between the quaternion subarrays. Comparing formulas (18) and (36) reveals that the quaternion method has evident advantages in calculating the rough estimation of the angle of arrival. The rotation invariant relationship between direct and delayed sampling data can be used to estimate the parameters for narrow-band far-field signals with different frequencies; in this scenario, no requirement is imposed for the array configuration. The translation invariant relationship between subarrays should be used for parameter estimation of the same-frequency signal; this approach is only suitable for uniform array and results in array aperture loss. For uniform L-array, the eigendecomposition should be performed using the translation-invariant relationship in the x- and y-axis directions to obtain the spatial steering vector in these directions.


## 5.2. Computational Complexity Analysis

The following is a comparison of the computational cost of the algorithms based on the two models. Given that the amount of algorithmic computation depends on hardware and software, we only focus on the computational complexity of the covariance matrix of the algorithm. Without loss of generality, an array ofM acoustic velocity vector sensors is considered. In data processing, N snapshots are used to estimate the data covariance matrix. The data covariance matrices based on the quaternion and plural model are expressed as follows:
(37)RS=1N∑n=1NZsnZsnH=1N∑n=1MRsn,RL=1N∑n=1NZLnZLnH=1N∑n=1NRLn.Each of the array elements is composed of three orthogonal velocity sensors along thex-, y-, and z-axis directions. Quaternion is Zsn∈HM, Rsn contains M2quaternion elements, storage matrix Rsn requires 4M2 machine memory units, and the multiplication of two quaternions is equivalent to 16 real multiplications and 12 real additions. The operation of the matrix Rsn requires 16M2 real multiplications and 12M2 real additions. Thus, 16M2N real multiplications, 4M2 real divisions, and 12M2N+N‐14M2=16M2N‐4M2 real additions are required for the operation of the matrix RS. The plural is ZLn∈C3M, RLn contains 9M2 complex units, storage matrix RLn requires 9M2 machine memory units, and the multiplication of two plurals is equivalent to 4 real multiplications and 2 real additions. The operation of the matrix RLn requires 36M2 real multiplications and 18M2 real additions. Thus, the operation of the matrix RL requires 36M2N real multiplications, 18M2real divisions, and 18M2N+N‐19M2=27M2N‐9M2 real additions. Under the same array conditions, the storage units needed for the calculation of the quaternion covariance matrix and the real number addition, multiplication, and division operations are less than those for the long-vector method.


## 6. Simulation Experiment and Analysis

The proposed algorithm in this paper is applicable to uniform linear array, uniform circular array, and so on. Without loss of generality, the simulation experiment takes uniform circular array as an example. Two far-field, narrow-band, non-Gaussian stationary sound signals are incident into the uniform circular acoustic vector sensor array shown in Figure1. The receiving array is composed of 9 array element s; the incident signal parameters are set as follows: θ1,ϕ1=72∘,85∘, θ2,ϕ2=30∘,40∘, digital frequencies f1′,f2′=0.75,0.45, the interelement d with ambiguity is 2.5λmin, 200 Monte Carlo experiments, and the noise is zero-mean additive white Gaussian noise.


### 6.1. Simulation Experiment

Figures2 and 3 are the scatter diagrams of the estimated results, the number of snapshots is 50, and the SNR is set to 20 dB. From the intuitive comparison of the two graphs, it can be seen that the elevation and azimuth angles of the quaternion deambiguity algorithm are more concentrated around the actual value, while the angle estimation of the long vector de-ambiguity algorithm is more scattered. It is shown that the quaternion de-ambiguity algorithm has a little improvement in parameter estimation accuracy compared with the long vector de-ambiguity algorithm, but the calculation amount and storage space of the quaternion de-ambiguity algorithm is much less than that of the long vector de-ambiguity algorithm.Figure 2
DOA scatter diagram using the quaternion method.Figure 3
DOA scatter diagram using the long vector method.Figures4 and 5 plot the root mean square error (RMSE) of elevation angle and azimuth angle, respectively estimated by long vector and quaternion de-ambiguity method, at different snapshot number levels. The SNR is set as 10 dB. It can be seen from Figures 4 and 5 that the RMSE of the elevation and azimuth angles of the two algorithms decreases with the increase of the number of snapshots. The RMSE performance of the quaternion-based de-ambiguity algorithm is better than that of the long vector de-ambiguity algorithm in the snapshot number range from 25 to 75.Figure 4
The RMSE of elevation versus the snapshot number.Figure 5
The RMSE of azimuth versus the snapshot number.Figures6 and 7 show the success probability of the elevation angle and azimuth angle, respectively, estimated by the long vector and quaternion de-ambiguity method, at different snapshot number levels. The SNR is set as 10 dB. It can be seen that the success probability of elevation and azimuth angles of the two algorithms increase with the increase of the number of snapshots. The success probability of the elevation angle and azimuth angle based on the quaternion-based de-ambiguity algorithm is better than that of the long vector de-ambiguity algorithm in the snapshot number range from 25 to 60 and from 25 to 55, respectively.Figure 6
The success probability of elevation.Figure 7
The success probability of azimuth.


### 6.2. Performance Analysis

The proposed method needs a quaternion eigendecomposition and a complex eigendecomposition. The eigenvalue of the plural eigendecomposition can be used to calculate the rough estimate of the direction of arrival of the signal, and the eigenvectors constitute the estimation of the spatial steering vector. The parameters are automatically paired, and no additional pairing operations are required. However, the long-vector method needs to divide the signal subspace into three subarrays corresponding to the steering vectors along thex, y, and z-axes. The complex eigendecomposition of the full data must be performed, and two additional complex eigendecompositions are needed. The eigenvectors of the two latter eigendecompositions are used to estimate the steering vector and the direction of arrival of the signal. The two eigendecompositions are performed separately, and additional pairing operations are required.


## 6.1. Simulation Experiment

Figures2 and 3 are the scatter diagrams of the estimated results, the number of snapshots is 50, and the SNR is set to 20 dB. From the intuitive comparison of the two graphs, it can be seen that the elevation and azimuth angles of the quaternion deambiguity algorithm are more concentrated around the actual value, while the angle estimation of the long vector de-ambiguity algorithm is more scattered. It is shown that the quaternion de-ambiguity algorithm has a little improvement in parameter estimation accuracy compared with the long vector de-ambiguity algorithm, but the calculation amount and storage space of the quaternion de-ambiguity algorithm is much less than that of the long vector de-ambiguity algorithm.Figure 2
DOA scatter diagram using the quaternion method.Figure 3
DOA scatter diagram using the long vector method.Figures4 and 5 plot the root mean square error (RMSE) of elevation angle and azimuth angle, respectively estimated by long vector and quaternion de-ambiguity method, at different snapshot number levels. The SNR is set as 10 dB. It can be seen from Figures 4 and 5 that the RMSE of the elevation and azimuth angles of the two algorithms decreases with the increase of the number of snapshots. The RMSE performance of the quaternion-based de-ambiguity algorithm is better than that of the long vector de-ambiguity algorithm in the snapshot number range from 25 to 75.Figure 4
The RMSE of elevation versus the snapshot number.Figure 5
The RMSE of azimuth versus the snapshot number.Figures6 and 7 show the success probability of the elevation angle and azimuth angle, respectively, estimated by the long vector and quaternion de-ambiguity method, at different snapshot number levels. The SNR is set as 10 dB. It can be seen that the success probability of elevation and azimuth angles of the two algorithms increase with the increase of the number of snapshots. The success probability of the elevation angle and azimuth angle based on the quaternion-based de-ambiguity algorithm is better than that of the long vector de-ambiguity algorithm in the snapshot number range from 25 to 60 and from 25 to 55, respectively.Figure 6
The success probability of elevation.Figure 7
The success probability of azimuth.


## 6.2. Performance Analysis

The proposed method needs a quaternion eigendecomposition and a complex eigendecomposition. The eigenvalue of the plural eigendecomposition can be used to calculate the rough estimate of the direction of arrival of the signal, and the eigenvectors constitute the estimation of the spatial steering vector. The parameters are automatically paired, and no additional pairing operations are required. However, the long-vector method needs to divide the signal subspace into three subarrays corresponding to the steering vectors along thex, y, and z-axes. The complex eigendecomposition of the full data must be performed, and two additional complex eigendecompositions are needed. The eigenvectors of the two latter eigendecompositions are used to estimate the steering vector and the direction of arrival of the signal. The two eigendecompositions are performed separately, and additional pairing operations are required.


## 7. Conclusion

The proposed phase ambiguity resolution method of the sparse array can effectively improve the estimation accuracy of parameters and the flexibility of array placement. A set of rough and unambiguous estimates of DOA can be given very conveniently through the ingenious combination of three acoustic velocity sensors into quaternion data. Thus, the phase ambiguity of the sparse array steering vector can be solved, and high-precision DOA estimates are herein obtained. The vector orthogonality between the quaternion data components is fully utilised, and the automatic matching of the eigenvalues and eigenvectors does not require additional matching operations. The quaternion method is less computationally complex than the long-vector method. ESPRIT of the long vector can be applied to arbitrarily distributed arrays in space and has a wide application range.



---

*Source: 2402489-2020-11-20.xml*