Skip to content

Commit 6362e1c

Browse files
robmaunderrobmaunder
committed
Added demodulation.
1 parent a50ebe1 commit 6362e1c

File tree

6 files changed

+224
-42
lines changed

6 files changed

+224
-42
lines changed

‎.gitignore‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,2 +1,3 @@
11
*.m~
2+
.DS_Store
23

‎jac.m‎

Lines changed: 4 additions & 41 deletions
Original file line numberDiff line numberDiff line change
@@ -18,47 +18,10 @@
1818

1919
function C = jac(A,B)
2020

21-
mode = 0; % Exact Jacobian logarithm
22-
% mode = 1; % Jacobian logarithm using a lookup table
23-
% mode = 2; % Approximate Jacobian logarithm
21+
C = max(A,B) + log(1+exp(-abs(A-B)));
22+
% C = max(A,B);
23+
24+
2425

25-
if(A == -inf && B == -inf)
26-
C = -inf;
27-
else
28-
if mode == 0
29-
30-
C = max(A,B) + log(1+exp(-abs(A-B)));
31-
32-
elseif mode == 1
33-
34-
difference = abs(A-B);
35-
if difference >= 4.5
36-
C = max(A,B);
37-
elseif difference >= 2.252
38-
C = max(A,B) + 0.05;
39-
elseif difference >= 1.508
40-
C = max(A,B) + 0.15;
41-
elseif difference >= 1.05
42-
C = max(A,B) + 0.25;
43-
elseif difference >= 0.71
44-
C = max(A,B) + 0.35;
45-
elseif difference >= 0.433
46-
C = max(A,B) + 0.45;
47-
elseif difference >= 0.196
48-
C = max(A,B) + 0.55;
49-
else % difference >= 0
50-
C = max(A,B) + 0.65;
51-
end
52-
53-
elseif mode == 2
54-
55-
C = max(A,B);
56-
57-
else
58-
59-
error('Invalid Jacobian mode');
60-
61-
end
62-
end
6326
end
6427

‎main_inner.m‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -22,7 +22,7 @@
2222
IA_count=11;
2323

2424
% Channel SNR in dB
25-
SNR = -6;
25+
SNR = -4;
2626

2727
% Noise variance
2828
N0 = 1/10^(SNR/10);

‎main_mod.m‎

Lines changed: 88 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,88 @@
1+
% EXIT function for a soft-input soft-output demodulator
2+
% Copyright (C) 2008 Robert G. Maunder
3+
4+
% This program is free software: you can redistribute it and/or modify it
5+
% under the terms of the GNU General Public License as published by the
6+
% Free Software Foundation, either version 3 of the License, or (at your
7+
% option) any later version.
8+
9+
% This program is distributed in the hope that it will be useful, but
10+
% WITHOUT ANY WARRANTY; without even the implied warranty of
11+
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
12+
% Public License for more details.
13+
14+
% The GNU General Public License can be seen at http://www.gnu.org/licenses/.
15+
16+
17+
18+
% Number of bits to encode
19+
bit_count=100000;
20+
21+
% Number of a priori mutual informations to consider
22+
IA_count=11;
23+
24+
% Channel SNR in dB
25+
SNR = 0;
26+
27+
% Noise variance
28+
N0 = 1/10^(SNR/10);
29+
30+
% A priori mutual informations to consider
31+
IA = 0.999*(0:1/(IA_count-1):1);
32+
33+
% Initialise results
34+
IE_av=zeros(1,IA_count);
35+
IE_hist=zeros(1,IA_count);
36+
area=0.0;
37+
38+
% Consider each a priori mutual information
39+
for IA_index = 1:IA_count
40+
41+
% Generate some random bits
42+
bits = round(rand(1,bit_count));
43+
44+
% Encode using a half-rate systematic recursive convolutional code having a single memory element
45+
tx = modulate(bits);
46+
47+
% Rayleigh fading
48+
h = sqrt(1/2)*(randn(size(tx))+1i*randn(size(tx)));
49+
50+
% Noise
51+
n = sqrt(N0/2)*(randn(size(tx))+1i*randn(size(tx)));
52+
53+
% Uncorrelated narrowband Rayleigh fading channel
54+
rx = h.*tx + n;
55+
56+
57+
% Generate the a priori LLRs having the a priori mutual information considered
58+
apriori_llrs = generate_llrs(bits, IA(IA_index));
59+
60+
% Do the BCJR
61+
extrinsic_llrs = soft_demodulate(rx, h, N0, apriori_llrs);
62+
63+
% Measure the mutual information of the extrinsic LLRs
64+
IE_hist(IA_index) = measure_mutual_information_histogram(extrinsic_llrs, bits);
65+
IE_av(IA_index) = measure_mutual_information_averaging(extrinsic_llrs);
66+
67+
% Update the area beneath the EXIT function
68+
if(IA_index > 1)
69+
area = area + (IE_av(IA_index)+IE_av(IA_index-1))*(IA(IA_index)-IA(IA_index-1))/2;
70+
end
71+
end
72+
73+
74+
% Plot EXIT function
75+
figure
76+
xlim([0 1]);
77+
ylim([0 1]);
78+
xlabel('Quality of input LLRs (a priori mutual information I_A)');
79+
ylabel('Quality of output LLRs (extrinsic mutual information I_E)');
80+
title(['EXIT function for SNR = ', num2str(SNR), ' dB']);
81+
hold on
82+
plot(IA,IE_hist,'r');
83+
plot(IA,IE_av,'b');
84+
legend({'True quality','Claimed quality'},'Location','northwest');
85+
86+
% Display the area beneath the EXIT function
87+
annotation('textbox','String',{['Area = ', num2str(area)]},'LineStyle','none','Position',[0.7 0.1 0.2 0.1]);
88+

‎modulate.m‎

Lines changed: 40 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,40 @@
1+
% QPSK modulator using natural mapping
2+
% Copyright (C) 2010 Robert G. Maunder
3+
4+
% This program is free software: you can redistribute it and/or modify it
5+
% under the terms of the GNU General Public License as published by the
6+
% Free Software Foundation, either version 3 of the License, or (at your
7+
% option) any later version.
8+
9+
% This program is distributed in the hope that it will be useful, but
10+
% WITHOUT ANY WARRANTY; without even the implied warranty of
11+
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
12+
% Public License for more details.
13+
14+
% The GNU General Public License can be seen at http://www.gnu.org/licenses/.
15+
16+
% bits is a k*n vector of bits
17+
% tx is a vector of complex symbols
18+
function tx = modulate(bits)
19+
20+
% Specify the constellation points and the bit mapping here
21+
constellation_points = [+1+1i; -1+1i; -1-1i; +1-1i]/sqrt(2);
22+
bit_labels = [0,0; 0,1; 1,1; 1,0];
23+
24+
% Determine the number of bits per symbol and the number of constellation points here
25+
k = size(bit_labels,2);
26+
M = 2^k;
27+
N = length(bits)/k;
28+
29+
30+
% Check that all the vectors and matrices have the correct dimensions
31+
if ~isequal(size(constellation_points),[M,1]) || ~isequal(size(bit_labels),[M,k])
32+
error('wrong dimensions');
33+
end
34+
35+
symbols = bin2dec(num2str(reshape(bits,[k,N])'))'+1;
36+
tx = constellation_points(symbols);
37+
38+
39+
end
40+

‎soft_demodulate.m‎

Lines changed: 90 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,90 @@
1+
% Soft QPSK demodulator using natural mapping
2+
% Copyright (C) 2010 Robert G. Maunder
3+
4+
% This program is free software: you can redistribute it and/or modify it
5+
% under the terms of the GNU General Public License as published by the
6+
% Free Software Foundation, either version 3 of the License, or (at your
7+
% option) any later version.
8+
9+
% This program is distributed in the hope that it will be useful, but
10+
% WITHOUT ANY WARRANTY; without even the implied warranty of
11+
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
12+
% Public License for more details.
13+
14+
% The GNU General Public License can be seen at http://www.gnu.org/licenses/.
15+
16+
% apriori_llrs is a 1xk vector of a priori LLRs
17+
% rx is a complex symbol
18+
% channel is a complex channel coefficient
19+
% N0 is the noise power spectral density
20+
% extrinsic_llrs is a 1xk vector of extrinsic LLRs
21+
function extrinsic_llrs = soft_demodulate(rx, channel, N0, apriori_llrs)
22+
23+
% Specify the constellation points and the bit mapping here
24+
constellation_points = [+1+1i; -1+1i; -1-1i; +1-1i]/sqrt(2);
25+
bit_labels = [0,0; 0,1; 1,0; 1,1];
26+
27+
% Determine the number of bits per symbol and the number of constellation points here
28+
k = size(bit_labels,2);
29+
M = 2^k;
30+
N = length(rx);
31+
32+
% Check that all the vectors and matrices have the correct dimensions
33+
if ~isequal(size(constellation_points),[M,1]) || ~isequal(size(bit_labels),[M,k])
34+
error('wrong dimensions');
35+
end
36+
37+
if length(channel) ~= length(rx) && length(channel) ~= 1
38+
error('wrong dimensions');
39+
end
40+
41+
42+
43+
aposteriori_symbol_LLRs = zeros(M,N);
44+
45+
% Put the influence of the received signals into the symbol LLRs
46+
for perm_index = 1:M
47+
aposteriori_symbol_LLRs(perm_index,:) = -abs(rx-channel*constellation_points(perm_index)).^2./N0;
48+
end
49+
50+
if exist('apriori_llrs','var')
51+
% Put the influence of the apriori LLRs into the symbol LLRs
52+
for bit_index = 1:k
53+
% aposteriori_symbol_LLRs(:,bit_permutations(bit_index,:) == 0) = aposteriori_symbol_LLRs(:,bit_permutations(bit_index,:) == 0) + repmat(apriori_llrs(bit_index:bits_per_symbol:end),[1,permutations/2]);
54+
55+
for perm_index = 1:M
56+
if bit_labels(perm_index,bit_index) == 0
57+
aposteriori_symbol_LLRs(perm_index,:) = aposteriori_symbol_LLRs(perm_index,:) + apriori_llrs(bit_index:k:end);
58+
end
59+
end
60+
end
61+
end
62+
63+
% Extract the aposteriori LLRs from the symbol LLRs
64+
aposteriori_llrs = zeros(1,N*k);
65+
for bit_index = 1:k
66+
p0 = -inf(1,N);
67+
p1 = -inf(1,N);
68+
69+
for perm_index = 1:M
70+
if bit_labels(perm_index,bit_index) == 0
71+
p0 = jac(p0, aposteriori_symbol_LLRs(perm_index,:));
72+
else
73+
p1 = jac(p1, aposteriori_symbol_LLRs(perm_index,:));
74+
end
75+
end
76+
77+
aposteriori_llrs(bit_index:k:end) = p0-p1;
78+
end
79+
80+
if exist('apriori_llrs','var')
81+
% Remove the apriori from the aposteriori to get the extrinsic
82+
extrinsic_llrs = aposteriori_llrs - apriori_llrs;
83+
else
84+
extrinsic_llrs = aposteriori_llrs;
85+
end
86+
87+
88+
89+
90+
end

0 commit comments

Comments
 (0)