From: drowe67 Date: Mon, 18 May 2015 07:04:21 +0000 (+0000) Subject: Octave and C agree with freq offset and freq drift and AWGN channel, yayyyyyyyy X-Git-Url: http://git.whiteaudio.com/gitweb/?a=commitdiff_plain;h=cbb870cb3cffb27bec228e085df1a966232031b2;p=freetel-svn-tracking.git Octave and C agree with freq offset and freq drift and AWGN channel, yayyyyyyyy git-svn-id: https://svn.code.sf.net/p/freetel/code@2134 01035d8c-6547-0410-b346-abe4f91aad63 --- diff --git a/codec2-dev/octave/cohpsk.m b/codec2-dev/octave/cohpsk.m index ba1c9503..31069972 100644 --- a/codec2-dev/octave/cohpsk.m +++ b/codec2-dev/octave/cohpsk.m @@ -638,12 +638,14 @@ function [ch_symb rx_timing rx_filt rx_baseband afdmdv f_est] = rate_Fs_rx_proce mod_strip += amod_strip; end %plot(mod_strip) - + %printf("modstrip: %f %f\n", real(mod_strip), imag(mod_strip)); + % loop filter made up of 1st order IIR plus integrator. Integerator % was found to be reqd afdmdv.filt = (1-beta)*afdmdv.filt + beta*angle(mod_strip); f_est += g*afdmdv.filt; + %printf("filt: %f angle: %f\n", afdmdv.filt, angle(mod_strip)); end end diff --git a/codec2-dev/octave/tcohpsk.m b/codec2-dev/octave/tcohpsk.m index 84c8da3a..36a1e370 100644 --- a/codec2-dev/octave/tcohpsk.m +++ b/codec2-dev/octave/tcohpsk.m @@ -10,6 +10,19 @@ % or (ii) can optionally be used to run an Octave version of the cohpsk % modem to tune and develop it. +% TODO: +% +% [ ] Test +% [X] AWGN channel +% [X] freq offset +% [X] fading channel +% [X] freq drift +% [ ] timing drift +% [X] tune perf/impl loss to get closer to ideal +% [X] linear interp of phase for better fading perf +% [X] freq offset/drift feedback loop +% [ ] PAPR measurement and reduction + graphics_toolkit ("gnuplot"); more off; @@ -30,8 +43,8 @@ test = 'compare to c'; if strcmp(test, 'compare to c') frames = 35; - foff = 0; - dfoff = 0; + foff = 55.5; + dfoff = -0.5/Fs; EsNodB = 8; fading_en = 0; hf_delay_ms = 2; @@ -111,7 +124,7 @@ afdmdv.Nfiltertiming = afdmdv.M + afdmdv.Nfilter + afdmdv.M; afdmdv.rx_filter_memory = zeros(afdmdv.Nc+1, afdmdv.Nfilter); afdmdv.filt = 0; -afdmdv.prev_rx_symb = zeros(1,afdmdv.Nc+1); +afdmdv.prev_rx_symb = ones(1,afdmdv.Nc+1); % COHPSK Init -------------------------------------------------------- @@ -156,7 +169,7 @@ ct_symb_ff_log = []; rx_timing_log = []; ratio_log = []; foff_log = []; -fest_log = []; +f_est_log = []; % Channel modeling and BER measurement ---------------------------------------- @@ -323,7 +336,7 @@ for f=1:frames % If in sync just do sample rate processing on latest frame if sync == 1 - [ch_symb rx_timing rx_filt rx_baseband afdmdv acohpsk.f_est] = rate_Fs_rx_processing(afdmdv, ch_fdm_frame, acohpsk.f_est, acohpsk.Nsymbrowpilot, nin, 0); + [ch_symb rx_timing rx_filt rx_baseband afdmdv acohpsk.f_est] = rate_Fs_rx_processing(afdmdv, ch_fdm_frame, acohpsk.f_est, acohpsk.Nsymbrowpilot, nin, 1); [next_sync acohpsk] = frame_sync_fine_freq_est(acohpsk, ch_symb, sync, next_sync); acohpsk.ct_symb_ff_buf(1:2,:) = acohpsk.ct_symb_ff_buf(acohpsk.Nsymbrowpilot+1:acohpsk.Nsymbrowpilot+2,:); @@ -333,7 +346,8 @@ for f=1:frames rx_filt_log = [rx_filt_log rx_filt]; ch_symb_log = [ch_symb_log; ch_symb]; rx_timing_log = [rx_timing_log rx_timing]; - fest_log = [fest_log acohpsk.f_est]; + f_est_log = [f_est_log acohpsk.f_est]; + %printf("%f\n", acohpsk.f_est); end % if we are in sync complete demodulation with symbol rate processing @@ -388,6 +402,10 @@ printf("\nOctave EsNodB: %4.1f ber..: %4.3f Nerrs..: %d Tbits..: %d\n", EsNodB, if compare_with_c + % Output vectors from C port --------------------------------------------------- + + load ../build_linux/unittest/tcohpsk_out.txt + % Determine bit error rate sz = length(tx_bits_log_c); @@ -396,10 +414,6 @@ if compare_with_c ber_c = Nerrs_c/Tbits_c; printf("C EsNodB.....: %4.1f ber_c: %4.3f Nerrs_c: %d Tbits_c: %d\n", EsNodB, ber_c, Nerrs_c, Tbits_c); - % Output vectors from C port --------------------------------------------------- - - load ../build_linux/unittest/tcohpsk_out.txt - stem_sig_and_error(1, 111, tx_bits_log_c, tx_bits_log - tx_bits_log_c, 'tx bits', [1 length(tx_bits) -1.5 1.5]) stem_sig_and_error(2, 211, real(tx_symb_log_c), real(tx_symb_log - tx_symb_log_c), 'tx symb re', [1 length(tx_symb_log_c) -1.5 1.5]) stem_sig_and_error(2, 212, imag(tx_symb_log_c), imag(tx_symb_log - tx_symb_log_c), 'tx symb im', [1 length(tx_symb_log_c) -1.5 1.5]) @@ -415,13 +429,15 @@ if compare_with_c [n m] = size(ch_symb_log); stem_sig_and_error(6, 211, real(ch_symb_log_c), real(ch_symb_log - ch_symb_log_c), 'ch symb re', [1 n -1.5 1.5]) stem_sig_and_error(6, 212, imag(ch_symb_log_c), imag(ch_symb_log - ch_symb_log_c), 'ch symb im', [1 n -1.5 1.5]) - %stem_sig_and_error(7, 211, real(ct_symb_ff_log_c), real(ct_symb_ff_log - ct_symb_ff_log_c), 'ct symb ff re', [1 n -1.5 1.5]) - %stem_sig_and_error(7, 212, imag(ct_symb_ff_log_c), imag(ct_symb_ff_log - ct_symb_ff_log_c), 'ct symb ff im', [1 n -1.5 1.5]) + + [n m] = size(rx_symb_log); stem_sig_and_error(8, 211, rx_amp_log_c, rx_amp_log - rx_amp_log_c, 'Amp Est', [1 n -1.5 1.5]) stem_sig_and_error(8, 212, rx_phi_log_c, rx_phi_log - rx_phi_log_c, 'Phase Est', [1 n -4 4]) stem_sig_and_error(9, 211, real(rx_symb_log_c), real(rx_symb_log - rx_symb_log_c), 'rx symb re', [1 n -1.5 1.5]) stem_sig_and_error(9, 212, imag(rx_symb_log_c), imag(rx_symb_log - rx_symb_log_c), 'rx symb im', [1 n -1.5 1.5]) - stem_sig_and_error(10, 111, rx_bits_log_c, rx_bits_log - rx_bits_log_c, 'rx bits', [1 n -1.5 1.5]) + + stem_sig_and_error(10, 111, rx_bits_log_c, rx_bits_log - rx_bits_log_c, 'rx bits', [1 length(rx_bits_log) -1.5 1.5]) + stem_sig_and_error(11, 111, f_est_log_c - Fcentre, f_est_log - f_est_log_c, 'f est', [1 length(f_est_log) foff-5 foff+5]) check(tx_bits_log, tx_bits_log_c, 'tx_bits'); check(tx_symb_log, tx_symb_log_c, 'tx_symb'); @@ -429,12 +445,13 @@ if compare_with_c check(ch_fdm_frame_log, ch_fdm_frame_log_c, 'ch_fdm_frame'); %check(rx_fdm_frame_bb_log, rx_fdm_frame_bb_log_c, 'rx_fdm_frame_bb', 0.01); - check(ch_symb_log, ch_symb_log_c, 'ch_symb',0.01); + check(ch_symb_log, ch_symb_log_c, 'ch_symb',0.05); %check(ct_symb_ff_log, ct_symb_ff_log_c, 'ct_symb_ff',0.01); check(rx_amp_log, rx_amp_log_c, 'rx_amp_log',0.01); - check(rx_phi_log, rx_phi_log_c, 'rx_phi_log',0.01); + check(rx_phi_log, rx_phi_log_c, 'rx_phi_log',0.05); check(rx_symb_log, rx_symb_log_c, 'rx_symb',0.01); check(rx_bits_log, rx_bits_log_c, 'rx_bits'); + check(f_est_log, f_est_log_c, 'f_est'); else @@ -466,13 +483,13 @@ else subplot(211) plot(foff_log,';freq offset;'); hold on; - plot(fest_log - Fcentre,'g;freq offset est;'); + plot(f_est_log - Fcentre,'g;freq offset est;'); hold off; title('freq offset'); legend("boxoff"); subplot(212) - plot(foff_log(1:length(fest_log)) - fest_log + Fcentre) + plot(foff_log(1:length(f_est_log)) - fest_log + Fcentre) title('freq offset estimation error'); end diff --git a/codec2-dev/src/cohpsk.c b/codec2-dev/src/cohpsk.c index 21967d97..f3b37e49 100644 --- a/codec2-dev/src/cohpsk.c +++ b/codec2-dev/src/cohpsk.c @@ -5,20 +5,7 @@ DATE CREATED: March 2015 Functions that implement a coherent PSK FDM modem. - - TODO: - - [ ] Code to plot EB/No v BER curves to char perf for - [ ] AWGN channel - [ ] freq offset - [ ] fading channel - [ ] freq drift - [ ] timing drift - [ ] tune perf/impl loss to get closer to ideal - [X] linear interp of phase for better fading perf - [ ] freq offset/drift feedback loop - [ ] smaller freq est block size to min ram req - + \*---------------------------------------------------------------------------*/ /* @@ -564,6 +551,8 @@ void rate_Fs_rx_processing(struct COHPSK *coh, COMP ch_symb[][COHPSK_NC*ND], COM mod_strip.real = 0.0; mod_strip.imag = 0.0; for(c=0; cNc+1; c++) { + //printf("rx_onesym[%d] %f %f prev_rx_symbols[%d] %f %f\n", c, rx_onesym[c].real, rx_onesym[c].imag, + // fdmdv->prev_rx_symbols[c].real, fdmdv->prev_rx_symbols[c].imag); adiff = cmult(rx_onesym[c], cconj(fdmdv->prev_rx_symbols[c])); fdmdv->prev_rx_symbols[c] = rx_onesym[c]; @@ -573,14 +562,16 @@ void rate_Fs_rx_processing(struct COHPSK *coh, COMP ch_symb[][COHPSK_NC*ND], COM amod_strip = cmult(adiff, adiff); amod_strip = cmult(amod_strip, amod_strip); - amod_strip.real = cabsolute(amod_strip); + amod_strip.real = fabsf(amod_strip.real); mod_strip = cadd(mod_strip, amod_strip); } - + //printf("modstrip: %f %f\n", mod_strip.real, mod_strip.imag); + /* loop filter made up of 1st order IIR plus integrator. Integerator was found to be reqd */ fdmdv->filt = (1.0-beta)*fdmdv->filt + beta*atan2(mod_strip.imag, mod_strip.real); + //printf("filt: %f angle: %f\n", fdmdv->filt, atan2(mod_strip.imag, mod_strip.real)); *f_est += g*fdmdv->filt; } @@ -663,7 +654,7 @@ void cohpsk_demod(struct COHPSK *coh, int rx_bits[], int *reliable_sync_bit, COM max_ratio = 0.0; for (coh->f_est = FDMDV_FCENTRE-40.0; coh->f_est <= FDMDV_FCENTRE+40.0; coh->f_est += 40.0) { - printf(" [%d] acohpsk.f_est: %f +/- 20\n", coh->frame, coh->f_est); + fprintf(stderr, " [%d] acohpsk.f_est: %f +/- 20\n", coh->frame, coh->f_est); /* we are out of sync so reset f_est and process two frames to clean out memories */ @@ -728,7 +719,7 @@ void cohpsk_demod(struct COHPSK *coh, int rx_bits[], int *reliable_sync_bit, COM /* If in sync just do sample rate processing on latest frame */ if (sync == 1) { - rate_Fs_rx_processing(coh, ch_symb, rx_fdm, &coh->f_est, NSYMROWPILOT, nin, 0); + rate_Fs_rx_processing(coh, ch_symb, rx_fdm, &coh->f_est, NSYMROWPILOT, nin, 1); frame_sync_fine_freq_est(coh, ch_symb, sync, &next_sync); for(r=0; r<2; r++) diff --git a/codec2-dev/unittest/tcohpsk.c b/codec2-dev/unittest/tcohpsk.c index c4990062..71d1109c 100644 --- a/codec2-dev/unittest/tcohpsk.c +++ b/codec2-dev/unittest/tcohpsk.c @@ -49,7 +49,8 @@ #define FRAMESL (SYNC_FRAMES*FRAMES) /* worst case is every frame is out of sync */ #define RS 50 -#define FOFF 0 +#define FOFF 55.5 +#define DFOFF (-0.5/(float)FS) #define ESNODB 8 extern float pilots_coh[][PILOTS_NC]; @@ -80,7 +81,8 @@ int main(int argc, char *argv[]) FILE *fout; int f, r, c, log_r, log_data_r, noise_r, ff_log_r; int *ptest_bits_coh, *ptest_bits_coh_end; - COMP phase_ch; + double foff; + COMP foff_rect, phase_ch; struct FDMDV *fdmdv; //COMP rx_filt[COHPSK_NC*ND][P+1]; @@ -91,6 +93,8 @@ int main(int argc, char *argv[]) //COMP rx_onesym[COHPSK_NC*ND]; //int rx_baseband_log_col_index = 0; //COMP rx_baseband_log[COHPSK_NC*ND][(M+M/P)*NSYMROWPILOT*FRAMES]; + float f_est_log[FRAMES]; + int f_est_samples; int log_bits; float EsNo, variance; @@ -114,12 +118,13 @@ int main(int argc, char *argv[]) /* init stuff */ - log_r = log_data_r = noise_r = log_bits = ff_log_r = 0; + log_r = log_data_r = noise_r = log_bits = ff_log_r = f_est_samples = 0; ptest_bits_coh = (int*)test_bits_coh; ptest_bits_coh_end = (int*)test_bits_coh + sizeof(test_bits_coh)/sizeof(int); memcpy(tx_bits, test_bits_coh, sizeof(int)*COHPSK_BITS_PER_FRAME); phase_ch.real = 1.0; phase_ch.imag = 0.0; + foff = FOFF; /* each carrier has power = 2, total power 2Nc, total symbol rate NcRs, noise BW B=Fs Es/No = (C/Rs)/(N/B), N = var = @@ -127,6 +132,7 @@ int main(int argc, char *argv[]) EsNo = pow(10.0, ESNODB/10.0); variance = 2.0*FS/(RS*EsNo); + //fprintf(stderr, "doff: %e\n", DFOFF); /* Main Loop ---------------------------------------------------------------------*/ @@ -156,8 +162,17 @@ int main(int argc, char *argv[]) Channel \*---------------------------------------------------------*/ - fdmdv_freq_shift(ch_fdm_frame, tx_fdm_frame, FOFF, &phase_ch, NSYMROWPILOT*M); + //fdmdv_freq_shift(ch_fdm_frame, tx_fdm_frame, FOFF, &phase_ch, NSYMROWPILOT*M); + for(r=0; rf_est; } assert(log_r <= NSYMROWPILOT*FRAMES); assert(noise_r <= NSYMROWPILOT*M*FRAMES); assert(log_data_r <= NSYMROW*FRAMES); - printf("\r[%d]", f+1); + printf("\r [%d]", f+1); } printf("\n"); @@ -232,6 +248,7 @@ int main(int argc, char *argv[]) octave_save_float(fout, "rx_phi_log_c", (float*)rx_phi_log, log_data_r, COHPSK_NC*ND, COHPSK_NC*ND); octave_save_complex(fout, "rx_symb_log_c", (COMP*)rx_symb_log, log_data_r, COHPSK_NC*ND, COHPSK_NC*ND); octave_save_int(fout, "rx_bits_log_c", rx_bits_log, 1, COHPSK_BITS_PER_FRAME*log_bits); + octave_save_float(fout, "f_est_log_c", &f_est_log[1], 1, f_est_samples-1, f_est_samples-1); fclose(fout); cohpsk_destroy(coh);