frames = nbits/(Nc*Nb);
prev_rx_symbols = ones(Nc+1,1);
- foff_phase = 1;
+ foff_phase_rect = 1;
% BER stats
rx_symbols_log = [];
rx_timing_log = [];
+ foff_coarse_log = [];
foff_log = [];
rx_fdm_log = [];
snr_est_log = [];
if sync == 0
foff = foff_coarse;
end
- foff_log = [ foff_log foff ];
+ foff_coarse_log = [foff_coarse_log foff_coarse];
+
foff_rect = exp(j*2*pi*foff/Fs);
for i=1:nin
- foff_phase *= foff_rect';
- rx_fdm(i) = rx_fdm(i)*foff_phase;
+ foff_phase_rect *= foff_rect';
+ rx_fdm_fcorr(i) = rx_fdm(i)*foff_phase_rect;
end
% baseband processing
- rx_baseband = fdm_downconvert(rx_fdm, nin);
- rx_filt = rx_filter(rx_baseband, nin);
-
+ rx_fdm_filter = rxdec_filter(rx_fdm_fcorr, nin);
+ rx_filt = down_convert_and_rx_filter(rx_fdm_filter, nin, M/Q);
[rx_symbols rx_timing] = rx_est_timing(rx_filt, nin);
-
rx_timing_log = [rx_timing_log rx_timing];
+
nin = M;
if rx_timing > 2*M/P
nin += M/P;
nin -= M/P;
end
- if strcmp(modulation,'dqpsk')
- rx_symbols_log = [rx_symbols_log rx_symbols.*conj(prev_rx_symbols./abs(prev_rx_symbols))*exp(j*pi/4)];
- else
- rx_symbols_log = [rx_symbols_log rx_symbols];
- endif
+ rx_symbols_log = [rx_symbols_log rx_symbols.*conj(prev_rx_symbols./abs(prev_rx_symbols))*exp(j*pi/4)];
[rx_bits sync_bit f_err pd] = psk_to_bits(prev_rx_symbols, rx_symbols, modulation);
% optionally save output symbols
% update some states
+ prev_rx_symbols = rx_symbols;
[sig_est noise_est] = snr_update(sig_est, noise_est, pd);
snr_est = calc_snr(sig_est, noise_est);
snr_est_log = [snr_est_log snr_est];
foff -= 0.5*f_err;
- prev_rx_symbols = rx_symbols;
+ foff_log = [foff_log foff];
% freq est state machine
% Version 2
%
-function fdmdv_demod_c(dumpfilename, bits)
- NumCarriers = 16;
+function fdmdv_demod_c(dumpfilename, bits, NumCarriers)
+
fdmdv; % include modem code
frames = bits/(Nc*Nb);
test_frame_sync_log = [test_frame_sync_log test_frame_sync_state];
end
+ ber = total_bit_errors / total_bits;
+ printf("%d bits %d errors BER: %1.4f\n",total_bits, total_bit_errors, ber);
+
% ---------------------------------------------------------------------
% Plots
% ---------------------------------------------------------------------
fbb_phase_rx /= mag;
[pilot prev_pilot pilot_lut_index prev_pilot_lut_index] = get_pilot(pilot_lut_index, prev_pilot_lut_index, nin);
-
[foff_coarse S1 S2] = rx_est_freq_offset(rx_fdm, pilot, prev_pilot, nin);
%sync = 0; % when debugging good idea to uncomment this to "open loop"
end
nin_log = [nin_log nin];
- [rx_bits sync_bit foff_fine pd] = psk_to_bits(prev_rx_symbols, rx_symbols, 'dqpsk');
+ [rx_bits sync_bit f_err pd] = psk_to_bits(prev_rx_symbols, rx_symbols, 'dqpsk');
phase_difference_log = [phase_difference_log pd];
[sig_est noise_est] = snr_update(sig_est, noise_est, pd);
rx_bits_log = [rx_bits_log rx_bits];
foff_fine_log = [foff_fine_log foff_fine];
sync_bit_log = [sync_bit_log sync_bit];
- foff -= 0.5*foff_fine;
+ foff -= 0.5*f_err;
foff_log = [foff_log foff];
% freq est state machine
add_executable(fdmdv_put_test_bits fdmdv_put_test_bits.c fdmdv.c kiss_fft.c)
target_link_libraries(fdmdv_put_test_bits ${CMAKE_REQUIRED_LIBRARIES})
+add_executable(fdmdv_channel fdmdv_channel.c fdmdv.c kiss_fft.c)
+target_link_libraries(fdmdv_channel ${CMAKE_REQUIRED_LIBRARIES})
+
add_executable(fdmdv_interleave fdmdv_interleave.c)
target_link_libraries(fdmdv_interleave ${CMAKE_REQUIRED_LIBRARIES})
fprintf(stderr, "\n\n");
}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: randn()
+ AUTHOR......: David Rowe
+ DATE CREATED: 2 August 2014
+
+ Simple approximation to normal (gaussian) random number generator
+ with 0 mean and unit variance.
+
+\*---------------------------------------------------------------------------*/
+
+#define RANDN_IT 12 /* This magic number of iterations gives us a
+ unit variance. I think beacuse var =
+ (b-a)^2/12 for one uniform random variable, so
+ for a sum of n random variables it's
+ n(b-a)^2/12, or for b=1, a = 0, n=12, we get
+ var = 12(1-0)^2/12 = 1 */
+
+static float randn() {
+ int i;
+ float rn = 0.0;
+
+ for(i=0; i<RANDN_IT; i++)
+ rn += (float)rand()/RAND_MAX;
+
+ rn -= (float)RANDN_IT/2.0;
+ return rn;
+}
+
+
/*---------------------------------------------------------------------------*\
FUNCTION....: fdmdv_simulate_channel()
void fdmdv_simulate_channel(struct FDMDV *f, COMP samples[], int nin, float target_snr)
{
- float sig_pwr, target_snr_linear, noise_pwr, noise_pwr_4k, var_uniform, noise_gain;
+ float sig_pwr, target_snr_linear, noise_pwr, noise_pwr_1Hz, noise_pwr_4000Hz, noise_gain;
int i;
/* estimate signal power */
for(i=0; i<nin; i++)
sig_pwr += samples[i].real*samples[i].real + samples[i].imag*samples[i].imag;
- f->sig_pwr_av = 0.9*f->sig_pwr_av + 0.1*(sig_pwr/nin);
+ sig_pwr /= nin;
+
+ f->sig_pwr_av = 0.9*f->sig_pwr_av + 0.1*sig_pwr;
/* det noise to meet target SNR */
target_snr_linear = powf(10.0, target_snr/10.0);
- noise_pwr = f->sig_pwr_av/target_snr_linear;
- noise_pwr_4k = 0.75*noise_pwr; /* this is the equivalent power scaled to a 4000 Hz BW */
-
- /* for convenience we are using a uniform random number generator */
-
- var_uniform = 1.0/12.0;
- noise_gain = sqrtf(0.5*noise_pwr_4k/var_uniform);
+ noise_pwr = f->sig_pwr_av/target_snr_linear; /* noise pwr in a 3000 Hz BW */
+ noise_pwr_1Hz = noise_pwr/3000.0; /* noise pwr in a 1 Hz bandwidth */
+ noise_pwr_4000Hz = noise_pwr_1Hz*4000.0; /* noise pwr in a 4000 Hz BW, which
+ due to fs=8000 Hz in our simulation noise BW */
+ noise_gain = sqrtf(0.5*noise_pwr_4000Hz); /* split noise pwr between real and imag sides */
+
for(i=0; i<nin; i++) {
- samples[i].real += noise_gain*(((float)rand()/RAND_MAX) - 0.5);
- samples[i].imag += noise_gain*(((float)rand()/RAND_MAX) - 0.5);
+ samples[i].real += noise_gain*randn();
+ samples[i].imag += noise_gain*randn();
}
-
- //printf("f->sig_pwr_av: %f target_snr_linear: %f noise_pwr_4k: %f noise_gain: %f\n",
- // f->sig_pwr_av, target_snr_linear, noise_pwr_4k, noise_gain);
+ /*
+ fprintf(stderr, "sig_pwr: %f f->sig_pwr_av: %e target_snr_linear: %f noise_pwr_4000Hz: %e noise_gain: %e\n",
+ sig_pwr, f->sig_pwr_av, target_snr_linear, noise_pwr_4000Hz, noise_gain);
+ */
}
--- /dev/null
+/*---------------------------------------------------------------------------*\
+
+ FILE........: fdmdv_channel.c
+ AUTHOR......: David Rowe
+ DATE CREATED: 2 August 2014
+
+ Given an input raw file (8kHz, 16 bit shorts) of FDMDV modem
+ samples, adds channel impairments and outputs to another raw file.
+
+\*---------------------------------------------------------------------------*/
+
+/*
+ Copyright (C) 2014 David Rowe
+
+ All rights reserved.
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU Lesser General Public License version 2, as
+ published by the Free Software Foundation. This program is
+ distributed in the hope that it will be useful, but WITHOUT ANY
+ WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
+ License for more details.
+
+ You should have received a copy of the GNU Lesser General Public License
+ along with this program; if not, see <http://www.gnu.org/licenses/>.
+*/
+
+#include <assert.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+#include <errno.h>
+
+#include "codec2_fdmdv.h"
+
+int main(int argc, char *argv[])
+{
+ FILE *fin, *fout;
+ short rx_fdm_buf[FDMDV_NOM_SAMPLES_PER_FRAME];
+ COMP rx_fdm[FDMDV_NOM_SAMPLES_PER_FRAME];
+ struct FDMDV *fdmdv;
+ float snrdB, sam;
+ int i;
+
+ if (argc < 3) {
+ printf("usage: %s InputModemRawFile OutputModemRawFile SNRdB\n", argv[0]);
+ printf("e.g %s test_in.raw test_out.raw 4\n", argv[0]);
+ exit(1);
+ }
+
+ if (strcmp(argv[1], "-") == 0) fin = stdin;
+ else if ( (fin = fopen(argv[1],"rb")) == NULL ) {
+ fprintf(stderr, "Error opening input modem sample file: %s: %s.\n",
+ argv[1], strerror(errno));
+ exit(1);
+ }
+
+ if (strcmp(argv[2], "-") == 0) fout = stdout;
+ else if ( (fout = fopen(argv[2],"wb")) == NULL ) {
+ fprintf(stderr, "Error opening output modem sample file: %s: %s.\n",
+ argv[2], strerror(errno));
+ exit(1);
+ }
+
+ snrdB = atof(argv[3]);
+ fdmdv = fdmdv_create(FDMDV_NC);
+
+ while(fread(rx_fdm_buf, sizeof(short), FDMDV_NOM_SAMPLES_PER_FRAME, fin) == FDMDV_NOM_SAMPLES_PER_FRAME) {
+
+ for(i=0; i<FDMDV_NOM_SAMPLES_PER_FRAME; i++) {
+ rx_fdm[i].real = (float)rx_fdm_buf[i]/FDMDV_SCALE;
+ rx_fdm[i].imag = 0.0;
+ }
+
+ /* real signal so we adjust SNR to suit. I think. I always get confused by this! */
+
+ fdmdv_simulate_channel(fdmdv, rx_fdm, FDMDV_NOM_SAMPLES_PER_FRAME, snrdB - 3.0);
+
+ for(i=0; i<FDMDV_NOM_SAMPLES_PER_FRAME; i++) {
+ sam = FDMDV_SCALE*rx_fdm[i].real;
+ if (sam > 32767.0) sam = 32767.0;
+ if (sam < -32767.0) sam = -32767.0;
+ rx_fdm_buf[i] = sam;
+ }
+
+ fwrite(rx_fdm_buf, sizeof(short), FDMDV_NOM_SAMPLES_PER_FRAME, fout);
+
+ /* if this is in a pipeline, we probably don't want the usual
+ buffering to occur */
+
+ if (fout == stdout) fflush(stdout);
+ if (fin == stdin) fflush(stdin);
+ }
+
+ fclose(fin);
+ fclose(fout);
+ fdmdv_destroy(fdmdv);
+
+ return 0;
+}
+
+++ /dev/null
-/*---------------------------------------------------------------------------*\
-
- FILE........: globals.c
- AUTHOR......: David Rowe
- DATE CREATED: 11/5/94
-
- Globals for sinusoidal speech coder.
-
-\*---------------------------------------------------------------------------*/
-
-/*
- Copyright (C) 2009 David Rowe
-
- All rights reserved.
-
- This program is free software; you can redistribute it and/or modify
- it under the terms of the GNU Lesser General Public License version 2.1, as
- published by the Free Software Foundation. This program is
- distributed in the hope that it will be useful, but WITHOUT ANY
- WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
- License for more details.
-
- You should have received a copy of the GNU Lesser General Public License
- along with this program; if not, see <http://www.gnu.org/licenses/>.
-*/
-
-#include "sine.h" /* global defines for coder */
-
-/* Globals used in encoder and decoder */
-
-int frames; /* number of frames processed so far */
-float Sn[M]; /* float input speech samples */
-MODEL model; /* model parameters for the current frame */
-int Nw; /* number of samples in analysis window */
-float sig; /* energy of current frame */
-
-/* Globals used in encoder */
-
-float w[M]; /* time domain hamming window */
-COMP W[FFT_ENC]; /* DFT of w[] */
-COMP Sw[FFT_ENC]; /* DFT of current frame */
-
-/* Globals used in decoder */
-
-COMP Sw_[FFT_ENC]; /* DFT of all voiced synthesised signal */
-float Sn_[AW_DEC]; /* synthesised speech */
-float Pn[AW_DEC]; /* time domain Parzen (trapezoidal) window */
-
+++ /dev/null
-/*---------------------------------------------------------------------------*\
-
- FILE........: globals.h
- AUTHOR......: David Rowe
- DATE CREATED: 1/11/94
-
- Globals for sinusoidal speech coder.
-
-\*---------------------------------------------------------------------------*/
-
-/*
- Copyright (C) 2009 David Rowe
-
- All rights reserved.
-
- This program is free software; you can redistribute it and/or modify
- it under the terms of the GNU Lesser General Public License version 2.1, as
- published by the Free Software Foundation. This program is
- distributed in the hope that it will be useful, but WITHOUT ANY
- WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
- License for more details.
-
- You should have received a copy of the GNU Lesser General Public License
- along with this program; if not, see <http://www.gnu.org/licenses/>.
-*/
-
-/* Globals used in encoder and decoder */
-
-extern int frames; /* number of frames processed so far */
-extern float Sn[]; /* float input speech samples */
-extern MODEL model; /* model parameters for the current frame */
-extern int Nw; /* number of samples in analysis window */
-extern float sig; /* energy of current frame */
-
-/* Globals used in encoder */
-
-extern float w[]; /* time domain hamming window */
-extern COMP W[]; /* frequency domain hamming window */
-extern COMP Sw[]; /* DFT of current frame */
-extern COMP Sw_[]; /* DFT of all voiced synthesised signal */
-
-/* Globals used in decoder */
-
-extern float Sn_[]; /* output synthesised speech samples */
-extern float Pn[]; /* time domain Parzen (trapezoidal) window */
-
+++ /dev/null
-const float glottal[]={
- 0.000000,
- -0.057687,
- -0.115338,
- -0.172917,
- -0.230385,
- -0.287707,
- -0.344845,
- -0.401762,
- -0.458419,
- -0.514781,
- -0.570809,
- -0.626467,
- -0.681721,
- -0.736537,
- -0.790884,
- -0.844733,
- -0.898057,
- -0.950834,
- -1.003044,
- -1.054670,
- -1.105700,
- -1.156124,
- -1.205936,
- -1.255132,
- -1.303711,
- -1.351675,
- -1.399026,
- -1.445769,
- -1.491908,
- -1.537448,
- -1.582393,
- -1.626747,
- -1.670514,
- -1.713693,
- -1.756285,
- -1.798288,
- -1.839697,
- -1.880507,
- -1.920712,
- -1.960302,
- -1.999269,
- -2.037603,
- -2.075295,
- -2.112335,
- -2.148716,
- -2.184430,
- -2.219472,
- -2.253839,
- -2.287531,
- -2.320550,
- -2.352900,
- -2.384588,
- -2.415623,
- -2.446019,
- -2.475788,
- -2.504946,
- -2.533512,
- -2.561501,
- -2.588934,
- -2.615827,
- -2.642198,
- -2.668064,
- -2.693439,
- -2.718337,
- -2.742767,
- -2.766738,
- -2.790256,
- -2.813322,
- -2.835936,
- -2.858094,
- -2.879790,
- -2.901016,
- -2.921759,
- -2.942008,
- -2.961747,
- -2.980961,
- -2.999632,
- -3.017745,
- -3.035282,
- -3.052228,
- -3.068567,
- -3.084285,
- -3.099371,
- -3.113813,
- -3.127605,
- -3.140738,
- 3.129975,
- 3.118167,
- 3.107022,
- 3.096537,
- 3.086709,
- 3.077531,
- 3.068996,
- 3.061096,
- 3.053821,
- 3.047159,
- 3.041102,
- 3.035636,
- 3.030753,
- 3.026441,
- 3.022690,
- 3.019491,
- 3.016836,
- 3.014718,
- 3.013132,
- 3.012072,
- 3.011535,
- 3.011521,
- 3.012028,
- 3.013057,
- 3.014612,
- 3.016695,
- 3.019310,
- 3.022463,
- 3.026160,
- 3.030407,
- 3.035212,
- 3.040580,
- 3.046520,
- 3.053038,
- 3.060141,
- 3.067836,
- 3.076128,
- 3.085023,
- 3.094525,
- 3.104639,
- 3.115367,
- 3.126712,
- 3.138674,
- -3.131930,
- -3.118731,
- -3.104915,
- -3.090485,
- -3.075444,
- -3.059795,
- -3.043543,
- -3.026695,
- -3.009254,
- -2.991229,
- -2.972625,
- -2.953449,
- -2.933710,
- -2.913414,
- -2.892567,
- -2.871176,
- -2.849248,
- -2.826787,
- -2.803798,
- -2.780284,
- -2.756247,
- -2.731689,
- -2.706609,
- -2.681005,
- -2.654875,
- -2.628213,
- -2.601015,
- -2.573272,
- -2.544977,
- -2.516121,
- -2.486694,
- -2.456686,
- -2.426084,
- -2.394879,
- -2.363060,
- -2.330616,
- -2.297538,
- -2.263816,
- -2.229444,
- -2.194416,
- -2.158727,
- -2.122375,
- -2.085359,
- -2.047682,
- -2.009347,
- -1.970361,
- -1.930732,
- -1.890470,
- -1.849587,
- -1.808098,
- -1.766017,
- -1.723360,
- -1.680145,
- -1.636388,
- -1.592105,
- -1.547313,
- -1.502025,
- -1.456256,
- -1.410016,
- -1.363314,
- -1.316157,
- -1.268547,
- -1.220486,
- -1.171971,
- -1.122997,
- -1.073555,
- -1.023636,
- -0.973227,
- -0.922312,
- -0.870875,
- -0.818899,
- -0.766366,
- -0.713257,
- -0.659554,
- -0.605242,
- -0.550303,
- -0.494723,
- -0.438492,
- -0.381598,
- -0.324036,
- -0.265800,
- -0.206889,
- -0.147303,
- -0.087046,
- -0.026121,
- 0.035463,
- 0.097698,
- 0.160576,
- 0.224087,
- 0.288221,
- 0.352969,
- 0.418323,
- 0.484276,
- 0.550822,
- 0.617958,
- 0.685681,
- 0.753991,
- 0.822889,
- 0.892378,
- 0.962462,
- 1.033144,
- 1.104430,
- 1.176325,
- 1.248833,
- 1.321956,
- 1.395696,
- 1.470051,
- 1.545019,
- 1.620593,
- 1.696763,
- 1.773516,
- 1.850837,
- 1.928705,
- 2.007097,
- 2.085987,
- 2.165347,
- 2.245145,
- 2.325347,
- 2.405919,
- 2.486824,
- 2.568025,
- 2.649485,
- 2.731167,
- 2.813033,
- 2.895045,
- 2.977167,
- 3.059362};
projects / freetel-svn-tracking.git / commitdiff