From: baobrien Date: Tue, 23 Jan 2018 19:39:19 +0000 (+0000) Subject: Finally have my head on straight about how to go about ofdm timing/freq est X-Git-Url: http://git.whiteaudio.com/gitweb/?a=commitdiff_plain;h=63e29d9b38aca0fb05eaecca2fb30f0e6887851c;p=freetel-svn-tracking.git Finally have my head on straight about how to go about ofdm timing/freq est git-svn-id: https://svn.code.sf.net/p/freetel/code@3392 01035d8c-6547-0410-b346-abe4f91aad63 --- diff --git a/codec2-dev/octave/ofdm_lib.m b/codec2-dev/octave/ofdm_lib.m index a618bda3..b04c9a30 100644 --- a/codec2-dev/octave/ofdm_lib.m +++ b/codec2-dev/octave/ofdm_lib.m @@ -168,6 +168,7 @@ function states = ofdm_init(bps, Rs, Tcp, Ns, Nc) states.sample_point = states.timing_est = 1; states.nin = states.Nsamperframe; + states.sync = 0; % generate OFDM pilot symbol, used for timing and freq offset est rate_fs_pilot_samples = states.pilots * W/states.M; @@ -430,162 +431,19 @@ endfunction function [rx_bits states aphase_est_pilot_log rx_np rx_amp] = ofdm_demod2(states, rxbuf_in) ofdm_load_const; - - % insert latest input samples into rxbuf - - rxbuf(1:Nrxbuf-states.nin) = rxbuf(states.nin+1:Nrxbuf); - rxbuf(Nrxbuf-states.nin+1:Nrxbuf) = rxbuf_in; - - % get latest freq offset estimate - - woff_est = 2*pi*foff_est_hz/Fs; - - % update timing estimate -------------------------------------------------- - - delta_t = 0; - if timing_en - % update timing at start of every frame - - st = M+Ncp + Nsamperframe + 1 - floor(ftwindow_width/2) + (timing_est-1); - en = st + Nsamperframe-1 + M+Ncp + ftwindow_width-1; - - ft_est = coarse_sync(states, rxbuf(st:en) .* exp(-j*woff_est*(st:en)), rate_fs_pilot_samples); - timing_est = timing_est + ft_est - ceil(ftwindow_width/2); - - if verbose > 1 - printf(" ft_est: %2d timing_est: %2d sample_point: %2d\n", ft_est, timing_est, sample_point); + [rx_bits states aphase_est_pilot_log rx_np rx_amp] =ofdm_demod(states,rxbuf_in); + + %If out of sync, try a coarse sync on RX buffer + if states.sync <= 0 + [t_est foff_est] = coarse_sync(states,states.rxbuf,rate_fs_pilot_samples); + %Correct to be within one frame and at head of frame instead of after first symbol + t_est = mod(t_est,Nsamperframe) - (M+Ncp); + if(t_est > 30) + states.nin = t_est; end - % Black magic to keep sample_point inside cyclic prefix. Or something like that. - - delta_t = ft_est - ceil(ftwindow_width/2); - sample_point = max(timing_est+Ncp/4, sample_point); - sample_point = min(timing_est+Ncp, sample_point); + printf("t_est is %d, foff_est is %d\n",t_est,foff_est) end - - % down convert at current timing instant---------------------------------- - - % todo: this cld be more efficent, as pilot r becomes r-Ns on next frame - - rx_sym = zeros(1+Ns+1+1, Nc+2); - - % previous pilot - - st = M+Ncp + Nsamperframe + (-Ns)*(M+Ncp) + 1 + sample_point; en = st + M - 1; - - for c=1:Nc+2 - acarrier = rxbuf(st:en) .* exp(-j*woff_est*(st:en)) .* conj(W(c,:)); - rx_sym(1,c) = sum(acarrier); - end - - % pilot - this frame - pilot - - for rr=1:Ns+1 - st = M+Ncp + Nsamperframe + (rr-1)*(M+Ncp) + 1 + sample_point; en = st + M - 1; - for c=1:Nc+2 - acarrier = rxbuf(st:en) .* exp(-j*woff_est*(st:en)) .* conj(W(c,:)); - rx_sym(rr+1,c) = sum(acarrier); - end - end - - % next pilot - - st = M+Ncp + Nsamperframe + (2*Ns)*(M+Ncp) + 1 + sample_point; en = st + M - 1; - for c=1:Nc+2 - acarrier = rxbuf(st:en) .* exp(-j*woff_est*(st:en)) .* conj(W(c,:)); - rx_sym(Ns+3,c) = sum(acarrier); - end - - % est freq err based on all carriers ------------------------------------ - - if foff_est_en - freq_err_rect = sum(rx_sym(2,:))' * sum(rx_sym(2+Ns,:)); - - % prevent instability in atan(im/re) when real part near 0 - - freq_err_rect += 1E-6; - - %printf("freq_err_rect: %f %f angle: %f\n", real(freq_err_rect), imag(freq_err_rect), angle(freq_err_rect)); - freq_err_hz = angle(freq_err_rect)*Rs/(2*pi*Ns); - foff_est_hz = foff_est_hz + foff_est_gain*freq_err_hz; - end - - % OK - now estimate and correct phase ---------------------------------- - - aphase_est_pilot = 10*ones(1,Nc+2); - aamp_est_pilot = zeros(1,Nc+2); - for c=2:Nc+1 - - % estimate phase using average of 6 pilots in a rect 2D window centred - % on this carrier - % PPP - % DDD - % DDD - % PPP - - cr = c-1:c+1; - aphase_est_pilot_rect = sum(rx_sym(2,cr)*pilots(cr)') + sum(rx_sym(2+Ns,cr)*pilots(cr)'); - - % use next step of pilots in past and future - - aphase_est_pilot_rect += sum(rx_sym(1,cr)*pilots(cr)'); - aphase_est_pilot_rect += sum(rx_sym(2+Ns+1,cr)*pilots(cr)'); - - aphase_est_pilot(c) = angle(aphase_est_pilot_rect); - aamp_est_pilot(c) = abs(aphase_est_pilot_rect/12); - end - - % correct phase offset using phase estimate, and demodulate - % bits, separate loop as it runs across cols (carriers) to get - % frame bit ordering correct - - aphase_est_pilot_log = []; - rx_bits = []; rx_np = []; rx_amp = []; - for rr=1:Ns-1 - for c=2:Nc+1 - if phase_est_en - rx_corr = rx_sym(rr+2,c) * exp(-j*aphase_est_pilot(c)); - else - rx_corr = rx_sym(rr+2,c); - end - rx_np = [rx_np rx_corr]; - rx_amp = [rx_amp aamp_est_pilot(c)]; - if bps == 1 - abit = real(rx_corr) > 0; - end - if bps == 2 - abit = qpsk_demod(rx_corr); - end - rx_bits = [rx_bits abit]; - end % c=2:Nc+1 - aphase_est_pilot_log = [aphase_est_pilot_log; aphase_est_pilot(2:Nc+1)]; - end - - % Adjust nin to take care of sample clock offset - - nin = Nsamperframe; - if timing_en - thresh = (M+Ncp)/8; - tshift = (M+Ncp)/4; - if timing_est > thresh - nin = Nsamperframe+tshift; - timing_est -= tshift; - sample_point -= tshift; - end - if timing_est < -thresh - nin = Nsamperframe-tshift; - timing_est += tshift; - sample_point += tshift; - end - end - - states.rx_sym = rx_sym; - states.rxbuf = rxbuf; - states.nin = nin; - states.timing_est = timing_est; - states.sample_point = sample_point; - states.delta_t = delta_t; - states.foff_est_hz = foff_est_hz; endfunction diff --git a/codec2-dev/octave/tofdm.m b/codec2-dev/octave/tofdm.m index 88a6dd66..c0177a06 100644 --- a/codec2-dev/octave/tofdm.m +++ b/codec2-dev/octave/tofdm.m @@ -10,8 +10,9 @@ function rx = ebno_awgn_channel(tx, Fs, Rb, EbNodB) EbNo = 10^(EbNodB/10); variance = Fs / (Rb*EbNo); noise = sqrt(variance)*randn(1,length(tx)); - avg = sum(abs(tx))/length(tx) + avg = sum(abs(tx))/length(tx); rx = noise.*avg + tx; + %rx = rx * .10; end function nums = im_re_interleave(nim) @@ -197,4 +198,4 @@ end % This works best on my machine -- Brady graphics_toolkit('fltk') -run_ofdm_test(60,100,.5,60 ,0) \ No newline at end of file +run_ofdm_test(60,20,0.5,40,100) \ No newline at end of file diff --git a/codec2-dev/src/fmfsk.c b/codec2-dev/src/fmfsk.c index ab56945f..2b523352 100644 --- a/codec2-dev/src/fmfsk.c +++ b/codec2-dev/src/fmfsk.c @@ -198,8 +198,8 @@ void fmfsk_demod(struct FMFSK *fmfsk, uint8_t rx_bits[],float fmfsk_in[]){ #endif /* Shift in nin samples */ - memcpy(&oldsamps[0] , &oldsamps[nmem-nold], sizeof(float)*nold); - memcpy(&oldsamps[nold], &fmfsk_in[0] , sizeof(float)*nin ); + memmove(&oldsamps[0] , &oldsamps[nmem-nold], sizeof(float)*nold); + memcpy (&oldsamps[nold], &fmfsk_in[0] , sizeof(float)*nin ); /* Allocate memory for filtering */ float *rx_filt = alloca(sizeof(float)*(nsym+1)*Ts); diff --git a/codec2-dev/src/ofdm.c b/codec2-dev/src/ofdm.c index bcca1eff..99734d27 100644 --- a/codec2-dev/src/ofdm.c +++ b/codec2-dev/src/ofdm.c @@ -192,7 +192,7 @@ static int coarser_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; @@ -278,7 +278,7 @@ static int coarser_sync(struct OFDM *ofdm, complex float *rx, int length, float } foff_est = (pmax > nmax) ? pmax_i : (nmax_i-Fs+1); - fprintf(stderr,"foff_est is %f, tpos is %d, ratio is %f\n",foff_est,t_est,max_mag/max_corr); + //fprintf(stderr,"foff_est is %f, tpos is %d, ratio is %f\n",foff_est,t_est,max_corr/max_mag); *foff_out = foff_est; @@ -388,7 +388,7 @@ static int coarse_sync(struct OFDM *ofdm, complex float *rx, int length, float * } foff_est = (pmax > nmax) ? pmax_i : (nmax_i-Fs+1); //foff_est = pmax_i; - //`fprintf(stderr,"foff_est is %f, tpos is %d\n",foff_est,t_est); + 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; @@ -778,6 +778,7 @@ void ofdm_mod(struct OFDM *ofdm, COMP result[], const int *tx_bits) { /* * Like ofdm_demod, but handles sync and large frame offset */ +//static int coarser_sync(struct OFDM *ofdm, complex float *rx, int length, float *foff_out, float *ratio_out) void ofdm_demod_coarse(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in){ const int SampsPerFrame = ofdm->config.SampsPerFrame; const int RxBufSize = ofdm->config.RxBufSize; @@ -789,16 +790,42 @@ void ofdm_demod_coarse(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in){ complex float * coarse_rxbuf = ofdm->coarse_rxbuf; complex float rxbuf_temp[RxBufSize]; + + /* Copy in nin samples to end of buffer */ + memmove(&coarse_rxbuf[0],&coarse_rxbuf[nin], (RxBufSize - nin) * sizeof(coarse_rxbuf[0])); + memcpy (&coarse_rxbuf[RxBufSize - nin],rxbuf_in, nin * sizeof(coarse_rxbuf[0])); + + int frame_pos_est; + float fshift; + float corr_ratio_max = 0; + float corr_ratio; + float foff_out; + float fshift_max; + int frame_pos_max = 0; /* NCO for freq. shifting */ complex float shift_nco = 1; complex float shift_nco_dph = 0; - /* Copy in nin samples to end of buffer */ - memmove(&coarse_rxbuf[0],&coarse_rxbuf[nin],nin); - memcpy (&coarse_rxbuf[RxBufSize - nin],rxbuf_in,nin); - - if(!sync) { + if(sync <= 0 ) { + /* Do coarse search over frequency range in 20hz slices*/ + for(fshift = -100; fshift < 100; fshift+=40){ + /* Shift block of samples by fshift */ + shift_nco_dph = cexpf(2*M_PI*(fshift/(float)Fs)*I); + shift_nco = 1; + for(int i=0; i corr_ratio_max){ + corr_ratio_max = corr_ratio; + fshift_max = -fshift + foff_out; + frame_pos_max = frame_pos_est; + } + } + fprintf(stderr,"Best ratio %f freq %f offset %d\n",corr_ratio_max,fshift_max,frame_pos_max); } } @@ -832,21 +859,21 @@ void ofdm_demod(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in) { complex float * rxbuf = ofdm->rxbuf; - //memmove(&rxbuf[0],&rxbuf[nin],nin); - //memcpy (&rxbuf[RxBufSize - nin],rxbuf_in,nin); + memmove(&rxbuf[0],&rxbuf[nin], (RxBufSize - nin) * sizeof(rxbuf[0])); + memcpy (&rxbuf[RxBufSize - nin],rxbuf_in, nin * sizeof(rxbuf[0])); /* shift the buffer left based on nin */ - for (i = 0, j = ofdm->nin; i < (RxBufSize - ofdm->nin); i++, j++) { - ofdm->rxbuf[i] = ofdm->rxbuf[j]; - } + // for (i = 0, j = ofdm->nin; i < (RxBufSize - ofdm->nin); i++, j++) { + // ofdm->rxbuf[i] = ofdm->rxbuf[j]; + // } /* insert latest input samples onto tail of rxbuf */ /* Note: COMP and complex float have the same memory layout */ /* see iso/iec 9899:1999 sec 6.2.5.13 */ - for (i = (RxBufSize - ofdm->nin), j = 0; i < RxBufSize; i++, j++) { - ofdm->rxbuf[i] = rxbuf_in[j].real + rxbuf_in[j].imag * I; - } + // for (i = (RxBufSize - ofdm->nin), j = 0; i < RxBufSize; i++, j++) { + // ofdm->rxbuf[i] = rxbuf_in[j].real + rxbuf_in[j].imag * I; + // } /* @@ -876,9 +903,12 @@ void ofdm_demod(struct OFDM *ofdm, int *rx_bits, COMP *rxbuf_in) { work[j] = ofdm->rxbuf[i] * cexpf(-I * woff_est * i); } + fprintf(stderr,"n>"); ft_est = coarse_sync(ofdm, work, (en - st), &freq_offset); - - coarser_sync(ofdm, work, (en - st), &freq_offset, NULL); + fprintf(stderr,"\nx>"); + coarse_sync(ofdm,ofdm->rxbuf,RxBufSize,&freq_offset); + fprintf(stderr,"\n"); + //coarser_sync(ofdm, work, (en - st), &freq_offset, NULL); //ofdm->foff_est_hz += (freq_offset/2); ofdm->timing_est += (ft_est - ceilf(FtWindowWidth / 2));