FreeDV
/
codec2
mirror of https://github.com/drowe67/codec2.git
1
0
Fork 0
Code Issues Releases Wiki Activity
codec2/octave/ofdm_helper.m

274 lines
8.6 KiB
Matlab
Raw Blame History

% ofdm_helper.m
%
% Misc functions that are used to support OFDM modem development, that
% aren't required for modem operation
1;
%------------------------------------------------------------------------------
% print_config - utility function to use ascii-art to describe the modem frame
%------------------------------------------------------------------------------
function print_config(states)
ofdm_load_const;
% ASCII-art packet visualisation
s=1; u=1; Nuwsyms=length(uw_ind_sym);
cr = 1:Nc+2;
for f=1:Np
for r=1:Ns
for c=cr
if r == 1
if (c==1) && states.edge_pilots
sym="P";
elseif (c==Nc+1) && states.edge_pilots
sym="P";
elseif c>1 && c <=(Nc+1)
sym="P";
else
sym=" ";
end
elseif c>1 && c <=(Nc+1)
sym=".";
if (u <= Nuwsyms) && (s == uw_ind_sym(u)) sym="U"; u++; end
s++;
else
sym=" ";
end
printf("%s",sym);
end
printf("\n");
end
end
printf("Nc=%d Ts=%4.3f Tcp=%4.3f Ns: %d Np: %d\n", Nc, 1/Rs, Tcp, Ns, Np);
printf("Nsymperframe: %d Nbitsperpacket: %d Nsamperframe: %d Ntxtbits: %d Nuwbits: %d Nuwframes: %d\n",
Ns*Nc, Nbitsperpacket, Nsamperframe, Ntxtbits, Nuwbits, Nuwframes);
printf("uncoded bits/s: %4.1f Duration (incl post/preamble): %4.2f s\n",
Nbitsperpacket*Fs/(Np*Nsamperframe), (Np+2)*Ns*(Tcp+1/Rs));
end
%-----------------------------------------------------------------------
% create_ldpc_test_frame - generate a test frame of bits
%-----------------------------------------------------------------------
function [tx_bits payload_data_bits codeword] = create_ldpc_test_frame(states, coded_frame=1)
ofdm_load_const;
ldpc;
gp_interleaver;
if coded_frame
% Set up LDPC code
mod_order = 4; bps = 2; modulation = 'QPSK'; mapping = 'gray';
init_cml(); % TODO: make this path sensible and portable
load HRA_112_112.txt
[code_param framesize rate] = ldpc_init_user(HRA_112_112, modulation, mod_order, mapping);
assert(Nbitsperframe == (code_param.coded_bits_per_frame + Nuwbits + Ntxtbits));
payload_data_bits = round(ofdm_rand(code_param.data_bits_per_frame)/32767);
codeword = LdpcEncode(payload_data_bits, code_param.H_rows, code_param.P_matrix);
Nsymbolsperframe = length(codeword)/bps;
% need all these steps to get actual raw codeword bits at demod ..
tx_symbols = [];
for s=1:Nsymbolsperframe
tx_symbols = [tx_symbols qpsk_mod( codeword(2*(s-1)+1:2*s) )];
end
tx_symbols = gp_interleave(tx_symbols);
codeword_raw = [];
for s=1:Nsymbolsperframe
codeword_raw = [codeword_raw qpsk_demod(tx_symbols(s))];
end
else
codeword_raw = round(ofdm_rand(Nbitsperpacket-(Nuwbits+Ntxtbits))/32767);
end
% insert UW and txt bits
tx_bits = assemble_modem_packet(states, codeword_raw, zeros(1,Ntxtbits));
assert(Nbitsperpacket == length(tx_bits));
endfunction
% automated test
function test_assemble_disassemble(states)
ofdm_load_const;
Nsymsperpacket = Nbitsperpacket/bps;
Ndatabitsperpacket = Nbitsperpacket-(Nuwbits+Ntxtbits);
Ndatasymsperpacket = Ndatabitsperpacket/bps;
codeword_bits = round(ofdm_rand(Ndatabitsperpacket)/32767);
tx_bits = assemble_modem_packet(states, codeword_bits, zeros(1,Ntxtbits));
tx_syms = zeros(1,Nsymsperpacket);
for s=1:Nsymsperpacket
if bps == 2
tx_syms(s) = qpsk_mod(tx_bits(bps*(s-1)+1:bps*s));
elseif bps == 4
tx_syms(s) = qam16_mod(states.qam16,tx_bits(bps*(s-1)+1:bps*s));
end
end
codeword_syms = zeros(1,Ndatasymsperpacket);
for s=1:Ndatasymsperpacket
if bps == 2
codeword_syms(s) = qpsk_mod(codeword_bits(bps*(s-1)+1:bps*s));
elseif bps == 4
codeword_syms(s) = qam16_mod(states.qam16,codeword_bits(bps*(s-1)+1:bps*s));
end
end
[rx_uw rx_codeword_syms payload_amps txt_bits] = disassemble_modem_packet(states, tx_syms, ones(1,Nsymsperpacket));
assert(rx_uw == states.tx_uw);
Ndatasymsperframe = (Nbitsperpacket-(Nuwbits+Ntxtbits))/bps;
assert(codeword_syms == rx_codeword_syms);
endfunction
% test function, kind of like a CRC for QPSK symbols, to compare two vectors
function acc = test_acc(v)
sre = 0; sim = 0;
for i=1:length(v)
x = v(i);
re = round(real(x)); im = round(imag(x));
sre += re; sim += im;
%printf("%d %10f %10f %10f %10f\n", i, re, im, sre, sim);
end
acc = sre + j*sim;
end
% Save test bits frame to a text file in the form of a C array
%
% usage:
% ofdm_lib; test_bits_ofdm_file
%
function test_bits_ofdm_file
Ts = 0.018; Tcp = 0.002; Rs = 1/Ts; bps = 2; Nc = 17; Ns = 8;
states = ofdm_init(bps, Rs, Tcp, Ns, Nc);
[test_bits_ofdm payload_data_bits codeword] = create_ldpc_test_frame(states);
printf("%d test bits\n", length(test_bits_ofdm));
f=fopen("../src/test_bits_ofdm.h","wt");
fprintf(f,"/* Generated by test_bits_ofdm_file() Octave function */\n\n");
fprintf(f,"const int test_bits_ofdm[]={\n");
for m=1:length(test_bits_ofdm)-1
fprintf(f," %d,\n",test_bits_ofdm(m));
endfor
fprintf(f," %d\n};\n",test_bits_ofdm(end));
fprintf(f,"\nconst int payload_data_bits[]={\n");
for m=1:length(payload_data_bits)-1
fprintf(f," %d,\n",payload_data_bits(m));
endfor
fprintf(f," %d\n};\n",payload_data_bits(end));
fprintf(f,"\nconst int test_codeword[]={\n");
for m=1:length(codeword)-1
fprintf(f," %d,\n",codeword(m));
endfor
fprintf(f," %d\n};\n",codeword(end));
fclose(f);
endfunction
% Get rid of nasty unfiltered stuff either side of OFDM signal
% This may need to be tweaked, or better yet made a function of Nc, if Nc changes
%
% usage:
% ofdm_lib; make_ofdm_bpf(1);
function bpf_coeff = make_ofdm_bpf(write_c_header_file)
filt_n = 100;
Fs = 8000;
bpf_coeff = fir2(filt_n,[0 900 1000 2000 2100 4000]/(Fs/2),[0.001 0.001 1 1 0.001 0.001]);
if write_c_header_file
figure(1)
clf;
h = freqz(bpf_coeff,1,Fs/2);
plot(20*log10(abs(h)))
grid minor
% save coeffs to a C header file
f=fopen("../src/ofdm_bpf_coeff.h","wt");
fprintf(f,"/* 1000 - 2000 Hz FIR filter coeffs */\n");
fprintf(f,"/* Generated by make_ofdm_bpf() in ofdm_lib.m */\n");
fprintf(f,"\n#define OFDM_BPF_N %d\n\n", filt_n);
fprintf(f,"float ofdm_bpf_coeff[]={\n");
for r=1:filt_n
if r < filt_n
fprintf(f, " %f,\n", bpf_coeff(r));
else
fprintf(f, " %f\n};", bpf_coeff(r));
end
end
fclose(f);
end
endfunction
% Helper function to help design UW error thresholds, in particular for raw
% data modes. See also https://www.rowetel.com/wordpress/?p=7467
function ofdm_determine_bad_uw_errors(Nuw)
figure(1); clf;
% Ideally the 10% and 50% BER curves are a long way apart
plot(0:Nuw, binocdf(0:Nuw,Nuw,0.1),';BER=0.1;'); hold on;
plot(binocdf(0:Nuw,Nuw,0.5),';BER=0.5;');
% Suggested threshold for raw data modes is the 5% probability
% level for the 50% BER curve. The pre/post-amble has a low chance
% of failure. If it does make an error, then we will have random
% bits presented as the UW (50% BER in UW). This threshold means
% there is only a 5% case of random bits being accepted as a valid UW
bad_uw_errors = max(find(binocdf(0:Nuw,Nuw,0.5) <= 0.05))+1;
plot([bad_uw_errors bad_uw_errors],[0 1],';bad uw errors;'); hold off; grid
xlabel('bits');
printf("for Nuw = %d, suggest bad_uw_errors = %d\n", Nuw, bad_uw_errors);
end
% Returns level threshold such that threshold_cdf of the tx magnitudes are
% beneath that level. Helper function that can be used to design
% the clipper level. See also https://www.rowetel.com/?p=7596
function threshold_level = ofdm_determine_clip_threshold(tx, threshold_cdf)
Nsteps = 25;
mx = max(abs(tx));
cdf = empirical_cdf(mx*(1:Nsteps)/Nsteps,abs(tx));
threshold_level = find(cdf >= threshold_cdf)(1)*mx/25;
printf("threshold_cdf: %f threshold_level: %f\n", threshold_cdf, threshold_level);
figure(1); clf; [hh nn] = hist(abs(tx),Nsteps,1);
plotyy(nn,hh,mx*(1:Nsteps)/Nsteps,cdf); title('PDF and CDF Estimates'); grid;
end
% helper function that adds channel simulation and ensures we don't saturate int16 output samples
function [rx_real rx] = ofdm_channel(states, tx, SNR3kdB, channel, freq_offset_Hz)
[rx_real rx sigma] = channel_simulate(states.Fs, SNR3kdB, freq_offset_Hz, channel, tx, states.verbose);
% multipath models can lead to clipping of int16 samples
num_clipped = length(find(abs(rx_real>32767)));
while num_clipped/length(rx_real) > 0.001
rx_real /= 2;
num_clipped = length(find(abs(rx_real>32767)));
printf("WARNING: output samples clipped, reducing level\n")
end
endfunction
Reference in New Issue View Git Blame Copy Permalink