Skip to content

Commit 9804b7e

Browse files
committed
doa-using-just-eigenvectors
1 parent 2402571 commit 9804b7e

1 file changed

Lines changed: 59 additions & 0 deletions

File tree

Lines changed: 59 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,59 @@
1+
# TLDR- works fine for 1 signal, but if there are more then the signals need to be at different amplitudes
2+
3+
import numpy as np
4+
import matplotlib.pyplot as plt
5+
6+
sample_rate = 1e6
7+
d = 0.5 # half wavelength spacing
8+
N = 10000 # number of samples to simulate
9+
t = np.arange(N)/sample_rate # time vector
10+
11+
Nr = 8 # elements
12+
theta1 = 15 / 180 * np.pi
13+
theta2 = 60 / 180 * np.pi
14+
theta3 = -50 / 180 * np.pi
15+
16+
tone1 = 1.0 * np.exp(2j*np.pi*0.0173e6*t).reshape(1,-1)
17+
tone2 = 0.5 * np.exp(2j*np.pi*0.0257e6*t).reshape(1,-1)
18+
tone3 = 0.25 * np.exp(2j*np.pi*0.0312e6*t).reshape(1,-1)
19+
20+
s1 = np.exp(-2j * np.pi * d * np.arange(Nr) * np.sin(theta1)).reshape(-1,1) # 8x1
21+
s2 = np.exp(-2j * np.pi * d * np.arange(Nr) * np.sin(theta2)).reshape(-1,1)
22+
s3 = np.exp(-2j * np.pi * d * np.arange(Nr) * np.sin(theta3)).reshape(-1,1)
23+
24+
# Simulate received signal
25+
r = s1 @ tone1 + s2 @ tone2 + s3 @ tone3
26+
n = np.random.randn(Nr, N) + 1j*np.random.randn(Nr, N)
27+
r = r + 0.00001*n # 8xN
28+
29+
# Compute covariance matrix
30+
R = r @ r.conj().T / N
31+
32+
# Eigenvalue decomposition
33+
w, V = np.linalg.eig(R)
34+
idx = np.argsort(np.abs(w))[::-1]
35+
w = w[idx]
36+
V = V[:, idx]
37+
38+
if False:
39+
# Plot eigenvalues
40+
plt.figure()
41+
plt.plot(np.sort(np.abs(w))[::-1], 'k*')
42+
plt.xlabel('Index')
43+
plt.ylabel('Eigenvalue')
44+
plt.show()
45+
46+
for i in range(3):
47+
phases = np.angle(V[:, i])
48+
# find phase between adjacent elements
49+
phase_diffs = []
50+
for i in range(len(phases)-1):
51+
phase_diffs.append(phases[i+1] - phases[i])
52+
phase_diffs = np.array(phase_diffs)
53+
phase_diffs = np.mod(phase_diffs + np.pi, 2*np.pi) - np.pi # make them all between -np.pi and np.pi
54+
print(phase_diffs)
55+
phase_diff = np.mean(phase_diffs)
56+
57+
# Convert to AoA
58+
result = np.arcsin(phase_diff / (-2*np.pi*d))
59+
print(np.rad2deg(result))

0 commit comments

Comments
 (0)