From: baobrien Date: Tue, 13 Feb 2018 08:05:15 +0000 (+0000) Subject: Testing and tuning C coarse ofdm acq X-Git-Url: http://git.whiteaudio.com/gitweb/?a=commitdiff_plain;h=6dbcc5293cf5fc3a1b0630d243056ec282ccd5f4;p=freetel-svn-tracking.git Testing and tuning C coarse ofdm acq git-svn-id: https://svn.code.sf.net/p/freetel/code@3400 01035d8c-6547-0410-b346-abe4f91aad63 --- diff --git a/codec2-dev/octave/tofdm.m b/codec2-dev/octave/tofdm.m index 2302958a..0d64e22a 100644 --- a/codec2-dev/octave/tofdm.m +++ b/codec2-dev/octave/tofdm.m @@ -4,9 +4,13 @@ % Octave script for comparing Octave and C versions of OFDZM modem % ------------------------------------------------------------------ -1; + +load rand_state_vec.mat +%rand("state",rand_state) function rx = ebno_awgn_channel(tx, Fs, Rb, EbNodB) + %rand("seed",2) + %randn("seed",3); EbNo = 10^(EbNodB/10); variance = Fs / (Rb*EbNo); noise = sqrt(variance)*randn(1,length(tx)); @@ -15,6 +19,23 @@ function rx = ebno_awgn_channel(tx, Fs, Rb, EbNodB) %rx = rx * .10; end +function rx = snr_awgn_channel(tx, Fs, BwS, snrdb) + snr = 10^(snrdb/10); + %signal power + psig = sum(abs(tx))^2 / length(tx); + %noise power in interesting band + pnoiseb = (1/snr) * psig; + %noise power in 1hz band + pnoise0 = pnoiseb/BwS; + %noise power over sample rate bandwidth + pnoiset = pnoise0 * (Fs/2); + + noise = sqrt(pnoiset)*randn(1,length(tx)); + + rx = noise + tx; + %rx = rx * .10; +end + function nums = im_re_interleave(nim) nums = zeros(1,length(nim)*2); nums(1:2:length(nums)) = real(nim); @@ -36,6 +57,8 @@ function [a b] = trim2(ain, bin) endfunction function pass = run_ofdm_test(Nframes,sample_clock_offset_ppm,foff_hz,EbNodB = 100,nheadsamp = 0) + rand("seed",2) + randn("seed",2) %Nframes = 30; %sample_clock_offset_ppm = 100; %foff_hz = 5; @@ -53,7 +76,6 @@ function pass = run_ofdm_test(Nframes,sample_clock_offset_ppm,foff_hz,EbNodB = 1 states = ofdm_init(bps, Rs, Tcp, Ns, Nc); ofdm_load_const; - rand('seed',2); tx_bits = round(rand(1,Nbitsperframe)); % Run tx loop @@ -121,7 +143,6 @@ function pass = run_ofdm_test(Nframes,sample_clock_offset_ppm,foff_hz,EbNodB = 1 end [rx_bits states aphase_est_pilot_log arx_np arx_amp] = ofdm_demod2(states, rxbuf_in); - % log some states for comparison to C rxbuf_in_log = [rxbuf_in_log rxbuf_in]; @@ -256,7 +277,189 @@ function pass = run_ofdm_test(Nframes,sample_clock_offset_ppm,foff_hz,EbNodB = 1 end -% This works best on my machine -- Brady -graphics_toolkit('fltk') +function [ber berc] = run_ber_test(Nframes,sample_clock_offset_ppm,foff_hz,EbNodB = 100,nheadsamp = 0) + rand("seed",2) + randn("seed",2) + %Nframes = 30; + %sample_clock_offset_ppm = 100; + %foff_hz = 5; + + more off; format; + ofdm_lib; + autotest; + + % --------------------------------------------------------------------- + % Run Octave version + % --------------------------------------------------------------------- + + Ts = 0.018; Tcp = 0.002; Rs = 1/Ts; bps = 2; Nc = 16; Ns = 8; + Rb_real = 1400.0; + states = ofdm_init(bps, Rs, Tcp, Ns, Nc); + ofdm_load_const; + + tx_bits = round(rand(1,Nbitsperframe)); + + % Run tx loop + + tx_bits_log = []; tx_log = []; + for f=1:Nframes + tx_bits_log = [tx_bits_log tx_bits]; + tx_log = [tx_log ofdm_mod(states, tx_bits)]; + end + + % Channel simulation ---------------------------------------------- + + rx_log = sample_clock_offset(tx_log, sample_clock_offset_ppm); + rx_log = freq_shift(rx_log, foff_hz, Fs); + #rx_log(1:nheadsamp) = zeros(nheadsamp,1); + rx_log_l = length(rx_log) + rx_log = [zeros(1,nheadsamp), rx_log]; + rx_log = rx_log(1:rx_log_l); + %rx_log = snr_awgn_channel(rx_log, Fs, Fs, EbNodB); + rx_log = ebno_awgn_channel(rx_log, Fs, Rb_real, EbNodB); + + rx_vec = fopen("tofdm_rx_vec","wb+"); + fwrite(rx_vec,im_re_interleave(rx_log),"single"); + fclose(rx_vec); + + % Rx --------------------------------------------------------------- + + % Init rx with ideal timing so we can test with timing estimation disabled + + Nsam = length(rx_log); + prx = 1; + nin = Nsamperframe+2*(M+Ncp); + %states.rxbuf(Nrxbuf-nin+1:Nrxbuf) = rx_log(prx:nin); + %prx += nin; + + rxbuf_log = []; rxbuf_in_log = []; rx_sym_log = []; foff_hz_log = []; + timing_est_log = []; sample_point_log = []; + phase_est_pilot_log = []; rx_amp_log = []; + rx_np_log = []; rx_bits_log = []; + + states.timing_en = 1; + states.foff_est_en = 1; + states.phase_est_en = 1; + + if states.timing_en == 0 + % manually set ideal timing instant + states.sample_point = Ncp; + end + + have_frame = true; + while have_frame + % insert samples at end of buffer, set to zero if no samples + % available to disable phase estimation on future pilots on last + % frame of simulation + + nin = states.nin; + lnew = min(Nsam-prx+1,nin); + rxbuf_in = zeros(1,nin); + %printf("nin: %d prx: %d lnew: %d\n", nin, prx, lnew); + if lnew + rxbuf_in(1:lnew) = rx_log(prx:prx+lnew-1); + end + prx += nin; + if prx > Nsam + have_frame = false; + end + + [rx_bits states aphase_est_pilot_log arx_np arx_amp] = ofdm_demod(states, rxbuf_in); + % log some states for comparison to C + + rxbuf_in_log = [rxbuf_in_log rxbuf_in]; + rxbuf_log = [rxbuf_log states.rxbuf]; + rx_sym_log = [rx_sym_log; states.rx_sym]; + phase_est_pilot_log = [phase_est_pilot_log; aphase_est_pilot_log]; + rx_amp_log = [rx_amp_log arx_amp]; + foff_hz_log = [foff_hz_log; states.foff_est_hz]; + timing_est_log = [timing_est_log; states.timing_est]; + sample_point_log = [sample_point_log; states.sample_point]; + rx_np_log = [rx_np_log arx_np]; + rx_bits_log = [rx_bits_log rx_bits]; + + end + + % --------------------------------------------------------------------- + % Run C version and plot Octave and C states and differences + % --------------------------------------------------------------------- + + % Override default path by setting path_to_tofdm = "/your/path/to/tofdm" + + if exist("path_to_tofdm", "var") == 0 + path_to_tofdm = "../build_linux/unittest/tofdm"; + end + system(path_to_tofdm); + + load tofdm_out.txt; + % Generated with modem probe thing + load ofdm_test.txt; + + front_cutoff = 10000; + rx_np_log = rx_np_log(front_cutoff/2:length(rx_np_log)); + rx_np_log_c = rx_np_log_c(front_cutoff/2:length(rx_np_log_c)); + + fg = 1; + figure(fg++); clf; plot(foff_hz_log_c); title('foff hz C'); + figure(fg++); clf; plot(rx_np_log,'+'); title('Octave Scatter Diagram'); axis([-1.5 1.5 -1.5 1.5]); + figure(fg++); clf; plot(rx_np_log_c,'+'); title('C Scatter Diagram'); axis([-1.5 1.5 -1.5 1.5]); + + + [rx_bits_log rx_bits_log_c] = trim1(rx_bits_log,rx_bits_log_c); + [tx_bits_log rx_bits_log] = trim1(tx_bits_log,rx_bits_log); + [tx_bits_log rx_bits_log_c] = trim1(tx_bits_log,rx_bits_log_c); + + bitcnt = length(rx_bits_log_c); + rx_bits = rx_bits_log_c(front_cutoff:length(rx_bits_log_c)); + tx_bits = tx_bits_log(front_cutoff:length(tx_bits_log)); + + ber = 1; + for offset = (1:100) + nerr = sum(xor(rx_bits(offset:length(rx_bits)),tx_bits(1:length(rx_bits)+1-offset))); + bern = nerr/(bitcnt-offset); + if(bern < ber) + ox = offset; + best_nerr = nerr; + end + ber = min([ber bern]); + end + offset = ox; + berc = ber; + + figure(fg++) + plot(xor(rx_bits(offset:length(rx_bits)),tx_bits(1:length(rx_bits)+1-offset))) + title('C bit errors') + + rx_bits = rx_bits_log(front_cutoff:length(rx_bits)); + ber = 1; + + for offset = (1:100) + nerr = sum(xor(rx_bits(offset:length(rx_bits)),tx_bits(1:length(rx_bits)+1-offset))); + bern = nerr/(bitcnt-offset); + if(bern < ber) + ox = offset; + best_nerr = nerr; + end + ber = min([ber bern]); + end + offset = ox; + ber = ber; + + figure(fg++) + plot(xor(rx_bits(offset:length(rx_bits)),tx_bits(1:length(rx_bits)+1-offset))) + title('Octave bit errors') + + + figure(fg++) + plot(20*log10(abs(fft(rx_log))), (1:length(rx_log)) .* (Fs/length(rx_log) )) + +endfunction + +% FLTK plots work better than gnuplot or the default on my machine -- Brady +try + graphics_toolkit('fltk') +catch err + warning(err.identifier, err.message); +end_try_catch -run_ofdm_test(60,100,5,8,10000) +%run_ofdm_test(60,100,60,5,10000) diff --git a/codec2-dev/src/ofdm.c b/codec2-dev/src/ofdm.c index f3bd188d..9229b176 100644 --- a/codec2-dev/src/ofdm.c +++ b/codec2-dev/src/ofdm.c @@ -307,7 +307,7 @@ static int coarse_sync(struct OFDM *ofdm, complex float *rx, int length, float * float max_mag = 0; float mag = 0.0f; int t_est = 0; - int fmax = 60; + int fmax = 20; int pmax_i,nmax_i; float pmax,nmax; float foff_est; @@ -333,6 +333,10 @@ static int coarse_sync(struct OFDM *ofdm, complex float *rx, int length, float * } } + if(ratio_out != NULL){ + *ratio_out = (max_corr + 1e-9)/(max_mag + 1e-9); + } + /* Skip coarse frequency est if we don't have anywhere to put result */ if(foff_out == NULL){ return t_est; @@ -391,10 +395,6 @@ static int coarse_sync(struct OFDM *ofdm, complex float *rx, int length, float * //fprintf(stderr,"foff_est is %f, tpos is %d, ratio is %f\n",foff_est,t_est%1280,max_corr/max_mag); *foff_out = foff_est; - if(ratio_out != NULL){ - *ratio_out = (max_corr + 1e-9)/(max_mag + 1e-9); - } - return t_est; } @@ -813,10 +813,10 @@ void ofdm_demod_coarse(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in){ complex float shift_nco = 1; complex float shift_nco_dph = 0; - if(sync <= 0 || 1 ) { + if(sync <= 0 ) { /* Do coarse search over frequency range in 20hz slices*/ - //for(fshift = -100; fshift < 100; fshift+=40){ - fshift = 0;{ + for(fshift = -103; fshift < 100; fshift+=40){ + //fshift = 0;{ /* Shift block of samples by fshift */ shift_nco_dph = cexpf(2*M_PI*(fshift/(float)Fs)*I); shift_nco = 1; @@ -826,21 +826,30 @@ void ofdm_demod_coarse(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in){ } /* Do a coarse search for the pilot */ - frame_pos_est = coarse_sync(ofdm,rxbuf_temp,RxBufSize,&foff_out,&corr_ratio); + frame_pos_est = coarse_sync(ofdm,rxbuf_temp,RxBufSize,NULL,&corr_ratio); if(corr_ratio > corr_ratio_max){ corr_ratio_max = corr_ratio; - fshift_max = -fshift + foff_out; + fshift_max = fshift; frame_pos_max = frame_pos_est; } } frame_pos_est = (frame_pos_max % SampsPerFrame) - (M+Ncp) + 1; if(corr_ratio_max > 0.5 && sync <=0){ - ofdm->sync_count = 3; + /* Now that we think we can sync, do a freq. domain search in the smaller range for our signal */ + shift_nco_dph = cexpf(2*M_PI*(fshift_max/(float)Fs)*I); + shift_nco = 1; + for(int i=0; isync_count = 2; ofdm->frame_point = frame_pos_est; - ofdm->foff_est_hz = fshift_max; - fprintf(stderr,"syncing\n"); + ofdm->foff_est_hz = -fshift_max + foff_out; + fprintf(stderr,"syncing, shift %f,fine foff %f\n",fshift_max,-fshift_max + foff_out); } - fprintf(stderr,"Best ratio %f freq %f offset %d\n",corr_ratio_max,fshift_max,frame_pos_est); + fprintf(stderr,"Best ratio %f freq %f offset %d fine %f\n",corr_ratio_max,fshift_max,frame_pos_est,ofdm->foff_est_hz); } } @@ -972,7 +981,9 @@ void ofdm_demod(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in) { ft_est = coarse_sync(ofdm, work, (en - st), &freq_offset,&ratio_out); fprintf(stderr,"good ratio %f\n",ratio_out); // fprintf(stderr,"\nx>"); - + if(ratio_out < .7){ + ofdm->sync_count--; + } //ofdm->foff_est_hz += (freq_offset/2); ofdm->timing_est += (ft_est - ceilf(FtWindowWidth / 2)); @@ -992,10 +1003,11 @@ void ofdm_demod(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in) { ft_est_other = (ft_est_other % SampsPerFrame) - (M+Ncp) + 1; // fprintf(stderr,"\n"); //coarser_sync(ofdm, work, (en - st), &freq_offset, NULL); - if(ratio_out<.5){ - ofdm->sync_count--; + + fprintf(stderr,"ratio is %f t_offset is %d f_off is %f freq_off is %f\n",ratio_out,ft_est_other,ofdm->foff_est_hz,freq_offset); + if(freq_offset > 1){ + //ofdm->foff_est_hz += freq_offset/4; } - fprintf(stderr,"ratio is %f t_offset is %d\n",ratio_out,ft_est_other); } /*