debugging C port of Octave modem, C version compiles OK
authordrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Fri, 15 May 2015 03:08:55 +0000 (03:08 +0000)
committerdrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Fri, 15 May 2015 03:08:55 +0000 (03:08 +0000)
git-svn-id: https://svn.code.sf.net/p/freetel/code@2130 01035d8c-6547-0410-b346-abe4f91aad63

codec2-dev/octave/cohpsk.m
codec2-dev/octave/tcohpsk.m
codec2-dev/src/cohpsk.c
codec2-dev/src/cohpsk_internal.h
codec2-dev/src/fdmdv_internal.h
codec2-dev/unittest/tcohpsk.c
codec2-dev/unittest/test_cohpsk_ch.c

index e742f652c75221d23a36489a8592d435b68cf465..78a1794afd174fd1f25007eab98dd7ab57bc8e38 100644 (file)
@@ -594,7 +594,7 @@ function [ch_symb rx_timing rx_filt rx_baseband afdmdv f_est] = rate_Fs_rx_proce
 
     % figure(10); clf; hold on;
     for r=1:nsymb
-      % downconvert entire signal to nominal baseband 
+      % shift signal to nominal baseband, this will put Nc/2 carriers either side of 0 Hz
 
       [rx_fdm_frame_bb afdmdv.fbb_phase_rx] = freq_shift(ch_fdm_frame(1+(r-1)*M:r*M), -f_est, afdmdv.Fs, afdmdv.fbb_phase_rx);
 
@@ -621,21 +621,21 @@ function [ch_symb rx_timing rx_filt rx_baseband afdmdv f_est] = rate_Fs_rx_proce
         for c=1:afdmdv.Nc+1
           adiff = rx_onesym(c) .* conj(afdmdv.prev_rx_symb(c));
           afdmdv.prev_rx_symb(c) = rx_onesym(c);
+
+          % 4th power strips QPSK modulation, by multiplying phase by 4
+          % Using the abs value of the real coord was found to help 
+          % non-linear issues when noise power was large
+
           amod_strip = adiff.^4;
           amod_strip = abs(real(amod_strip)) + j*imag(amod_strip);
           mod_strip += amod_strip;
         end
         %plot(mod_strip)
-
-        % 4th power strips QPSK modulation, by multiplying phase by 4
-        % Using the abs value of the real coord was found to help 
-        % non-linear issues when noise power was large
-
         
         % 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));
+        afdmdv.filt = (1-beta)*afdmdv.filt + beta*angle(mod_strip);
         f_est += g*afdmdv.filt;
 
       end
@@ -643,7 +643,7 @@ function [ch_symb rx_timing rx_filt rx_baseband afdmdv f_est] = rate_Fs_rx_proce
 endfunction
 
 
-function ct_symb_buf = update_ct_symb_buf(ct_symb_buf, ch_symb, Nct_sym_buf, Nsymbrowpilot)
+void update_ct_symb_buf(COMP ct_symb_buf, ch_symb, Nct_sym_buf, Nsymbrowpilot)
 
   % update memory in symbol buffer
 
@@ -714,79 +714,20 @@ function [next_sync cohpsk] = frame_sync_fine_freq_est(cohpsk, ch_symb, sync, ne
     end
     cohpsk.ratio = abs(max_corr/max_mag);
   end
-
-  if sync == 1
-
-    % we are in sync so just sample correlation over 1D grid of fine freq points
-
-    max_corr = 0;
-    st = cohpsk.f_fine_est - 1;
-    en = cohpsk.f_fine_est + 1;
-    for f_fine = st:0.25:en
-        f_fine_rect = exp(-j*f_fine*2*pi*sampling_points/Rs)';
-        corr = 0; mag = 0;
-        for c=1:Nc*Nd
-          f_corr_vec = f_fine_rect .* ct_symb_buf(cohpsk.ct+sampling_points,c);
-          for p=1:length(sampling_points)
-            corr += pilot2(p,c-Nc*floor((c-1)/Nc)) * f_corr_vec(p);
-            mag  += abs(f_corr_vec(p));
-          end
-        end
-        if corr > max_corr
-          max_corr = corr;
-          max_mag = mag;
-          f_fine_est = f_fine;
-        end
-    end
-
-    %cohpsk.f_est -= 0.5*f_fine_est;
-    %printf("  coarse: %f  fine: %f\n", cohpsk.f_est, f_fine_est);
-    cohpsk.ratio = abs(max_corr/max_mag);
-  end
   
+  % single point correlation just to see if we are still in sync
 
-endfunction
-
-
-% fine freq correction
-
-function acohpsk = fine_freq_correct(acohpsk, sync, next_sync);
-  ct_symb_ff_buf = acohpsk.ct_symb_ff_buf;
-
-  % We can decode first frame that we achieve sync.  Need to fine freq
-  % correct all of it's symbols, including pilots.  From then on, just
-  % correct new symbols into frame.  make copy, so if we lose sync we
-  % havent fine freq corrected ct_symb_buf if next_sync == 4 correct
-  % all 8 if sync == 2 correct latest 6
-
-
-  if (next_sync == 4) && (sync == 2)
-      
-      % first frame, we've just gotten sync so fine freq correct all Nsymbrowpilot+2 samples
-
-      ct_symb_ff_buf = acohpsk.ct_symb_buf(acohpsk.ct+1:acohpsk.ct+acohpsk.Nsymbrowpilot+2,:);
-      for r=1:acohpsk.Nsymbrowpilot+2
-        acohpsk.ff_phase *= acohpsk.ff_rect';
-        ct_symb_ff_buf(r,:) *= acohpsk.ff_phase;
-      end
-  end
-
-  if sync == 4
-      % second and subsequent frames, just fine freq correct the latest Nsymbrowpilot
-
-      ct_symb_ff_buf(1:2,:) = ct_symb_ff_buf(acohpsk.Nsymbrowpilot+1:acohpsk.Nsymbrowpilot+2,:);
-      ct_symb_ff_buf(3:acohpsk.Nsymbrowpilot+2,:) = acohpsk.ct_symb_buf(acohpsk.ct+3:acohpsk.ct+acohpsk.Nsymbrowpilot+2,:);
-      for r=3:acohpsk.Nsymbrowpilot+2
-        acohpsk.ff_phase *= acohpsk.ff_rect';
-        ct_symb_ff_buf(r,:) *= acohpsk.ff_phase;
+  if sync == 1
+    corr = 0; mag = 0;
+    for c=1:Nc*Nd
+      for p=1:length(sampling_points)
+        corr += pilot2(p, c-Nc*floor((c-1)/Nc)) * ct_symb_buf(cohpsk.ct + sampling_points,c);
+        mag  += abs(f_corr_vec(p));
       end
+    end
+    cohpsk.ratio = abs(corr)/mag;
   end
 
-  mag = abs(acohpsk.ff_phase);
-  acohpsk.ff_phase /= mag;
-
-  acohpsk.ct_symb_ff_buf = ct_symb_ff_buf;
-
 endfunction
 
 
index f313a93a68b2c085a3d11cb29c025e79fe4bdb66..42954845570ed0e1faf6a53a124eeb89ea50087d 100644 (file)
@@ -225,7 +225,7 @@ for f=1:frames
   next_sync = sync;
   nin = M;
 
-  % if out of sync do Initial Freq offset estimation over last two frames
+  % if out of sync do Initial Freq offset estimation over NSW frames to flush out memories
 
   if (sync == 0) && (f>1)
 
index deca9c49cc8bc06043ea277fda082542b2ea016f..e3e88bf89d511142f30012e7b70f9ccfb8e38f59 100644 (file)
@@ -69,6 +69,7 @@ static COMP qpsk_mod[] = {
 static int sampling_points[] = {0, 1, 6, 7};
 
 void corr_with_pilots(COMP *corr_out, float *mag_out, struct COHPSK *coh, int t, float f_fine);
+void update_ct_symb_buf(COMP ct_symb_buf[][COHPSK_NC*ND], COMP ch_symb[][COHPSK_NC*ND]);
 
 /*---------------------------------------------------------------------------*\
                                                                              
@@ -131,7 +132,16 @@ struct COHPSK *cohpsk_create(void)
     }
 
     coh->ff_phase.real = 1.0; coh->ff_phase.imag = 0.0;
-    coh->sync = 0;
+    coh->sync  = 0;
+    coh->frame = 0;
+    coh->ratio = 0.0;
+
+    /* clear sync window buffer */
+
+    for (i=0; i<NSW*NSYMROWPILOT*M; i++) {
+        coh->ch_fdm_frame_buf[i].real = 0.0;
+        coh->ch_fdm_frame_buf[i].imag = 0.0;
+    }
 
     /* set up fdmdv states so we can use those modem functions */
 
@@ -153,6 +163,15 @@ struct COHPSK *cohpsk_create(void)
     }
     coh->fdmdv = fdmdv;
 
+    /* disable optional logging by default */
+
+    coh->rx_baseband_log = NULL;
+    coh->rx_baseband_log_col_index = 0;
+    coh->rx_filt_log = NULL;
+    coh->rx_filt_log_col_index = 0;
+    coh->ch_symb_log = NULL;
+    coh->ch_symb_log_r = 0;
+
     return coh;
 }
 
@@ -260,8 +279,8 @@ void qpsk_symbols_to_bits(struct COHPSK *coh, int rx_bits[], COMP ct_symb_buf[][
     float x[NPILOTSFRAME+2], x1;
     COMP  y[NPILOTSFRAME+2], yfit;
     COMP  m, b;
-    COMP  corr, rot, pi_on_4, phi_rect, div_symb;
-    float mag, phi_, amp_;
+    COMP   __attribute__((unused)) corr, rot, pi_on_4, phi_rect, div_symb;
+    float mag,  __attribute__((unused)) phi_,  __attribute__((unused)) amp_;
 
     pi_on_4.real = cosf(M_PI/4); pi_on_4.imag = sinf(M_PI/4);
    
@@ -345,63 +364,6 @@ void qpsk_symbols_to_bits(struct COHPSK *coh, int rx_bits[], COMP ct_symb_buf[][
 }
 
 
-/*---------------------------------------------------------------------------*\
-                                                       
-  FUNCTION....: coarse_freq_offset_est()            
-  AUTHOR......: David Rowe                           
-  DATE CREATED: March 2015
-
-  Determines an estimate of frequency offset, advances to next sync state.
-
-  TODO: This is currently very stack heavy for an embedded uC.  Tweak
-        algorithm in Octave and C to use a smaller DFT and test.
-
-\*---------------------------------------------------------------------------*/
-
-void coarse_freq_offset_est(struct COHPSK *coh, struct FDMDV *fdmdv, COMP ch_fdm_frame[], int sync, int *next_sync)
-{
-    float f_start, f_stop, sc, h, magsq, num, den, bin_est;
-    COMP  s[COARSE_FEST_NDFT], S[COARSE_FEST_NDFT];
-    int   bin_start, bin_stop, i;
-
-    if (sync == 0) {
-        f_start = FDMDV_FCENTRE - (((float)(COHPSK_NC*ND-1)/2)+2)*FSEP; 
-        f_stop = FDMDV_FCENTRE + (((float)(COHPSK_NC*ND-1)/2)+2)*FSEP;
-        sc = (float)COARSE_FEST_NDFT/FS;
-        bin_start = floorf(f_start*sc+0.5);
-        bin_stop = floorf(f_stop*sc+0.5);
-        //printf("f_start: %f f_stop: %f sc: %f bin_start: %d bin_stop: %d\n",
-        //        f_start, f_stop, sc, bin_start, bin_stop);
-
-        for(i=0; i<NSYMROWPILOT*M; i++) {
-            h = 0.5 - 0.5*cosf(2*M_PI*i/(NSYMROWPILOT*M-1));
-            s[i] = fcmult(h, ch_fdm_frame[i]); 
-        }
-        
-        for (; i<COARSE_FEST_NDFT; i++) {
-            s[i].real = 0.0;
-            s[i].imag = 0.0;
-        }
-
-        kiss_fft(coh->fft_coarse_fest, (kiss_fft_cpx *)s, (kiss_fft_cpx *)S);
-        
-        /* find centroid of signal energy inside search window */
-
-        num = den = 0.0;
-        for(i=bin_start; i<=bin_stop; i++) {
-            magsq = S[i].real*S[i].real + S[i].imag*S[i].imag;
-            num += (float)(i)*magsq;
-            den += magsq;
-        }
-        bin_est = num/den;
-        coh->f_est = floor(bin_est/sc+0.5);
-
-        fprintf(stderr, "  coarse freq est: %f\n", coh->f_est);
-        
-        *next_sync = 1;
-    }
-}
-
 
 void corr_with_pilots(COMP *corr_out, float *mag_out, struct COHPSK *coh, int t, float f_fine) 
 {
@@ -436,40 +398,24 @@ void corr_with_pilots(COMP *corr_out, float *mag_out, struct COHPSK *coh, int t,
   frequency offset, advances to next sync state if we have a reliable
   match for frame sync.
 
-  TODO: This is very stack heavy for an embedded uC.  Tweak algorthim to use
-        a smaller DFT and test.
-
 \*---------------------------------------------------------------------------*/
 
 void frame_sync_fine_freq_est(struct COHPSK *coh, COMP ch_symb[][COHPSK_NC*ND], int sync, int *next_sync)
 {
-    int   r,c,i,t;
+    int   t;
     float f_fine, mag, max_corr, max_mag;
     COMP  corr;
 
-    /* update memory in symbol buffer */
-
-    for (r=0; r<NCT_SYMB_BUF-NSYMROWPILOT; r++) {
-        for(c=0; c<COHPSK_NC*ND; c++) {
-            coh->ct_symb_buf[r][c] = coh->ct_symb_buf[r+NSYMROWPILOT][c];
-        }
-    }
-
-    for (r=NCT_SYMB_BUF-NSYMROWPILOT,i=0; r<NCT_SYMB_BUF; r++,i++) {
-       for(c=0; c<COHPSK_NC*ND; c++) {
-           coh->ct_symb_buf[r][c] = ch_symb[i][c];
-           //printf("%d %d %f %f\n", i,c,ch_symb[i][c].real, ch_symb[i][c].imag);
-       }
-    }
+    update_ct_symb_buf(coh->ct_symb_buf, ch_symb);
 
     /* sample pilots at start of this frame and start of next frame */
 
-    if (sync == 2) {
+    if (sync == 0) {
 
         /* sample correlation over 2D grid of time and fine freq points */
 
         max_corr = 0;
-        for (f_fine=-20; f_fine<=20; f_fine+=1.0) {
+        for (f_fine=-20; f_fine<=20; f_fine+=0.25) {
             for (t=0; t<NSYMROWPILOT; t++) {
                 corr_with_pilots(&corr, &mag, coh, t, f_fine);
                 //printf("  f: %f  t: %d corr: %f %f\n", f_fine, t, corr.real, corr.imag);
@@ -487,81 +433,35 @@ void frame_sync_fine_freq_est(struct COHPSK *coh, COMP ch_symb[][COHPSK_NC*ND],
         coh->ff_rect.imag = -sinf(coh->f_fine_est*2.0*M_PI/RS);
         fprintf(stderr, "  fine freq f: %f max_corr: %f max_mag: %f ct: %d\n", coh->f_fine_est, max_corr, max_mag, coh->ct);
  
-        if (max_corr/max_mag > 0.9) {
-            fprintf(stderr, "  in sync!\n");
-            *next_sync = 4;
+        if (max_corr/max_mag > 0.7) {
+            fprintf(stderr, "  [%d] encouraging sync word!\n", coh->frame);
             coh->sync_timer = 0;
+            *next_sync = 1;
         }
         else {
             *next_sync = 0;
-            fprintf(stderr, "  back to coarse freq offset est...\n");
+            fprintf(stderr, "  [%d]  back to coarse freq offset est...\n", coh->frame);
         }
-        
+        coh->ratio = max_corr/max_mag;        
     }
 }
 
 
-/*---------------------------------------------------------------------------*\
-                                                       
-  FUNCTION....: fine_freq_correct()         
-  AUTHOR......: David Rowe                           
-  DATE CREATED: April 2015
-
-  Fine frequency correction of symbols at rate Rs.
-
-\*---------------------------------------------------------------------------*/
-
-void fine_freq_correct(struct COHPSK *coh, int sync, int next_sync) {
-    int   r,c;
-    float mag;
-
-  /*
-    We can decode first frame that we achieve sync.  Need to fine freq
-    correct all of it's symbols, including pilots.  From then on, just
-    correct new symbols into frame.  make copy, so if we lose sync we
-    havent fine freq corrected ct_symb_buf if next_sync == 4 correct
-    all 8 if sync == 2 correct latest 6.
-  */
-
-  if ((next_sync == 4) || (sync == 4)) {
-
-      if ((next_sync == 4) && (sync == 2)) {
-          
-          /* first frame, we've just gotten sync so fine freq correct all Nsymbrowpilot+2 samples */
-
-          for(r=0; r<NSYMROWPILOT+2; r++) {
-              coh->ff_phase = cmult(coh->ff_phase, cconj(coh->ff_rect));
-              for(c=0; c<COHPSK_NC*ND; c++) {
-                  coh->ct_symb_ff_buf[r][c] = coh->ct_symb_buf[coh->ct+r][c];
-                  coh->ct_symb_ff_buf[r][c] = cmult(coh->ct_symb_ff_buf[r][c], coh->ff_phase);
-              }
-          }
-      }
-      else {
-          
-          /* second and subsequent frames, just fine freq correct the latest Nsymbrowpilot */
-
-          for(r=0; r<2; r++) {
-              for(c=0; c<COHPSK_NC*ND; c++) {
-                  coh->ct_symb_ff_buf[r][c] = coh->ct_symb_ff_buf[r+NSYMROWPILOT][c];
-              }
-          }
-          
-          for(; r<NSYMROWPILOT+2; r++) {
-              coh->ff_phase = cmult(coh->ff_phase, cconj(coh->ff_rect));
-              for(c=0; c<COHPSK_NC*ND; c++) {
-                  coh->ct_symb_ff_buf[r][c] = coh->ct_symb_buf[coh->ct+r][c];
-                  //printf("%d %d %f %f\n", r,c,coh->ct_symb_ff_buf[r][c].real, coh->ct_symb_ff_buf[r][c].imag);
-                  coh->ct_symb_ff_buf[r][c] = cmult(coh->ct_symb_ff_buf[r][c], coh->ff_phase);
-              }
-          }
+void update_ct_symb_buf(COMP ct_symb_buf[][COHPSK_NC*ND], COMP ch_symb[][COHPSK_NC*ND])
+{
+    int r, c, i;
 
-      }
+    /* update memory in symbol buffer */
 
-      mag = cabsolute(coh->ff_phase);
-      coh->ff_phase.real /= mag;
-      coh->ff_phase.imag /= mag;
-  }
+    for(r=0; r<NCT_SYMB_BUF-NSYMROWPILOT; r++) {
+        for(c=0; c<COHPSK_NC*ND; c++)
+            ct_symb_buf[r][c] = ct_symb_buf[r+NSYMROWPILOT][c];
+    }
+  
+    for(r=NCT_SYMB_BUF-NSYMROWPILOT, i=0; r<NSYMROWPILOT; r++, i++) {
+        for(c=0; c<COHPSK_NC*ND; c++)
+            ct_symb_buf[r][c] = ch_symb[i][c];
+    }
 }
 
 
@@ -570,14 +470,11 @@ int sync_state_machine(struct COHPSK *coh, int sync, int next_sync)
     COMP  corr;
     float mag;
 
-    if (sync == 1)
-        next_sync = 2;
+    if (sync == 1) {
 
-    if (sync == 4) {
+        /* check that sync is still good, fall out of sync on consecutive bad frames */
 
-        /* check that sync is still good, fall out of sync on 3 consecutive bad frames */
-
-        corr_with_pilots(&corr, &mag, coh, coh->ct, coh->f_fine_est);
+        corr_with_pilots(&corr, &mag, coh, coh->ct, 0.0);
         // printf("%f\n", cabsolute(corr)/mag);
 
         if (cabsolute(corr)/mag < 0.5) 
@@ -585,8 +482,8 @@ int sync_state_machine(struct COHPSK *coh, int sync, int next_sync)
         else
             coh->sync_timer = 0;            
 
-        if (coh->sync_timer == 3) {
-            fprintf(stderr,"  lost sync ....\n");
+        if (coh->sync_timer == 10) {
+            fprintf(stderr,"  [%d] lost sync ....\n", coh->frame);
             next_sync = 0;
         }
     }
@@ -630,6 +527,94 @@ void cohpsk_mod(struct COHPSK *coh, COMP tx_fdm[], int tx_bits[])
 }
 
 
+void rate_Fs_rx_processing(struct COHPSK *coh, COMP ch_symb[NSYMROWPILOT][COHPSK_NC*ND], COMP ch_fdm_frame[], float *f_est, int nsymb, int nin, int freq_track)
+{
+    struct FDMDV *fdmdv = coh->fdmdv;
+    int   r, c, i;
+    COMP  rx_fdm_frame_bb[M];
+    COMP  rx_baseband[COHPSK_NC*ND][M+M/P];
+    COMP  rx_filt[COHPSK_NC*ND][P+1];
+    float env[NT*P], __attribute__((unused)) rx_timing;
+    COMP  rx_onesym[COHPSK_NC*ND];
+    float beta, g;
+    COMP  adiff, amod_strip, mod_strip;
+
+    assert(nin < M*NSYMROWPILOT);
+
+    for (r=0; r<nsymb; r++) {
+        fdmdv_freq_shift(rx_fdm_frame_bb, &ch_fdm_frame[r*M], -(*f_est), &fdmdv->fbb_phase_rx, nin);
+        fdm_downconvert(rx_baseband, fdmdv->Nc, rx_fdm_frame_bb, fdmdv->phase_rx, fdmdv->freq, nin);
+        rx_filter(rx_filt, fdmdv->Nc, rx_baseband, coh->rx_filter_memory, nin);
+        rx_timing = rx_est_timing(rx_onesym, fdmdv->Nc, rx_filt, fdmdv->rx_filter_mem_timing, env, nin);
+          
+        for(c=0; c<COHPSK_NC*ND; c++) {
+            ch_symb[r][c] = rx_onesym[c];
+        }
+
+        /* freq tracking, see test_ftrack.m for unit test.  Placed in
+           this function as it needs to work on a symbol by symbol
+           abasis rather than frame by frame.  This means the control
+           loop operates at a sample rate of Rs = 50Hz for say 1 Hz/s
+           drift. */
+
+        if (freq_track) {
+            beta = 0.005;
+            g = 0.2;
+
+            /* combine difference on phase from last symbol over Nc carriers */
+
+            mod_strip.real = 0.0; mod_strip.imag = 0.0;
+            for(c=0; c<fdmdv->Nc+1; c++) {
+                adiff = cmult(rx_onesym[c], cconj(fdmdv->prev_rx_symbols[c]));
+                fdmdv->prev_rx_symbols[c] = rx_onesym[c];
+
+                /* 4th power strips QPSK modulation, by multiplying phase by 4
+                   Using the abs value of the real coord was found to help 
+                   non-linear issues when noise power was large. */
+
+                amod_strip = cmult(adiff, adiff);
+                amod_strip = cmult(amod_strip, amod_strip);
+                amod_strip.real = cabsolute(amod_strip);
+                mod_strip = cadd(mod_strip, amod_strip);
+            }
+        
+            /* 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);
+            *f_est += g*fdmdv->filt;
+        }    
+
+        /* Optional logging used for testing against Octave version */
+
+        if (coh->rx_baseband_log) {
+            for(c=0; c<COHPSK_NC*ND; c++) {       
+            for(i=0; i<nin; i++) {
+              coh->rx_baseband_log[c*(M+M/P)*LOG_FRAMES*NSYMROWPILOT + coh->rx_baseband_log_col_index + i] = rx_baseband[c][i]; 
+            }
+          }
+          coh->rx_baseband_log_col_index += nin;        
+        }
+
+        if (coh->rx_filt_log) {
+         for(c=0; c<COHPSK_NC*ND; c++) {       
+            for(i=0; i<P; i++) {
+              coh->rx_filt_log[c*(P+1)*LOG_FRAMES*NSYMROWPILOT + coh->rx_filt_log_col_index + i] = rx_filt[c][i]; 
+            }
+         }
+         coh->rx_filt_log_col_index += P;        
+        }
+    }
+
+    if (coh->ch_symb_log) {
+       for(r=0; r<NSYMROWPILOT; r++, coh->ch_symb_log_r++) {
+            for(c=0; c<COHPSK_NC*ND; c++) {
+               coh->ch_symb_log[coh->ch_symb_log_r*COHPSK_NC*ND* + c] = ch_symb[r][c]; 
+            }
+        }
+    }
+}
+
 /*---------------------------------------------------------------------------*\
                                                        
   FUNCTION....: cohpsk_demod()      
@@ -653,45 +638,103 @@ void cohpsk_mod(struct COHPSK *coh, COMP tx_fdm[], int tx_bits[])
 
 void cohpsk_demod(struct COHPSK *coh, int rx_bits[], int *reliable_sync_bit, COMP rx_fdm[])
 {
-    struct FDMDV *fdmdv = coh->fdmdv;
-    COMP  rx_fdm_frame_bb[M*NSYMROWPILOT];
-    COMP  rx_baseband[COHPSK_NC*ND][M+M/P];
-    COMP  rx_filt[COHPSK_NC*ND][P+1];
-    float env[NT*P], __attribute__((unused)) rx_timing;
-    COMP  ch_symb[NSYMROWPILOT][COHPSK_NC*ND];
-    COMP  rx_onesym[COHPSK_NC*ND];
-    int   sync, next_sync, nin, r, c;
+    COMP  ch_symb[NSW*NSYMROWPILOT][COHPSK_NC*ND];
+    int   i, j, sync, anext_sync, next_sync, nin, r, c;
+    float max_ratio, f_est;
+
+    /* store two frames of received samples so we can rewind if we get a good candidate */
+
+    for (i=0; i<(NSW-1)*NSYMROWPILOT*M; i++)
+        coh->ch_fdm_frame_buf[i] = coh->ch_fdm_frame_buf[i+NSYMROWPILOT*M];
+    for (j=0; i<NSYMROWPILOT*M; i++,j++)
+        coh->ch_fdm_frame_buf[i] = rx_fdm[j];
 
     next_sync = sync = coh->sync;
+    nin = M;
 
-    coarse_freq_offset_est(coh, fdmdv, rx_fdm, sync, &next_sync);
+    /* if out of sync do Initial Freq offset estimation using NSW frames to flush out filter memories */
 
-    /* sample rate demod processing */
+    if (sync == 0) {
 
-    nin = M;
-    for (r=0; r<NSYMROWPILOT; r++) {
-        fdmdv_freq_shift(&rx_fdm_frame_bb[r*M], &rx_fdm[r*M], -coh->f_est, &fdmdv->fbb_phase_rx, nin);
-        fdm_downconvert(rx_baseband, fdmdv->Nc, &rx_fdm_frame_bb[r*M], fdmdv->phase_rx, fdmdv->freq, nin);
-        rx_filter(rx_filt, fdmdv->Nc, rx_baseband, coh->rx_filter_memory, nin);
-        rx_timing = rx_est_timing(rx_onesym, fdmdv->Nc, rx_filt, fdmdv->rx_filter_mem_timing, env, nin);
-          
-        for(c=0; c<COHPSK_NC*ND; c++) {
-            ch_symb[r][c] = rx_onesym[c];
+        /* we can test +/- 20Hz, so we break this up into 3 tests to cover +/- 60Hz */
+
+        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);
+
+            /* we are out of sync so reset f_est and process two frames to clean out memories */
+
+            rate_Fs_rx_processing(coh, ch_symb, coh->ch_fdm_frame_buf, &coh->f_est, NSW*NSYMROWPILOT, nin, 0);
+            for (i=0; i<NSW-1; i++) {
+                update_ct_symb_buf(coh->ct_symb_buf, &ch_symb[i*NSYMROWPILOT]);
+            }
+            frame_sync_fine_freq_est(coh, &ch_symb[(NSW-1)*NSYMROWPILOT], sync, &anext_sync);
+
+            if (anext_sync == 1) {
+                //printf("  [%d] acohpsk.ratio: %f\n", f, coh->ratio);
+                if (coh->ratio > max_ratio) {
+                    max_ratio   = coh->ratio;
+                    f_est       = coh->f_est - coh->f_fine_est;
+                    next_sync   = anext_sync;
+                }
+            }
         }
+
+        if (next_sync == 1) {
+
+            /* we've found a sync candidate!
+               re-process last NSW frames with adjusted f_est then check again */
+
+            coh->f_est = f_est;
+
+            fprintf(stderr, "  [%d] trying sync and f_est: %f\n", coh->frame, coh->f_est);
+
+            rate_Fs_rx_processing(coh, ch_symb, coh->ch_fdm_frame_buf, &coh->f_est, NSW*NSYMROWPILOT, nin, 0);
+            for (i=0; i<NSW-1; i++) {
+                update_ct_symb_buf(coh->ct_symb_buf, &ch_symb[i*NSYMROWPILOT]);
+            }
+            frame_sync_fine_freq_est(coh, &ch_symb[(NSW-1)*NSYMROWPILOT], sync, &next_sync);
+
+            if (fabs(coh->f_fine_est) > 2.0) {
+                fprintf(stderr, "  [%d] Hmm %f is a bit big so back to coarse est ...\n", coh->frame, coh->f_fine_est);
+                next_sync = 0;
+            }
+        }
+
+        if (next_sync == 1) {
+            /* OK we are in sync!
+               demodulate first frame (demod completed below) */
+
+            fprintf(stderr, "  [%d] in sync!\n", coh->frame);
+            for(r=0; r<NSYMROWPILOT+2; r++)
+                for(c=0; c<COHPSK_NC*ND; c++)
+                    coh->ct_symb_ff_buf[r][c] = coh->ct_symb_buf[coh->ct+r][c];
+        }
+    }  
+
+    /* 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, 1);
+        frame_sync_fine_freq_est(coh, &ch_symb[(NSW-1)*NSYMROWPILOT], sync, &next_sync);
+        for(r=0; r<2; r++)
+            for(c=0; c<COHPSK_NC*ND; c++)
+                coh->ct_symb_ff_buf[r][c] = coh->ct_symb_ff_buf[r+NSYMROWPILOT][c];
+        for(; r<NSYMROWPILOT+2; r++)
+            for(c=0; c<COHPSK_NC*ND; c++)
+                coh->ct_symb_ff_buf[r][c] = coh->ct_symb_buf[coh->ct+r][c];
     }
-     
-    /* coarse timing (frame sync) and initial fine freq est */
+
+    /* if we are in sync complete demodulation with symbol rate processing */
   
-    frame_sync_fine_freq_est(coh, ch_symb, sync, &next_sync);
-    fine_freq_correct(coh, sync, next_sync);
-        
     *reliable_sync_bit = 0;
-    if ((sync == 4) || (next_sync == 4)) {
+    if ((next_sync == 1) || (sync == 1)) {
         qpsk_symbols_to_bits(coh, rx_bits, coh->ct_symb_ff_buf);
         *reliable_sync_bit = 1;
     }
 
     sync = sync_state_machine(coh, sync, next_sync);        
-
     coh->sync = sync;
 }
index 1fee91249458ff089ec380565138b47f32db8a1a..76d81525dd3db1f8e2cda44745ef2c37acd8f6d8 100644 (file)
 #define NCT_SYMB_BUF     (2*NSYMROWPILOT+2)
 #define ND               2                           /* diversity factor ND 1 is no diveristy, ND we have orginal plus 
                                                         one copy */
+#define NSW              3                           /* number of sync window frames */
+
+#define LOG_FRAMES       35
 
 #include "fdmdv_internal.h"
 #include "kiss_fft.h"
 
 struct COHPSK {
+    COMP         ch_fdm_frame_buf[NSW*NSYMROWPILOT*M];  /* buffer of several frames of symbols from channel      */
     float        pilot2[2*NPILOTSFRAME][COHPSK_NC];    
     float        phi_[NSYMROW][COHPSK_NC*ND];           /* phase estimates for this frame of rx data symbols     */
     float        amp_[NSYMROW][COHPSK_NC*ND];           /* amplitude estimates for this frame of rx data symbols */
@@ -53,7 +57,19 @@ struct COHPSK {
     int          sync;
     int          sync_timer;
 
+    int          frame;
+    float        ratio;
+
     struct FDMDV *fdmdv;
+    
+    /* optional log variables used for tetsing Octave to C port */
+
+    COMP           *rx_baseband_log;
+    int            rx_baseband_log_col_index;
+    COMP           *rx_filt_log;
+    int            rx_filt_log_col_index;
+    COMP           *ch_symb_log;
+    int            ch_symb_log_r;
 };
 
 
index 4476aea0c1830b3d7e7da8c7107a63ff7f15f674..ba27be342b94bd08cdaa2798033fbc140e5777b3 100644 (file)
@@ -124,7 +124,8 @@ struct FDMDV {
 
     float foff;
     COMP foff_phase_rect;
-    
+    float filt;
+
     /* Demodulator */
 
     COMP  rxdec_lpf_mem[NRXDEC-1+M]; 
index d48730b7a3f90a5d44a3fe421a284b9351002e56..74371dc402868a6a58d52f6a442af8ac50f40d0b 100644 (file)
@@ -44,7 +44,7 @@
 #include "comp_prim.h"
 #include "noise_samples.h"
 
-#define FRAMES 35
+#define FRAMES LOG_FRAMES  /* #defined in cohpsk_internal.h */
 #define RS     50
 #define FOFF   10.5
 #define ESNODB 8.0
@@ -58,16 +58,16 @@ int main(int argc, char *argv[])
     COMP           tx_symb[NSYMROWPILOT][COHPSK_NC*ND];
     COMP           tx_fdm_frame[M*NSYMROWPILOT];
     COMP           ch_fdm_frame[M*NSYMROWPILOT];
-    COMP           rx_fdm_frame_bb[M*NSYMROWPILOT];
-    COMP           ch_symb[NSYMROWPILOT][COHPSK_NC*ND];
+    //COMP           rx_fdm_frame_bb[M*NSYMROWPILOT];
+    //COMP           ch_symb[NSYMROWPILOT][COHPSK_NC*ND];
     int            rx_bits[COHPSK_BITS_PER_FRAME];
     
     int            tx_bits_log[COHPSK_BITS_PER_FRAME*FRAMES];
     COMP           tx_symb_log[NSYMROWPILOT*FRAMES][COHPSK_NC*ND];
     COMP           tx_fdm_frame_log[M*NSYMROWPILOT*FRAMES];
     COMP           ch_fdm_frame_log[M*NSYMROWPILOT*FRAMES];
-    COMP           rx_fdm_frame_bb_log[M*NSYMROWPILOT*FRAMES];
-    COMP           ch_symb_log[NSYMROWPILOT*FRAMES][COHPSK_NC*ND];
+    //COMP           rx_fdm_frame_bb_log[M*NSYMROWPILOT*FRAMES];
+    //COMP           ch_symb_log[NSYMROWPILOT*FRAMES][COHPSK_NC*ND];
     COMP           ct_symb_ff_log[NSYMROWPILOT*FRAMES][COHPSK_NC*ND];
     float          rx_amp_log[NSYMROW*FRAMES][COHPSK_NC*ND];
     float          rx_phi_log[NSYMROW*FRAMES][COHPSK_NC*ND];
@@ -75,40 +75,40 @@ int main(int argc, char *argv[])
     int            rx_bits_log[COHPSK_BITS_PER_FRAME*FRAMES];
                                           
     FILE          *fout;
-    int            f, r, c, log_r, log_data_r, noise_r, i;
+    int            f, r, c, log_r, log_data_r, noise_r;
     int           *ptest_bits_coh, *ptest_bits_coh_end;
     COMP           phase_ch;
 
     struct FDMDV  *fdmdv;
-    COMP           rx_baseband[COHPSK_NC*ND][M+M/P];
-    int            nin;
-    COMP           rx_filt[COHPSK_NC*ND][P+1];
-    COMP           rx_filt_log[COHPSK_NC*ND][(P+1)*NSYMROWPILOT*FRAMES];
-    int            rx_filt_log_col_index = 0;
-    float          env[NT*P];
-    float           __attribute__((unused)) rx_timing;
+    //COMP           rx_filt[COHPSK_NC*ND][P+1];
+    //int            rx_filt_log_col_index = 0;
+    //float          env[NT*P];
+    //float           __attribute__((unused)) rx_timing;
     COMP           tx_onesym[COHPSK_NC*ND];
-    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];
+    //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];
 
-    int            sync, next_sync, log_bits;
+    int            log_bits;
     float          EsNo, variance;
     COMP           scaled_noise;
+    int            reliable_sync_bit;
 
     coh = cohpsk_create();
     fdmdv = coh->fdmdv;
     assert(coh != NULL);
 
+    coh->rx_filt_log = (COMP *)malloc(sizeof(COMP)*COHPSK_NC*ND*(P+1)*NSYMROWPILOT*FRAMES);
+    coh->rx_baseband_log = (COMP *)malloc(sizeof(COMP)*COHPSK_NC*ND*(M+M/P)*NSYMROWPILOT*FRAMES);
+    coh->ch_symb_log = (COMP *)malloc(sizeof(COMP)*NSYMROWPILOT*FRAMES*COHPSK_NC*ND);
+
     log_r = log_data_r = noise_r = log_bits = 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; 
-    sync = 0;
-    
+     
     /*  each carrier has power = 2, total power 2Nc, total symbol rate
         NcRs, noise BW B=Fs Es/No = (C/Rs)/(N/B), N = var =
         2NcFs/NcRs(Es/No) = 2Fs/Rs(Es/No) */
@@ -154,69 +154,25 @@ int main(int argc, char *argv[])
                                  Demod
        \*---------------------------------------------------------*/
 
-        next_sync = sync;
-
-        coarse_freq_offset_est(coh, fdmdv, ch_fdm_frame, sync, &next_sync);
-
-        /* sample rate demod processing */
-
-        nin = M;
-        for (r=0; r<NSYMROWPILOT; r++) {
-          fdmdv_freq_shift(&rx_fdm_frame_bb[r*M], &ch_fdm_frame[r*M], -coh->f_est, &fdmdv->fbb_phase_rx, nin);
-          fdm_downconvert(rx_baseband, fdmdv->Nc, &rx_fdm_frame_bb[r*M], fdmdv->phase_rx, fdmdv->freq, nin);
-          rx_filter(rx_filt, fdmdv->Nc, rx_baseband, coh->rx_filter_memory, nin);
-         rx_timing = rx_est_timing(rx_onesym, fdmdv->Nc, rx_filt, fdmdv->rx_filter_mem_timing, env, nin);
-          
-          for(c=0; c<COHPSK_NC*ND; c++) {
-             ch_symb[r][c] = rx_onesym[c];
-          }
-          
-         for(c=0; c<COHPSK_NC*ND; c++) {       
-            for(i=0; i<nin; i++) {
-              rx_baseband_log[c][rx_baseband_log_col_index + i] = rx_baseband[c][i]; 
-            }
-         }
-         rx_baseband_log_col_index += nin;        
+        cohpsk_demod(coh, rx_bits, &reliable_sync_bit, ch_fdm_frame);
 
-         for(c=0; c<COHPSK_NC*ND; c++) {       
-            for(i=0; i<P; i++) {
-              rx_filt_log[c][rx_filt_log_col_index + i] = rx_filt[c][i]; 
-            }
-         }
-         rx_filt_log_col_index += P;        
-
-        }
-
-        /* coarse timing (frame sync) and initial fine freq est */
-  
-        frame_sync_fine_freq_est(coh, ch_symb, sync, &next_sync);
-        fine_freq_correct(coh, sync, next_sync);
-        
-        if ((sync == 4) || (next_sync == 4)) {
-           qpsk_symbols_to_bits(coh, rx_bits, coh->ct_symb_ff_buf);
-        }
-
-        //printf("f: %d sync: %d next_sync: %d\n", f, sync, next_sync);
-        sync = sync_state_machine(coh, sync, next_sync);
-        
-       /* --------------------------------------------------------*\
+       /* --------------------------------------------------------*\
                               Log each vector 
        \*---------------------------------------------------------*/
 
        memcpy(&tx_bits_log[COHPSK_BITS_PER_FRAME*f], tx_bits, sizeof(int)*COHPSK_BITS_PER_FRAME);
        memcpy(&tx_fdm_frame_log[M*NSYMROWPILOT*f], tx_fdm_frame, sizeof(COMP)*M*NSYMROWPILOT);
        memcpy(&ch_fdm_frame_log[M*NSYMROWPILOT*f], ch_fdm_frame, sizeof(COMP)*M*NSYMROWPILOT);
-               memcpy(&rx_fdm_frame_bb_log[M*NSYMROWPILOT*f], rx_fdm_frame_bb, sizeof(COMP)*M*NSYMROWPILOT);
+               //memcpy(&rx_fdm_frame_bb_log[M*NSYMROWPILOT*f], rx_fdm_frame_bb, sizeof(COMP)*M*NSYMROWPILOT);
 
        for(r=0; r<NSYMROWPILOT; r++, log_r++) {
             for(c=0; c<COHPSK_NC*ND; c++) {
                tx_symb_log[log_r][c] = tx_symb[r][c]; 
-               ch_symb_log[log_r][c] = ch_symb[r][c]; 
                ct_symb_ff_log[log_r][c] = coh->ct_symb_ff_buf[r][c]; 
             }
         }
 
-        if ((sync == 4) || (next_sync == 4)) {
+        if (coh->sync == 1) {
 
            for(r=0; r<NSYMROW; r++, log_data_r++) {
                 for(c=0; c<COHPSK_NC*ND; c++) {
@@ -246,10 +202,10 @@ int main(int argc, char *argv[])
     octave_save_complex(fout, "tx_symb_log_c", (COMP*)tx_symb_log, NSYMROWPILOT*FRAMES, COHPSK_NC*ND, COHPSK_NC*ND);  
     octave_save_complex(fout, "tx_fdm_frame_log_c", (COMP*)tx_fdm_frame_log, 1, M*NSYMROWPILOT*FRAMES, M*NSYMROWPILOT*FRAMES);  
     octave_save_complex(fout, "ch_fdm_frame_log_c", (COMP*)ch_fdm_frame_log, 1, M*NSYMROWPILOT*FRAMES, M*NSYMROWPILOT*FRAMES);  
-    octave_save_complex(fout, "rx_fdm_frame_bb_log_c", (COMP*)rx_fdm_frame_bb_log, 1, M*NSYMROWPILOT*FRAMES, M*NSYMROWPILOT*FRAMES);  
-    octave_save_complex(fout, "rx_baseband_log_c", (COMP*)rx_baseband_log, COHPSK_NC*ND, rx_baseband_log_col_index, (M+M/P)*FRAMES*NSYMROWPILOT);  
-    octave_save_complex(fout, "rx_filt_log_c", (COMP*)rx_filt_log, COHPSK_NC*ND, rx_filt_log_col_index, (P+1)*FRAMES*NSYMROWPILOT);  
-    octave_save_complex(fout, "ch_symb_log_c", (COMP*)ch_symb_log, NSYMROWPILOT*FRAMES, COHPSK_NC*ND, COHPSK_NC*ND);  
+    //octave_save_complex(fout, "rx_fdm_frame_bb_log_c", (COMP*)rx_fdm_frame_bb_log, 1, M*NSYMROWPILOT*FRAMES, M*NSYMROWPILOT*FRAMES);  
+    octave_save_complex(fout, "rx_baseband_log_c", (COMP*)coh->rx_baseband_log, COHPSK_NC*ND, coh->rx_baseband_log_col_index, (M+M/P)*FRAMES*NSYMROWPILOT);  
+    octave_save_complex(fout, "rx_filt_log_c", (COMP*)coh->rx_filt_log, COHPSK_NC*ND, coh->rx_filt_log_col_index, (P+1)*FRAMES*NSYMROWPILOT);  
+    octave_save_complex(fout, "ch_symb_log_c", (COMP*)coh->ch_symb_log, NSYMROWPILOT*FRAMES, COHPSK_NC*ND, COHPSK_NC*ND);  
     octave_save_complex(fout, "ct_symb_ff_log_c", (COMP*)ct_symb_ff_log, NSYMROWPILOT*FRAMES, COHPSK_NC*ND, COHPSK_NC*ND);  
     octave_save_float(fout, "rx_amp_log_c", (float*)rx_amp_log, log_data_r, COHPSK_NC*ND, COHPSK_NC*ND);  
     octave_save_float(fout, "rx_phi_log_c", (float*)rx_phi_log, log_data_r, COHPSK_NC*ND, COHPSK_NC*ND);  
index bed611a23d50d4329c0e8562b3467e5802e02dc9..2ef6f3c5019ed37251712bae7753ac8581d17aaa 100644 (file)
 #include "comp_prim.h"
 #include "noise_samples.h"
 
-#define FRAMES   350
-#define FOFF_HZ  10.5
-#define ES_NO_DB  12.0
+#define FS          8000
+#define FRAMES      100
+#define FOFF_HZ     10.5
+#define ES_NO_DB     8.0
+#define HF_DELAY_MS  2.0
+#define HF_GAIN      1.423599
+
+#define FADING_FILE_NAME "../../raw/fading_samples.float"
 
 int main(int argc, char *argv[])
 {
@@ -54,19 +59,23 @@ int main(int argc, char *argv[])
     int           *ptest_bits_coh, *ptest_bits_coh_end, *ptest_bits_coh_rx;
     COMP           phase_ch;
     int            noise_r, noise_end;
-    float          corr;
+    int            errors;
     int            state, next_state, nerrors, nbits, reliable_sync_bit;
     float          EsNo, variance;
     COMP           scaled_noise;
     float          EsNodB, foff_hz;
+    int            fading_en, nhfdelay, ret;
+    COMP          *ch_fdm_delay, aspread, aspread_2ms, delayed, direct;
+    FILE          *ffading;
 
     EsNodB = ES_NO_DB;
     foff_hz =  FOFF_HZ;
-    if (argc == 3) {
+    if (argc == 4) {
         EsNodB = atof(argv[1]);
         foff_hz = atof(argv[2]);
+        fading_en = atoi(argv[3]);
     }
-    fprintf(stderr, "EsNodB: %4.2f foff: %4.2f Hz\n", EsNodB, foff_hz);
+    fprintf(stderr, "EsNodB: %4.2f foff: %4.2f Hz fading: %d\n", EsNodB, foff_hz, fading_en);
 
     coh = cohpsk_create();
     assert(coh != NULL);
@@ -84,6 +93,23 @@ int main(int argc, char *argv[])
     EsNo = pow(10.0, EsNodB/10.0);
     variance = 2.0*COHPSK_FS/(COHPSK_RS*EsNo);
 
+    /* init HF fading model */
+
+    if (fading_en) {
+        ffading = fopen(FADING_FILE_NAME, "rb");
+        if (ffading == NULL) {
+            printf("Can't find fading file: %s\n", FADING_FILE_NAME);
+            exit(1);
+        }
+        nhfdelay = floor(HF_DELAY_MS*FS/1000);
+        ch_fdm_delay = (COMP*)malloc((nhfdelay+COHPSK_SAMPLES_PER_FRAME)*sizeof(COMP));
+        assert(ch_fdm_delay != NULL);
+        for(i=0; i<nhfdelay+COHPSK_SAMPLES_PER_FRAME; i++) {
+            ch_fdm_delay[i].real = 0.0;
+            ch_fdm_delay[i].imag = 0.0;
+        }
+    }
+
     /* Main Loop ---------------------------------------------------------------------*/
 
     for(f=0; f<FRAMES; f++) {
@@ -106,6 +132,35 @@ int main(int argc, char *argv[])
 
         fdmdv_freq_shift(ch_fdm, tx_fdm, foff_hz, &phase_ch, COHPSK_SAMPLES_PER_FRAME);
 
+        /* optional HF fading -------------------------------------*/
+
+        if (fading_en) {
+
+            /* update delayed signal buffer */
+
+            for(i=0; i<nhfdelay; i++)
+                ch_fdm_delay[i] = ch_fdm_delay[i+COHPSK_SAMPLES_PER_FRAME];
+            for(; i<COHPSK_SAMPLES_PER_FRAME; i++)
+                ch_fdm_delay[i] = ch_fdm[i];
+
+            /* combine direct and delayed paths, both multiplied by
+               "spreading" (doppler) functions */
+
+            for(i=0; i<COHPSK_SAMPLES_PER_FRAME; i++) {
+                ret = fread(&aspread, sizeof(COMP), 1, ffading);
+                assert(ret == 1);
+                ret = fread(&aspread_2ms, sizeof(COMP), 1, ffading);
+                assert(ret == 1);
+                //printf("%f %f %f %f\n", aspread.real, aspread.imag, aspread_2ms.real, aspread_2ms.imag);
+                
+                direct    = cmult(aspread, ch_fdm[i]);
+                delayed   = cmult(aspread_2ms, ch_fdm_delay[i]);
+                ch_fdm[i] = fcmult(HF_GAIN, cadd(direct, delayed));
+            }
+        }
+
+        /* AWGN noise ------------------------------------------*/
+
         for(r=0; r<COHPSK_SAMPLES_PER_FRAME; r++) {
             scaled_noise = fcmult(sqrt(variance), noise[noise_r]);
             ch_fdm[r] = cadd(ch_fdm[r], scaled_noise);
@@ -120,9 +175,9 @@ int main(int argc, char *argv[])
 
        cohpsk_demod(coh, rx_bits, &reliable_sync_bit, ch_fdm);
 
-        corr = 0.0;
+        errors = 0;
         for(i=0; i<COHPSK_BITS_PER_FRAME; i++) {
-            corr += (1.0 - 2.0*(rx_bits[i] & 0x1)) * (1.0 - 2.0*ptest_bits_coh_rx[i]);
+            errors += (rx_bits[i] & 0x1) ^ ptest_bits_coh_rx[i];
         }
 
         /* state logic to sync up to test data */
@@ -130,10 +185,10 @@ int main(int argc, char *argv[])
         next_state = state;
 
         if (state == 0) {
-            if (reliable_sync_bit && (corr == COHPSK_BITS_PER_FRAME)) {
+            if (reliable_sync_bit && (errors == 0)) {
                 next_state = 1;
                 ptest_bits_coh_rx += COHPSK_BITS_PER_FRAME;
-                nerrors = COHPSK_BITS_PER_FRAME - corr;
+                nerrors = errors;
                 nbits = COHPSK_BITS_PER_FRAME;
                 fprintf(stderr, "  test data sync\n");            
 
@@ -141,7 +196,7 @@ int main(int argc, char *argv[])
         }
 
         if (state == 1) {
-            nerrors += COHPSK_BITS_PER_FRAME - corr;
+            nerrors += errors;
             nbits   += COHPSK_BITS_PER_FRAME;
             ptest_bits_coh_rx += COHPSK_BITS_PER_FRAME;
             if (ptest_bits_coh_rx >= ptest_bits_coh_end) {
@@ -154,6 +209,10 @@ int main(int argc, char *argv[])
     
     printf("%4.3f %d %d\n", (float)nerrors/nbits, nbits, nerrors);
 
+    if (fading_en) {
+        free(ch_fdm_delay);
+        fclose(ffading);
+    }
     cohpsk_destroy(coh);
 
     return 0;