From 0b0b29bdb264cfba7b8031a23f646fa5381075f4 Mon Sep 17 00:00:00 2001 From: drowe67 Date: Sun, 30 Mar 2014 04:59:21 +0000 Subject: [PATCH] first pass at a Octave HF channelsimulator git-svn-id: https://svn.code.sf.net/p/freetel/code@1484 01035d8c-6547-0410-b346-abe4f91aad63 --- codec2-dev/octave/hf_sim.m | 74 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 74 insertions(+) create mode 100644 codec2-dev/octave/hf_sim.m diff --git a/codec2-dev/octave/hf_sim.m b/codec2-dev/octave/hf_sim.m new file mode 100644 index 00000000..90458283 --- /dev/null +++ b/codec2-dev/octave/hf_sim.m @@ -0,0 +1,74 @@ +% hf_sim.m +% David Rowe March 2014 +% +% HF channel simulation. + +function sim_out = hf_sim(sim_in, snr3kHz) + + % Init HF channel model from stored sample files of spreading signal ---------------------------------- + + % convert "spreading" samples from 1kHz carrier at Fs to complex + % baseband, generated by passing a 1kHz sine wave through PathSim + % with the ccir-poor model, enabling one path at a time. + + Fc = 1000; Fs=8000; + fspread = fopen("../raw/sine1k_2Hz_spread.raw","rb"); + spread1k = fread(fspread, "int16")/10000; + fclose(fspread); + fspread = fopen("../raw/sine1k_2ms_delay_2Hz_spread.raw","rb"); + spread1k_2ms = fread(fspread, "int16")/10000; + fclose(fspread); + + % down convert to complex baseband + spreadbb = spread1k.*exp(-j*(2*pi*Fc/Fs)*(1:length(spread1k))'); + spreadbb_2ms = spread1k_2ms.*exp(-j*(2*pi*Fc/Fs)*(1:length(spread1k_2ms))'); + + % remove -2000 Hz image + b = fir1(50, 5/Fs); + spread = filter(b,1,spreadbb); + spread_2ms = filter(b,1,spreadbb_2ms); + + % discard first 1000 samples as these were near 0, probably as + % PathSim states were ramping up + + spread = spread(1000:length(spread)); + spread_2ms = spread_2ms(1000:length(spread_2ms)); + + hf_gain = 1.0/sqrt(var(spread)+var(spread_2ms)); + + % 300 - 3000 Hz filter + + b = fir1(100,[300/4000, 3000/4000], 'pass'); + + % det power of unit variance noise passed through this filter + + filter_var = var(filter(b,1,randn(1000,1))); + + % Start simulation + + s = hilbert(filter(b,1,sim_in)); + n1 = length(s); n2 = length(spread); + n = min(n1,n2); + path1 = s(1:n) .* spread(1:n); + path2 = s(1:n) .* spread_2ms(1:n); + delay = floor(0.002*Fs); + + combined = path1(delay+1:n) + path2(1:n-delay); + + snr = 10 .^ (snr3kHz/10); + variance = (combined'*combined)/(snr*n); + noise = sqrt(variance*0.5/filter_var)*(randn(n-delay,1) + j*randn(n-delay,1)); + filtered_noise = filter(b,1,noise); + + sim_out = real(combined+filtered_noise); + printf("measured SNR: %3.2fdB\n", 10*log10(var(real(combined))/var(real(filtered_noise)))); + + figure(1); + plot(s); + figure(2) + plot(real(combined)) + figure(2) + plot(sim_out) + +endfunction + -- 2.25.1