next_sync = 1;
else
next_sync = 0;
- printf(" [%d] back to coarse freq offset est...\n", cohpsk.frame) ;
+ %printf(" [%d] back to coarse freq offset est...\n", cohpsk.frame) ;
end
cohpsk.ratio = abs(max_corr/max_mag);
end
% test parameters ----------------------------------------------------------
-frames = 200;
-foff = -20;
-dfoff = 0/Fs;
+% TODO
+% [ ] set up various tests we use to characterise modem for easy
+% repeating when we change modem
+
+frames = 100;
+foff = -40;
+dfoff = -0.5/Fs;
EsNodB = 12;
fading_en = 1;
hf_delay_ms = 2;
if (sync == 0) && (f>1)
- % we are out of sync so reset f_est and process two frames to clean out memories
+ % we can test +/- 20Hz, so we break this up into 3 tests to cover +/- 60Hz
+
+ max_ratio = 0;
+ for acohpsk.f_est = Fcentre-40:40:Fcentre+40
+
+ printf(" [%d] acohpsk.f_est: %f +/- 20\n", f, acohpsk.f_est);
- acohpsk.f_est = Fcentre;
- %[rx_fdm_frame_bb afdmdv.fbb_phase_rx] = freq_shift(ch_fdm_frame_buf, -acohpsk.f_est, Fs, afdmdv.fbb_phase_rx);
- [ch_symb rx_timing rx_filt rx_baseband afdmdv f_est] = rate_Fs_rx_processing(afdmdv, ch_fdm_frame_buf, acohpsk.f_est, Nsw*acohpsk.Nsymbrowpilot, nin, 0);
- for i=1:Nsw-1
- acohpsk.ct_symb_buf = update_ct_symb_buf(acohpsk.ct_symb_buf, ch_symb((i-1)*acohpsk.Nsymbrowpilot+1:i*acohpsk.Nsymbrowpilot,:), acohpsk.Nct_sym_buf, acohpsk.Nsymbrowpilot);
+ % we are out of sync so reset f_est and process two frames to clean out memories
+
+ [ch_symb rx_timing rx_filt rx_baseband afdmdv acohpsk.f_est] = rate_Fs_rx_processing(afdmdv, ch_fdm_frame_buf, acohpsk.f_est, Nsw*acohpsk.Nsymbrowpilot, nin, 0);
+ for i=1:Nsw-1
+ acohpsk.ct_symb_buf = update_ct_symb_buf(acohpsk.ct_symb_buf, ch_symb((i-1)*acohpsk.Nsymbrowpilot+1:i*acohpsk.Nsymbrowpilot,:), acohpsk.Nct_sym_buf, acohpsk.Nsymbrowpilot);
+ end
+ [anext_sync acohpsk] = frame_sync_fine_freq_est(acohpsk, ch_symb((Nsw-1)*acohpsk.Nsymbrowpilot+1:Nsw*acohpsk.Nsymbrowpilot,:), sync, next_sync);
+
+ if anext_sync == 1
+ %printf(" [%d] acohpsk.ratio: %f\n", f, acohpsk.ratio);
+ if acohpsk.ratio > max_ratio
+ max_ratio = acohpsk.ratio;
+ f_est = acohpsk.f_est - acohpsk.f_fine_est;
+ next_sync = anext_sync;
+ end
+ end
end
- [next_sync acohpsk] = frame_sync_fine_freq_est(acohpsk, ch_symb((Nsw-1)*acohpsk.Nsymbrowpilot+1:Nsw*acohpsk.Nsymbrowpilot,:), sync, next_sync);
if next_sync == 1
% we've found a sync candidate!
% re-process last two frames with adjusted f_est then check again
- acohpsk.f_est -= acohpsk.f_fine_est;
+ acohpsk.f_est = f_est;
printf(" [%d] trying sync and f_est: %f\n", f, acohpsk.f_est);
- %[rx_fdm_frame_bb afdmdv.fbb_phase_rx] = freq_shift(ch_fdm_frame_buf, -acohpsk.f_est, Fs, afdmdv.fbb_phase_rx);
[ch_symb rx_timing rx_filt rx_baseband afdmdv f_est] = rate_Fs_rx_processing(afdmdv, ch_fdm_frame_buf, acohpsk.f_est, Nsw*acohpsk.Nsymbrowpilot, nin, 0);
for i=1:Nsw-1
acohpsk.ct_symb_buf = update_ct_symb_buf(acohpsk.ct_symb_buf, ch_symb((i-1)*acohpsk.Nsymbrowpilot+1:i*acohpsk.Nsymbrowpilot,:), acohpsk.Nct_sym_buf, acohpsk.Nsymbrowpilot);
% some other useful plots
figure(1)
+ clf;
plot(rx_symb_log*exp(j*pi/4),'+')
title('Scatter');
figure(2)
+ clf;
plot(rx_timing_log)
title('rx timing');
figure(3)
+ clf;
stem(nerr_log)
title('Bit Errors');
figure(4)
+ clf;
stem(ratio_log)
title('Sync ratio');
figure(5);
+ clf;
subplot(211)
- plot(foff_log);
+ plot(foff_log,';freq offset;');
hold on;
- plot(fest_log - Fcentre,'g');
+ plot(fest_log - Fcentre,'g;freq offset est;');
hold off;
-
- title('freq offset error');
+ title('freq offset');
+ legend("boxoff");
subplot(212)
plot(foff_log(1:length(fest_log)) - fest_log + Fcentre)
\r
next_sync = sync;\r
sync = 0;\r
+\r
if sync == 0\r
- acohpsk.f_est = Fcentre;\r
- end\r
\r
- [rx_fdm_frame_bb afdmdv.fbb_phase_rx] = freq_shift(ch_fdm_frame_buf, -acohpsk.f_est, Fs, afdmdv.fbb_phase_rx);\r
- [ch_symb rx_timing rx_filt rx_baseband afdmdv] = rate_Fs_rx_processing(afdmdv, rx_fdm_frame_bb, Nsw*acohpsk.Nsymbrowpilot, nin);\r
- ch_symb_log = [ch_symb_log; ch_symb((Nsw-1)*acohpsk.Nsymbrowpilot+1:Nsw*acohpsk.Nsymbrowpilot,:)];\r
+ % we can test +/- 20Hz, so we break this up into 3 tests to cover +/- 60Hz\r
\r
- % coarse timing (frame sync) and initial fine freq est ---------------------------------------------\r
+ max_ratio = 0;\r
+ for acohpsk.f_est = Fcentre-40:40:Fcentre+40\r
+ \r
+ printf(" [%d] acohpsk.f_est: %f +/- 20\n", f, acohpsk.f_est);\r
+ [ch_symb rx_timing rx_filt rx_baseband afdmdv acohpsk.f_est ] = rate_Fs_rx_processing(afdmdv, ch_fdm_frame_buf, acohpsk.f_est, Nsw*acohpsk.Nsymbrowpilot, nin, 0);\r
+ ch_symb_log = [ch_symb_log; ch_symb((Nsw-1)*acohpsk.Nsymbrowpilot+1:Nsw*acohpsk.Nsymbrowpilot,:)];\r
+\r
+ % coarse timing (frame sync) and initial fine freq est ---------------------------------------------\r
\r
- for i=1:Nsw-1\r
- acohpsk.ct_symb_buf = update_ct_symb_buf(acohpsk.ct_symb_buf, ch_symb((i-1)*acohpsk.Nsymbrowpilot+1:i*acohpsk.Nsymbrowpilot,:), acohpsk.Nct_sym_buf, acohpsk.Nsymbrowpilot);\r
+ for i=1:Nsw-1\r
+ acohpsk.ct_symb_buf = update_ct_symb_buf(acohpsk.ct_symb_buf, ch_symb((i-1)*acohpsk.Nsymbrowpilot+1:i*acohpsk.Nsymbrowpilot,:), acohpsk.Nct_sym_buf, acohpsk.Nsymbrowpilot);\r
+ end\r
+ [anext_sync acohpsk] = frame_sync_fine_freq_est(acohpsk, ch_symb((Nsw-1)*acohpsk.Nsymbrowpilot+1:Nsw*acohpsk.Nsymbrowpilot,:), sync, next_sync);\r
+\r
+ if anext_sync == 1\r
+ %printf(" [%d] acohpsk.ratio: %f\n", f, acohpsk.ratio);\r
+ if acohpsk.ratio > max_ratio\r
+ max_ratio = acohpsk.ratio;\r
+ f_est = acohpsk.f_est - acohpsk.f_fine_est;\r
+ next_sync = anext_sync;\r
+ end\r
+ end\r
+ end\r
end\r
- [next_sync acohpsk] = frame_sync_fine_freq_est(acohpsk, ch_symb((Nsw-1)*acohpsk.Nsymbrowpilot+1:Nsw*acohpsk.Nsymbrowpilot,:), sync, next_sync);\r
\r
% we've found a sync candidate\r
\r
\r
% rewind and re-process last few frames with f_est\r
\r
- acohpsk.f_est -= acohpsk.f_fine_est;\r
+ acohpsk.f_est = f_est;\r
printf(" [%d] trying sync and f_est: %f\n", f, acohpsk.f_est);\r
- [rx_fdm_frame_bb afdmdv.fbb_phase_rx] = freq_shift(ch_fdm_frame_buf, -acohpsk.f_est, Fs, afdmdv.fbb_phase_rx);\r
- [ch_symb rx_timing rx_filt rx_baseband afdmdv] = rate_Fs_rx_processing(afdmdv, rx_fdm_frame_bb, Nsw*acohpsk.Nsymbrowpilot, nin);\r
+ [ch_symb rx_timing rx_filt rx_baseband afdmdv acohpsk.f_est] = rate_Fs_rx_processing(afdmdv, ch_fdm_frame_buf, acohpsk.f_est, Nsw*acohpsk.Nsymbrowpilot, nin, 0);\r
for i=1:Nsw-1\r
acohpsk.ct_symb_buf = update_ct_symb_buf(acohpsk.ct_symb_buf, ch_symb((i-1)*acohpsk.Nsymbrowpilot+1:i*acohpsk.Nsymbrowpilot,:), acohpsk.Nct_sym_buf, acohpsk.Nsymbrowpilot);\r
end\r
function freq_off_est_test_single\r
sim_in.frames = 100;\r
sim_in.EsNodB = 12;\r
- sim_in.foff = -20;\r
+ sim_in.foff = -59;\r
sim_in.dfoff = 0;\r
sim_in.fading_en = 1;\r
\r
EsNodBSet = [20 12 8];\r
\r
sim_in.frames = 100;\r
- sim_in.foff = -20;\r
+ sim_in.foff = -40;\r
sim_in.dfoff = 0;\r
sim_in.fading_en = 1;\r
freq_off_log = 1E6*ones(sim_in.frames, length(EsNodBSet) );\r
endfunction\r
\r
\r
-%freq_off_est_test_single;\r
-freq_off_est_test_curves;\r
+freq_off_est_test_single;\r
+%freq_off_est_test_curves;\r
\r
% 1. start with +/- 20Hz offset\r
% 2. Measure frames to sync. How to define sync? Foff to withn 1 Hz. Sync state\r
% Octave script that to test and develop symbol rate frequency offset
% tracking
%
-% [ ] add Es/No noise and measure BER
-% [ ] QPSK mod of random symbols
-% [ ] tracking 1 Hz/s freq drift with negl impl loss
-% [ ] simulate delay through filter, or implement R/N filter
+% [X] add Es/No noise and measure BER
+% [X] QPSK mod of random symbols
+% [X] tracking 1 Hz/s freq drift shiel staying in lock
+% [X] simulate R/N filter delay through filter
rand('state',1);
randn('state',1);
graphics_toolkit ("gnuplot");
-Fs = 8000;
-Rs = 50;
-Nbits = 5000;
-Nsymb = Nbits/2;
-foff = 1;
-dfoff = -1/Fs;
-M = Fs/Rs;
+Fs = 8000; % sample rate
+Rs = 50; % QPSK modem symbol rate
+Nbits = 5000; % number of bits to simulute over
+Nsymb = Nbits/2;
+foff = 1; % initial freq offset
+dfoff = -1/Fs; % rate of change of frequency offset in Hz/sample
+M = Fs/Rs; % oversampling rate
-phase_ch = 1;
-prev_rx_symb = 1;
-
-EsNodB = 8;
+EsNodB = 8; % noise level
EsNo = 10^(EsNodB/10);
variance = 1/EsNo;
-Nsym = 5;
+Nsym = 5; % simulated delay of root nyquist filter in symbols
+
+% control system constants, these can be tweaked
beta = 0.005;
-g = 0.2;
+g = 0.2;
+
+% init state variables
+
+phase_ch = 1;
+prev_rx_symb = 1;
ftrack = 0;
filt = 0;
-% Second order IIR filter coeffs, derived with pencil and paper
+% QPSK symbols ------> x -------> x ----------G(z) ------------\
+% | | |
+% exp(j*w_off*n) | |
+% | |
+% exp(-j*w_track*n) |
+% \____________________________/
+%
+% for small phase/freq errors, where phase/freq detector is linear:
+%
+% -1 beta 1 1
+% G(z) = 1 - z g --------------- -------- -------
+% -1 -1 -1
+% 1 - (1-beta)z 1 - z 1 - z
+%
+% (1) (2) (3) (4)
+%
+% (1) - differentiator that gives us an error based on frequency difference
+% (2) - first order IIR loop filter
+% (3) - integrator I added so loop could track ramp inputs (i.e. freq drift)
+% (4) - frequency to phase integrator implicit in VCO
+%
+% Note (1) and (3) cancel out, H(z) = 1, so we get a 2nd order loop. We need to keep (1)
+% in the time domain implementation to make the modulation stripper work, it likes angles
+% near 0.
+%
+% Transfer function, derived with pencil and paper:
+%
+% G(z) g*beta
+% ------------ = ---------- * ------------------------------------
+% 1 + G(z)H(z) 1 + g*beta 2 - beta -1 1 - beta -2
+% 1 - ---------- z + ----------- z
+% 1 + g*beta 1 + g*beta
+
+% Second order IIR filter coeffs
a1 = (2-beta)/(1+g*beta);
a2 = (1-beta)/(1+g*beta);
printf("2nd order loop w: %f gamma: %f\n", w, gamma);
+
function symbol = qpsk_mod(two_bits)
two_bits_decimal = sum(two_bits .* [2 1]);
switch(two_bits_decimal)
endswitch
endfunction
-% Gray coded QPSK demodulation function
-
-function two_bits = qpsk_demod(symbol)
- if isscalar(symbol) == 0
- printf("only works with scalars\n");
- return;
- end
- bit0 = real(symbol*exp(j*pi/4)) < 0;
- bit1 = imag(symbol*exp(j*pi/4)) < 0;
- two_bits = [bit1 bit0];
-endfunction
% modulator
tx_symb(i) = qpsk_mod(tx_bits(2*(i-1)+1:2*i));
end
-% run symbol by symbol
+% run rest of simulation symbol by symbol
rx_symb = zeros(1, Nsymb);
rx_filt = zeros(1,Nsym);
mod_strip = diff.^4;
mod_strip = abs(real(mod_strip)) + j*imag(mod_strip);
- % loop filter made up of 1st order IIR plus integrator. Integerator
- % was found to be reqd
+ % loop filter made up of 1st order IIR plus integrator plus gain
+ % term. Integrator was found to be reqd to track ramp inputs
+ % (i.e. freq drift)
+
filt = ((1-beta)*filt + beta*angle(mod_strip));
ftrack += g*filt;
projects / freetel-svn-tracking.git / commitdiff