From: drowe67 Date: Sun, 20 Apr 2014 08:28:35 +0000 (+0000) Subject: variable power simulation for energy modulation works quite well, some quality improv... X-Git-Url: http://git.whiteaudio.com/gitweb/?a=commitdiff_plain;h=4984339c13c1cc5a5e3f1ad2f4438ba60a0b8c00;p=freetel-svn-tracking.git variable power simulation for energy modulation works quite well, some quality improvement git-svn-id: https://svn.code.sf.net/p/freetel/code@1513 01035d8c-6547-0410-b346-abe4f91aad63 --- diff --git a/codec2-dev/octave/fuzzy_gray.m b/codec2-dev/octave/fuzzy_gray.m index abb9e12e..7fcca72c 100644 --- a/codec2-dev/octave/fuzzy_gray.m +++ b/codec2-dev/octave/fuzzy_gray.m @@ -2,8 +2,9 @@ % David Rowe % 10 April 2014 % -% Experiments in fuzzy gray codes. Idea is that with one bit error -% in the codeword only changes the encoded value by at most 1. +% Various experiments in fuzzy gray codes and quantising and +% transmitting scalars. Idea is that with one bit error in the +% codeword only changes the encoded value by at most 1. 1; @@ -69,11 +70,12 @@ function natural = gray_to_binary(gray) end endfunction -function sim_out = test_baseline_uncoded(Ebvec, Nbits, Ntrials) +function sim_out = test_baseline_uncoded(Ebvec, Nbits, Ntrials, enable_error_log) Nlevels = 2.^ Nbits; powersOfTwo = 2 .^ fliplr(0:(Nbits-1)); Nsymb = Nbits; sim_out.qnoise_log = zeros(length(Ebvec),Ntrials); + sim_out.error_log = []; for ne = 1:length(Ebvec) EbNodB = Ebvec(ne); @@ -103,6 +105,65 @@ function sim_out = test_baseline_uncoded(Ebvec, Nbits, Ntrials) Nerrs = sum(error_positions); Terrs += Nerrs; Tbits += length(tx_bits); + if enable_error_log + sim_out.error_log = [sim_out.error_log error_positions]; + end + + rx_index = (powersOfTwo * rx_bits'); + rx_value = unquantise_value(rx_index, 0, 1, Nlevels); + + qsignal += tx_value*tx_value; + qnoise += (tx_value - rx_value) .^ 2; + sim_out.qnoise_log(ne,nn) = tx_value - rx_value; + end + + sim_out.TERvec(ne) = Terrs; + sim_out.BERvec(ne) = Terrs/Tbits; + sim_out.QSNRvec(ne) = 10*log10(qsignal/qnoise); + printf("EbNo (dB): %3.2f Terrs: %6d BER %1.4f QSNR (dB): %3.2f\n", + EbNodB, Terrs, Terrs/Tbits, 10*log10(qsignal/qnoise)); + end + +endfunction + +function sim_out = test_varpower(Ebvec, Nbits, Ntrials, amps, enable_error_log) + Nlevels = 2.^ Nbits; powersOfTwo = 2 .^ fliplr(0:(Nbits-1)); + Nsymb = Nbits; + + sim_out.qnoise_log = zeros(length(Ebvec), Ntrials); + sim_out.error_log = []; + + for ne = 1:length(Ebvec) + EbNodB = Ebvec(ne); + EbNo = 10^(EbNodB/10); + + variance = 1/EbNo; + + Terrs = 0; Tbits = 0; + qsignal = qnoise = 0; + + for nn = 1:Ntrials + + tx_value = rand(1,1); + tx_index = quantise_value(tx_value, 0, 1, Nlevels); + tx_bits = dec2bin(tx_index, Nbits) - '0'; + tx_symbols = (-1 + 2*tx_bits) .* amps; + + % AWGN noise and phase/freq offset channel simulation + % 0.5 factor ensures var(noise) == variance , i.e. splits power between Re & Im + + noise = sqrt(variance*0.5)*(randn(1,Nsymb) + j*randn(1,Nsymb)); + rx_symbols = tx_symbols + noise; + + rx_bits = rx_symbols > 0; + + error_positions = xor(rx_bits, tx_bits); + if enable_error_log + sim_out.error_log = [sim_out.error_log error_positions]; + end + Nerrs = sum(error_positions); + Terrs += Nerrs; + Tbits += length(tx_bits); rx_index = (powersOfTwo * rx_bits'); rx_value = unquantise_value(rx_index, 0, 1, Nlevels); @@ -132,7 +193,7 @@ function valid_codewords = fuzzy_code_create(ndata,nparity) bad_distance = 0; for i=1:Nvalid for k=i+1:Nvalid - distance = sum(bitxor(valid_codewords(i,:), valid_codewords(k,:))); + distance = sum(bitxor(valid_codewords(i,:), valid_codewords(k,:))) if distance < 2 bad_distance++; end @@ -162,33 +223,47 @@ function value = fuzzy_code_decode(codewords, rx_symbols) end endfunction -function valid_codewords = fuzzy_code_create(ndata, nparity) +function set_codewords = fuzzy_code_create(ndata, nparity) Nbits = ndata + nparity; Nvalid = 2 .^ ndata; - codewords = binary_to_gray(0:(2 .^ Nbits)-1); - valid_codewords = dec2bin(codewords(1:2:(2 .^ Nbits)), Nbits) - '0'; + M = 2 .^ Nbits; - % check all valid codewords have a hamming distance of at least 2^nparity + codewords = dec2bin(0:(M-1), Nbits) - '0'; - bad_distance = 0; - for i=1:Nvalid - for k=i+1:Nvalid - distance = sum(bitxor(valid_codewords(i,:), valid_codewords(k,:))); - if distance < 2 - bad_distance++; + set_codewords(1,:) = codewords(1,:); + + for i=1:Nvalid-1 + + % measure distance from current set to all codewords + + for m=1:i + printf("%d: ", m); + for k=1:length(codewords); + d(m,k) = sum(bitxor(set_codewords(m,:), codewords(k,:))); + printf("%d ", d(m,k)); end + printf("\n"); end - end - if bad_distance != 0 - printf("Error: Nvalid: %d bad_distance: %d\n", Nvalid, bad_distance); - return; - end + % choose new codeword for set that has maximum distance from current set + + dmax = 0; best_cw = 1; + for k=1:length(codewords); + dmin = min(d(:,k)); + if dmin > dmax + dmax = dmin; + best_cw = k; + end + end + + printf("[%d %d]\n", dmax, best_cw); + set_codewords(i+1,:) = codewords(best_cw,:) + end + endfunction -function sim_out = test_fuzzy_code(Ebvec, Ndata, Ntrials) - Nparity = 1; +function sim_out = test_fuzzy_code(Ebvec, Ndata, Nparity, Ntrials) Nbits = Ndata + Nparity; Nlevels = 2 .^ Ndata; Nsymb = Nbits; @@ -232,38 +307,129 @@ function sim_out = test_fuzzy_code(Ebvec, Ndata, Ntrials) end endfunction -Ebvec = 0:7; -Ntrials = 1000; -Nbits = 5; - -baseline = test_baseline_uncoded(Ebvec, Nbits, Ntrials); -fuzzy = test_fuzzy_code(Ebvec, Nbits, Ntrials); - -figure(1); -clf; -semilogy(Ebvec, baseline.BERvec) -xlabel('Eb/N0') -ylabel('BER') -grid("minor") - -figure(2); -clf; -plot(Ebvec, baseline.QSNRvec,'b;baseline;') -hold on; -plot(Ebvec, fuzzy.QSNRvec,'r;fuzzy;') -hold off; -xlabel('Eb/N0') -ylabel('SNR') -grid("minor") - -figure(3); -subplot(211) -hist(baseline.qnoise_log(4,:),50); -subplot(212) -hist(fuzzy.qnoise_log(4,:),50); - -figure(4) -plot(baseline.qnoise_log(4,1:250),'b;baseline;') -hold on; -plot(fuzzy.qnoise_log(4,1:250),'r;fuzzy;') -hold off; +function compare_baseline_fuzzy + Ebvec = 0:7; + Ntrials = 1000; + Nbits = 3; Nparity = 1; + + baseline = test_baseline_uncoded(Ebvec, Nbits, Ntrials); + fuzzy = test_fuzzy_code(Ebvec, Nbits, Nparity, Ntrials); + + figure(1); + clf; + semilogy(Ebvec, baseline.BERvec) + xlabel('Eb/N0') + ylabel('BER') + grid("minor") + + figure(2); + clf; + plot(Ebvec, baseline.QSNRvec,'b;baseline;') + hold on; + plot(Ebvec, fuzzy.QSNRvec,'r;fuzzy;') + hold off; + xlabel('Eb/N0') + ylabel('SNR') + grid("minor") + + figure(3); + subplot(211) + hist(baseline.qnoise_log(4,:),50); + subplot(212) + hist(fuzzy.qnoise_log(4,:),50); + + figure(4) + plot(baseline.qnoise_log(1,1:250),'b;baseline;') + hold on; + plot(fuzzy.qnoise_log(1,1:250),'r;fuzzy;') + hold off; +endfunction + +function compare_baseline_varpower_plot + Ebvec = -2:5; + Ntrials = 5000; + Nbits = 5; + + baseline = test_baseline_uncoded(Ebvec, Nbits, Ntrials, 0); + amps = [2 1.5 1.0 0.5 0.5]; + av_pwr = (amps*amps')/length(amps); + amps_norm = amps/sqrt(av_pwr) + varpower = test_varpower(Ebvec, Nbits, Ntrials, amps_norm, 0); + + figure(1); + clf; + semilogy(Ebvec, baseline.BERvec) + xlabel('Eb/No (dB)') + ylabel('BER') + grid("minor") + title('BER versus Eb/No') + + figure(2); + clf; + plot(Ebvec, baseline.QSNRvec,'b;baseline;') + hold on; + plot(Ebvec, varpower.QSNRvec,'r;varpower;') + hold off; + xlabel('Eb/No (dB)') + ylabel('SNR (dB)') + grid("minor") + title('Quantiser SNR versus Eb/No') + + figure(3); + subplot(211) + hist(baseline.qnoise_log(1,:),50); + title('Baseline and Variable Power Error Histograms') + subplot(212) + hist(varpower.qnoise_log(1,:),50); + + figure(4) + subplot(211) + plot(baseline.qnoise_log(1,1:250),'b;baseline;') + title('Baseline and Variable Power Error plots for Eb/No = -2dB') + subplot(212) + plot(varpower.qnoise_log(1,1:250),'r;varpower;') +endfunction + +function compare_baseline_varpower_error_files + Ebvec = -2; + Fs = 25; % number of samples per second + Nsec = 15; % seconds to simulate + Ntrials = Fs*Nsec; + Nbits = 5; + bits_per_frame = 52; + bits_per_frame_rounded = ceil(bits_per_frame/8)*8; % c2enc uses integer number of bytes/frame + start_bit = 12; % first energy bit (after 4 voicing, 7 Wo bits) + + baseline = test_baseline_uncoded(Ebvec, Nbits, Ntrials, 1); + amps = [2 1.5 1.0 0.5 0.5]; + av_pwr = (amps*amps')/length(amps); + amps_norm = amps/sqrt(av_pwr); + varpower = test_varpower(Ebvec, Nbits, Ntrials, amps_norm, 1); + + % construct error patterns to apply to c2enc bit stream + + baseline_errors = []; + for i=1:Ntrials + error_positions = baseline.error_log(Nbits*(i-1)+1:Nbits*i); + baseline_errors = [baseline_errors zeros(1,start_bit-1) error_positions ... + zeros(1, bits_per_frame_rounded - Nbits - (start_bit-1))]; + end + + varpower_errors = []; + for i=1:Ntrials + error_positions = varpower.error_log(Nbits*(i-1)+1:Nbits*i); + varpower_errors = [varpower_errors zeros(1,start_bit-1) error_positions ... + zeros(1, bits_per_frame_rounded - Nbits - (start_bit-1))]; + end + + % save error patterns + + fep=fopen("energy_errors_baseline.bin","wb"); fwrite(fep, baseline_errors, "short"); fclose(fep); + fep=fopen("energy_errors_varpower.bin","wb"); fwrite(fep, varpower_errors, "short"); fclose(fep); +endfunction + +%compare_baseline_varpower_plot +%compare_baseline_varpower_error_files + +%compare_baseline_fuzzy +%fuzzy_code_create(3,1)