-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathfunc_theta.m
122 lines (96 loc) · 4.15 KB
/
func_theta.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
function [Theta, V] = func_theta(hRI_norm,hIT_norm,NG)
% Compute the scattering matrix Theta given the normalized channels
% hRI_norm and hIT_norm, and the group size NG
% Inputs: hRI_norm: normalized channel hRI
% hIT_norm: normalized channel hIT
% NG: group size (NG=0 means fully connected)
% Outputs: Theta: scattering matrix
NI = size(hIT_norm,1); % Number of RIS elements
if NG == 0 % Fully connected
NG = NI;
end
G = NI/NG; % Number of groups
Theta = [];
%XI = [];
if NG == 1 % Single connected
% Compute Theta
theta = - angle(hRI_norm) - angle(hIT_norm.');
Theta = diag(exp(1i * theta));
% Compute XI
%lambda = 50 * cot(theta / 2);
%XI = diag(lambda);
else % Group or fully connected
for g = 1:G
% Truncated channels
hRI_g = hRI_norm(NG*(g-1)+1:NG*g);
hIT_g = hIT_norm(NG*(g-1)+1:NG*g);
hRI_g_norm = hRI_g / norm(hRI_g);
hIT_g_norm = hIT_g / norm(hIT_g);
% Matrix A
RRI = hRI_g_norm' * hRI_g_norm;
RIT = hIT_g_norm * hIT_g_norm';
ARI = (RRI + RRI.') / 2;
AIT = (RIT + RIT.') / 2;
A = ARI - AIT;
% Alternatively:
%A = zeros(NI);
%for i = 1:NI
% for j = i:NI
% A(i,j) = abs(hIT_norm(i)) * abs(hIT_norm(j)) * cos(angle(hIT_norm(i)) - angle(hIT_norm(j))) ...
% - abs(hRI_norm(i)) * abs(hRI_norm(j)) * cos(angle(hRI_norm(i)) - angle(hRI_norm(j)));
% A(j,i) = A(i,j);
% end
%end
% Eigenvalue decomposition of A
[U,Delta] = eig(A);
delta = flip(diag(Delta)); % Order delta in decreasing order
% Compute matrix T distinguishing three cases
T = zeros(NG);
if NG == 2
T = [sqrt(1/2), sqrt(1/2);
sqrt(1/2), -sqrt(1/2)];
elseif NG == 3
T = [sqrt(-delta(3)/(delta(1)-delta(3))), sqrt(delta(1)/(2*(delta(1)-delta(3)))), -sqrt(delta(1)/(2*(delta(1)-delta(3))));
0, sqrt(1/2), sqrt(1/2);
sqrt(delta(1)/(delta(1)-delta(3))), -sqrt(-delta(3)/(2*(delta(1)-delta(3)))), sqrt(-delta(3)/(2*(delta(1)-delta(3))))];
else
T(1,1) = sqrt(-delta(NG-1) / (delta(1) - delta(NG-1)));
T(NG-1,1) = sqrt(delta(1) / (delta(1) - delta(NG-1)));
T(2,2) = sqrt(-delta(NG) / (delta(2) - delta(NG)));
T(NG,2) = sqrt(delta(2) / (delta(2) - delta(NG)));
T(1,3) = sqrt(1/2) * T(NG-1,1);
T(2,3) = sqrt(1/2) * T(NG,2);
T(NG-1,3) = -sqrt(1/2) * T(1,1);
T(NG,3) = -sqrt(1/2) * T(2,2);
T(1,4) = sqrt(1/2) * T(NG-1,1);
T(2,4) = -sqrt(1/2) * T(NG,2);
T(NG-1,4) = -sqrt(1/2) * T(1,1);
T(NG,4) = sqrt(1/2) * T(2,2);
% Alternatively
%T(1,1) = sqrt(abs(delta(NG-1)) / (abs(delta(1)) + abs(delta(NG-1))));
%T(NG-1,1) = sqrt(abs(delta(1)) / (abs(delta(1)) + abs(delta(NG-1))));
%T(2,2) = sqrt(abs(delta(NG)) / (abs(delta(2)) + abs(delta(NG))));
%T(NG,2) = sqrt(abs(delta(2)) / (abs(delta(2)) + abs(delta(NG))));
%B = null(T(:,1:2)');
%T(:,3) = sqrt(1/2) * B(:,NG-3) + sqrt(1/2) * B(:,NG-2);
%T(:,4) = sqrt(1/2) * B(:,NG-3) - sqrt(1/2) * B(:,NG-2);
T(3:NG-2,5:NG) = eye(NG-4);
end
% Compute matrices V, D, and Theta
V = flip(U,2) * T; % Order U according to delta
hRI_g_bar = hRI_g_norm * V;
hIT_g_bar = V' * hIT_g_norm;
theta = - angle(hRI_g_bar) - angle(hIT_g_bar.');
d = exp(1i * theta);
D = diag(d);
Theta_tmp = V * D * V';
% Compute XI
%lambda = 50 * cot(theta / 2);
%Lambda = diag(lambda);
%XI_tmp = V * Lambda * V';
Theta = blkdiag(Theta,Theta_tmp); % The scattering matrix is block diagonal
%XI = blkdiag(XI,XI_tmp); % The reactance matrix is block diagonal
end
end
%Theta2 = (1i*XI + 50*eye(NI)) \ (1i*XI - 50*eye(NI)); % Sanity check
end