% This function converts the covariace matrices in a Gaussian MAC to corresponding
% covariance matrices in its Dual BC. These covariance matrices achieve the
% same set of rates in the dual BC by dirty paper coding scheme. The
% encoding order in the BC is the reverse of the decoding order in the MAC.
% The total number of users is denoted by K. The decoding order in the MAC
% is given by K by 1 vector ind. ind(k) is the user that is decoded kth in
% the successive decoding. H is a n by K by N matrix containing all the
% channel matrices of the MAC. H(:,k,i) is the channel matrix for user k on tone i in
% the MAC. n is the number transmit antennas in the BC (receiver antennas
% in the MAC). S is a n by n by K by N matrix containing the covariance matrices of the
% each user of BC on each tone. Gtot is a n by n by K matrix containing the
% total covariance matrices of the BC on each tone. This function is
% written for just one receive antenna per user for the BC and is for
% parallel BCs and MACs (multi tone).
% the inputs are:
% 1) H, an n by K by N channel matrix. H(:,k, n) is the channel
% vector for user k on tone n. This code assumes each user has just one
% transmit antenna. n is the number of receive antennas and n=K in DSL.
% 2) P, a K by N vector of powers for each user on each tone.
% 3) ind, a K by 1 vector that specified the decoding order in the MAC.
% ind(k) is the user that is decoded kth in the successive decoding.
function [S, Gtot, pBC] = mac2BcMultiTone(H, P, ind)
% The code is written for decoding order [1 2 3 ... K]
[n, K, N] = size(H);
H = H(:,ind,:);
P = P(ind,:);
Gtot = zeros(n,n,N);
S = zeros(n,n,K,N);
% after re-ordering the matrices, the decoding order is 1,2,3...,K
% Gtot is the transmit covariance matrix of the BC on each tone
pBC = zeros(K,N);
for i = 1:N
r = zeros(K,1);
for k = 1:K
r(k,1) = log2(det(eye(n) + H(:,k:K,i) * diag(P(k:K, i)) * H(:,k:K,i)') ...
/ det(eye(n) + H(:,k+1:K,i) * diag(P(k+1:K, i)) * H(:,k+1:K,i)'));
end
g = 2.^r - 1;
w = zeros(n,K);
L = zeros(K,K);
for k = 1:K
w(:,k) = inv(eye(n) + H(:,k+1:K,i) * diag(P(k+1:K, i)) * H(:,k+1:K,i)') * H(:,k,i);
w(:,k) = w(:,k) / norm(w(:,k));
L(k,k) = abs(w(:,k)' * H(:,k,i))^2/g(k);
for j = 1:k-1
L(k,j) = -abs(w(:,j)' * H(:,k,i))^2;
end
end
pBC(:,i) = L \ ones(K,1);
for k = 1:K
S(:,:,k,i) = w(:,k) * w(:,k)' * pBC(k,i);
Gtot(:,:,i) = Gtot(:,:,i) + S(:,:,k,i);
end
end
pBC(ind, :) = real(pBC);
S(:,:,ind,:) = S;