-% ldpc.m\r
+% ldpc_short.m\r
%\r
% David Rowe Mar 2017\r
% Based on Bill Cowley's LDPC simulations\r
\r
function sim_out = run_sim(sim_in, HRA, Ntrials)\r
\r
+ % Note this is effective Eb/No of payload data bits, sorta thing we\r
+ % plot on BER versus Eb/No graphs of decoded data. So if we have a\r
+ % rate 1/2 code, each codeword bit will have Eb/No - 3dB.\r
+\r
+ EbNovec = sim_in.EbNovec;\r
+\r
genie_Es = sim_in.genie_Es;\r
packet_size = sim_in.packet_size;\r
code = sim_in.code;\r
hf_en = sim_in.hf_en;\r
- Esvec = sim_in.Esvec;\r
+ verbose = sim_in.verbose;\r
\r
if strcmp(code, 'golay')\r
rate = 0.5;\r
\r
if hf_en\r
\r
- % some typical values\r
+ % some typical values, assume symbols spread acroos Nc OFDM carriers\r
\r
- Rs = 50; dopplerSpreadHz = 1.0; path_delay = 1E-3*Rs;\r
+ Nc = 28;\r
+ Rs = 50; dopplerSpreadHz = 0.5; path_delay = 2E-3*Rs;\r
\r
- nsymb = Ntrials*framesize*bps;\r
- spread1 = doppler_spread(dopplerSpreadHz, Rs, nsymb);\r
- spread2 = doppler_spread(dopplerSpreadHz, Rs, nsymb);\r
+ nsymb = Ntrials*framesize;\r
+ Tp = (framesize/Nc)/Rs;\r
+ spread1 = doppler_spread(dopplerSpreadHz, Rs, 1.1*ceil(nsymb/Nc));\r
+ spread2 = doppler_spread(dopplerSpreadHz, Rs, 1.1*ceil(nsymb/Nc));\r
hf_gain = 1.0/sqrt(var(spread1)+var(spread2));\r
% printf("nsymb: %d lspread1: %d\n", nsymb, length(spread1));\r
- end\r
+ hf_model = zeros(ceil(nsymb/Nc),Nc);\r
+ end\r
\r
- % loop around each Esvec point\r
+ % loop around each EbNovec point\r
\r
- for ne = 1:length(Esvec)\r
- Es = Esvec(ne);\r
- EsNo = 10^(Es/10);\r
- EbNodB = Es - 10*log10(code_param.bits_per_symbol * rate);\r
+ for ne = 1:length(EbNovec)\r
+\r
+ % Given Eb/No of payload data bits, work out Es/No we need to\r
+ % apply to each codeword symbol, e.g. for rate 1/2 code:\r
+ %\r
+ % Es/No = Eb/No - 3 dB\r
\r
- Terrs = 0; Tbits = 0; Ferrs = 0; \r
+ EbNodB = EbNovec(ne);\r
+ EsNodB = EbNodB + 10*log10(code_param.bits_per_symbol*rate);\r
+ EsNo = 10^(EsNodB/10);\r
+ hf_r = 1;\r
+\r
+ Terrs = 0; Tbits = 0; Ferrs = 0; Terrs_raw = Tbits_raw = 0;\r
tx_bits = rx_bits = [];\r
- hfi = 1;\r
+ r_log = []; Nerrs_log = []; Nerrs_raw_log = [];\r
\r
for nn = 1: Ntrials \r
data = round( rand( 1, code_param.data_bits_per_frame ) );\r
% modulate\r
\r
s = 1 - 2 * codeword; \r
- code_param.symbols_per_frame = length( s );\r
+ Ns = code_param.symbols_per_frame = length(s);\r
\r
if hf_en\r
\r
- % simplified rate Rs simulation model that doesn't include\r
- % ISI, just freq filtering. We assume perfect phase estimation\r
- % so it's just amplitude distortion.\r
-\r
- for i=1:length(s)\r
-\r
- % 1.5*Rs carrier spacing, symbols mapped to 14 carriers\r
- % OK we'd probably use QPSK in practice but meh a few approximations....\r
-\r
- w = 1.5*mod(i,14)*2*pi; \r
- hf_model(i) = hf_gain*(spread1(hfi) + exp(-j*w*path_delay)*spread2(hfi));\r
- s(i) *= abs(hf_model(i));\r
- hfi++;\r
+ % Simplified rate Rs OFDM simulation model that doesn't\r
+ % include ISI, just freq filtering. We assume perfect phase\r
+ % estimation so it's just amplitude distortion. We distribute\r
+ % symbols across Nc carriers\r
+\r
+ Nr = ceil(length(s)/Nc);\r
+ w = (1:Nc)*2*pi; \r
+ s = [s zeros(1,Nr*Nc-Ns)]; % pad out to integer number of rows\r
+\r
+ for r=1:Nr\r
+ for c=1:Nc\r
+ hf_model(hf_r,c) = hf_gain*(spread1(hf_r) + exp(-j*w(c)*path_delay)*spread2(hf_r));\r
+ s(Nc*(r-1)+c) *= abs(hf_model(hf_r,c));\r
+ end\r
+ hf_r++;\r
end\r
+ s = s(1:Ns); % remove padding\r
end\r
-\r
- variance = 1/(2*EsNo);\r
- noise = sqrt(variance)* randn(1,code_param.symbols_per_frame); \r
+ \r
+ variance = 1/(EsNo);\r
+ noise = sqrt(variance/2)*(randn(1,Ns) + j*randn(1,Ns));\r
r = s + noise;\r
\r
Nr = length(r); \r
+ r_log = [r_log r];\r
+\r
+ % raw BER\r
+\r
+ detected_data = real(r) < 0;\r
+ error_positions = xor(detected_data, codeword);\r
+\r
+ Nerrs_raw = sum(error_positions);\r
+ Terrs_raw += Nerrs_raw;\r
+ Tbits_raw += length(codeword);\r
+ Nerrs_raw_log = [Nerrs_raw_log Nerrs_raw];\r
+\r
+ % other decoders\r
\r
if strcmp(code, 'golay')\r
- detected_data = egolaydec(r < 0);\r
+ detected_data = egolaydec(real(r) < 0);\r
detected_data = detected_data(code_param.data_bits_per_frame+1:framesize);\r
end\r
\r
end\r
\r
if strcmp(code, 'diversity')\r
- detected_data = (r(1:7) + r(8:14)) < 0;\r
+ detected_data = real(r(1:7) + r(8:14)) < 0;\r
end\r
\r
rx_bits = [rx_bits detected_data];\r
\r
error_positions = xor( detected_data, data );\r
Nerrs = sum(error_positions);\r
- \r
+ Nerrs_log = [Nerrs_log Nerrs];\r
+ \r
if Nerrs>0, Ferrs = Ferrs +1; end\r
Terrs = Terrs + Nerrs;\r
Tbits = Tbits + code_param.data_bits_per_frame;\r
-\r
end\r
\r
% count packet errors using supplied packet size. Operating one\r
Tpackets++;\r
end\r
\r
- printf("EbNo: %3.1f dB BER: %5.4f PER: %5.4f Nbits: %4d Nerrs: %4d Tpackets: %4d Perr: %4d\n", EbNodB, Terrs/Tbits, Perrs/ Tpackets, Tbits, Terrs, Tpackets, Perrs);\r
+ printf("Coded EbNo: %3.1f dB BER: %5.4f PER: %5.4f Nbits: %4d Nerrs: %4d Tpackets: %4d Perr: %4d\n", \r
+ EbNodB, Terrs/Tbits, Perrs/ Tpackets, Tbits, Terrs, Tpackets, Perrs);\r
+ EbNodB_raw = EsNodB - 10*log10(code_param.bits_per_symbol);\r
+ printf("Raw EbNo..: %3.1f dB BER: %5.4f Nbits: %4d Nerrs: %4d\n", EbNodB_raw, \r
+ Terrs_raw/Tbits_raw, Tbits_raw, Terrs_raw);\r
\r
TERvec(ne) = Terrs;\r
FERvec(ne) = Ferrs;\r
BERvec(ne) = Terrs/ Tbits;\r
PERvec(ne) = Perrs/ Tpackets;\r
- Ebvec = Esvec - 10*log10(code_param.bits_per_symbol * rate);\r
- \r
+\r
sim_out.BERvec = BERvec;\r
sim_out.PERvec = PERvec;\r
- sim_out.Ebvec = Ebvec;\r
sim_out.FERvec = FERvec;\r
sim_out.TERvec = TERvec;\r
sim_out.framesize = framesize;\r
sim_out.error_positions = error_positions;\r
+\r
+ if verbose\r
+ figure(3); clf;\r
+ plot(real(r_log),imag(r_log),'+')\r
+ axis([-2 2 -2 2])\r
+ title('Scatter');\r
+\r
+ figure(4)\r
+\r
+ subplot(211);\r
+ stem((1:Ntrials)*Tp, Nerrs_raw_log);\r
+ subplot(212);\r
+ stem((1:Ntrials)*Tp, Nerrs_log);\r
+ %axis([1 Ntrials*Tp -0.2 1.2])\r
+\r
+ if hf_en\r
+ figure(5); clf;\r
+\r
+ % limit mesh plot to Np points to plot quickly\r
+ \r
+ Np = 500;\r
+ step = ceil(hf_r/Np);\r
+ mesh(1:Nc, (1:step:hf_r-1)/Rs, abs(hf_model(1:step:hf_r-1,:)))\r
+ title('HF channel amplitude');\r
+ xlabel('Carrier');\r
+ ylabel('Time (s)');\r
+ end\r
+\r
+ end\r
end\r
endfunction\r
\r
sim_in.packet_size = 28;\r
sim_in.code = 'ldpc';\r
sim_in.hf_en = hf_en;\r
+ sim_in.verbose = 0;\r
\r
if hf_en\r
- Esvec = 0:0.5:9; \r
+ EbNovec = 2:0.5:10; \r
else\r
- Esvec = -3:0.5:3; \r
+ EbNovec = 0:0.5:8; \r
end\r
- sim_in.Esvec = Esvec;\r
+ sim_in.EbNovec = EbNovec;\r
\r
load HRA_112_112.txt\r
load HRA_112_56.txt\r
load HRA_56_56.txt\r
load HRA_56_28.txt\r
\r
+ printf("HRA_112_112-------------\n");\r
sim_out1 = run_sim(sim_in, HRA_112_112, Ntrials);\r
+ printf("HRA_112_56-------------\n");\r
sim_out2 = run_sim(sim_in, HRA_112_56 , Ntrials);\r
+#{\r
+ printf("HRA_56_56-------------\n");\r
sim_out3 = run_sim(sim_in, HRA_56_56 , Ntrials*2);\r
+ printf("HRA_56_28-------------\n");\r
sim_out4 = run_sim(sim_in, HRA_56_28 , Ntrials*2);\r
- sim_in.code = 'golay';\r
- sim_out5 = run_sim(sim_in, [], Ntrials*10);\r
- sim_in.code = 'diversity';\r
- sim_out6 = run_sim(sim_in, [], Ntrials*10);\r
+ %printf("Golay -------------\n");\r
+ %sim_in.code = 'golay';\r
+ %sim_out5 = run_sim(sim_in, [], Ntrials*10);\r
+ %printf("Diversity -------------\n");\r
+ %sim_in.code = 'diversity';\r
+ %sim_out6 = run_sim(sim_in, [], Ntrials*10);\r
+#}\r
\r
if hf_en\r
Ebvec_theory = 2:0.5:12;\r
figure(1); clf;\r
semilogy(Ebvec_theory, uncoded_BER_theory, 'b+-;BPSK theory;','markersize', 10, 'linewidth', 2)\r
hold on;\r
- semilogy(sim_out1.Ebvec, sim_out1.BERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out2.Ebvec, sim_out2.BERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out3.Ebvec, sim_out3.BERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out4.Ebvec, sim_out4.BERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out5.Ebvec, sim_out5.BERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out6.Ebvec, sim_out6.BERvec, 'go-;rate 1/2 Diversity;','markersize', 10, 'linewidth', 2)\r
+ semilogy(EbNovec, sim_out1.BERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
+ semilogy(EbNovec, sim_out2.BERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out3.BERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out4.BERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out5.BERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out6.BERvec, 'go-;rate 1/2 Diversity;','markersize', 10, 'linewidth', 2)\r
hold off;\r
xlabel('Eb/No')\r
ylabel('BER')\r
figure(2); clf;\r
semilogy(Ebvec_theory, uncoded_PER_theory, 'b+-;BPSK theory;','markersize', 10, 'linewidth', 2)\r
hold on;\r
- semilogy(sim_out1.Ebvec, sim_out1.PERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out2.Ebvec, sim_out2.PERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out3.Ebvec, sim_out3.PERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out4.Ebvec, sim_out4.PERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out5.Ebvec, sim_out5.PERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
- semilogy(sim_out6.Ebvec, sim_out6.PERvec, 'go-;rate 1/2 Diversity;','markersize', 10, 'linewidth', 2)\r
+ semilogy(EbNovec, sim_out1.PERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
+ semilogy(EbNovec, sim_out2.PERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out3.PERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out4.PERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out5.PERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
+ %semilogy(EbNovec, sim_out6.PERvec, 'go-;rate 1/2 Diversity;','markersize', 10, 'linewidth', 2)\r
hold off;\r
xlabel('Eb/No')\r
ylabel('PER')\r
\r
if strcmp(channel, 'awgn')\r
sim_in.hf_en = 0;\r
- sim_in.Esvec = EbNodB + 10*log10(rate);\r
+ sim_in.EbNovec = EbNodB;\r
else\r
sim_in.hf_en = 1;\r
- sim_in.Esvec = EbNodB + 10*log10(rate);\r
+ sim_in.EbNovec = EbNodB;\r
end\r
\r
+ sim_in.verbose = 1;\r
sim_out = run_sim(sim_in, HRA_112_112, Ntrials);\r
\r
if nargin == 5\r
fclose(fep);\r
end\r
\r
- figure\r
- plot(sim_out.error_positions);\r
- axis([1 length(sim_out.error_positions) -0.2 1.2])\r
endfunction\r
\r
\r
Ntrials = bits/data_bits_per_frame;\r
sim_in.genie_Es = 1;\r
sim_in.packet_size = 28;\r
- EbNodB = 2;\r
sim_in.hf_en = 0;\r
- sim_in.Esvec = EbNodB + 10*log10(rate);\r
+ sim_in.Esvec = 2;\r
sim_in.c_include_file = "../src/HRA_112_112.h";\r
\r
sim_out = run_sim(sim_in, HRA_112_112, Ntrials);\r
randn('seed',1);\r
more off;\r
format;\r
-close all;\r
+\r
+% simple single point test\r
+\r
+run_single(700*120, 'ldpc', 'hf', 12)\r
\r
% plotting curves (may take a while)\r
\r
\r
% generating error files \r
\r
+#{\r
run_single(700*10, 'ldpc', 'awgn', 2, 'awgn_2dB_ldpc.err')\r
run_single(700*10, 'diversity', 'awgn', 2, 'awgn_2dB_diversity.err')\r
run_single(700*10, 'ldpc', 'hf', 6, 'hf_6dB_ldpc.err')\r
run_single(700*10, 'diversity', 'hf', 6, 'hf_6dB_diversity.err')\r
+#}\r
\r
% generate C header for C port of code\r
\r
projects / freetel-svn-tracking.git / commitdiff