diff options
| author | Marin Ivanov <[email protected]> | 2025-07-25 10:17:14 +0300 |
|---|---|---|
| committer | Marin Ivanov <[email protected]> | 2026-01-18 20:09:26 +0200 |
| commit | 0168586485e6310c598713c911b1dec5618d61a1 (patch) | |
| tree | 6aabc2a12ef8fef70683f5389bea00f948015f77 /unittest | |
* codec2 cut-down version 1.2.0
* Remove codebook and generation of sources
* remove c2dec c2enc binaries
* prepare for emscripten
Diffstat (limited to 'unittest')
43 files changed, 5953 insertions, 0 deletions
diff --git a/unittest/CMakeLists.txt b/unittest/CMakeLists.txt new file mode 100644 index 0000000..656afd0 --- /dev/null +++ b/unittest/CMakeLists.txt @@ -0,0 +1,123 @@ +add_definitions(-DFLOATING_POINT -DVAR_ARRAYS) +include_directories(../src) + +add_executable(tfdmdv tfdmdv.c ../src/fdmdv.c ../src/kiss_fft.c ../src/octave.c) +target_link_libraries(tfdmdv codec2) + +add_executable(tcohpsk tcohpsk.c ../src/cohpsk.c ../src/octave.c) +target_link_libraries(tcohpsk codec2) + +add_executable(tfsk tfsk.c ../src/kiss_fft.c ../src/kiss_fftr.c ../src/octave.c ../src/modem_probe.c) +target_link_libraries(tfsk m) + +add_executable(tfreedv_data_channel tfreedv_data_channel.c ../src/freedv_data_channel.c) + +add_executable(tfmfsk tfmfsk.c ../src/octave.c ../src/modem_probe.c) +target_link_libraries(tfmfsk m) + +add_definitions(-DMODEMPROBE_ENABLE -DXXXXX) + +add_executable(tofdm tofdm.c ../src/octave.c) +target_link_libraries(tofdm m codec2) + +add_executable(tofdm_acq tofdm_acq.c ../src/octave.c) +target_link_libraries(tofdm_acq m codec2) + +if(UNIX) # Uses pthreads +add_executable(tfifo tfifo.c ../src/codec2_fifo.c) +target_link_libraries(tfifo codec2 ${CMAKE_THREAD_LIBS_INIT}) +endif() + +add_definitions(-D__UNITTEST__) + +add_executable(tnewamp1 tnewamp1.c ../src/quantise.c ../src/newamp1.c ../src/mbest.c ../src/kiss_fft.c ../src/sine.c ../src/nlp.c ../src/dump.c ../src/octave.c ${CODEBOOKS}) +target_link_libraries(tnewamp1 codec2) + +add_executable(compare_ints compare_ints.c) +add_executable(compare_floats compare_floats.c) +add_executable(tvq_mbest tvq_mbest.c) + +add_executable(tfreedv_800XA_rawdata tfreedv_800XA_rawdata.c) +target_link_libraries(tfreedv_800XA_rawdata codec2) + +add_executable(tfreedv_2400A_rawdata tfreedv_2400A_rawdata.c) +target_link_libraries(tfreedv_2400A_rawdata codec2) + +add_executable(tfreedv_2400B_rawdata tfreedv_2400B_rawdata.c) +target_link_libraries(tfreedv_2400B_rawdata codec2) + +add_executable(tfsk_llr tfsk_llr.c) +target_link_libraries(tfsk_llr codec2 m) + +add_executable(thash thash.c) +target_link_libraries(thash codec2 m) + +add_executable(tqam16 tqam16.c) +target_link_libraries(tqam16 codec2 m) + +add_executable(t16_8 t16_8.c ../src/fdmdv.c ../src/kiss_fft.c) +target_link_libraries(t16_8 codec2) + +add_executable(t16_8_short t16_8_short.c ../src/fdmdv.c ../src/kiss_fft.c) +target_link_libraries(t16_8_short codec2) + +add_executable(t48_8 t48_8.c ../src/fdmdv.c ../src/kiss_fft.c) +target_link_libraries(t48_8 codec2) + +add_executable(t48_8_short t48_8_short.c ../src/fdmdv.c ../src/kiss_fft.c) +target_link_libraries(t48_8_short codec2) + +add_executable(tquisk_filter tquisk_filter.c) +target_link_libraries(tquisk_filter codec2) + +# Build CML as part of unit test setup +find_program(OCTAVE_CMD octave-cli REQUIRED) +message("Octave command: ${OCTAVE_CMD}") + +include(ExternalProject) +set(CML_PATH ${CMAKE_CURRENT_BINARY_DIR}/../cml) +ExternalProject_Add(cml + GIT_REPOSITORY https://github.com/drowe67/cml.git + SOURCE_DIR ${CML_PATH} + BUILD_IN_SOURCE 1 + CONFIGURE_COMMAND true # No configuration required + BUILD_COMMAND cd ${CMAKE_CURRENT_BINARY_DIR}/../cml && make + INSTALL_COMMAND true # No installation required +) + +# Create fading files (used for channel simulation) as part of unit test setup +add_custom_target(fading_files ALL + DEPENDS ${CMAKE_CURRENT_BINARY_DIR}/fast_fading_samples.float + DEPENDS ${CMAKE_CURRENT_BINARY_DIR}/faster_fading_samples.float + ) +add_custom_command( + OUTPUT ${CMAKE_CURRENT_BINARY_DIR}/fast_fading_samples.float + OUTPUT ${CMAKE_CURRENT_BINARY_DIR}/faster_fading_samples.float + COMMAND cd ${CMAKE_CURRENT_SOURCE_DIR} && ./fading_files.sh ${CMAKE_CURRENT_BINARY_DIR} +) + +add_executable(freedv_700d_comptx freedv_700d_comptx.c) +add_executable(freedv_700d_comprx freedv_700d_comprx.c) + +if(LPCNET AND lpcnetfreedv_FOUND) + target_link_libraries(freedv_700d_comptx m codec2 lpcnetfreedv) + target_link_libraries(freedv_700d_comprx m codec2 lpcnetfreedv) +else() + target_link_libraries(freedv_700d_comptx m codec2) + target_link_libraries(freedv_700d_comprx m codec2) +endif() + +add_executable(golay23 ../src/golay23.c) +target_compile_options(golay23 PUBLIC -DGOLAY23_UNITTEST) + +add_executable(golay23_runtime_tables ../src/golay23.c) +target_compile_options(golay23_runtime_tables PUBLIC -DGOLAY23_UNITTEST -DRUN_TIME_TABLES) + +add_executable(mksine mksine.c) +target_link_libraries(mksine m) + +add_executable(vq_mbest vq_mbest.c) +target_link_libraries(vq_mbest codec2) + +add_executable(tesno_est tesno_est.c) +target_link_libraries(tesno_est m codec2)
\ No newline at end of file diff --git a/unittest/check_comp.sh b/unittest/check_comp.sh new file mode 100755 index 0000000..b559a2e --- /dev/null +++ b/unittest/check_comp.sh @@ -0,0 +1,29 @@ +#!/usr/bin/env bash +# +# Check Octave and C raw data mode waveforms have about the same +# compression - sanity check for C port of raw data modes. +# +# For manual run outside of ctest: +# cd codec/build_linux +# ../unittest/check_comp.sh ${CODEC2} ${PATH}:${CODEC2}/build_linux/src + +CODEC2=$1 +PATH=$2:$PATH +set -x +octave_log=$(mktemp) +ch_log=$(mktemp) + +echo "warning ('off', 'Octave:data-file-in-path'); + ofdm_ldpc_tx('test_datac0.raw','datac0',1,100,'awgn','bursts',10,'txclip'); + quit" | DISPLAY="" octave-cli -p ${CODEC2}/octave 1>${octave_log} +oct_rms=$(cat ${octave_log} | tr -s ' ' | grep 'RMS:' | cut -d' ' -f4) +oct_cpapr=$(cat ${octave_log} | grep 'RMS:' | tr -s ' ' | cut -d' ' -f6) + +freedv_data_raw_tx datac0 /dev/zero - --delay 1000 --testframes 10 --bursts 10 --clip 1 --txbpf 1 | \ +ch - /dev/null 2>${ch_log} +ch_rms=$(cat ${ch_log} | grep RMS | tr -s ' ' | cut -d' ' -f5) +ch_cpapr=$(cat ${ch_log} | grep RMS | tr -s ' ' | cut -d' ' -f7) + +# Allow 5% difference +python3 -c "import sys; sys.exit(0) if abs((${oct_rms} - ${ch_rms})/${oct_rms}) < 0.05 else sys.exit(1)" +python3 -c "import sys; sys.exit(0) if abs((${oct_cpapr} - ${ch_cpapr})/${oct_cpapr}) < 0.05 else sys.exit(1)" diff --git a/unittest/check_peak.sh b/unittest/check_peak.sh new file mode 100755 index 0000000..b357534 --- /dev/null +++ b/unittest/check_peak.sh @@ -0,0 +1,58 @@ +#!/usr/bin/env bash +# +# Check peak level of each FreeDV waveform is about the same to present +# consistent drive to transmitters. +# +# For manual run outside of ctest: +# cd codec/build_linux +# PATH=${PATH}:${HOME}/codec2/build_linux/src +# ./unittest/check_peak.sh +# OR: +# + +voice_test() { + mode=$1 + echo -n "$mode " + f=$(mktemp) + freedv_tx $mode ../raw/ve9qrp_10s.raw $f --clip 1 + octave_cmd="cd ../octave; + t=load_raw('${f}'); + mx=max(t); printf('%d ',max(t)); + if (mx > 16000) && (mx < 17000) printf('PASS\n') else printf('FAIL\n') end" + octave-cli -qf --eval "$octave_cmd" +} + +data_test() { + mode=$1 + echo -n "$mode " + f=$(mktemp) + freedv_data_raw_tx --framesperburst 2 --bursts 3 --testframes 6 $mode /dev/zero $f 2>/dev/null + octave_cmd="cd ../octave; + t=load_raw('${f}'); + mx=max(t); printf('%d ',max(t)); + if (mx > 16000) && (mx < 17000) printf('PASS\n') else printf('FAIL\n') end" + octave-cli -qf --eval "$octave_cmd" +} + +if [ "$1" == "LPCNet" ]; then + # these don't get run unless we build with LPCNet + voice_test "2020" + voice_test "2020B" + else + voice_test "1600" + voice_test "700C" + voice_test "700D" + voice_test "700E" + voice_test "800XA" + voice_test "2400A" + voice_test "2400B" + data_test "datac0" + data_test "datac1" + data_test "datac3" + data_test "datac4" + data_test "datac13" + data_test "datac14" +fi + +exit 0 + diff --git a/unittest/check_real_comp.sh b/unittest/check_real_comp.sh new file mode 100755 index 0000000..a7b57db --- /dev/null +++ b/unittest/check_real_comp.sh @@ -0,0 +1,15 @@ +#!/usr/bin/env bash +# check_real_comp.sh +# Check the output of freedv_tx() and the real part of freedv_comptx() match, +# as they use different code paths. Run from codec2/unittest, set path to +# include codec2/build/misc and codec2/build/unittest + +set -x +cat ../raw/ve9qrp_10s.raw | freedv_700d_tx > tx_700d.int16 +cat ../raw/ve9qrp_10s.raw | freedv_700d_comptx > tx_700d.iq16 + +echo "tx_real=load_raw('tx_700d.int16'); tx_comp=load_raw('tx_700d.iq16'); \ + tx_comp=tx_comp(1:2:end)+j*tx_comp(2:2:end); \ + diff = sum(real(tx_comp)-tx_real); printf('diff: %f\n', diff); \ + if diff < 1, quit(0), end; \ + quit(1)" | octave-cli -p ../octave -qf diff --git a/unittest/compare_floats.c b/unittest/compare_floats.c new file mode 100644 index 0000000..1d98968 --- /dev/null +++ b/unittest/compare_floats.c @@ -0,0 +1,87 @@ +/* compare floats - a test utility */ + +#include <errno.h> +#include <getopt.h> +#include <math.h> +#include <stdint.h> +#include <stdio.h> +#include <stdlib.h> + +/* Declarations */ + +/* Globals */ + +/* Main */ + +int main(int argc, char *argv[]) { + char usage[] = "Usage: %s [-t tolerance] file1 file2\n"; + + float tol = .001; + + int opt; + while ((opt = getopt(argc, argv, "t:")) != -1) { + switch (opt) { + case 't': + tol = atof(optarg); + break; + default: + fprintf(stderr, usage, argv[0]); + exit(1); + } + } + + if ((optind + 2) > argc) { + fprintf(stderr, usage, argv[0]); + exit(1); + } + char *fname1 = argv[optind++]; + char *fname2 = argv[optind++]; + + FILE *f1 = fopen(fname1, "rb"); + if (f1 == NULL) { + fprintf(stderr, "Error opening file1 \"%s\": ", fname1); + perror(NULL); + exit(1); + } + + FILE *f2 = fopen(fname2, "rb"); + if (f2 == NULL) { + fprintf(stderr, "Error opening file2 \"%s\": ", fname2); + perror(NULL); + exit(1); + } + + float data1, data2; + int count = 0; + int errors = 0; + double rms_sum = 0; + + while (fread(&data1, sizeof(float), 1, f1)) { + if (!fread(&data2, sizeof(float), 1, f2)) { + fprintf(stderr, "Error: file2 is shorter!"); + exit(1); + } + float err = fabsf((data1 - data2) / data1); + if (err > tol) { + errors++; + printf("%d %g %g %g\n", count, data1, data2, err); + } + rms_sum += (err * err); + count++; + } + if (fread(&data2, sizeof(float), 1, f2)) { + fprintf(stderr, "Error: file1 is shorter\n"); + exit(1); + } + + if (errors) { + printf("Fail: %d errors\n", errors); + printf(" rms error = %g\n", ((double)rms_sum / count)); + exit(1); + } else + printf("Pass\n"); + exit(0); + +} // main + +/* vi:set ts=4 et sts=4: */ diff --git a/unittest/compare_ints.c b/unittest/compare_ints.c new file mode 100644 index 0000000..85ef2f3 --- /dev/null +++ b/unittest/compare_ints.c @@ -0,0 +1,160 @@ +/* compare ints - a test utility */ + +#include <errno.h> +#include <getopt.h> +#include <inttypes.h> +#include <stdint.h> +#include <stdio.h> +#include <stdlib.h> + +/* Declarations */ + +/* Globals */ + +/* Functions */ +int get_data(FILE *f, int64_t *dd, int signed_flag, int bytes) { + int res; + int8_t d_8; + int16_t d_16; + uint8_t d_u8; + uint16_t d_u16; + // TODO Loop on reads until, but catch EOF!! + if (signed_flag) { + switch (bytes) { + case 1: + res = fread(&d_8, bytes, 1, f); + *dd = d_8; + break; + case 2: + res = fread(&d_16, bytes, 1, f); + *dd = d_16; + break; + default: + fprintf(stderr, "Error: unsupported size %d bytes\n", bytes); + exit(1); + } + } else { // unsigned + switch (bytes) { + case 1: + res = fread(&d_u8, bytes, 1, f); + *dd = d_u8; + break; + case 2: + res = fread(&d_u16, bytes, 1, f); + *dd = d_u16; + break; + default: + fprintf(stderr, "Error: unsupported size %d bytes\n", bytes); + exit(1); + } + } + + if (res != 1) + return (0); + else + return (1); +} + +/* Main */ + +int main(int argc, char *argv[]) { + char usage[] = + "Usage: %s [-b size_in_bytes] [-c] [-s] [-t tolerance] [-n " + "numerrorstoexit] file1 file2\n"; + + int bytes = 1; + int count_errors = 0; + int signed_flag = 0; + int tol = 1; + int numerrorstoexit = -1; + + int opt; + while ((opt = getopt(argc, argv, "b:cst:n:")) != -1) { + switch (opt) { + case 'b': + bytes = atoi(optarg); + break; + case 'c': + count_errors = 1; + break; + case 's': + signed_flag = 1; + break; + case 'n': + numerrorstoexit = atoi(optarg); + break; + case 't': + tol = atof(optarg); + break; + default: + fprintf(stderr, usage, argv[0]); + exit(1); + } + } + + if ((optind + 2) > argc) { + fprintf(stderr, usage, argv[0]); + exit(1); + } + char *fname1 = argv[optind++]; + char *fname2 = argv[optind++]; + + FILE *f1 = fopen(fname1, "rb"); + if (f1 == NULL) { + fprintf(stderr, "Error opening file1 \"%s\": ", fname1); + perror(NULL); + exit(1); + } + + FILE *f2 = fopen(fname2, "rb"); + if (f2 == NULL) { + fprintf(stderr, "Error opening file2 \"%s\": ", fname2); + perror(NULL); + exit(1); + } + + // Convert inputs to SIGNED long values + int64_t data1, data2; + + int count = 0; + int errors = 0; + int rms_sum = 0; + + while (get_data(f1, &data1, signed_flag, bytes)) { + if (!get_data(f2, &data2, signed_flag, bytes)) { + fprintf(stderr, "Error: file2 is shorter\n"); + exit(1); + } + uint64_t err = llabs(data1 - data2); + if (err > tol) { + errors++; + printf("%d %" PRId64 " %" PRId64 "\n", count, data1, data2); + if (numerrorstoexit != -1) + if (errors > numerrorstoexit) { + printf("reached errors: %d, bailing!", numerrorstoexit); + exit(1); + } + } + rms_sum += (err * err); + count++; + } + if (get_data(f2, &data2, signed_flag, bytes)) { + fprintf(stderr, "Error: file1 is shorter\n"); + exit(1); + } + + if (count_errors) + exit(errors); + else { + if (errors) { + printf("Fail: %d errors\n", errors); + printf(" rms error = %f\n", ((double)rms_sum / count)); + exit(1); + } else + printf("Pass\n"); + exit(0); + } + +} // main + +/* vi:set ts=4 et sts=4: */ diff --git a/unittest/fading_files.sh b/unittest/fading_files.sh new file mode 100755 index 0000000..9fe17f3 --- /dev/null +++ b/unittest/fading_files.sh @@ -0,0 +1,14 @@ +#!/usr/bin/env bash +# +# Generate fading files used for channel simulation + +output_path=$1 +echo "Generating fading files ......" +cmd='cd ../octave; pkg load signal; ch_fading("'${output_path}'/fast_fading_samples.float", 8000, 1.0, 8000*60)' +octave --no-gui -qf --eval "$cmd" +[ ! $? -eq 0 ] && { echo "octave failed to run correctly .... exiting"; exit 1; } +cmd='cd ../octave; pkg load signal; ch_fading("'${output_path}'/faster_fading_samples.float", 8000, 2.0, 8000*60)' +octave --no-gui -qf --eval "$cmd" +[ ! $? -eq 0 ] && { echo "octave failed to run correctly .... exiting"; exit 1; } +exit 0 + diff --git a/unittest/freedv_700d_comprx.c b/unittest/freedv_700d_comprx.c new file mode 100644 index 0000000..0cf8407 --- /dev/null +++ b/unittest/freedv_700d_comprx.c @@ -0,0 +1,140 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: freedv_700d_comprx.c + AUTHOR......: David Rowe + DATE CREATED: July 2022 + + Complex valued rx to support ctests. Includes a few operations that will + only work if complex Tx and Rx signals are being handled correctly. + +\*---------------------------------------------------------------------------*/ + +#include <assert.h> +#include <stdio.h> +#include <stdlib.h> + +#include "codec2_cohpsk.h" +#include "comp_prim.h" +#include "freedv_api.h" +#include "freedv_api_internal.h" +#include "ofdm_internal.h" + +int main(int argc, char *argv[]) { + /* with no arguments then run with no test code */ + int test_num = 0; + if (argc == 2) { + if (strcmp(argv[1], "tx") == 0) { + test_num = 1; + } + if (strcmp(argv[1], "rx") == 0) { + test_num = 2; + } + } + fprintf(stderr, "%d\n", test_num); + + struct freedv *freedv; + freedv = freedv_open(FREEDV_MODE_700D); + assert(freedv != NULL); + + /* note API functions to tell us how big our buffers need to be */ + short speech_out[freedv_get_n_max_speech_samples(freedv)]; + short demod_in[2 * freedv_get_n_max_modem_samples(freedv)]; + COMP demod_in_comp[2 * freedv_get_n_max_modem_samples(freedv)]; + + /* set up small freq offset */ + float foff_hz = 25; + COMP phase_ch; + phase_ch.real = 1.0; + phase_ch.imag = 0.0; + + /* set complex sine wave interferer at -fc */ + COMP interferer_phase = {1.0, 0.0}; + COMP interferer_freq; + interferer_freq.real = + cos(2.0 * M_PI * freedv->ofdm->tx_centre / FREEDV_FS_8000); + interferer_freq.imag = + sin(2.0 * M_PI * freedv->ofdm->tx_centre / FREEDV_FS_8000); + interferer_freq = cconj(interferer_freq); + + /* log a file of demod input samples for plotting in Octave */ + FILE *fdemod = fopen("demod.f32", "wb"); + assert(fdemod != NULL); + + /* measure demod input power, interferer input power */ + float power_d = 0.0; + float power_interferer = 0.0; + + int frames = 0, sum_sync = 0, frames_snr = 0; + float sum_snr = 0.0; + size_t nin, nout; + nin = freedv_nin(freedv); + + while (fread(demod_in, sizeof(short), 2 * nin, stdin) == 2 * nin) { + for (int i = 0; i < nin; i++) { + demod_in_comp[i].real = (float)demod_in[2 * i]; + demod_in_comp[i].imag = (float)demod_in[2 * i + 1]; + // demod_in_comp[i].imag = 0; + } + + if (test_num == 1) { + /* So Tx is a complex OFDM signal centered at +fc. A small + shift fd followed by Re{} will only work if Tx is complex. + If Tx is real, neg freq components at -fc+fd will be + aliased on top of fc+fd wanted signal by Re{} operation. + This can be tested by setting demod_in_comp[i].imag = 0 + above */ + fdmdv_freq_shift_coh(demod_in_comp, demod_in_comp, foff_hz, + FREEDV_FS_8000, &phase_ch, nin); + for (int i = 0; i < nin; i++) demod_in_comp[i].imag = 0.0; + } + + if (test_num == 2) { + /* a complex sinewave (carrier) at -fc will only be ignored if + Rx is treating signal as complex, otherwise if real a +fc + alias will appear in the middle of our wanted signal at + +fc, this can be tested by setting demod_in_comp[i].imag = + 0 below */ + for (int i = 0; i < nin; i++) { + COMP a = fcmult(2E4, interferer_phase); + interferer_phase = cmult(interferer_phase, interferer_freq); + power_interferer += a.real * a.real + a.imag * a.imag; + COMP d = demod_in_comp[i]; + power_d += d.real * d.real + d.imag * d.imag; + demod_in_comp[i] = cadd(d, a); + // demod_in_comp[i].imag = 0; + } + } + + /* useful to take a look at this with Octave */ + fwrite(demod_in_comp, sizeof(COMP), nin, fdemod); + + nout = freedv_comprx(freedv, speech_out, demod_in_comp); + nin = freedv_nin(freedv); /* call me on every loop! */ + fwrite(speech_out, sizeof(short), nout, stdout); + int sync; + float snr_est; + freedv_get_modem_stats(freedv, &sync, &snr_est); + fprintf(stderr, "sync: %d snr_est: %f\n", sync, snr_est); + frames++; + sum_sync += sync; + if (sync) { + sum_snr += snr_est; + frames_snr++; + } + } + + fclose(fdemod); + freedv_close(freedv); + + if (test_num == 2) + fprintf(stderr, "Demod/Interferer power ratio: %3.2f dB\n", + 10 * log10(power_d / power_interferer)); + float snr_av = sum_snr / frames_snr; + fprintf(stderr, "frames: %d sum_sync: %d snr_av: %3.2f dB\n", frames, + sum_sync, snr_av); + + if (snr_av > 8.0) + return 0; + else + return 1; +} diff --git a/unittest/freedv_700d_comptx.c b/unittest/freedv_700d_comptx.c new file mode 100644 index 0000000..9f66c76 --- /dev/null +++ b/unittest/freedv_700d_comptx.c @@ -0,0 +1,44 @@ +/*---------------------------------------------------------------------------*\ + + freedv_comptx.c + + Complex valued Tx to support ctests. + +\*---------------------------------------------------------------------------*/ + +#include <assert.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "freedv_api.h" + +int main(int argc, char *argv[]) { + struct freedv *freedv; + + freedv = freedv_open(FREEDV_MODE_700D); + assert(freedv != NULL); + + /* handy functions to set buffer sizes */ + int n_speech_samples = freedv_get_n_speech_samples(freedv); + short speech_in[n_speech_samples]; + int n_nom_modem_samples = freedv_get_n_nom_modem_samples(freedv); + COMP mod_out[n_nom_modem_samples]; + short mod_out_short[2 * n_nom_modem_samples]; + + /* OK main loop --------------------------------------- */ + + while (fread(speech_in, sizeof(short), n_speech_samples, stdin) == + n_speech_samples) { + freedv_comptx(freedv, mod_out, speech_in); + for (int i = 0; i < n_nom_modem_samples; i++) { + mod_out_short[2 * i] = mod_out[i].real; + mod_out_short[2 * i + 1] = mod_out[i].imag; + } + fwrite(mod_out_short, sizeof(short), 2 * n_nom_modem_samples, stdout); + } + + freedv_close(freedv); + + return 0; +} diff --git a/unittest/hts1a.h b/unittest/hts1a.h new file mode 100644 index 0000000..fa721ff --- /dev/null +++ b/unittest/hts1a.h @@ -0,0 +1,890 @@ +short hts1a_raw[] = { + -14, -14, -8, -7, -11, -15, -14, -16, -24, + -26, -25, -26, -22, -22, -24, -19, -19, -19, + -26, -28, -28, -21, -16, -14, -19, -19, -18, + -18, -16, -18, -26, -28, -35, -28, -19, -12, + -12, -14, -15, -21, -16, -12, -9, -11, -5, + -8, -7, -5, -8, -8, -7, 3, 3, -1, + -2, -5, -1, -7, -5, -4, -4, -7, -5, + -9, -8, -12, -21, -21, -28, -28, -24, -25, + -29, -29, -31, -32, -28, -31, -35, -26, -35, + -31, -28, -32, -26, -21, -22, -16, -15, -14, + -18, -12, -19, -12, -12, -16, -15, -16, -16, + -16, -12, -15, -12, -18, -18, -15, -19, -18, + -16, -14, -15, -16, -16, -16, -14, -16, -11, + -4, -8, -8, -9, -8, -15, -12, -11, -12, + -9, -11, -8, -11, -14, -11, -18, -16, -14, + -14, -9, -5, -8, -15, -16, -14, -16, -18, + -15, -31, -32, -19, -15, -12, -16, -15, -18, + -14, -12, -12, -16, -24, -25, -19, -18, -22, + -21, -19, -16, -16, -14, -16, -24, -24, -19, + -24, -24, -19, -21, -24, -25, -28, -25, -25, + -26, -18, -12, -22, -25, -26, -25, -24, -24, + -22, -21, -19, -16, -15, -12, -12, -14, -8, + -12, -11, -5, -1, 0, 0, -1, -4, -4, + -5, -7, -7, -11, -8, -11, -5, -4, -2, + -8, -5, -12, -14, -14, -14, -12, -12, -7, + -16, -21, -22, -22, -25, -28, -24, -31, -32, + -33, -33, -35, -31, -29, -32, -36, -36, -35, + -35, -32, -26, -35, -29, -22, -18, -11, -16, + -14, -14, -11, -8, -8, -12, -11, -24, -25, + -12, -8, -7, -9, -5, -9, -8, -4, -7, + -5, -11, -11, -4, -5, -5, -8, -4, -8, + -4, -8, -14, -8, -9, -12, -11, -15, -22, + -21, -19, -22, -15, -22, -18, -15, -16, -18, + -12, -14, -21, -19, -16, -19, -21, -22, -21, + -25, -19, -26, -21, -19, -21, -19, -16, -15, + -18, -9, -8, -9, -8, -14, -19, -15, -16, + -16, -16, -12, -11, -12, -11, -11, -7, -19, + -18, -14, -28, -28, -26, -28, -31, -31, -24, + -25, -26, -29, -28, -31, -28, -24, -21, -19, + -21, -16, -24, -24, -18, -16, -19, -12, -9, + -12, -14, -12, -7, -9, -9, -5, -5, -9, + -7, -5, -5, -4, -21, -25, -12, -5, 5, + 10, 10, -1, -2, 3, -4, -7, -8, -5, + -11, -12, -25, -26, -24, -33, -32, -29, -35, + -36, -33, -38, -42, -38, -38, -36, -36, -33, + -32, -38, -32, -28, -22, -18, -14, -9, -15, + -9, -8, -11, -2, -2, -5, -2, 3, -2, + -1, -4, -7, -12, -12, -12, -16, -15, -11, + -11, -8, -7, -5, -9, -12, -12, -18, -19, + -26, -26, -22, -22, -24, -21, -15, -12, -14, + -18, -16, -15, -26, -26, -28, -26, -26, -25, + -28, -25, -14, -12, -14, -18, -24, -14, -9, + -5, -7, -9, -7, -8, -14, -11, -8, -15, + -8, -7, -5, -2, -5, -5, -8, -15, -15, + -16, -33, -35, -25, -15, -14, -18, -22, -18, + -22, -24, -29, -31, -32, -33, -31, -36, -31, + -25, -31, -33, -28, -26, -22, -25, -25, -16, + -18, -16, -15, -15, -14, -11, -7, -5, -4, + -8, -4, -5, -8, -4, 0, -7, -1, -1, + -9, -11, -12, -14, -15, -9, -8, -7, -5, + -11, -12, -12, -19, -15, -16, -16, -12, -18, + -19, -19, -18, -15, -19, -24, -19, -25, -28, + -26, -26, -26, -28, -26, -22, -21, -25, -19, + -16, -16, -9, -9, -12, -7, -11, -7, -12, + -15, -18, -24, -16, -16, -12, -19, -18, -16, + -25, -24, -22, -25, -21, -25, -24, -24, -19, + -15, -26, -19, -14, -19, -16, -18, -26, -28, + -16, -4, 10, 15, 12, 13, 10, 8, 8, + -7, -11, -19, -29, -31, -33, -29, -26, -26, + -26, -22, -15, -16, -15, -19, -18, -19, -12, + -15, -21, -21, -24, -22, -19, -19, -15, -16, + -16, -24, -18, -25, -24, -24, -26, -21, -22, + -22, -18, -16, -19, -22, -22, -19, -24, -19, + -16, -19, -15, -16, -12, -4, -9, -19, -12, + -15, -19, -16, -16, -15, -14, -12, -12, -11, + -8, -9, -9, -12, -11, -11, -8, -9, -5, + 5, -1, -1, -4, -8, -8, -8, -9, -8, + -7, -11, -19, -19, -18, -26, -21, -24, -24, + -26, -35, -32, -36, -31, -26, -28, -26, -25, + -22, -26, -35, -35, -36, -45, -45, -33, -28, + -19, -16, -9, -14, -12, -5, -11, -8, -7, + -5, -5, 2, -2, -5, 0, -7, -5, -11, + -14, -14, -8, -8, -7, -9, -11, -11, -15, + -15, -14, -12, -14, -18, -16, -9, -12, -12, + -14, -16, -22, -25, -26, -31, -29, -26, -29, + -25, -22, -19, -18, -24, -21, -24, -19, -12, + -15, -15, -16, -15, -14, -16, -15, -16, -24, + -19, -22, -26, -24, -19, -18, -19, -15, -11, + -5, -2, -2, -4, -7, -4, -8, -8, -11, + -15, -11, -9, -7, -9, -4, 2, -8, -1, + -2, -12, -9, -15, -21, -31, -38, -32, -32, + -35, -31, -28, -33, -32, -35, -33, -33, -36, + -36, -42, -45, -43, -33, -29, -25, -21, -14, + -12, -11, -16, -15, -12, -14, -15, -15, -9, + -2, -2, 0, -5, -4, -2, 0, 8, 9, + 10, 12, 3, 6, 5, 9, 3, -7, -9, + -25, -32, -25, -11, -1, 2, -1, -9, -15, + -12, -15, -21, -29, -35, -39, -39, -31, -33, + -33, -26, -28, -29, -31, -33, -26, -24, -22, + -24, -21, -18, -15, -18, -26, -25, -22, -18, + -21, -24, -26, -35, -28, -26, -26, -24, -22, + -16, -18, -22, -15, -22, -24, -16, -14, -11, + -4, 3, 5, 3, 8, 8, 6, 0, 6, + 3, -5, 0, 6, -5, -5, -8, -11, -14, + -19, -21, -24, -25, -28, -28, -22, -28, -38, + -38, -26, -26, -22, -32, -31, -26, -18, -12, + -11, -9, -16, -21, -19, -16, -16, -18, -12, + -12, -8, -7, -15, -16, -16, -19, -21, -22, + -22, -22, -25, -32, -29, -35, -32, -33, -33, + -29, -25, -18, -15, -16, -12, -12, -8, -2, + 0, 6, 2, 0, 5, 2, -5, -1, -5, + -8, -8, -8, -9, -15, -12, -11, -5, -9, + -18, -19, -22, -16, -14, -11, -2, -2, -8, + -16, -21, -22, -19, -24, -24, -16, -16, -12, + -7, -8, -5, -8, -16, -24, -29, -35, -36, + -36, -38, -41, -38, -43, -41, -39, -35, -32, + -31, -32, -32, -28, -24, -24, -19, -16, -18, + -12, -2, -1, 5, 5, -9, -25, -21, -15, + -14, -15, -14, -18, -16, -12, -11, -5, -2, + -7, -11, -7, -18, -11, -5, -2, 0, -2, + 0, 0, -5, -11, -14, -12, -14, -14, -16, + -18, -18, -29, -38, -38, -42, -46, -38, -35, + -32, -31, -16, -5, 0, 12, 19, 20, 17, + 20, 20, 16, 16, 12, 5, -7, -15, -14, + -22, -25, -26, -26, -32, -38, -43, -45, -49, + -55, -56, -52, -56, -48, -39, -33, -28, -31, + -24, -29, -32, -26, -16, -7, -14, -11, -11, + -16, -18, -29, -28, -24, -14, -11, -15, -4, + -11, -12, -5, -2, 3, 5, 2, 10, 5, + 6, 2, -1, -7, -7, -12, -14, -16, -22, + -21, -18, -19, -22, -11, -1, 3, 2, 2, + 3, 6, -2, -7, -15, -18, -24, -26, -31, + -38, -41, -39, -36, -39, -33, -26, -24, -18, + -19, -21, -19, -18, -16, -21, -21, -15, -14, + -18, -24, -25, -31, -38, -43, -45, -46, -43, + -39, -33, -28, -19, -11, -8, -4, 5, 12, + 12, 17, 16, 9, 9, 10, 6, 8, -4, + 3, 0, -5, -11, -14, -28, -26, -15, -24, + -32, -32, -28, -32, -28, -18, -22, -22, -15, + -15, -24, -25, -26, -25, -16, -16, -18, -22, + -21, -26, -29, -25, -22, -19, -16, -9, -4, + 2, 6, 10, 3, 2, 0, -7, -7, -14, + -16, -15, -22, -26, -29, -25, -25, -22, -29, + -35, -25, -19, -14, -15, -12, -15, -26, -24, + -29, -28, -26, -26, -29, -32, -38, -42, -38, + -33, -29, -25, -25, -21, -14, -5, 5, 6, + 8, 6, 2, 0, -1, -9, -16, -18, -19, + -22, -22, -21, -15, -18, -22, -11, -8, -11, + -7, 5, 2, -2, -4, -5, -7, -5, 6, + -7, -9, -8, -19, -22, -24, -26, -29, -33, + -29, -25, -24, -21, -21, -24, -29, -28, -26, + -25, -21, -26, -26, -25, -31, -31, -38, -39, + -38, -33, -21, -9, -5, -4, 5, 3, -2, + -7, -8, -9, -15, -16, -16, -22, -21, -24, + -24, -16, -22, -21, -24, -26, -22, -16, -9, + -11, -2, 6, 2, -4, -9, -16, -21, -21, + -21, -24, -22, -22, -19, -18, -16, -12, -14, + -12, -4, -2, -4, -8, -5, -5, -11, 0, + 2, -7, -8, -12, -14, -19, -24, -25, -28, + -9, -5, -29, -33, -22, -22, -21, -15, -18, + -24, -21, -26, -29, -25, -33, -29, -29, -29, + -31, -28, -28, -25, -24, -21, -25, -14, -7, + -5, -16, -19, -21, -28, -33, -38, -36, -26, + -25, -22, -8, -5, 0, 5, 10, 16, 13, + 10, 8, 5, -4, -1, -7, -11, -18, -28, + -31, -42, -43, -38, -38, -22, -11, -8, -7, + 6, 6, 3, 13, 8, -7, -2, -9, -16, + -11, -15, -18, -24, -28, -24, -25, -22, -25, + -28, -25, -38, -39, -35, -36, -14, -25, -42, + -16, -9, -29, -28, -26, -31, -29, -39, -29, + -12, -7, -1, -2, -4, 0, 2, -5, -15, + -21, -35, -32, -22, -19, -18, -19, -5, 6, + 0, -14, -26, -11, 0, -19, -24, -24, -25, + -25, -31, -35, -32, -19, -5, -7, -1, 0, + 6, 8, 0, -8, -5, 9, 9, 6, 2, + 9, 13, -7, -26, -36, -35, -42, -56, -49, + -42, -42, -36, -28, -12, -12, -21, -18, -24, + -19, -22, -25, -24, -21, -18, -15, -15, -8, + -7, 3, -4, -11, -22, -22, -16, -24, -21, + -7, -22, -31, -16, -21, -11, -12, -21, -26, + -28, -19, -28, -31, -25, -38, -38, -29, -33, + -38, -33, -9, 10, 19, 5, -4, -4, -1, + -12, -21, -18, -16, -16, -19, -8, -5, -7, + -2, 0, 12, 16, 15, 2, -14, -12, -22, + -29, -42, -36, -25, -16, -18, -19, -12, 6, + 2, 2, 9, -4, -11, -19, -25, -24, -28, + -35, -43, -35, -25, -41, -45, -42, -39, -35, + -41, -33, -29, -18, -2, -19, -32, -12, 3, + -8, -11, -26, -35, -29, -29, -24, -15, -9, + -8, -1, 2, 0, -1, -2, 2, -1, -8, + -9, -18, -21, -26, -32, -35, -26, -15, -9, + -7, -2, -2, 2, 13, 12, -2, -11, -15, + -18, -28, -29, -24, -28, -32, -31, -28, -15, + -9, -21, -4, 10, -7, -5, -14, -18, -9, + -21, -25, -29, -32, -36, -31, -16, -24, -18, + -18, -22, -25, -32, -31, -12, -8, -8, 2, + -2, 15, 9, -14, -9, -9, 0, -1, -22, + -24, -18, -32, -29, -29, -39, -33, -24, -41, + -33, -18, -33, -28, -32, -28, -24, -35, -32, + -24, -21, -14, -21, -25, -18, -9, -8, -11, + 0, 3, 5, -5, -14, -19, -15, -9, 2, + -2, 0, 12, 10, 19, 22, 8, -2, -5, + -9, -22, -16, -5, -16, -21, -19, -16, -16, + -26, -29, -28, -25, -1, -2, -15, -19, -29, + -35, -35, -38, -36, -38, -29, -46, -56, -38, + -45, -62, -55, -48, -33, -42, -79, -93, 33, + 207, 203, 112, 30, -29, -28, -70, -148, -176, + -96, -69, -87, -32, -45, -24, 46, 64, 76, + 46, 43, 60, 39, 20, 12, 6, -35, -26, + -28, -67, -8, -18, -55, 6, -38, -75, -24, + -5, -5, -25, -25, -35, 22, 93, 23, -28, + -55, -83, -42, -83, -103, -56, -22, 44, 3, + -33, 6, 8, 12, 17, -5, -131, -189, 100, + 265, 40, 49, 135, -59, 12, 51, -123, -87, + -182, -298, -226, -192, -205, -198, -49, 30, 39, + 183, 238, 183, 141, 187, 132, 83, 176, 16, + -116, -90, -118, -138, -189, -123, -137, -147, -9, + -28, 49, 125, 66, 29, 43, 46, -70, -100, + -75, -121, -117, -109, -58, -28, 29, 77, 74, + 128, 118, 63, 12, 56, 158, 90, 0, -29, + -127, -103, -33, -137, -140, -18, -43, -66, -7, + -45, -73, 6, 42, -35, -111, -104, -92, -147, + -526, -1264, -1575, -765, 438, 997, 1207, 1339, 1320, + 1748, 2310, 2055, 1176, 299, -434, -1009, -1257, -1632, + -2439, -2960, -2928, -2615, -2054, -1308, -736, -303, 445, + 1377, 2112, 2587, 2728, 2464, 2061, 1838, 1525, 888, + 182, -511, -1203, -1516, -1425, -1447, -1556, -1386, -1076, + -642, -16, 437, 547, 662, 915, 1092, 1055, 883, + 519, 53, -259, -472, -771, -1037, -1070, -981, -814, + -511, -203, 2, 251, 519, 616, 631, 648, 582, + 476, 391, 271, -18, -305, -348, -399, -478, -389, + -348, -351, -222, -121, -386, -1735, -3442, -2533, 573, + 2049, 2601, 3434, 2708, 2736, 4175, 4209, 2551, 296, + -1119, -2477, -3494, -3361, -4179, -5485, -5636, -4580, -3316, + -2031, -22, 1181, 1763, 3271, 4874, 5510, 5224, 4612, + 3186, 1425, 652, -147, -1598, -2834, -3484, -3818, -3613, + -2537, -1594, -1278, -574, 734, 1997, 2818, 3281, 3233, + 2488, 1960, 1561, 520, -672, -1616, -2288, -2651, -2450, + -2115, -2137, -1772, -998, -341, 312, 891, 1208, 1322, + 1574, 1687, 1285, 810, 427, -18, -404, -664, -1030, + -1306, -1064, -656, -397, -182, -16, 170, 514, 922, + 1023, 956, 915, 626, -42, -1701, -4508, -5074, -1033, + 2559, 2488, 3256, 3525, 2940, 5224, 6208, 3992, -50, + -2480, -3034, -4628, -4714, -5427, -7560, -7667, -5645, -2855, + -1414, 621, 2572, 3066, 5038, 7469, 7951, 6155, 4438, + 2841, 544, -373, -1233, -3376, -5177, -5271, -4414, -3653, + -2121, -788, -392, 720, 2701, 3972, 4022, 3914, 3451, + 2206, 1426, 796, -632, -2163, -3381, -4159, -4113, -3045, + -1776, -1312, -743, 275, 1098, 1851, 2335, 2131, 1527, + 1259, 1200, 713, 73, -516, -1176, -1677, -1755, -1672, + -1670, -1336, -625, 37, 718, 1435, 1620, 1412, 1445, + 1137, 438, 388, 645, 505, 194, -907, -3965, -7195, + -5524, 1111, 4509, 3145, 3535, 3294, 3846, 7424, 7315, + 3111, -1601, -3095, -3366, -4972, -5281, -7332, -9647, -8219, + -4862, -1956, -808, 1125, 2677, 3750, 7120, 9202, 7949, + 5515, 3979, 2538, 716, -4, -1881, -5029, -6072, -5264, + -4448, -3573, -2254, -1472, -802, 1513, 3981, 4536, 4289, + 4083, 3336, 2484, 2141, 963, -1417, -3088, -3937, -4523, + -3991, -2694, -2023, -1694, -593, 580, 1431, 2310, 2417, + 1777, 1433, 1598, 1241, 374, -69, -753, -1490, -1449, + -1357, -1625, -1628, -1176, -845, -356, 723, 1418, 1370, + 1530, 1663, 1222, 1057, 1140, 345, -576, -628, -1435, + -4343, -6010, -2139, 2712, 2651, 2832, 3744, 2781, 5433, + 7295, 4184, 117, -2377, -2990, -4553, -4853, -5022, -7143, + -6439, -4074, -2314, -1340, -260, 822, 720, 2570, 4741, + 4967, 4932, 4619, 4026, 2856, 2222, 1221, -911, -1898, + -2510, -2962, -2593, -1870, -1447, -1553, -1085, -648, -342, + 240, 418, 505, 623, 1054, 1416, 1337, 1193, 544, + -101, -348, -407, -426, -604, -573, -519, -297, 9, + -235, -589, -927, -952, -563, -310, 9, 168, 219, + 369, 245, 40, -235, -451, -287, 46, 415, 683, + 836, 842, 670, 489, 316, 166, 54, -36, -65, + 36, -96, -883, -1693, -1570, -756, -117, 340, 975, + 1918, 2596, 2086, 1068, 199, -509, -625, -662, -529, + 71, -33, -597, -1067, -1512, -1670, -1870, -1901, -1547, + -974, -269, 142, 260, 272, 56, -164, -24, 509, + 1019, 1292, 1496, 1616, 1660, 1613, 1357, 966, 660, + 605, 597, 480, 275, -134, -671, -1134, -1376, -1473, + -1534, -1564, -1496, -1237, -812, -402, -128, 64, 267, + 431, 522, 479, 332, 320, 468, 565, 548, 407, + 180, -63, -264, -455, -638, -685, -574, -358, -32, + 238, 328, 329, 371, 427, 438, 448, 424, 383, + 441, 482, 29, -1189, -2337, -1898, -396, 261, 599, + 1708, 2529, 2347, 1500, 312, -430, -148, 278, 227, + 432, 471, -532, -1710, -2265, -2405, -2323, -2014, -1597, + -1182, -690, -406, -567, -617, -358, -89, 343, 1088, + 1780, 1952, 1833, 1659, 1418, 1391, 1493, 1446, 1302, + 1180, 911, 415, -45, -492, -962, -1186, -1142, -1095, + -1100, -1196, -1422, -1556, -1353, -995, -586, -59, 257, + 287, 287, 197, 3, -36, 156, 418, 737, 932, + 759, 427, 165, -84, -252, -225, -148, -84, 34, + 39, -96, -169, -168, -72, 63, 179, 360, 539, + 485, -67, -1439, -2707, -1793, 379, 1014, 1030, 2144, + 2711, 2106, 1211, -4, -525, 364, 955, 628, 584, + 170, -1432, -2769, -2892, -2525, -2075, -1551, -1390, -1363, + -1025, -944, -1087, -615, 145, 594, 1099, 1758, 1969, + 1712, 1453, 1289, 1384, 1843, 2083, 1800, 1387, 975, + 400, -86, -260, -446, -681, -781, -918, -1183, -1481, + -1768, -1877, -1558, -1020, -515, -28, 192, 9, -261, + -322, -178, 98, 444, 706, 805, 761, 492, 145, + -5, 6, 20, 53, 81, 27, -103, -249, -325, + -219, -9, 139, 291, 473, 547, 400, 42, -904, + -2367, -2330, -46, 1336, 924, 1518, 2371, 1895, 1372, + 564, -271, 454, 1397, 1020, 578, 328, -986, -2542, + -2821, -2398, -1994, -1564, -1564, -1956, -1922, -1598, -1434, + -879, 175, 771, 983, 1429, 1726, 1615, 1534, 1616, + 1721, 2049, 2369, 2073, 1426, 888, 304, -192, -219, + -94, -284, -630, -961, -1359, -1618, -1696, -1748, -1454, + -778, -324, -128, -18, -195, -423, -382, -226, 70, + 483, 701, 628, 428, 210, 10, -22, 114, 227, + 272, 257, 169, 25, -75, -138, -118, 102, 328, + 287, 159, 135, -135, -1224, -2588, -2139, 243, 1339, + 861, 1548, 2238, 1746, 1333, 575, 20, 917, 1538, + 977, 662, 405, -964, -2463, -2640, -2228, -1890, -1544, + -1616, -1970, -1980, -1750, -1567, -867, 258, 771, 834, + 1153, 1445, 1448, 1465, 1600, 1732, 2025, 2283, 1981, + 1387, 883, 305, -131, -60, 112, -138, -472, -713, + -1204, -1632, -1686, -1616, -1274, -676, -513, -574, -314, + -246, -516, -545, -349, -28, 512, 853, 680, 427, + 319, 192, 159, 289, 343, 264, 254, 193, -29, + -143, -148, -107, 152, 371, 271, -62, -913, -2303, + -2531, -572, 992, 788, 1179, 2054, 1888, 1608, 1079, + 285, 625, 1391, 1211, 911, 815, -239, -1829, -2473, + -2303, -1946, -1557, -1548, -1949, -2143, -1976, -1775, -1202, + -178, 458, 614, 910, 1272, 1374, 1395, 1466, 1518, + 1745, 2114, 2090, 1680, 1224, 657, 53, -175, -52, + -45, -203, -454, -892, -1329, -1492, -1447, -1204, -814, + -655, -631, -406, -247, -334, -468, -526, -358, 83, + 466, 547, 455, 343, 228, 230, 309, 335, 371, + 394, 352, 231, 64, -60, -48, 104, 202, -11, + -747, -2004, -2576, -1166, 550, 643, 737, 1650, 1868, + 1642, 1375, 674, 546, 1071, 1077, 943, 1115, 595, + -841, -1963, -2249, -2154, -1810, -1444, -1523, -1822, -1890, + -1813, -1512, -726, 32, 333, 650, 1095, 1306, 1343, + 1429, 1450, 1429, 1627, 1802, 1692, 1440, 1016, 360, + -106, -130, -84, -157, -290, -645, -1068, -1270, -1301, + -1210, -984, -715, -532, -436, -353, -370, -471, -464, + -363, -165, 206, 496, 517, 400, 306, 192, 166, + 284, 373, 383, 386, 383, 319, 224, 135, 127, + 115, -491, -1799, -2388, -1161, 363, 522, 420, 1024, + 1442, 1559, 1517, 1026, 750, 863, 628, 452, 800, + 684, -348, -1405, -1942, -2143, -2044, -1704, -1476, -1461, + -1526, -1621, -1422, -716, 5, 325, 585, 1013, 1288, + 1380, 1476, 1487, 1409, 1382, 1354, 1255, 1174, 968, + 519, 129, -52, -240, -438, -441, -399, -604, -941, + -1258, -1284, -315, 394, -628, -1085, -532, -675, -126, + -277, -1573, -522, 598, 237, 407, 805, 1075, 653, + 318, 362, -59, -56, -240, -233, 312, -489, -2514, + -4669, -4202, -712, 2020, 3805, 4136, 2529, 3440, 5718, + 5515, 3867, 1626, -647, -2193, -2176, -1803, -3276, -5465, + -7012, -7661, -6402, -3818, -2200, -1338, 105, 1535, 3292, + 5837, 7231, 6633, 5443, 4310, 3120, 2167, 863, -1674, + -4124, -5139, -5305, -4852, -3473, -1936, -1080, -52, 1583, + 2812, 3833, 4309, 3278, 2191, 1920, 1330, 168, -1132, + -2559, -3548, -3413, -3045, -3000, -2440, -1105, 424, 1562, + 2003, 1718, 1409, 1794, 1922, 959, -191, -649, -329, + -22, -2232, -7587, -10133, -6123, 3815, 11010, 10048, 7655, + 4465, 7798, 10768, 5766, -1615, -9642, -12287, -10166, -6666, + -6229, -9606, -9109, -5832, -1326, 4475, 7162, 6468, 5722, + 6975, 8255, 7576, 5333, 919, -3576, -5125, -4751, -4305, + -4140, -4257, -4346, -2541, 1238, 4948, 7386, 6832, 4496, + 3434, 2708, 1527, -206, -3454, -6091, -5990, -4237, -2710, + -1963, -1253, -324, 1191, 2740, 2787, 2061, 1796, 1911, + 1408, -65, -1632, -2642, -2612, -2111, -2132, -2119, -703, + 1824, 2113, -3052, -10217, -9285, 2474, 15080, 14696, 9168, + 3569, 1871, 8840, 7733, -863, -10981, -16373, -12835, -7172, + -2606, -4832, -7242, -3728, 1232, 7185, 9862, 7265, 2730, + 1101, 4013, 4560, 2412, -866, -5642, -6610, -3630, -702, + 380, 20, -726, -825, 2400, 6645, 7202, 4332, 728, + -1582, -1799, -535, -1312, -4539, -5701, -4256, -1013, 2246, + 2495, 1214, 859, 1630, 1894, 788, -715, -1510, -1008, + -797, -1376, -2139, -2347, -1166, 64, 874, 1826, 2927, + 1167, -7529, -13865, -7119, 8693, 17405, 11376, 5806, -706, + 3910, 12384, 4898, -5841, -15567, -15268, -9199, -2915, -2, + -6828, -7298, -645, 4931, 9593, 8510, 3404, -1029, 1545, + 5927, 3179, 704, -3058, -7267, -5092, -1456, 523, 578, + 803, 476, 128, 4503, 7132, 4843, 2484, -314, -2661, + -1629, -312, -2667, -4845, -4645, -3556, -106, 3033, 2307, + 815, 1327, 2202, 1608, 238, -1347, -2299, -1676, -1349, + -2026, -2123, -993, 185, 601, 1363, 2536, 1356, -6385, + -14216, -8118, 9341, 18113, 10028, 5714, -662, 2201, 13868, + 5555, -6998, -16122, -13932, -7484, -3936, -511, -7031, -7875, + 1569, 7356, 9852, 7233, 2382, -812, 2168, 7258, 2096, + -2579, -3745, -6101, -3546, -123, 70, -1707, -134, 2587, + 1895, 4233, 5028, 2090, 1861, 1389, -720, -1444, -1278, + -3112, -4277, -3164, -2898, -1520, 861, 1640, 1528, 1967, + 2726, 1695, 329, -492, -1697, -1969, -1929, -2346, -2041, + -600, 480, 449, 1232, 2001, -2082, -11654, -11999, 3894, + 17853, 11683, 5279, 1994, -1245, 12112, 9934, -5184, -14616, + -13918, -5837, -3437, -2238, -7345, -9704, 1457, 9286, 8837, + 4980, 1010, 1222, 4271, 7312, 1819, -5336, -4002, -3308, + -1866, -985, -3300, -3466, -199, 5187, 4104, 1751, 2997, + 2126, 3600, 3670, -678, -3106, -1891, -1261, -2368, -2264, + -3919, -3476, 176, 1251, 1354, 400, 691, 2138, 1896, + 1695, -373, -2231, -1989, -1816, -1622, -1393, -557, -261, + 839, 1234, -6099, -13830, -6704, 12207, 17650, 5219, 4873, + 2652, 5438, 14840, 2038, -11781, -14722, -7010, -2482, -5104, + -8182, -11393, -3796, 9197, 8387, 2298, 752, 2771, 7659, + 8564, 3986, -5080, -5073, 893, -749, -2984, -6554, -6884, + -950, 4891, 4912, -866, 561, 5333, 7261, 5833, -1037, + -4326, -1758, 1554, -393, -5394, -5454, -3454, 179, 2665, + -1061, -2457, 677, 3635, 4058, 1347, -451, -846, -346, + -403, -2221, -2535, -1778, -671, 454, -2159, -10071, -12945, + -179, 17116, 12061, 403, 7550, 7610, 11287, 11279, -5843, + -14361, -7034, -1442, -6719, -11613, -12495, -8241, 4100, 9619, + -682, -1891, 5146, 9786, 11057, 5758, -974, -3229, 3852, + 2712, -6304, -8209, -8290, -4368, 1412, 1119, -2455, -1243, + 5528, 8521, 7421, 3537, -1505, 772, 3590, -341, -4257, + -4940, -4165, -2414, -583, -3011, -4147, 39, 2313, 2528, + 2055, 1207, 1947, 2715, 1007, -2408, -2377, -414, -1230, + -1867, -3633, -9724, -12863, -3041, 14045, 10335, -2909, 9310, + 13964, 12254, 12191, -3180, -10105, -345, -192, -11763, -16662, + -13357, -7749, 648, 2417, -6327, -892, 9776, 11498, 8643, + 4306, 4143, 5477, 7138, 935, -8335, -6130, -4015, -5890, + -6088, -6166, -3580, 1133, 4766, 3210, 4129, 8211, 5810, + 3999, 3155, 502, -186, -1683, -4503, -5169, -3667, -3607, + -4424, -3364, -1626, 496, 2096, 1967, 2167, 3503, 3642, + 1586, -283, -332, -32, -1544, -4229, -7899, -12648, -8414, + 7338, 10181, -4368, 4177, 17458, 13517, 13236, 4243, -5162, + 1964, 2892, -11114, -17844, -12323, -7832, -6496, -4604, -7565, + -3331, 7284, 7605, 3611, 5653, 10209, 9934, 7845, 4885, + -679, -230, -329, -7610, -9785, -6780, -5146, -4804, -3114, + -808, 1753, 6140, 6461, 4152, 5323, 6386, 4387, 857, + -710, -1095, -2162, -3546, -6282, -6035, -3066, -2787, -2956, + -1881, 510, 2940, 3258, 2654, 2028, 2429, 2474, 1033, + -233, -2762, -7034, -11056, -7771, 4646, 5433, -5715, 3990, + 16374, 11856, 10648, 5494, -18, 5462, 1544, -10835, -12677, + -6704, -7160, -10246, -8169, -7359, -3156, 2089, -420, 260, + 6647, 8892, 6911, 7202, 8552, 5729, 4090, 1688, -3204, + -3246, -3532, -7259, -7538, -4539, -3069, -2792, -884, 932, + 3421, 5464, 3526, 3536, 5345, 4133, 2168, 584, -281, + -1044, -2500, -4101, -4951, -3804, -3109, -3549, -2503, -910, + 410, 1643, 2037, 1937, 2188, 2788, 2315, -572, -6503, + -9657, -1302, 7699, -1541, -6680, 9781, 13307, 6279, 7430, + 3430, 3365, 5749, -4342, -10037, -4328, -4039, -10940, -9452, + -4717, -5445, -3381, -2392, -2448, 2723, 4173, 2037, 4710, + 9437, 7904, 4223, 5554, 4404, 1385, -43, -3303, -4253, + -3654, -5194, -5691, -3290, -1456, -2081, -498, 1909, 2226, + 3394, 3325, 2560, 4033, 3472, 1280, 1159, 354, -1428, + -2003, -2912, -3889, -3756, -3522, -3245, -2163, -801, -58, + 556, 1671, 2327, 2137, 1292, -3214, -7879, -365, 8792, + -2058, -6971, 10764, 11376, 2080, 7475, 5194, 2926, 4953, + -5039, -7315, -715, -5646, -12680, -6937, -1755, -6613, -5901, + -2174, -1343, 2238, 575, -682, 5661, 9086, 4503, 2478, + 7910, 6936, 1933, 1266, -19, -835, -2541, -5805, -4897, + -2034, -3317, -5455, -1761, 1654, -167, -301, 2165, 3287, + 3091, 2318, 1799, 3063, 3081, -334, -481, 970, -1720, + -3691, -3239, -2820, -2991, -3206, -2470, -811, 529, 173, + 716, 2623, 1564, -3845, -6581, 2529, 8041, -5370, -4499, + 13582, 7503, -109, 8820, 5684, 4438, 3673, -6724, -2925, + 1272, -9126, -12206, -2956, -2350, -9067, -5194, -1653, -1288, + 390, -2551, 530, 6921, 4657, 798, 5477, 10107, 4739, + 2049, 4361, 2961, 561, -2843, -4325, -1612, -2273, -6381, + -5260, -109, -1162, -3773, -522, 1928, 1069, 328, 1702, + 3025, 2946, 2046, 955, 2895, 2338, -1672, -1136, -243, + -2378, -3644, -2902, -1987, -2112, -1495, -920, 270, 1683, + 595, -992, -4240, -4096, 5163, 3365, -8693, 2644, 13871, + 721, 2637, 10991, 4279, 5176, 1215, -4254, 1188, -1917, + -10538, -7706, -284, -6409, -9875, -1970, -2217, -3228, -2424, + -2575, 1869, 3423, 1156, 1937, 7458, 7611, 2627, 5239, + 6736, 2968, 897, -186, -48, -1270, -3511, -4529, -2943, + -1507, -4284, -3616, -532, -870, -1332, -877, 1520, 2161, + 605, 2045, 3158, 2852, 1636, 554, 1527, 479, -1720, + -2091, -1287, -1420, -2853, -1994, -642, -702, -651, -305, + 809, -1228, -6310, -3284, 6444, 57, -9694, 6692, 12322, + -2264, 4655, 10553, 3995, 5217, -66, -1294, 2509, -3061, + -9247, -4914, 871, -7968, -9271, -1035, -3443, -5114, -3791, + -1539, 1162, 829, 660, 2385, 7135, 5285, 1660, 6675, + 7080, 2594, 1368, 2909, 2358, -1581, -2121, -1548, -1810, + -3079, -4944, -2309, -1422, -3854, -2656, -14, 115, -634, + 735, 2188, 2383, 1821, 1231, 2298, 1998, 6, -378, + 210, -244, -1924, -1735, -317, -1118, -1969, -886, 5, + -368, -410, -812, -3948, -5540, 1489, 4735, -7002, -3001, + 13628, 3594, -3100, 10227, 7590, 2375, 2916, 1211, 1479, + -631, -3134, -6350, -1802, -978, -10207, -6177, -1046, -4499, + -6925, -3423, 1443, -1421, -584, 2444, 2882, 5009, 3294, + 3499, 5559, 5411, 3077, 1096, 4451, 2750, -1703, -631, + 383, -1183, -4035, -2959, -1350, -3069, -3804, -2694, -236, + -904, -1864, 513, 1545, 1463, 455, 1019, 2475, 970, + 230, 551, 798, 211, -916, -468, -216, -352, -1036, + -1298, -138, -232, -750, -580, 40, -1001, -5521, -4688, + 5193, 3492, -9220, 1094, 13917, 1137, -3168, 7172, 9020, + 2400, -2739, 3367, 3104, -2930, -2350, -4632, -2215, -2432, + -6993, -6593, -3793, 228, -5689, -6176, 2457, 1840, -1003, + -1275, 4305, 5439, 333, 2841, 4258, 4558, 3512, 1238, + 3254, 2559, 1884, -421, -1639, 1150, -1170, -3158, -3000, + -1500, -1025, -3848, -2200, -277, -804, -1328, -773, 1228, + 162, -96, 884, 812, 1234, 386, 645, 708, 381, + 633, -624, -126, 146, -872, -957, -586, 223, -698, + -781, 544, 261, 200, -328, -2742, -4169, 1605, 6106, + -3718, -5298, 7649, 6832, -1261, -1390, 4067, 6767, -470, + -1764, 1864, -38, 257, -2297, -3900, -2558, -1958, -1561, + -6583, -4334, 516, -2266, -3226, -2200, 2206, 1794, -1454, + 1915, 2658, 3157, 2576, 1949, 3410, 1896, 3450, 2310, + 117, 1588, 1002, 917, -1493, -1639, 129, -2106, -2237, + -2421, -1808, -1254, -2307, -727, -923, -651, 449, 46, + 471, 268, 1098, 710, -117, 1031, 362, -192, -90, + 398, 471, -775, -349, 271, 125, -277, -403, -31, + 162, 682, -104, -477, 441, 534, 682, 136, -219, + 114, -447, -2159, -2436, 1852, 1779, -2721, -567, 1960, + 1799, 1184, -179, 1457, 1266, 629, 1852, -82, -866, + -345, 61, -696, -3194, -1676, -726, -1898, -1857, -1936, + -516, -494, -253, 808, -301, 554, 1218, 1193, 1286, + 398, 1408, 1364, 843, 1215, 506, 400, 94, 496, + 632, -447, -192, -387, -386, -505, -1087, -390, -591, + -458, -196, -703, -257, -225, 87, 90, -470, 210, + 258, -121, -270, -150, 267, -203, -43, 302, -101, + 47, 244, 209, 63, 142, 353, 64, -158, -79, + 312, -121, -985, -253, 224, -389, -603, -351, -28, + 56, 376, 220, -222, 353, 635, 478, 441, 85, + -219, -768, -570, 876, 771, -478, -189, 482, 529, + 329, 211, 618, 210, -94, 575, 101, -511, -508, + -454, -437, -1064, -996, -845, -1071, -617, -593, -579, + -382, -58, 471, 253, 135, 289, 926, 1221, 374, + 432, 626, 645, 834, 325, 203, 84, -69, -36, + -344, -247, -287, -308, -257, -498, -158, -109, -198, + 51, -60, -96, -138, -56, 98, -11, 54, -73, + -188, -216, -209, -5, -219, -244, -28, -33, 68, + -84, -84, 80, 73, 119, -16, 66, 277, 234, + 102, -55, -16, 47, 2, -93, -193, -113, -55, + 50, -53, -307, -24, 213, 98, 37, -31, -7, + 158, 335, 243, 22, -138, -114, 100, -46, -12, + 59, -209, -89, -387, -869, -338, 349, 214, -232, + 57, 417, 335, 312, 414, 604, 151, 63, 478, + 22, -329, -402, -315, -355, -754, -553, -526, -538, + -332, -424, -199, -28, 5, -42, 94, 415, 271, + 316, 136, 90, 485, 233, 202, 189, -7, 119, + 61, 238, 200, -237, -254, -49, 136, 64, -79, + 20, 129, 93, 9, 163, 95, -203, -49, -24, + -295, -223, -110, -223, -396, -499, -400, -240, -33, + 98, -1, 54, 179, 134, 97, 78, 77, 83, + 197, 296, 85, 70, 180, -124, -99, 192, -128, + -218, 57, -193, -181, 74, -113, -140, -164, -338, + -179, 243, 366, -42, -277, -104, 326, 618, -14, + -451, -82, 427, 529, -158, -205, 119, 80, 131, + -334, -213, 29, -440, -203, -304, -836, -890, -169, + 646, -67, -447, 172, 636, 871, 527, 643, 608, + 262, 642, 401, -199, -569, -523, 23, -393, -1005, + -848, -468, -303, -376, -237, -406, -287, 127, -24, + 32, 107, 294, 410, 155, 369, 735, 1157, 894, + 228, 386, 77, -355, -358, -547, -555, -441, -222, + -250, -179, 108, 136, 345, 287, 162, -33, -328, + -152, -376, -598, -557, -479, -116, -18, -29, -179, + -58, 306, 129, 206, 138, 67, 337, 284, 462, + 255, 127, 401, 97, 33, -124, -90, 194, -120, + -174, -185, -107, -46, -546, -587, -11, 101, -488, + -414, -205, -223, 604, 237, -778, -128, 197, -56, + 346, -42, -471, 543, 625, 42, 398, 270, 352, + 612, 226, -97, -545, -986, -945, 180, 328, -916, + -487, 557, 1254, 725, -351, 751, 1211, 679, 478, + -325, -431, -253, -48, -569, -1721, -1248, -523, -457, + -1384, -1587, 255, 138, -553, -114, -66, 592, 657, + 497, 609, 570, 1263, 1060, 490, 405, 734, 1057, + -32, -239, 210, 68, -104, -819, -495, -219, -666, + -557, -737, -389, -352, -536, -222, -397, 50, 74, + -90, 59, -155, 548, 417, -70, 325, 411, 649, + 248, 204, 645, 186, 219, 112, -70, -109, -446, + -278, -542, -702, -586, -596, -352, -617, -462, -216, + -222, 25, -247, -18, 456, 500, 543, 379, 531, + 935, 1053, 335, 25, 672, 318, -28, -816, -2456, + -1516, 306, -121, -1684, -1292, 1694, 2342, 236, 60, + 2252, 3274, 1020, -188, 1232, 1227, -252, -1636, -1516, + -904, -1924, -2779, -2848, -2019, -1439, -2150, -1854, -683, + 309, 118, -339, 1142, 2092, 1617, 1234, 1834, 2791, + 2007, 1289, 1677, 1789, 1272, 289, 217, 270, -317, + -828, -1357, -1162, -1102, -1626, -1775, -1529, -903, -1022, + -1267, -600, -113, -67, -274, 162, 766, 571, 349, + 359, 1031, 1385, 534, 396, 1299, 1418, 144, -242, + 669, 265, -870, -965, -404, -662, -1479, -1228, -777, + -889, -1008, -862, -344, 42, 110, 54, 404, 1167, + 1098, 571, 942, 1528, 1010, 363, 766, 584, -897, + -1931, -1206, -25, -637, -1718, -533, 1422, 900, -67, + 1286, 2403, 1732, 659, 1009, 1558, 228, -802, -877, + -925, -1503, -2632, -2514, -2009, -2026, -2343, -2281, -906, + -525, -892, -334, 718, 1351, 805, 1249, 2375, 2283, + 1903, 1872, 2457, 2290, 1439, 1279, 1129, 860, 43, + -607, -550, -942, -1455, -1851, -1708, -1513, -1901, -1768, + -1313, -951, -920, -816, -144, 112, 93, 227, 495, + 816, 757, 554, 582, 1324, 1474, 458, 711, 1433, + 798, -86, 145, 626, -513, -1102, -339, -651, -1347, + -1325, -907, -1033, -1323, -828, -642, -444, -113, -33, + 411, 789, 969, 827, 980, 1337, 1112, 1075, 727, + 374, 27, -933, -1200, -696, -79, -913, -1343, 905, + 997, -110, 963, 1732, 1617, 723, 866, 1238, 71, + -341, -839, -1040, -1214, -2238, -2186, -1978, -1808, -2140, + -2121, -819, -845, -894, -189, 489, 830, 605, 1397, + 1894, 1789, 1928, 1903, 2286, 2038, 1530, 1394, 1204, + 891, -12, -249, -227, -934, -1386, -1439, -1292, -1607, + -1796, -1363, -1213, -1135, -1015, -678, -307, -260, -90, + 165, 466, 520, 431, 665, 667, 502, 727, 1115, + 806, 367, 1021, 1166, 173, 83, 558, -42, -816, + -562, -464, -1091, -1190, -917, -975, -1025, -924, -692, + -450, -291, -193, 142, 595, 626, 686, 970, 1082, + 990, 805, 867, 843, 319, 17, -195, -859, -952, + -270, -365, -1005, -184, 815, 199, 328, 1211, 1085, + 894, 789, 633, 384, -89, -576, -976, -924, -1469, + -2074, -1662, -1558, -1801, -1713, -1090, -713, -758, -223, + 244, 524, 791, 986, 1443, 1626, 1612, 1636, 1802, + 1817, 1365, 1225, 1156, 713, 267, 23, -169, -603, + -894, -988, -1146, -1236, -1325, -1250, -1073, -1022, -940, + -717, -448, -366, -295, -28, 160, 192, 305, 522, + 591, 553, 554, 502, 492, 703, 751, 373, 354, + 735, 388, -164, 102, 87, -564, -675, -444, -706, + -992, -816, -678, -750, -689, -505, -259, -150, -116, + 176, 438, 465, 560, 815, 880, 638, 615, 851, + 466, -58, 207, 19, -777, -959, -577, -66, -379, + -678, 444, 939, 363, 565, 1241, 1237, 565, 407, + 598, 155, -564, -972, -849, -1078, -1833, -1883, -1422, + -1428, -1737, -1373, -620, -553, -545, 46, 621, 725, + 734, 1186, 1600, 1530, 1380, 1579, 1777, 1428, 1053, + 1111, 992, 452, 61, 49, -175, -716, -967, -907, + -1025, -1335, -1370, -1115, -1063, -1161, -979, -622, -477, + -475, -219, 156, 240, 243, 503, 785, 734, 650, + 795, 826, 557, 315, 452, 605, 60, -329, 190, + 153, -550, -494, -135, -389, -811, -593, -314, -597, + -669, -386, -167, -210, -266, 78, 299, 187, 206, + 471, 554, 353, 455, 582, 393, 261, 217, 119, + -53, -186, -361, -647, -732, -604, -236, -86, -325, + 388, 1004, 599, 871, 1273, 1200, 945, 645, 653, + 224, -325, -740, -1025, -1149, -1776, -1976, -1733, -1752, + -1781, -1622, -985, -698, -570, 27, 500, 829, 990, + 1336, 1726, 1749, 1748, 1770, 1867, 1667, 1258, 1130, + 910, 447, -38, -286, -509, -991, -1326, -1383, -1461, + -1615, -1615, -1384, -1196, -1080, -838, -482, -175, -14, + 216, 546, 721, 771, 874, 1017, 977, 800, 710, + 674, 476, 139, 50, 37, -290, -491, -485, -547, + -501, -539, -647, -450, -303, -443, -355, 5, 6, + -121, 158, 337, 213, 237, 445, 488, 373, 371, + 309, 216, 257, 166, -93, -67, -25, -487, -577, + -237, -757, -1373, -832, -334, -368, -174, 226, 1135, + 1470, 1102, 1724, 2099, 1709, 1214, 900, 893, -62, + -961, -1193, -1544, -2098, -2837, -2650, -2327, -2490, -2272, + -1748, -890, -508, -179, 742, 1418, 1785, 1945, 2451, + 2852, 2549, 2342, 2263, 2071, 1439, 723, 422, -70, + -780, -1397, -1652, -1795, -2228, -2373, -2126, -1836, -1672, + -1380, -736, -179, 146, 496, 999, 1399, 1449, 1470, + 1602, 1565, 1249, 881, 663, 364, -144, -546, -764, + -1006, -1224, -1274, -1323, -1236, -962, -746, -407, -120, + 80, 463, 721, 747, 843, 993, 900, 669, 618, + 533, 236, -48, -107, -201, -458, -556, -542, -529, + -465, -437, -307, -188, -438, -662, -305, -39, -849, + -1609, -263, 1428, 228, -375, 2323, 3098, 1664, 1745, + 2553, 2535, 1112, 77, -77, -756, -1720, -3222, -3587, + -2878, -3439, -4100, -3459, -1755, -1227, -1445, 83, 1728, + 2314, 2379, 2965, 3941, 3821, 3266, 2732, 2490, 2154, + 755, -389, -760, -1169, -2230, -3181, -2864, -2619, -2872, + -2717, -1946, -927, -441, 68, 922, 1777, 2357, 2318, + 2405, 2664, 2408, 1643, 895, 553, -104, -1131, -1816, + -2166, -2399, -2609, -2615, -2416, -1697, -736, -598, 40, + 1632, 2103, 1854, 2137, 2549, 2345, 1456, 902, 564, + -128, -821, -1510, -1762, -1629, -1827, -1897, -1407, -639, + -298, -256, 342, 1237, 1516, 948, 881, 1300, 285, + -601, -106, -1005, -2152, -765, 366, -225, 105, 1800, + 2667, 2263, 2358, 2609, 1902, 1309, 383, -1178, -2028, + -2431, -3415, -4550, -4205, -3344, -3426, -2845, -1330, 27, + 949, 1939, 3073, 3622, 4112, 4210, 3486, 2987, 2535, + 1477, 2, -942, -1447, -2432, -3129, -3225, -3139, -2816, + -2288, -1557, -826, 190, 1229, 1697, 2171, 2713, 2848, + 2494, 1998, 1574, 810, -38, -777, -1573, -2133, -2414, + -2554, -2640, -2453, -1750, -916, -403, 165, 1489, 2218, + 1984, 2950, 3358, 1847, 1268, 1173, 3, -1163, -1755, + -2176, -2690, -2518, -1972, -1908, -1192, 22, 618, 1019, + 1695, 2107, 1666, 1636, 2010, 1092, 94, -263, -1376, + -1768, -457, -1731, -5151, -3981, 955, 2337, -400, 750, + 5303, 6160, 5173, 3503, 1711, 2310, 1724, -1812, -6147, + -6323, -3861, -5816, -7464, -5990, -3493, -1001, 97, 1540, + 2770, 5367, 7786, 5898, 4386, 4606, 4302, 1937, -1267, + -2217, -3112, -3790, -4215, -5298, -4526, -2551, -773, -546, + -32, 2599, 3965, 3955, 3518, 3064, 2937, 2131, 893, + -1169, -2477, -2268, -2844, -3763, -3940, -2899, -1616, -1023, + -181, 607, 1586, 2617, 3023, 3070, 2914, 2991, 2172, + 623, -310, -1374, -2183, -2752, -3338, -3214, -2667, -1629, + -816, -278, 799, 1765, 2420, 2415, 2095, 1811, 1020, + 509, 175, -542, -1258, -1676, -1622, -1721, -1740, -1068, + -852, -2054, -2617, 929, 5290, 4343, 2568, 4316, 5511, + 5146, 3002, -412, -2276, -2717, -3204, -6159, -8896, -7253, + -4846, -3602, -2969, -1694, 1159, 3995, 6267, 6021, 4739, + 5447, 5599, 3839, 720, -1621, -2390, -3349, -4039, -4975, + -5170, -3395, -1315, 139, 701, 2083, 4056, 4738, 4561, + 3559, 2429, 1602, 507, -1149, -3232, -4042, -3810, -3787, + -3580, -3000, -1677, -99, 1293, 2238, 2416, 2580, 2706, + 3382, 4128, 2082, -743, -502, -239, -2169, -3426, -3827, + -3585, -2181, -999, -751, -658, 1255, 3304, 2794, 1816, + 1619, 1680, 1241, 8, -1322, -2153, -1180, -321, -1526, + -2302, -1478, -72, 363, -230, -426, -82, -103, -77, + 2541, 4581, 3326, 3631, 4197, 3161, 2512, 785, -1567, + -3641, -4417, -4491, -6115, -6640, -5357, -3878, -2181, -630, + 1004, 2187, 3811, 5725, 5413, 4276, 3600, 2869, 1625, + -298, -1819, -3206, -3817, -3270, -3022, -2734, -1844, -175, + 1408, 2225, 2944, 3158, 3145, 3144, 2440, 1057, -421, + -1236, -1857, -2701, -3165, -3225, -2780, -1898, -937, -298, + -206, 87, 1419, 3719, 4816, 3523, 2068, 1482, 881, + -351, -2295, -3896, -4219, -3195, -2016, -1842, -1410, 9, + 1300, 1732, 1636, 1591, 1442, 1292, 1060, 905, 951, + 40, -693, -286, -671, -1437, -1264, -383, 371, 148, + 623, 1538, -1352, -5744, -4090, 1312, 2921, 1799, 1479, + 2488, 5347, 7009, 4782, 210, -1976, -191, -607, -4050, + -6603, -6961, -5451, -3585, -2599, -3395, -3123, 320, 3569, + 4271, 3628, 3813, 4892, 5331, 4548, 2148, -569, -1301, + -1040, -1894, -3494, -4162, -3511, -2160, -805, -235, -93, + 737, 2313, 3162, 2711, 2027, 1466, 942, 786, 684, + -712, -2331, -1677, -1330, -2395, -1867, -539, -52, -62, + 434, 820, 233, 735, 1092, -206, -628, -8, 124, + -546, -681, -390, -814, -535, 90, -5, -80, 381, + 1054, 919, 810, 759, 221, 373, 445, -358, -1278, + -1213, 15, 517, 240, 469, 618, 117, -366, -1994, + -3771, -679, 3229, 1528, -90, 1057, 1545, 2429, 2774, + 1183, -1064, -1624, 94, -409, -2670, -2848, -2292, -1846, + -1057, -812, -1316, -1070, 454, 1003, 468, 747, 1275, + 1268, 1598, 1792, 975, 808, 942, 410, 441, -32, + -622, -261, -436, -603, -474, -709, -641, -368, -305, + -310, -242, -165, -29, 114, 107, 193, 187, 158, + 287, 194, 100, 51, -28, -31, -70, -19, -33, + -143, -165, -225, -240, -244, -273, -283, -269, -118, + -7, 34, 68, 94, 196, 267, 258, 202, 136, + 102, 29, -72, -179, -281, -298, -274, -227, -161, + -15, 162, 227, 236, 268, 332, 335, 237, 132, + 8, -77, -130, -179, -220, -267, -270, -218, -159, + -65, 46, 111, 144, 187, 226, 197, 129, 56, + -14, -93, -179, -259, -280, -256, -236, -198, -168, + -134, -41, 51, 66, 42, 50, 76, 107, 156, + 168, 144, 124, 110, 76, -5, -50, -72, -130, + -219, -284, -253, -181, -118, -75, -56, -9, 77, + 144, 175, 176, 206, 247, 240, 187, 85, -15, + -79, -158, -240, -274, -283, -257, -178, -96, -36, + 44, 131, 196, 206, 180, 175, 134, 61, -26, + -116, -179, -195, -188, -215, -250, -209, -137, -97, + -62, -22, -5, 32, 84, 114, 93, 100, 142, + 135, 84, 33, 3, -8, -9, -22, -60, -79, + -50, -5, 17, 2, -1, 6, 23, 40, 43, + 25, 26, 42, 29, -15, -55, -104, -134, -147, + -203, -244, -220, -175, -118, -58, -5, 44, 95, + 136, 135, 134, 132, 115, 60, -4, -46, -72, + -110, -124, -135, -159, -145, -93, -46, -21, 44, + 107, 141, 151, 132, 111, 110, 100, 44, -15, + -43, -67, -92, -107, -110, -96, -80, -58, -15, + 13, 44, 67, 68, 42, 23, 19, -15, -49, + -86, -101, -121, -131, -141, -159, -165, -152, -120, + -92, -52, 13, 67, 107, 139, 151, 141, 117, + 76, 36, 10, -16, -35, -73, -75, -56, -46, + -49, -58, -46, -15, 17, 40, 43, 53, 59, + 49, 34, 12, -9, -19, -43, -79, -116, -130, + -134, -123, -100, -96, -77, -48, -12, 30, 44, + 59, 60, 47, 40, 25, 20, 22, 15, 6, + -14, -24, -29, -35, -36, -43, -45, -41, -33, + -18, -9, 3, 6, -2, -4, 3, 9, 17, + 20, 17, 10, 3, 5, 3, -12, -24, -31, + -36, -45, -52, -48, -55, -59, -60, -55, -53, + -53, -63, -65, -60, -50, -43, -48, -35, -24, + -25, -35, -41, -39, -33, -16, -2, 12, 27, + 42, 57, 68, 66, 61, 51, 44, 25, 0, + -9, -31, -49, -53, -52, -36, -39, -41, -7, + 27, 25, 15, 27, 32, 23, 20, 5, -18, + -26, -36, -39, -50, -62, -63, -62, -55, -50, + -41, -33, -25, -14, -1, 5, 9, 15, 5, + -7, -18, -26, -33, -46, -56, -73, -82, -66, + -67, -50, -38, -22, -7, 17, 30, 37, 36, + 29, 19, 13, 12, -7, -25, -41, -42, -38, + -49, -49, -35, -26, -8, 6, 0, 0, -2, + -5, -1, -4, -2, -12, -12, -15, -16, -19, + -24, -35, -39, -48, -46, -35, -25, -4, 13, + 25, 27, 3, 2, 30, 27, -1, -1, -16, + -25, -18, -36, -42, -48, -62, -56, -48, -53, + -46, -36, -26, -15, -8, 12, 23, 29, 25, + 13, 5, 5, 0, -16, -12, -19, -25, -24, + -29, -29, -39, -33, -25, -33, -19, -5, -11, + -22, -24, -18, -25, -32, -39, -39, -32, -41, + -36, -43, -39, -18, -11, -7, -12, -8, -4, + -7, 0, 5, -4, 2, 6, -2, -5, -19, + -25, -29, -33, -39, -36, -25, -24, -12, -7, + 5, 13, 16, 34, 29, 3, -1, -1, -14, + -16, -18, -19, -19, -26, -22, -32, -26, -26, + -28, -14, -11, -5, -2, -4, 0, -9, -15, + -16, -28, -28, -32, -48, -56, -52, -55, -56, + -56, -48, -41, -31, -22, -8, -2, -1, 9, + 10, -1, 3, -7, -1, 0, 2, -1, -11, + -16, -28, -38, -45, -41, -42, -24, -15, 5, + 9, -9, -19, -62, -192, -328, 125, 1057, 514, + -579, 277, -138, -1499, -18, 595, -276, 67, -332, + -853, 761, 1310, -301, -237, 513, 762, 540, -164, + 57, -1481, -1500, 396, -1350, -365, 1445, -933, 136, + 1200, -686, 876, 1901, 166, -349, -244, -198, -1073, + -1085, 187, -481, 294, 893, -877, -4, 659, -7, + -7, -178, 73, -62, 463, 223, -523, 500, 209, + 318, 199, -818, 398, -28, -508, 320, -635, -127, + 226, -688, 83, -65, -426, 214, -318, -297, 425, + 134, 29, 311, -21, 115, 316, -417, -84, 127, + -294, 407, 219, -470, -110, 221, -148, -176, 231, + -189, 180, 194, -562, 172, 71, 110, 90, -645, + 107, -58, 267, 646, -467, 349, 628, -366, -205, + -16, -235, -386, -322, -31, 22, 54, 316, -226, + 219, 689, -186, 272, -288, -679, 37, -351, 17, + -284, -130, 95, -312, 347, -15, 563, 439, -696, + 54, -113, -233, -212, -246, 190, -33, 91, -76, + -106, 192, 88, 117, 196, 432, 56, -100, 279, + 0, -65, -243, -412, -49, -193, -361, 121, 104, + -280, 102, -89, -120, 488, -7, -145, 194, 261, + -270, -314, 396, -185, -93, 271, 81, 134, -317, + 13, -158, -549, 445, 388, -53, -2, -550, -332, + 88, -237, -131, -45, 36, 386, -128, -243, 90, + 145, 686, -259, -511, 217, -213, 422, -318, -277, + 809, -339, 131, -73, -577, 354, 296, 527, -158, + -488, -89, -196, 203, -252, -86, 146, -242, 264, + -39, 56, 287, -470, -624, -666, 71, 632, 221, + 292, 189, 390, 39, -600, -80, -182, -79, 182, + 44, -33, -700, 101, 319, -417, 895, 459, -155, + 550, -325, -645, -325, -94, -491, -460, 284, -193, + 182, 669, 325, -58, 122, 520, -504, -423, 169, + 169, 360, -165, -195, -404, -161, 42, -1156, -581, + 449, 434, 301, -41, 405, 544, 298, -87, -716, + 46, 322, -460, -590, -164, -28, -530, -140, 141, + 138, 776, 514, -121, -461, 380, 1170, 257, -32, + 10, -174, -114, -369, -771, -925, 148, 571, -413, + -89, 539, 17, -87, 6, -716, -392, 683, 124, + -436, 233, 473, 66, -198, -329, -478, -26, 588, + 0, -475, 362, 483, -189, -145, -310, -516, -4, + 63, -461, -269, 413, 360, 507, 805, 291, 447, + 580, -383, -723, -312, -141, -270, -638, -824, -416, + -89, -277, -126, 476, 1006, 1010, 898, 1040, 320, + -185, -148, -794, -989, -1019, -972, -450, -210, 17, + -31, 83, 90, -138, 415, 219, -162, 597, 1156, + 1003, 108, -376, -202, -291, -216, -644, -1095, -404, + 282, 248, -150, 145, 371, -104, 305, 261, -512, + -267, -101, -549, -553, 59, 260, 350, 557, -114, + -223, 495, 388, 278, 414, 177, 22, 23, 29, + -12, -65, 90, -250, -756, -60, 533, 142, -76, + -157, -181, 193, -77, -1185, -1054, 54, 333, -124, + -750, -431, 700, 1449, 1017, -877, -1236, 415, 277, + -586, -126, 373, 326, 85, 410, 132, -361, 185, + -293, -440, 915, 1344, 155, -698, 262, 519, -270, + -431, -678, -693, -556, -294, 78, -182, -325, -478, + -675, -55, 248, -404, -395, 415, 907, 614, -206, + -14, 699, 387, -135, -276, 192, 759, 456, -263, + -420, 340, 585, -128, -208, -114, 111, 781, 93, + -675, -99, -7, -560, -835, -647, -638, -199, 853, + 291, -198, 666, 439, -127, -266, -457, -229, -134, + -103, -678, -1277, -120, 782, 461, -101, 124, 985, + 421, -76, -128, -329, 493, 806, 653, 363, 349, + 687, -593, -945, -67, -495, -1020, -1158, -869, -436, + 359, 757, -753, -1234, 42, 680, 217, -504, -137, + 609, 934, 333, -750, -252, 867, 632, -327, -104, + 951, 1184, 548, -529, -618, 565, 969, -179, -967, + -184, 680, 735, -87, -1227, -985, 125, -7, -1102, + -794, 676, 669, -28, -332, -1333, -1237, 513, 1204, + 335, 318, 1758, 420, -2489, -1776, 76, -359, -1987, + -1489, 925, 2204, 2995, 2332, 306, 1254, 3049, 1171, + -1367, -1189, -468, -1397, -2357, -2116, -1536, -1060, -750, + -790, -788, -75, 638, 683, 699, 383, 91, 642, + 911, 258, -376, -362, -271, -167, 241, 296, 308, + 650, 945, 999, 883, 953, 946, 418, -453, -808, + -461, -706, -1258, -988, -434, -52, 66, -123, 151, + 679, 428, -481, -509, 771, 714, -512, -1035, -1486, + -726, -31, -1281, -1212, 1309, 2854, -717, -5541, -2528, + 4411, 4364, -965, -2135, 1605, 6128, 6542, 901, -2670, + 955, 3056, -1970, -6739, -4304, -1117, -2714, -3398, -2501, + -995, 1508, 2242, 360, -608, 1734, 2558, 452, -246, + 119, -35, -26, 448, 260, 60, 1385, 1602, 292, + 541, 1344, 673, -386, -482, -556, -744, -250, -352, + -906, -236, 1023, 238, -1847, -1223, 598, 466, 107, + 156, 44, 587, 1174, 431, -1296, -1827, -647, -169, + -780, -1298, -1319, -43, 2335, 3246, 1074, -2402, -6047, + -7079, -1059, 4745, 1480, -734, 4822, 9405, 9224, 5118, + -67, -1939, -822, -1629, -6357, -8638, -5218, -2163, -2139, + -2847, -1910, 727, 2801, 3025, 1442, 682, 1915, 2640, + 1630, -328, -1203, -831, -587, -253, 102, 145, 388, + 1416, 1997, 1157, 750, 902, 285, -634, -968, -935, + -1489, -1441, -559, -307, 98, 255, -79, 393, 757, + 383, -92, 224, 977, 612, -250, -768, -1037, -1196, + -1434, -1257, -528, 53, -83, 128, 1382, 1759, 1323, + -341, -8125, -13054, -157, 15976, 6087, -6341, 5951, 14103, + 11032, 4159, -6743, -10446, -3027, 1775, -9133, -16626, -6476, + 3326, 1674, -3781, -1209, 4156, 6635, 7114, 2682, -1791, + 1685, 5147, 898, -4427, -3907, -1369, -1268, -375, -530, + -1009, 1666, 3818, 2716, 1439, 2971, 2807, -334, -836, + 548, -1022, -3429, -2739, -1867, -1887, -494, 346, 44, + 972, 2167, 1438, -239, -396, -202, -971, -1043, -826, + -696, 13, 401, -335, -805, 129, 226, -515, 410, + 1208, 717, 88, -5350, -13627, -8527, 10641, 17424, -276, + -2734, 12987, 14423, 8109, -1717, -11484, -8209, -128, -4410, + -16396, -14255, -1867, 5154, 1136, -2316, 3394, 8584, 10376, + 6841, -603, -2497, 2138, 1684, -4646, -7354, -4663, -2188, + -675, 945, 795, 1892, 5103, 6024, 3291, 1198, 1634, + 500, -1851, -2889, -3691, -3739, -1530, 228, -833, -215, + 2780, 3205, 1946, 1246, 642, -103, -849, -1823, -3171, + -2957, -1587, -1554, -1490, -322, 710, 1122, 1530, 1998, + 1508, 689, 476, 669, 754, -720, -4025, -6837, -7325, + -6367, -2698, 6174, 12964, 9699, 7442, 11835, 10703, 3935, + -2704, -7836, -10023, -10369, -10922, -12144, -10017, -3153, 3001, + 3870, 4798, 9021, 10308, 7988, 4810, 1285, -1649}; diff --git a/unittest/mksine.c b/unittest/mksine.c new file mode 100644 index 0000000..a77830c --- /dev/null +++ b/unittest/mksine.c @@ -0,0 +1,54 @@ +/* + mksine.c + David Rowe + 10 Nov 2010 + + Creates a file of sine wave samples. +*/ + +#include <assert.h> +#include <errno.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#define TWO_PI 6.283185307 +#define FS 8000.0 + +int main(int argc, char *argv[]) { + FILE *f; + int i, n; + float freq, length; + short *buf; + float amp = 1E4; + + if (argc < 4) { + printf("usage: %s outputFile frequencyHz lengthSecs [PeakAmp]\n", argv[0]); + exit(1); + } + + if (strcmp(argv[1], "-") == 0) { + f = stdout; + } else if ((f = fopen(argv[1], "wb")) == NULL) { + fprintf(stderr, "Error opening output file: %s: %s.\n", argv[3], + strerror(errno)); + exit(1); + } + freq = atof(argv[2]); + length = atof(argv[3]); + if (argc == 5) amp = atof(argv[4]); + + n = length * FS; + buf = (short *)malloc(sizeof(short) * n); + assert(buf != NULL); + + for (i = 0; i < n; i++) buf[i] = amp * cos(freq * i * (TWO_PI / FS)); + + fwrite(buf, sizeof(short), n, f); + + fclose(f); + free(buf); + + return 0; +} diff --git a/unittest/ofdm_check b/unittest/ofdm_check new file mode 100755 index 0000000..a4b80ed --- /dev/null +++ b/unittest/ofdm_check @@ -0,0 +1,68 @@ +#! /bin/bash + +# ofdm_check +# +# A series of checks of the ofdm functions, mostly decode. +# +# This uses ofdm_mod to supply test data to ofdm_demod and mostly +# assumes that the encode function is correct. + +# Define macros to (later) allow testing alternate versions. +alias OFDM_MOD=ofdm_mod +alias OFDM_DEMOD=ofdm_demod +shopt -s expand_aliases + +# PATH +PATH=$PATH:../src + +PASS=1 + +############################### +echo +echo "Simple test, plain, ideal" +OFDM_MOD --in /dev/zero --testframes 100 | + OFDM_DEMOD --out /dev/null --testframes --verbose 1 2> tmp +cat tmp +p1=$(grep '^BER\.*: 0.000' tmp | wc -l) +p2=$(grep '^BER2\.*: 0.000' tmp | wc -l) +if [[ $p1 -eq 1 && $p2 -eq 1 ]]; then echo "OK"; else echo "BAD"; PASS=0; fi + +############################### +echo +echo "Simple test, plain, AWGN" +OFDM_MOD --in /dev/zero --testframes 100 | + cohpsk_ch - - -20 -Fs 8000 -f -5 | + OFDM_DEMOD --out /dev/null --testframes --verbose 1 2>tmp +cat tmp +n=$(grep '^BER\.*:' tmp | cut -d ' ' -f 2) +p1=$(echo $n '<=' 0.10 | bc) +n=$(grep '^BER2\.*:' tmp | cut -d ' ' -f 2) +p2=$(echo $n '<=' 0.10 | bc) +if [[ $p1 -eq 1 && $p2 -eq 1 ]]; then echo "OK"; else echo "BAD"; PASS=0; fi + +############################### +echo +echo "Simple test, LDPC, ideal" +OFDM_MOD --in /dev/zero --ldpc --testframes 100 | + OFDM_DEMOD --out /dev/null --ldpc --testframes --verbose 1 2>tmp +cat tmp +p1=$(grep '^BER\.*: 0.000' tmp | wc -l) +p2=$(grep '^Coded BER: 0.000' tmp | wc -l) +if [[ $p1 -eq 1 && $p2 -eq 1 ]]; then echo "OK"; else echo "BAD"; PASS=0; fi + +############################### +echo +echo "Simple test, LDPC, AWGN" +OFDM_MOD --in /dev/zero --ldpc --testframes 100 | + cohpsk_ch - - -20 -Fs 8000 -f -5 --fading_dir ../../build_linux/unittest | + OFDM_DEMOD --out /dev/null --ldpc --testframes --verbose 1 2>tmp +cat tmp +n=$(grep '^BER\.*:' tmp | cut -d ' ' -f 2) +p1=$(echo $n '<=' 0.10 | bc) +n=$(grep '^Coded.*BER\.*:' tmp | cut -d ' ' -f 3) +p2=$(echo $n '<=' 0.01 | bc) +if [[ $p1 -eq 1 && $p2 -eq 1 ]]; then echo "OK"; else echo "BAD"; PASS=0; fi + + +echo +if [[ $PASS == 1 ]]; then echo "PASSED"; else echo "FAILED"; fi diff --git a/unittest/ofdm_fade.sh b/unittest/ofdm_fade.sh new file mode 100755 index 0000000..03f290e --- /dev/null +++ b/unittest/ofdm_fade.sh @@ -0,0 +1,12 @@ +#!/usr/bin/env bash +# +# David June 2019 +# Tests 700D OFDM modem fading channel performance, using a simulated channel + +results=$(mktemp) +fading_dir=$1 +# BER should be around 4% for this test (it's better for larger interleavers but no one uses interleaving in practice) +ofdm_mod --in /dev/zero --ldpc 1 --testframes 60 --txbpf | ch - - --No -24 -f -10 --mpp --fading_dir $fading_dir | ofdm_demod --out /dev/null --testframes --verbose 2 --ldpc 1 2> $results +cat $results +cber=$(cat $results | sed -n "s/^Coded BER.* \([0-9..]*\) Tbits.*/\1/p") +python3 -c "import sys; sys.exit(0) if $cber<=0.05 else sys.exit(1)" diff --git a/unittest/ofdm_phase_est_bw.sh b/unittest/ofdm_phase_est_bw.sh new file mode 100755 index 0000000..7eb826a --- /dev/null +++ b/unittest/ofdm_phase_est_bw.sh @@ -0,0 +1,24 @@ +#!/usr/bin/env bash +# +# ofdm_phase_est_bw.sh +# David August 2019 + +# Tests 2020 OFDM modem phase est bandwidth mode option. Locking the +# phase est bandwidth to "high" is useful for High SNR channels with +# fast fading or high phase noise. In this test we show that with +# high bandwidth phase est mode, the BER is < 5% for the "--mpp" (1 +# Hz fading) channel model on a fairly high SNR channel. +# +# To run manually outside of ctest: +# $ cd codec2/unittest +# $ PATH=$PATH:../build_linux/src ./ofdm_phase_est_bw.sh + +results=$(mktemp) +fading_dir=$1 +# BER should be < 5% for this test +ofdm_mod --in /dev/zero --testframes 300 --mode 2020 --ldpc --verbose 0 | \ +ch - - --No -40 -f 10 --ssbfilt 1 --mpp --fading_dir $fading_dir | \ +ofdm_demod --out /dev/null --testframes --mode 2020 --verbose 2 --ldpc --bandwidth 1 2> $results +cat $results +cber=$(cat $results | sed -n "s/^Coded BER.* \([0-9..]*\) Tbits.*/\1/p") +python3 -c "import sys; sys.exit(0) if $cber<=0.05 else sys.exit(1)" diff --git a/unittest/ofdm_time_sync.sh b/unittest/ofdm_time_sync.sh new file mode 100755 index 0000000..fec51c9 --- /dev/null +++ b/unittest/ofdm_time_sync.sh @@ -0,0 +1,30 @@ +#!/usr/bin/env bash +# Shell script version of ofdm_time_sync() +# David June 2019 +# Tests ofdm modem sync time, using real, off air files + +onerun=$(mktemp) +results=$(mktemp) + +# bunch of runs at different time offsets into a sample off air file with low SNR and semi-staionary fading +for start_secs in `seq 0 29`; +do + if [ "$1" = "700D" ]; then + ofdm_demod --in ../wav/vk2tpm_004.wav --out /dev/null --verbose 2 --ldpc \ + --start_secs $start_secs --len_secs 5 2>/dev/null > $onerun + fi; + if [ "$1" = "2020" ]; then + ofdm_demod --mode 2020 --in ../wav/david4.wav --out /dev/null --verbose 2 --ldpc \ + --start_secs $start_secs --len_secs 5 2>/dev/null > $onerun + fi; + [ ! $? -eq 0 ] && { echo "error running ofdm_demod"; exit 1; } + cat $onerun | sed -n "s/time_to_sync: \([0-9..]*\)/\1/p" >> $results +done +# a pass is we never take longer than 5 secs to sync (mean is much smaller) +python3 -c " +import sys; import numpy as np +x=np.loadtxt(\"$results\") +fails=sum(x == -1) +print(\"fails: %d mean: %5.2f var: %5.2f \" % (fails, np.mean(x), np.var(x))) +sys.exit(0) if fails==0 else sys.exit(1) +" diff --git a/unittest/raw_data_curves/Makefile b/unittest/raw_data_curves/Makefile new file mode 100644 index 0000000..09aa7ca --- /dev/null +++ b/unittest/raw_data_curves/Makefile @@ -0,0 +1,149 @@ +# Makefile +# Dec 2022 +# +# Automates PER/BER curve generation for raw data mode: +# +# 1. Compare "ch" noise injection/SNR measurement against reference Octave Tx +# 2. Compare C Tx against Octave Tx (with and without compression) +# 3. Plot curves for SNR estimated from C Rx against actual SNR +# 4. Plot AWGN/PER C Tx curves for end user documentation +# +# usage: +# +# 1. Create 20 minutes of MPP fading samples: +# $ cd codec2/octave/ +# $ octave-cli +# octave:24> pkg load signal +# octave:24> time_secs=60*20 +# octave:26> ch_fading("~/codec2/build_linux/unittest/fast_fading_samples.float", 8000, 1.0, 8000*time_secs) +# +# 2. Run scripts: +# +# $ make + +SHELL := /bin/bash +CODEC2 := $(HOME)/codec2 + +all: test \ + octave_ch_noise_awgn.png octave_c_tx_awgn.png octave_c_tx_comp_awgn.png \ + octave_ch_noise_mpp.png octave_c_tx_mpp.png octave_c_tx_comp_mpp.png \ + snrest_snr_ctx.png snrest_snr_ctxc.png \ + c_tx_comp.png c_tx_comp_thruput.png + +clean: + rm -f *.txt *.png *.raw + +# run this first, traps common CML setup error +test: + source snr_curves.sh; test_ldpc + +# subset of files generated, but enough to set up Makefile dependencies +snr_oct = snr_oct_datac0_awgn.txt snr_oct_datac1_awgn.txt snr_oct_datac3_awgn.txt +snr_ch = snr_ch_datac0_awgn.txt snr_ch_datac1_awgn.txt snr_ch_datac3_awgn.txt +snr_ctx = snr_ctx_datac0_awgn.txt snr_ctx_datac1_awgn.txt snr_ctx_datac3_awgn.txt +snr_ctxc = snr_ctxc_datac0_awgn.txt snr_ctxc_datac3_awgn.txt + +snr_oct_mpp = snr_oct_datac0_mpp.txt snr_oct_datac1_mpp.txt snr_oct_datac3_mpp.txt +snr_ch_mpp = snr_ch_datac0_mpp.txt snr_ch_datac1_mpp.txt snr_ch_datac3_mpp.txt +snr_ctx_mpp = snr_ctx_datac0_mpp.txt snr_ctx_datac1_mpp.txt snr_ctx_datac3_mpp.txt +snr_ctxc_mpp = snr_ctxc_datac0_mpp.txt snr_ctxc_datac3_mpp.txt + +$(snr_oct): + source snr_curves.sh; generate_octave_tx_data datac0 awgn + source snr_curves.sh; generate_octave_tx_data datac1 awgn + source snr_curves.sh; generate_octave_tx_data datac3 awgn + +$(snr_oct_mpp): + source snr_curves.sh; generate_octave_tx_data datac0 mpp + source snr_curves.sh; generate_octave_tx_data datac1 mpp + source snr_curves.sh; generate_octave_tx_data datac3 mpp + +$(snr_ch): + source snr_curves.sh; generate_ch_data datac0 awgn + source snr_curves.sh; generate_ch_data datac1 awgn + source snr_curves.sh; generate_ch_data datac3 awgn + +$(snr_ch_mpp): + source snr_curves.sh; generate_ch_data datac0 mpp + source snr_curves.sh; generate_ch_data datac1 mpp + source snr_curves.sh; generate_ch_data datac3 mpp + +# C without compression + +$(snr_ctx): + source snr_curves.sh; generate_snrest_v_snr_data datac0 awgn + source snr_curves.sh; generate_snrest_v_snr_data datac1 awgn + source snr_curves.sh; generate_snrest_v_snr_data datac3 awgn + source snr_curves.sh; generate_snrest_v_snr_data datac4 awgn + source snr_curves.sh; generate_snrest_v_snr_data datac13 awgn + source snr_curves.sh; generate_snrest_v_snr_data datac14 awgn + +$(snr_ctx_mpp): + source snr_curves.sh; generate_snrest_v_snr_data datac0 mpp + source snr_curves.sh; generate_snrest_v_snr_data datac1 mpp + source snr_curves.sh; generate_snrest_v_snr_data datac3 mpp + source snr_curves.sh; generate_snrest_v_snr_data datac4 mpp + source snr_curves.sh; generate_snrest_v_snr_data datac13 mpp + source snr_curves.sh; generate_snrest_v_snr_data datac14 mpp + +# C with compression + +$(snr_ctxc): + source snr_curves.sh; generate_snrest_v_snr_data datac0 awgn 1 + source snr_curves.sh; generate_snrest_v_snr_data datac1 awgn 1 + source snr_curves.sh; generate_snrest_v_snr_data datac3 awgn 1 + source snr_curves.sh; generate_snrest_v_snr_data datac4 awgn 1 + source snr_curves.sh; generate_snrest_v_snr_data datac13 awgn 1 + source snr_curves.sh; generate_snrest_v_snr_data datac14 awgn 1 + +$(snr_ctxc_mpp): + source snr_curves.sh; generate_snrest_v_snr_data datac0 mpp 1 + source snr_curves.sh; generate_snrest_v_snr_data datac1 mpp 1 + source snr_curves.sh; generate_snrest_v_snr_data datac3 mpp 1 + source snr_curves.sh; generate_snrest_v_snr_data datac4 mpp 1 + source snr_curves.sh; generate_snrest_v_snr_data datac13 mpp 1 + source snr_curves.sh; generate_snrest_v_snr_data datac14 mpp 1 + +# Octave and C curves should be on top of each other, indicating Octave +# and ch noise injection/SNR measurement are equivalent +octave_ch_noise_awgn.png: $(snr_oct) $(snr_ch) + echo "snr_curves_plot; octave_ch_noise_print('awgn'); quit" | \ + octave-cli -p $(CODEC2)/octave +octave_ch_noise_mpp.png: $(snr_oct_mpp) $(snr_ch_mpp) + echo "snr_curves_plot; octave_ch_noise_print('mpp'); quit" | \ + octave-cli -p $(CODEC2)/octave + +# Octave Tx and C Tx curves should be on top of each other +octave_c_tx_awgn.png: $(snr_oct) $(snr_ctx) + echo "snr_curves_plot; octave_c_tx_print('awgn'); quit" | \ + octave-cli -p $(CODEC2)/octave +octave_c_tx_mpp.png: $(snr_oct_mpp) $(snr_ctx_mpp) + echo "snr_curves_plot; octave_c_tx_print('mpp'); quit" | \ + octave-cli -p $(CODEC2)/octave + +# Octave Tx and C Tx (compressed) curves should be close, but C may be 1dB +# poorer +octave_c_tx_comp_awgn.png: $(snr_oc) $(snr_ctxc) + echo "snr_curves_plot; octave_c_tx_comp_print('awgn'); quit" | \ + octave-cli -p $(CODEC2)/octave +octave_c_tx_comp_mpp.png: $(snr_oct_mpp) $(snr_ctxc_mpp) + echo "snr_curves_plot; octave_c_tx_comp_print('mpp'); quit" | \ + octave-cli -p $(CODEC2)/octave + +# combined AWGN and MPP from C Tx (compressed) - what end users would run +c_tx_comp.png: $(snr_ctxc) $(snr_ctxc_mpp) + echo "snr_curves_plot; c_tx_comp_print; quit" | \ + octave-cli -p $(CODEC2)/octave + +# Curves of SNR estimates from C Rx compared to actual SNR, useful for "gear shifting" +snrest_snr_ctx.png: $(snr_ctx) + echo "snr_curves_plot; snrest_snr_print('ctx', 'awgn'); quit" | \ + octave-cli -p $(CODEC2)/octave +snrest_snr_ctxc.png: $(snr_ctxc) + echo "snr_curves_plot; snrest_snr_print('ctxc', 'awgn'); quit" | \ + octave-cli -p $(CODEC2)/octave + +# Throughput of payload data in bits/s of modes against SNR +c_tx_comp_thruput.png: $(snr_ctxc) $(snr_ctxc_mpp) + echo "snr_curves_plot; c_tx_comp_thruput_print; quit" | \ + octave-cli -p $(CODEC2)/octave diff --git a/unittest/raw_data_curves/snr_curves.sh b/unittest/raw_data_curves/snr_curves.sh new file mode 100755 index 0000000..62f0a52 --- /dev/null +++ b/unittest/raw_data_curves/snr_curves.sh @@ -0,0 +1,191 @@ +# snr_curves.sh +# +# Library of bash functions to generate data for SNR curves. +# +# testing a function example: +# $ bash -c "source ./snr_curves.sh; generate_octave_tx_data datac0 awgn" + +set -x + +PATH=${PATH}:${HOME}/codec2/build_linux/src +CODEC2=${HOME}/codec2 +FADING_DIR=${CODEC2}/build_linux/unittest + +snr_list='-5 -4 -3 -2 0 1 2 4' +No_list='-13 -14 -15 -16 -18 -20 -22 -24 -26' +Nbursts_awgn=20 +Nbursts_mpp=100 + +# Octave Tx injects noise and is source of truth for SNR, measure BER/PER v SNR +function generate_octave_tx_data { + mode=$1 + channel=$2 + + Nbursts=$Nbursts_awgn + snr_nudge=0 + if [ "$channel" == "mpp" ]; then + Nbursts=$Nbursts_mpp + snr_nudge=4 + fi + + rx_log=$(mktemp) + + i=1 + rm -f snr_oct_${mode}_${channel}*.txt + rm -f ber_oct_${mode}_${channel}*.txt + rm -f per_oct_${mode}_${channel}*.txt + for snr in $snr_list + do + snr_adj=$((${snr}+${snr_nudge})) + echo "warning ('off', 'Octave:data-file-in-path'); + ofdm_ldpc_tx('test_${mode}.raw','${mode}',1,${snr_adj},'${channel}','bursts',${Nbursts},'crc'); + quit" | DISPLAY="" octave-cli -p ${CODEC2}/octave + freedv_data_raw_rx --testframes $mode test_${mode}.raw /dev/null 2>${rx_log} -v + BERmeas=$(cat ${rx_log} | grep 'BER......:' | cut -d' ' -f2) + PERmeas=$(cat ${rx_log} | grep 'Coded FER' | cut -d' ' -f3) + + echo ${snr_adj} >> snr_oct_${mode}_${channel}.txt + echo ${BERmeas} >> ber_oct_${mode}_${channel}.txt + echo ${PERmeas} >> per_oct_${mode}_${channel}.txt + i=$((i+1)) + done + echo 0 > offset_oct_${mode}_${channel}.txt +} + +# ch injects noise and is source of truth for SNR, measure BER/PER v SNR +# Octave Tx +function generate_ch_data { + mode=$1 + channel=$2 + + ch_multipath='' + Nbursts=$Nbursts_awgn + snr_nudge=0 + if [ "$channel" == "mpp" ]; then + ch_multipath='--mpp' + Nbursts=$Nbursts_mpp + snr_nudge=4 + fi + + octave_log=$(mktemp) + ch_log=$(mktemp) + rx_log=$(mktemp) + + i=1 + rm -f snr_ch_${mode}_${channel}*.txt + rm -f ber_ch_${mode}_${channel}*.txt + rm -f per_ch_${mode}_${channel}*.txt + for No in $No_list + do + No_adj=$((${No}-${snr_nudge})) + echo "warning ('off', 'Octave:data-file-in-path'); + ofdm_ldpc_tx('test_${mode}.raw','${mode}',1,100,'awgn','bursts',${Nbursts},'crc'); + quit" | DISPLAY="" octave-cli -p ${CODEC2}/octave 1>${octave_log} + SNRoffset=$(cat ${octave_log} | grep 'Burst offset:' | cut -d' ' -f5) + + ch test_${mode}.raw - --No $No_adj ${ch_multipath} --fading_dir ${FADING_DIR} 2>>${ch_log} | \ + freedv_data_raw_rx --testframes $mode - /dev/null -v 2>${rx_log} + BERmeas=$(cat ${rx_log} | grep 'BER......:' | cut -d' ' -f2) + PERmeas=$(cat ${rx_log} | grep 'Coded FER' | cut -d' ' -f3) + + echo ${BERmeas} >> ber_ch_${mode}_${channel}.txt + echo ${PERmeas} >> per_ch_${mode}_${channel}.txt + i=$((i+1)) + + # trap not enough fading file samples (with mpp) + grep "Fading file finished" ${ch_log} + if [ $? -eq 0 ]; then + cat ${ch_log} + exit 1 + fi + done + + echo ${SNRoffset} > offset_ch_${mode}_${channel}.txt + SNRch=$(cat ${ch_log} | grep SNR3k | tr -s ' ' | cut -d' ' -f3) + echo ${SNRch} > snr_ch_${mode}_${channel}.txt +} + +# ch injects noise and is source of truth for SNR, measure BER/PER v SNR and +# SNR estimates v SNR from rx, C Tx +function generate_snrest_v_snr_data { + mode=$1 + channel=$2 + + snr_nudge=0 + aNo_list=$No_list + + # nudge SNR test range to get meaningful results for these tests + if [ "$mode" == "datac1" ]; then + snr_nudge=4 + fi + if [[ "$mode" == "datac4" || "$mode" == "datac13" || "$mode" == "datac14" ]]; then + snr_nudge=-6 + fi + + ch_multipath='' + Nbursts=$Nbursts_awgn + if [ "$channel" == "mpp" ]; then + ch_multipath='--mpp' + Nbursts=$Nbursts_mpp + snr_nudge=$((${snr_nudge}+4)) + fi + + clip=0 + id='ctx' + if [ "$#" -eq 3 ]; then + clip=$3 + id='ctxc' + snr_nudge=$((${snr_nudge}-4)) + fi + + tx_log=$(mktemp) + ch_log=$(mktemp) + rx_log=$(mktemp) + + i=1 + rm -f snrest_${id}_${mode}_${channel}*.txt + rm -f ber_${id}_${mode}_${channel}*.txt + rm -f per_${id}_${mode}_${channel}*.txt + for No in $aNo_list + do + No_adj=$((${No}-${snr_nudge})) + freedv_data_raw_tx --clip ${clip} --delay 1000 --txbpf ${clip} --bursts $Nbursts --testframes $Nbursts $mode /dev/zero - 2>${tx_log} | \ + ch - - --No $No_adj ${ch_multipath} --fading_dir ${FADING_DIR} 2>>${ch_log} | \ + freedv_data_raw_rx --testframes $mode - /dev/null 2>${rx_log} -v + SNRoffset=$(cat ${tx_log} | grep "mark:space" | tr -s ' ' | cut -d' ' -f 5) + + SNRest=$(cat ${rx_log} | grep '\-BS\-' | tr -s ' ' | cut -d' ' -f17) + if [ ! -z "$SNRest" ]; then + echo ${SNRest} > snrest_${id}_${mode}_${channel}_${i}.txt + fi + BERmeas=$(cat ${rx_log} | grep 'BER......:' | cut -d' ' -f2) + PERmeas=$(cat ${rx_log} | grep 'Coded FER' | cut -d' ' -f3) + echo ${BERmeas} >> ber_${id}_${mode}_${channel}.txt + echo ${PERmeas} >> per_${id}_${mode}_${channel}.txt + i=$((i+1)) + done + + echo ${SNRoffset} > offset_${id}_${mode}_${channel}.txt + + # trap not enough fading file samples (with mpp) + grep "Fading file finished" ${ch_log} + if [ $? -eq 0 ]; then + cat ${ch_log} + exit 1 + fi + SNRch=$(cat ${ch_log} | grep SNR3k | tr -s ' ' | cut -d' ' -f3) + echo ${SNRch} > snr_${id}_${mode}_${channel}.txt +} + +# Sanity check to make sure Octave/CML is set up OK +function test_ldpc { + echo "ldpcut; quit" | DISPLAY="" octave-cli -p ${CODEC2}/octave + if [ "$?" -ne 0 ]; then + echo "basic octave test failed, you may need to" + echo "(a) run ctests to create build_xxx/cml" + echo "(b) set up ~/.octaverc as per octave/ldpc.m" + exit 1 + else + echo "OK" + fi +} diff --git a/unittest/reliable_text_fade.sh b/unittest/reliable_text_fade.sh new file mode 100755 index 0000000..8e86c2a --- /dev/null +++ b/unittest/reliable_text_fade.sh @@ -0,0 +1,27 @@ +#!/usr/bin/env bash +# +# Tests reliable_text fading channel performance, using a simulated channel + +results=$(mktemp -d) +mode=$1 +snr=$2 +min_text_packets=$3 +clip=$4 +build_folder=$5 +fading_dir=${build_folder}/../unittest +rx=$build_folder/freedv_rx +tx=$build_folder/freedv_tx + +if [ $clip -eq 1 ]; then + clip_args="--txbpf 1 --clip 1" +else + clip_args= +fi + +$tx $mode ../raw/ve9qrp.raw - --reliabletext AB1CDEF $clip_args | $build_folder/ch - - --No $snr --mpp -f -5 --fading_dir $fading_dir > $results/reliable_fade.raw +$rx $mode $results/reliable_fade.raw /dev/null --txtrx $results/reliable_fade.txt --reliabletext +if [ `cat $results/reliable_fade.txt | wc -l` -ge $min_text_packets ]; then + exit 0 +else + exit -1 +fi diff --git a/unittest/sum_debug_alloc b/unittest/sum_debug_alloc new file mode 100755 index 0000000..21e3137 --- /dev/null +++ b/unittest/sum_debug_alloc @@ -0,0 +1,79 @@ +#!/usr/bin/env python3 + +""" sum_debug_alloc + + Sum lines reported from codec2-dev/src/debug_alloc.h and report + + Lines like: + + CALLOC: run_ldpc_decoder(112, 32) + + + MALLOC: other_func(238) + + """ + +import fileinput + +by_func = {} +by_addr = {} + +def new_func(): + rec = {} + rec['cur'] = 0 + rec['max'] = 0 + return(rec) + +def new_addr(): + rec = {} + rec['func'] = 0 # Where allocated, may not be same as where free'd! + rec['size'] = 0 + return(rec) + +for line in fileinput.input(): + + if ((line.startswith("MALLOC:")) or (line.startswith("CALLOC:"))): + if (line.startswith("MALLOC:")): + words = line.translate( str.maketrans("()", " ") ).strip().split() + func = words[1] + addr = words[2] + size = int(words[3]) + + elif (line.startswith("CALLOC:")): + words = line.translate( str.maketrans("(,)", " ") ).strip().split() + func = words[1] + addr = words[2] + size = int(words[3]) * int(words[4]) + + #print("Alloc: {} to {} in {}".format(size, addr, func)) + if (not (func in by_func)): by_func[func] = new_func() + data = by_func[func] + data['cur'] += size + if (data['cur'] > data['max']): + data['max'] = data['cur'] + if (addr in by_addr): + print("Error: duplicate allocation to {} in {}".format(addr, func)) + else: + by_addr[addr] = new_addr() + by_addr[addr]['func'] = func + by_addr[addr]['size'] = size + + elif (line.startswith("FREE:")): + words = line.translate( str.maketrans("(,)", " ") ).strip().split() + func = words[1] + addr = words[2] + if (addr in by_addr): + func_a = by_addr[addr]['func'] + by_func[func_a]['cur'] -= by_addr[addr]['size'] + del(by_addr[addr]) + else: + print("Error: free of unallocated location {} in {}".format(addr, func)) + #print("Free: {}:{} in {} to {}".format(func_a, addr, func, by_func[func_a]['cur'])) + +##### +total = 0 +for func, sum in by_func.items(): + max = by_func[func]['max'] + print("{} = {}".format(func, max)) + total += max +print("Total = {}".format(total)) diff --git a/unittest/t16_8.c b/unittest/t16_8.c new file mode 100644 index 0000000..1463b30 --- /dev/null +++ b/unittest/t16_8.c @@ -0,0 +1,99 @@ +/* + t16_8.c + David Rowe + May 10 2012 + + Unit test for 16 <-> 8 kHz sample rate conversion functions. I + + Evaluated output by plotting using Octave and looking for jaggies: + + pl("../unittest/out16.raw",1,3000) + pl("../unittest/out8.raw",1,3000) + + Listening to it also shows up anything nasty: + + $ play -s -2 -r 16000 out16.raw + $ play -s -2 -r 8000 out8.raw + + */ + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> + +#include "codec2_fdmdv.h" + +#define N8 159 /* processing buffer size at 8 kHz (odd number deliberate) */ +#define N16 (N8 * FDMDV_OS) +#define FRAMES 50 +#define TWO_PI 6.283185307 +#define FS 16000 + +#define SINE + +int main() { + float in8k[FDMDV_OS_TAPS_8K + N8]; + short in8k_short[N8]; + float out16k[N16]; + short out16k_short[N16]; + FILE *f16; + + float in16k[FDMDV_OS_TAPS_16K + N16]; + float out8k[N16]; + short out8k_short[N8]; + FILE *f8, *f8in; + + int i, f, t, t1; + float freq = 800.0; + + f16 = fopen("out16.raw", "wb"); + assert(f16 != NULL); + f8 = fopen("out8.raw", "wb"); + assert(f8 != NULL); + f8in = fopen("in8.raw", "wb"); + assert(f8in != NULL); + + /* clear filter memories */ + for (i = 0; i < FDMDV_OS_TAPS_8K; i++) in8k[i] = 0.0; + for (i = 0; i < FDMDV_OS_TAPS_16K; i++) in16k[i] = 0.0; + + t = t1 = 0; + for (f = 0; f < FRAMES; f++) { +#ifdef DC + for (i = 0; i < N8; i++) in8k[FDMDV_OS_TAPS_8K + i] = 16000.0; +#endif +#ifdef SINE + for (i = 0; i < N8; i++, t++) + in8k[FDMDV_OS_TAPS_8K + i] = + 16000.0 * cos(TWO_PI * t * freq / (FS / FDMDV_OS)); +#endif + for (i = 0; i < N8; i++) in8k_short[i] = (short)in8k[i]; + fwrite(in8k_short, sizeof(short), N8, f8in); + + /* upsample */ + fdmdv_8_to_16(out16k, &in8k[FDMDV_OS_TAPS_8K], N8); + + /* save 16k to disk for plotting and check out */ + for (i = 0; i < N16; i++) out16k_short[i] = (short)out16k[i]; + fwrite(out16k_short, sizeof(short), N16, f16); + + /* add a 6 kHz spurious signal, down sampler should + knock this out */ + for (i = 0; i < N16; i++, t1++) + in16k[i + FDMDV_OS_TAPS_16K] = + out16k[i] + 16000.0 * cos(TWO_PI * t1 * 6000.0 / FS); + + /* downsample */ + fdmdv_16_to_8(out8k, &in16k[FDMDV_OS_TAPS_16K], N8); + + /* save 8k to disk for plotting and check out */ + for (i = 0; i < N8; i++) out8k_short[i] = (short)out8k[i]; + fwrite(out8k_short, sizeof(short), N8, f8); + } + + fclose(f16); + fclose(f8); + fclose(f8in); + return 0; +} diff --git a/unittest/t16_8_short.c b/unittest/t16_8_short.c new file mode 100644 index 0000000..040d9a5 --- /dev/null +++ b/unittest/t16_8_short.c @@ -0,0 +1,90 @@ +/* + t16_8_short.c + David Rowe + 19 August 2014 + + Unit test for 16 <-> 8 kHz sample rate conversion functions. I + evaluated output by plotting using Octave and looking for jaggies: + + pl("../unittest/out16_short.raw",1,3000) + pl("../unittest/out8.raw",1,3000) + + Listening to it also shows up anything nasty: + + $ play -s -2 -r 16000 out16_short.raw + $ play -s -2 -r 8000 out8.raw + + */ + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> + +#include "codec2_fdmdv.h" + +#define N8 159 /* procssing buffer size at 8 kHz */ +#define N16 (N8 * FDMDV_OS) +#define FRAMES 100 +#define TWO_PI 6.283185307 +#define FS 16000 + +#define SINE + +int main() { + short in8k_short[FDMDV_OS_TAPS_8K + N8]; + short out16k_short[N16]; + FILE *f16; + + short in16k_short[FDMDV_OS_TAPS_16K + N16]; + short out8k_short[N16]; + FILE *f8, *f8in; + + int i, f, t, t1; + float freq = 800.0; + + f16 = fopen("out16_short.raw", "wb"); + assert(f16 != NULL); + f8 = fopen("out8_short.raw", "wb"); + assert(f8 != NULL); + f8in = fopen("in8_short.raw", "wb"); + assert(f8in != NULL); + + /* clear filter memories */ + for (i = 0; i < FDMDV_OS_TAPS_8K; i++) in8k_short[i] = 0; + for (i = 0; i < FDMDV_OS_TAPS_16K; i++) in16k_short[i] = 0; + + t = t1 = 0; + for (f = 0; f < FRAMES; f++) { +#ifdef DC + for (i = 0; i < N8; i++) in8k_short[FDMDV_OS_TAPS_8K + i] = 16000.0; +#endif +#ifdef SINE + for (i = 0; i < N8; i++, t++) + in8k_short[FDMDV_OS_TAPS_8K + i] = 8000.0 * cos(TWO_PI * t * freq / FS); +#endif + fwrite(in8k_short, sizeof(short), N8, f8in); + + /* upsample */ + fdmdv_8_to_16_short(out16k_short, &in8k_short[FDMDV_OS_TAPS_8K], N8); + + fwrite(out16k_short, sizeof(short), N16, f16); + + /* add a 6 kHz spurious signal for fun, we want down sampler to + knock this out */ + for (i = 0; i < N16; i++, t1++) + in16k_short[i + FDMDV_OS_TAPS_16K] = + out16k_short[i] + 8000.0 * cos(TWO_PI * t1 * 6000.0 / FS); + + /* downsample */ + fdmdv_16_to_8_short(out8k_short, &in16k_short[FDMDV_OS_TAPS_16K], N8); + + /* save 8k to disk for plotting and check out */ + fwrite(out8k_short, sizeof(short), N8, f8); + } + + fclose(f16); + fclose(f8); + fclose(f8in); + return 0; +} diff --git a/unittest/t48_8.c b/unittest/t48_8.c new file mode 100644 index 0000000..88aa929 --- /dev/null +++ b/unittest/t48_8.c @@ -0,0 +1,104 @@ +/* + t48_8.c + David Rowe + May 10 2012 + + Unit test for 48 to 8 kHz sample rate conversion functions. I + evaluated output by plotting using Octave and looking for jaggies: + + pl("../unittest/out48.raw",1,3000) + pl("../unittest/out8.raw",1,3000) + + Listening to it also shows up anything nasty: + + $ play -s -2 -r 48000 out48.raw + $ play -s -2 -r 8000 out8.raw + + */ + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> + +#include "codec2_fdmdv.h" + +#define N8 180 /* processing buffer size at 8 kHz */ +#define N48 (N8 * FDMDV_OS_48) +#define MEM8 FDMDV_OS_TAPS_48_8K +#define FRAMES 50 +#define TWO_PI 6.283185307 +#define FS 48000 + +#define SINE + +int main() { + float in8k[MEM8 + N8]; + short in8k_short[N8]; + float out48k[N48]; + short out48k_short[N48]; + FILE *f48; + + float in48k[FDMDV_OS_TAPS_48K + N48]; + float out8k[N48]; + short out8k_short[N8]; + FILE *f8, *f8in; + + int i, f, t, t1; + float freq = 800.0; + + f48 = fopen("out48.raw", "wb"); + assert(f48 != NULL); + f8 = fopen("out8.raw", "wb"); + assert(f8 != NULL); + f8in = fopen("in8.raw", "wb"); + assert(f8in != NULL); + + /* clear filter memories */ + + for (i = 0; i < MEM8; i++) in8k[i] = 0.0; + for (i = 0; i < FDMDV_OS_TAPS_48K; i++) in48k[i] = 0.0; + + t = t1 = 0; + for (f = 0; f < FRAMES; f++) { +#ifdef DC + for (i = 0; i < N8; i++) in8k[MEM8 + i] = 16000.0; +#endif +#ifdef SINE + for (i = 0; i < N8; i++, t++) + in8k[MEM8 + i] = 16000.0 * cos(TWO_PI * t * freq / (FS / FDMDV_OS_48)); +#endif + for (i = 0; i < N8; i++) in8k_short[i] = (short)in8k[i]; + fwrite(in8k_short, sizeof(short), N8, f8in); + + /* upsample */ + + fdmdv_8_to_48(out48k, &in8k[MEM8], N8); + + /* save 48k to disk for plotting and check out */ + + for (i = 0; i < N48; i++) out48k_short[i] = (short)out48k[i]; + fwrite(out48k_short, sizeof(short), N48, f48); + + /* add a 10 kHz spurious signal for fun, we want down sampler to + knock this out */ + + for (i = 0; i < N48; i++, t1++) + in48k[i + FDMDV_OS_TAPS_48K] = + out48k[i] + 16000.0 * cos(TWO_PI * t1 * 1E4 / FS); + + /* downsample */ + + fdmdv_48_to_8(out8k, &in48k[FDMDV_OS_TAPS_48K], N8); + + /* save 8k to disk for plotting and check out */ + + for (i = 0; i < N8; i++) out8k_short[i] = (short)out8k[i]; + fwrite(out8k_short, sizeof(short), N8, f8); + } + + fclose(f48); + fclose(f8); + fclose(f8in); + return 0; +} diff --git a/unittest/t48_8_short.c b/unittest/t48_8_short.c new file mode 100644 index 0000000..9469655 --- /dev/null +++ b/unittest/t48_8_short.c @@ -0,0 +1,81 @@ +/* + t48_8_short.c + David Rowe + Dec 2021 +*/ + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> + +#include "codec2_fdmdv.h" + +#define N8 180 /* processing buffer size at 8 kHz */ +#define N48 (N8 * FDMDV_OS_48) +#define MEM8 FDMDV_OS_TAPS_48_8K +#define FRAMES 50 +#define TWO_PI 6.283185307 +#define FS 48000 + +#define SINE + +int main() { + short in8k[MEM8 + N8]; + short out48k[N48]; + FILE *f48; + + short in48k[FDMDV_OS_TAPS_48K + N48]; + short out8k[N48]; + FILE *f8, *f8in; + + int i, f, t, t1; + float freq = 800.0; + + f48 = fopen("out48.raw", "wb"); + assert(f48 != NULL); + f8 = fopen("out8.raw", "wb"); + assert(f8 != NULL); + f8in = fopen("in8.raw", "wb"); + assert(f8in != NULL); + + /* clear filter memories */ + + for (i = 0; i < MEM8; i++) in8k[i] = 0.0; + for (i = 0; i < FDMDV_OS_TAPS_48K; i++) in48k[i] = 0.0; + + t = t1 = 0; + for (f = 0; f < FRAMES; f++) { +#ifdef DC + for (i = 0; i < N8; i++) in8k[MEM8 + i] = 16000.0; +#endif +#ifdef SINE + for (i = 0; i < N8; i++, t++) + in8k[MEM8 + i] = 16000.0 * cos(TWO_PI * t * freq / (FS / FDMDV_OS_48)); +#endif + fwrite(&in8k[MEM8], sizeof(short), N8, f8in); + + /* upsample */ + fdmdv_8_to_48_short(out48k, &in8k[MEM8], N8); + + /* save 48k to disk for plotting and check out */ + fwrite(out48k, sizeof(short), N48, f48); + + /* add a 10 kHz spurious signal for fun, we want down sampler to + knock this out */ + for (i = 0; i < N48; i++, t1++) + in48k[i + FDMDV_OS_TAPS_48K] = + out48k[i] + 16000.0 * cos(TWO_PI * t1 * 1E4 / FS); + + /* downsample */ + fdmdv_48_to_8_short(out8k, &in48k[FDMDV_OS_TAPS_48K], N8); + + /* save 8k to disk for plotting and check out */ + fwrite(out8k, sizeof(short), N8, f8); + } + + fclose(f48); + fclose(f8); + fclose(f8in); + return 0; +} diff --git a/unittest/tcohpsk.c b/unittest/tcohpsk.c new file mode 100644 index 0000000..1519d7d --- /dev/null +++ b/unittest/tcohpsk.c @@ -0,0 +1,324 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tcohpsk.c + AUTHOR......: David Rowe + DATE CREATED: March 2015 + + Tests for the C version of the coherent PSK FDM modem. This program + outputs a file of Octave vectors that are loaded and automatically + tested against the Octave version of the modem by the Octave script + tcohpsk.m + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2015 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 <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "codec2_cohpsk.h" +#include "codec2_fdmdv.h" +#include "cohpsk_defs.h" +#include "cohpsk_internal.h" +#include "comp_prim.h" +#include "fdmdv_internal.h" +#include "noise_samples.h" +#include "octave.h" + +#define FRAMES \ + 30 /* LOG_FRAMES is #defined in cohpsk_internal.h */ +#define SYNC_FRAMES \ + 12 /* sync state uses up extra log storage as we reprocess several times */ +#define FRAMESL \ + (SYNC_FRAMES * FRAMES) /* worst case is every frame is out of sync */ + +#define FOFF 58.7 +#define DFOFF (-0.5 / (float)COHPSK_FS) +#define ESNODB 8 +#define PPM -1500 + +extern float pilots_coh[][PILOTS_NC]; + +int main(int argc, char *argv[]) { + struct COHPSK *coh; + int tx_bits[COHPSK_BITS_PER_FRAME]; + COMP tx_symb[NSYMROWPILOT][COHPSK_NC * COHPSK_ND]; + COMP tx_fdm_frame[COHPSK_M * NSYMROWPILOT]; + COMP ch_fdm_frame[COHPSK_M * NSYMROWPILOT]; + // COMP rx_fdm_frame_bb[M*NSYMROWPILOT]; + // COMP ch_symb[NSYMROWPILOT][COHPSK_NC*COHPSK_ND]; + float rx_bits_sd[COHPSK_BITS_PER_FRAME]; + int rx_bits[COHPSK_BITS_PER_FRAME]; + + int tx_bits_log[COHPSK_BITS_PER_FRAME * FRAMES]; + COMP tx_symb_log[NSYMROWPILOT * FRAMES][COHPSK_NC * COHPSK_ND]; + COMP tx_fdm_frame_log[COHPSK_M * NSYMROWPILOT * FRAMES]; + COMP ch_fdm_frame_log[COHPSK_M * NSYMROWPILOT * FRAMES]; + COMP ch_fdm_frame_log_out[(COHPSK_M * NSYMROWPILOT + 1) * FRAMES]; + // COMP rx_fdm_frame_bb_log[M*NSYMROWPILOT*FRAMES]; + // COMP ch_symb_log[NSYMROWPILOT*FRAMES][COHPSK_NC*COHPSK_ND]; + COMP ct_symb_ff_log[NSYMROWPILOT * FRAMES][COHPSK_NC * COHPSK_ND]; + float rx_amp_log[NSYMROW * FRAMES][COHPSK_NC * COHPSK_ND]; + float rx_phi_log[NSYMROW * FRAMES][COHPSK_NC * COHPSK_ND]; + COMP rx_symb_log[NSYMROW * FRAMES][COHPSK_NC * COHPSK_ND]; + int rx_bits_log[COHPSK_BITS_PER_FRAME * FRAMES]; + + FILE *fout; + int f, r, c, log_r, log_data_r, noise_r, ff_log_r, i; + double foff; + COMP foff_rect, phase_ch; + + struct FDMDV *fdmdv; + // COMP rx_filt[COHPSK_NC*COHPSK_ND][P+1]; + // int rx_filt_log_col_index = 0; + // float env[NT*P]; + // float __attribute__((unused)) rx_timing; + COMP tx_onesym[COHPSK_NC * COHPSK_ND]; + // COMP rx_onesym[COHPSK_NC*COHPSK_ND]; + // int rx_baseband_log_col_index = 0; + // COMP rx_baseband_log[COHPSK_NC*COHPSK_ND][(M+M/P)*NSYMROWPILOT*FRAMES]; + float f_est_log[FRAMES], sig_rms_log[FRAMES], noise_rms_log[FRAMES]; + int f_est_samples; + + int log_bits; + float EsNo, variance; + COMP scaled_noise; + int reliable_sync_bit; + int ch_fdm_frame_log_index, nin_frame, tmp, nout; + + coh = cohpsk_create(); + fdmdv = coh->fdmdv; + assert(coh != NULL); + cohpsk_set_verbose(coh, 1); + + /* these puppies are used for logging data in the bowels on the modem */ + + coh->rx_baseband_log_col_sz = + (COHPSK_M + COHPSK_M / P) * NSYMROWPILOT * FRAMESL; + coh->rx_baseband_log = (COMP *)malloc(sizeof(COMP) * COHPSK_NC * COHPSK_ND * + coh->rx_baseband_log_col_sz); + + coh->rx_filt_log_col_sz = (P + 1) * NSYMROWPILOT * FRAMESL; + coh->rx_filt_log = (COMP *)malloc(sizeof(COMP) * COHPSK_NC * COHPSK_ND * + coh->rx_filt_log_col_sz); + + coh->ch_symb_log_col_sz = COHPSK_NC * COHPSK_ND; + coh->ch_symb_log = (COMP *)malloc(sizeof(COMP) * NSYMROWPILOT * FRAMESL * + coh->ch_symb_log_col_sz); + + coh->rx_timing_log = (float *)malloc(sizeof(float) * NSYMROWPILOT * FRAMESL); + + /* init stuff */ + + log_r = log_data_r = noise_r = log_bits = ff_log_r = f_est_samples = 0; + 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 = + 2NcFs/NcRs(Es/No) = 2Fs/Rs(Es/No) */ + + EsNo = pow(10.0, ESNODB / 10.0); + variance = 2.0 * COHPSK_FS / (COHPSK_RS * EsNo); + // fprintf(stderr, "doff: %e\n", DFOFF); + + /* Main Loop + * ---------------------------------------------------------------------*/ + + for (f = 0; f < FRAMES; f++) { + /* --------------------------------------------------------*\ + Mod + \*---------------------------------------------------------*/ + + cohpsk_get_test_bits(coh, tx_bits); + bits_to_qpsk_symbols(tx_symb, (int *)tx_bits, COHPSK_BITS_PER_FRAME); + + for (r = 0; r < NSYMROWPILOT; r++) { + for (c = 0; c < COHPSK_NC * COHPSK_ND; c++) tx_onesym[c] = tx_symb[r][c]; + tx_filter_and_upconvert_coh( + &tx_fdm_frame[r * COHPSK_M], COHPSK_NC * COHPSK_ND, tx_onesym, + fdmdv->tx_filter_memory, fdmdv->phase_tx, fdmdv->freq, + &fdmdv->fbb_phase_tx, fdmdv->fbb_rect); + } + cohpsk_clip(tx_fdm_frame, COHPSK_CLIP, NSYMROWPILOT * COHPSK_M); + + /* --------------------------------------------------------*\ + Channel + \*---------------------------------------------------------*/ + + for (r = 0; r < NSYMROWPILOT * COHPSK_M; r++) { + foff_rect.real = cos(2.0 * M_PI * foff / COHPSK_FS); + foff_rect.imag = sin(2.0 * M_PI * foff / COHPSK_FS); + foff += DFOFF; + phase_ch = cmult(phase_ch, foff_rect); + ch_fdm_frame[r] = cmult(tx_fdm_frame[r], phase_ch); + } + phase_ch.real /= cabsolute(phase_ch); + phase_ch.imag /= cabsolute(phase_ch); + // fprintf(stderr, "%f\n", foff); + for (r = 0; r < NSYMROWPILOT * COHPSK_M; r++, noise_r++) { + scaled_noise = fcmult(sqrt(variance), noise[noise_r]); + ch_fdm_frame[r] = cadd(ch_fdm_frame[r], scaled_noise); + } + + /* optional gain to test level sensitivity */ + + for (r = 0; r < NSYMROWPILOT * COHPSK_M; r++) { + ch_fdm_frame[r] = fcmult(1.0, ch_fdm_frame[r]); + } + + /* tx vector logging */ + + memcpy(&tx_bits_log[COHPSK_BITS_PER_FRAME * f], tx_bits, + sizeof(int) * COHPSK_BITS_PER_FRAME); + memcpy(&tx_fdm_frame_log[COHPSK_M * NSYMROWPILOT * f], tx_fdm_frame, + sizeof(COMP) * COHPSK_M * NSYMROWPILOT); + memcpy(&ch_fdm_frame_log[COHPSK_M * NSYMROWPILOT * f], ch_fdm_frame, + sizeof(COMP) * COHPSK_M * NSYMROWPILOT); + + for (r = 0; r < NSYMROWPILOT; r++, log_r++) { + for (c = 0; c < COHPSK_NC * COHPSK_ND; c++) { + tx_symb_log[log_r][c] = tx_symb[r][c]; + } + } + } + + /* Fs offset simulation */ + + nout = cohpsk_fs_offset(ch_fdm_frame_log_out, ch_fdm_frame_log, + COHPSK_M * NSYMROWPILOT * FRAMES, PPM); + assert(nout < (COHPSK_M * NSYMROWPILOT + 1) * FRAMES); + + nin_frame = COHPSK_NOM_SAMPLES_PER_FRAME; + ch_fdm_frame_log_index = 0; + + /* --------------------------------------------------------*\ + Demod + \*---------------------------------------------------------*/ + + for (f = 0; f < FRAMES; f++) { + coh->frame = f; + + // printf("nin_frame: %d\n", nin_frame); + + assert(ch_fdm_frame_log_index < COHPSK_M * NSYMROWPILOT * FRAMES); + tmp = nin_frame; + cohpsk_demod(coh, rx_bits_sd, &reliable_sync_bit, + &ch_fdm_frame_log_out[ch_fdm_frame_log_index], &nin_frame); + for (i = 0; i < COHPSK_BITS_PER_FRAME; i++) + rx_bits[i] = rx_bits_sd[i] < 0.0; + + ch_fdm_frame_log_index += tmp; + + /* --------------------------------------------------------*\ + Log each vector + \*---------------------------------------------------------*/ + + if (coh->sync == 1) { + for (r = 0; r < NSYMROWPILOT; r++, ff_log_r++) { + for (c = 0; c < COHPSK_NC * COHPSK_ND; c++) { + ct_symb_ff_log[ff_log_r][c] = coh->ct_symb_ff_buf[r][c]; + } + } + + for (r = 0; r < NSYMROW; r++, log_data_r++) { + for (c = 0; c < COHPSK_NC * COHPSK_ND; c++) { + rx_amp_log[log_data_r][c] = coh->amp_[r][c]; + rx_phi_log[log_data_r][c] = coh->phi_[r][c]; + rx_symb_log[log_data_r][c] = coh->rx_symb[r][c]; + } + } + memcpy(&rx_bits_log[COHPSK_BITS_PER_FRAME * log_bits], rx_bits, + sizeof(int) * COHPSK_BITS_PER_FRAME); + log_bits++; + f_est_log[f_est_samples] = coh->f_est; + sig_rms_log[f_est_samples] = coh->sig_rms; + noise_rms_log[f_est_samples] = coh->noise_rms; + f_est_samples++; + ; + } + + assert(log_r <= NSYMROWPILOT * FRAMES); + assert(noise_r <= NSYMROWPILOT * COHPSK_M * FRAMES); + assert(log_data_r <= NSYMROW * FRAMES); + + printf("\r [%d]", f + 1); + } + printf("\n"); + + /*---------------------------------------------------------*\ + Dump logs to Octave file for evaluation + by tcohpsk.m Octave script + \*---------------------------------------------------------*/ + + fout = fopen("tcohpsk_out.txt", "wt"); + assert(fout != NULL); + fprintf(fout, "# Created by tcohpsk.c\n"); + octave_save_int(fout, "tx_bits_log_c", tx_bits_log, 1, + COHPSK_BITS_PER_FRAME * FRAMES); + octave_save_complex(fout, "tx_symb_log_c", (COMP *)tx_symb_log, + NSYMROWPILOT * FRAMES, COHPSK_NC * COHPSK_ND, + COHPSK_NC * COHPSK_ND); + octave_save_complex(fout, "tx_fdm_frame_log_c", (COMP *)tx_fdm_frame_log, 1, + COHPSK_M * NSYMROWPILOT * FRAMES, + COHPSK_M * NSYMROWPILOT * FRAMES); + octave_save_complex(fout, "ch_fdm_frame_log_c", (COMP *)ch_fdm_frame_log_out, + 1, nout - 1, nout - 1); + // 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 * COHPSK_ND, coh->rx_baseband_log_col_index, + coh->rx_baseband_log_col_sz); + octave_save_complex(fout, "rx_filt_log_c", (COMP *)coh->rx_filt_log, + COHPSK_NC * COHPSK_ND, coh->rx_filt_log_col_index, + coh->rx_filt_log_col_sz); + octave_save_complex(fout, "ch_symb_log_c", (COMP *)coh->ch_symb_log, + coh->ch_symb_log_r, COHPSK_NC * COHPSK_ND, + COHPSK_NC * COHPSK_ND); + octave_save_float(fout, "rx_timing_log_c", (float *)coh->rx_timing_log, 1, + coh->rx_timing_log_index, coh->rx_timing_log_index); + octave_save_complex(fout, "ct_symb_ff_log_c", (COMP *)ct_symb_ff_log, + NSYMROWPILOT * FRAMES, COHPSK_NC * COHPSK_ND, + COHPSK_NC * COHPSK_ND); + octave_save_float(fout, "rx_amp_log_c", (float *)rx_amp_log, log_data_r, + COHPSK_NC * COHPSK_ND, COHPSK_NC * COHPSK_ND); + octave_save_float(fout, "rx_phi_log_c", (float *)rx_phi_log, log_data_r, + COHPSK_NC * COHPSK_ND, COHPSK_NC * COHPSK_ND); + octave_save_complex(fout, "rx_symb_log_c", (COMP *)rx_symb_log, log_data_r, + COHPSK_NC * COHPSK_ND, COHPSK_NC * COHPSK_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); + octave_save_float(fout, "sig_rms_log_c", sig_rms_log, 1, f_est_samples, + f_est_samples - 1); + octave_save_float(fout, "noise_rms_log_c", noise_rms_log, 1, f_est_samples, + f_est_samples); +#ifdef XX +#endif + fclose(fout); + + cohpsk_destroy(coh); + + return 0; +} diff --git a/unittest/tesno_est.c b/unittest/tesno_est.c new file mode 100644 index 0000000..8937b1b --- /dev/null +++ b/unittest/tesno_est.c @@ -0,0 +1,32 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tesno_est.c + AUTHORS.....: David Rowe + DATE CREATED: Mar 2021 + + Test for C port of Es/No estimator. + +\*---------------------------------------------------------------------------*/ + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> + +#include "ofdm_internal.h" + +int main(int argc, char *argv[]) { + FILE *fin = fopen(argv[1], "rb"); + assert(fin != NULL); + size_t nsym = atoi(argv[2]); + assert(nsym >= 0); + complex float rx_sym[nsym]; + size_t nread = fread(rx_sym, sizeof(complex float), nsym, fin); + assert(nread == nsym); + fclose(fin); + + float EsNodB = ofdm_esno_est_calc(rx_sym, nsym); + printf("%f\n", EsNodB); + + return 0; +} diff --git a/unittest/test_700c_eq.sh b/unittest/test_700c_eq.sh new file mode 100755 index 0000000..1a106f0 --- /dev/null +++ b/unittest/test_700c_eq.sh @@ -0,0 +1,12 @@ +#!/usr/bin/env bash +# test_700c_eq.sh +# make sure 700C EQ is reducing VQ distortion + +results=$(mktemp) + +c2enc 700C ../raw/kristoff.raw /dev/null --var 2> $results +var=$(cat $results | sed -n "s/.*var: \([0-9..]*\) .*/\1/p") +c2enc 700C ../raw/kristoff.raw /dev/null --var --eq 2> $results +var_eq=$(cat $results | sed -n "s/.*var: \([0-9..]*\) .*/\1/p") +printf "var: %5.2f var_eq: %5.2f\n" $var $var_eq +python3 -c "import sys; sys.exit(0) if $var_eq<=$var else sys.exit(1)" diff --git a/unittest/tfdmdv.c b/unittest/tfdmdv.c new file mode 100644 index 0000000..e7fbca8 --- /dev/null +++ b/unittest/tfdmdv.c @@ -0,0 +1,319 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfdmdv.c + AUTHOR......: David Rowe + DATE CREATED: April 16 2012 + + Tests for the C version of the FDMDV modem. This program outputs a + file of Octave vectors that are loaded and automatically tested + against the Octave version of the modem by the Octave script + tfmddv.m + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2012 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 <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "codec2_fdmdv.h" +#include "fdmdv_internal.h" +#include "octave.h" + +#define FRAMES 35 +#define CHANNEL_BUF_SIZE (10 * M_FAC) + +extern float pilot_coeff[]; + +int main(int argc, char *argv[]) { + struct FDMDV *fdmdv; + int tx_bits[FDMDV_BITS_PER_FRAME]; + COMP tx_symbols[FDMDV_NC + 1]; + COMP tx_fdm[M_FAC]; + float channel[CHANNEL_BUF_SIZE]; + int channel_count; + COMP rx_fdm[M_FAC + M_FAC / P]; + float foff_coarse; + int nin, next_nin; + COMP rx_fdm_fcorr[M_FAC + M_FAC / P]; + COMP rx_fdm_filter[M_FAC + M_FAC / P]; + COMP rx_filt[NC + 1][P + 1]; + float rx_timing; + float env[NT * P]; + COMP rx_symbols[FDMDV_NC + 1]; + int rx_bits[FDMDV_BITS_PER_FRAME]; + float foff_fine; + int sync_bit, reliable_sync_bit; + + int tx_bits_log[FDMDV_BITS_PER_FRAME * FRAMES]; + COMP tx_symbols_log[(FDMDV_NC + 1) * FRAMES]; + COMP tx_fdm_log[M_FAC * FRAMES]; + COMP pilot_baseband1_log[NPILOTBASEBAND * FRAMES]; + COMP pilot_baseband2_log[NPILOTBASEBAND * FRAMES]; + COMP pilot_lpf1_log[NPILOTLPF * FRAMES]; + COMP pilot_lpf2_log[NPILOTLPF * FRAMES]; + COMP S1_log[MPILOTFFT * FRAMES]; + COMP S2_log[MPILOTFFT * FRAMES]; + float foff_coarse_log[FRAMES]; + float foff_log[FRAMES]; + COMP rx_fdm_filter_log[(M_FAC + M_FAC / P) * FRAMES]; + int rx_fdm_filter_log_index; + COMP rx_filt_log[NC + 1][(P + 1) * FRAMES]; + int rx_filt_log_col_index; + float env_log[NT * P * FRAMES]; + float rx_timing_log[FRAMES]; + COMP rx_symbols_log[FDMDV_NC + 1][FRAMES]; + COMP phase_difference_log[FDMDV_NC + 1][FRAMES]; + float sig_est_log[FDMDV_NC + 1][FRAMES]; + float noise_est_log[FDMDV_NC + 1][FRAMES]; + int rx_bits_log[FDMDV_BITS_PER_FRAME * FRAMES]; + float foff_fine_log[FRAMES]; + int sync_bit_log[FRAMES]; + int sync_log[FRAMES]; + int nin_log[FRAMES]; + + FILE *fout; + int f, c, i, j; + + fdmdv = fdmdv_create(FDMDV_NC); + next_nin = M_FAC; + channel_count = 0; + + rx_fdm_filter_log_index = 0; + rx_filt_log_col_index = 0; + + printf("sizeof FDMDV states: %zd bytes\n", sizeof(struct FDMDV)); + + for (f = 0; f < FRAMES; f++) { + /* --------------------------------------------------------*\ + Modulator + \*---------------------------------------------------------*/ + + fdmdv_get_test_bits(fdmdv, tx_bits); + bits_to_dqpsk_symbols(tx_symbols, FDMDV_NC, fdmdv->prev_tx_symbols, tx_bits, + &fdmdv->tx_pilot_bit, 0); + memcpy(fdmdv->prev_tx_symbols, tx_symbols, sizeof(COMP) * (FDMDV_NC + 1)); + tx_filter_and_upconvert(tx_fdm, FDMDV_NC, tx_symbols, + fdmdv->tx_filter_memory, fdmdv->phase_tx, + fdmdv->freq, &fdmdv->fbb_phase_tx, fdmdv->fbb_rect); + + /* --------------------------------------------------------*\ + Channel + \*---------------------------------------------------------*/ + + nin = next_nin; + + // nin = M_FAC; // when debugging good idea to uncomment this to "open + // loop" + + /* add M_FAC tx samples to end of buffer */ + + assert((channel_count + M_FAC) < CHANNEL_BUF_SIZE); + for (i = 0; i < M_FAC; i++) channel[channel_count + i] = tx_fdm[i].real; + channel_count += M_FAC; + + /* take nin samples from start of buffer */ + + for (i = 0; i < nin; i++) { + rx_fdm[i].real = channel[i]; + rx_fdm[i].imag = 0; + } + + /* shift buffer back */ + + for (i = 0, j = nin; j < channel_count; i++, j++) channel[i] = channel[j]; + channel_count -= nin; + + /* --------------------------------------------------------*\ + Demodulator + \*---------------------------------------------------------*/ + + /* shift down to complex baseband */ + + fdmdv_freq_shift(rx_fdm, rx_fdm, -FDMDV_FCENTRE, &fdmdv->fbb_phase_rx, nin); + + /* freq offset estimation and correction */ + + // fdmdv->sync = 0; // when debugging good idea to uncomment this to "open + // loop" + + foff_coarse = rx_est_freq_offset(fdmdv, rx_fdm, nin, !fdmdv->sync); + + if (fdmdv->sync == 0) fdmdv->foff = foff_coarse; + fdmdv_freq_shift(rx_fdm_fcorr, rx_fdm, -fdmdv->foff, + &fdmdv->foff_phase_rect, nin); + + /* baseband processing */ + + rxdec_filter(rx_fdm_filter, rx_fdm_fcorr, fdmdv->rxdec_lpf_mem, nin); + down_convert_and_rx_filter(rx_filt, fdmdv->Nc, rx_fdm_filter, + fdmdv->rx_fdm_mem, fdmdv->phase_rx, fdmdv->freq, + fdmdv->freq_pol, nin, M_FAC / Q); + rx_timing = rx_est_timing(rx_symbols, FDMDV_NC, rx_filt, + fdmdv->rx_filter_mem_timing, env, nin, M_FAC); + foff_fine = + qpsk_to_bits(rx_bits, &sync_bit, FDMDV_NC, fdmdv->phase_difference, + fdmdv->prev_rx_symbols, rx_symbols, 0); + + // for(i=0; i<FDMDV_NC;i++) + // printf("rx_symbols: %f %f prev_rx_symbols: %f %f phase_difference: %f + // %f\n", rx_symbols[i].real, rx_symbols[i].imag, + // fdmdv->prev_rx_symbols[i].real, fdmdv->prev_rx_symbols[i].imag, + // fdmdv->phase_difference[i].real, + // fdmdv->phase_difference[i].imag); + // if (f==1) + // exit(0); + + snr_update(fdmdv->sig_est, fdmdv->noise_est, FDMDV_NC, + fdmdv->phase_difference); + memcpy(fdmdv->prev_rx_symbols, rx_symbols, sizeof(COMP) * (FDMDV_NC + 1)); + + next_nin = M_FAC; + + if (rx_timing > 2 * M_FAC / P) next_nin += M_FAC / P; + + if (rx_timing < 0) next_nin -= M_FAC / P; + + fdmdv->sync = freq_state(&reliable_sync_bit, sync_bit, &fdmdv->fest_state, + &fdmdv->timer, fdmdv->sync_mem); + fdmdv->foff -= TRACK_COEFF * foff_fine; + + /* --------------------------------------------------------*\ + Log each vector + \*---------------------------------------------------------*/ + + memcpy(&tx_bits_log[FDMDV_BITS_PER_FRAME * f], tx_bits, + sizeof(int) * FDMDV_BITS_PER_FRAME); + memcpy(&tx_symbols_log[(FDMDV_NC + 1) * f], tx_symbols, + sizeof(COMP) * (FDMDV_NC + 1)); + memcpy(&tx_fdm_log[M_FAC * f], tx_fdm, sizeof(COMP) * M_FAC); + + memcpy(&pilot_baseband1_log[f * NPILOTBASEBAND], fdmdv->pilot_baseband1, + sizeof(COMP) * NPILOTBASEBAND); + memcpy(&pilot_baseband2_log[f * NPILOTBASEBAND], fdmdv->pilot_baseband2, + sizeof(COMP) * NPILOTBASEBAND); + memcpy(&pilot_lpf1_log[f * NPILOTLPF], fdmdv->pilot_lpf1, + sizeof(COMP) * NPILOTLPF); + memcpy(&pilot_lpf2_log[f * NPILOTLPF], fdmdv->pilot_lpf2, + sizeof(COMP) * NPILOTLPF); + memcpy(&S1_log[f * MPILOTFFT], fdmdv->S1, sizeof(COMP) * MPILOTFFT); + memcpy(&S2_log[f * MPILOTFFT], fdmdv->S2, sizeof(COMP) * MPILOTFFT); + foff_coarse_log[f] = foff_coarse; + foff_log[f] = fdmdv->foff; + + /* rx filtering */ + + for (i = 0; i < nin; i++) + rx_fdm_filter_log[rx_fdm_filter_log_index + i] = rx_fdm_filter[i]; + rx_fdm_filter_log_index += nin; + + for (c = 0; c < NC + 1; c++) { + for (i = 0; i < (P * nin) / M_FAC; i++) + rx_filt_log[c][rx_filt_log_col_index + i] = rx_filt[c][i]; + } + rx_filt_log_col_index += (P * nin) / M_FAC; + + /* timing estimation */ + + memcpy(&env_log[NT * P * f], env, sizeof(float) * NT * P); + rx_timing_log[f] = rx_timing; + nin_log[f] = nin; + + for (c = 0; c < FDMDV_NC + 1; c++) { + rx_symbols_log[c][f] = rx_symbols[c]; + phase_difference_log[c][f] = fdmdv->phase_difference[c]; + } + + /* qpsk_to_bits() */ + + memcpy(&rx_bits_log[FDMDV_BITS_PER_FRAME * f], rx_bits, + sizeof(int) * FDMDV_BITS_PER_FRAME); + for (c = 0; c < FDMDV_NC + 1; c++) { + sig_est_log[c][f] = fdmdv->sig_est[c]; + noise_est_log[c][f] = fdmdv->noise_est[c]; + } + foff_fine_log[f] = foff_fine; + sync_bit_log[f] = sync_bit; + + sync_log[f] = fdmdv->sync; + } + + /*---------------------------------------------------------*\ + Dump logs to Octave file for evaluation + by tfdmdv.m Octave script + \*---------------------------------------------------------*/ + + fout = fopen("tfdmdv_out.txt", "wt"); + assert(fout != NULL); + fprintf(fout, "# Created by tfdmdv.c\n"); + octave_save_int(fout, "tx_bits_log_c", tx_bits_log, 1, + FDMDV_BITS_PER_FRAME * FRAMES); + octave_save_complex(fout, "tx_symbols_log_c", tx_symbols_log, 1, + (FDMDV_NC + 1) * FRAMES, (FDMDV_NC + 1) * FRAMES); + octave_save_complex(fout, "tx_fdm_log_c", (COMP *)tx_fdm_log, 1, + M_FAC * FRAMES, M_FAC * FRAMES); + octave_save_complex(fout, "pilot_lut_c", (COMP *)fdmdv->pilot_lut, 1, + NPILOT_LUT, NPILOT_LUT); + octave_save_complex(fout, "pilot_baseband1_log_c", pilot_baseband1_log, 1, + NPILOTBASEBAND * FRAMES, NPILOTBASEBAND * FRAMES); + octave_save_complex(fout, "pilot_baseband2_log_c", pilot_baseband2_log, 1, + NPILOTBASEBAND * FRAMES, NPILOTBASEBAND * FRAMES); + octave_save_float(fout, "pilot_coeff_c", pilot_coeff, 1, NPILOTCOEFF, + NPILOTCOEFF); + octave_save_complex(fout, "pilot_lpf1_log_c", pilot_lpf1_log, 1, + NPILOTLPF * FRAMES, NPILOTLPF * FRAMES); + octave_save_complex(fout, "pilot_lpf2_log_c", pilot_lpf2_log, 1, + NPILOTLPF * FRAMES, NPILOTLPF * FRAMES); + octave_save_complex(fout, "S1_log_c", S1_log, 1, MPILOTFFT * FRAMES, + MPILOTFFT * FRAMES); + octave_save_complex(fout, "S2_log_c", S2_log, 1, MPILOTFFT * FRAMES, + MPILOTFFT * FRAMES); + octave_save_float(fout, "foff_log_c", foff_log, 1, FRAMES, FRAMES); + octave_save_float(fout, "foff_coarse_log_c", foff_coarse_log, 1, FRAMES, + FRAMES); + octave_save_complex(fout, "rx_fdm_filter_log_c", (COMP *)rx_fdm_filter_log, 1, + rx_fdm_filter_log_index, rx_fdm_filter_log_index); + octave_save_complex(fout, "rx_filt_log_c", (COMP *)rx_filt_log, + (FDMDV_NC + 1), rx_filt_log_col_index, (P + 1) * FRAMES); + octave_save_float(fout, "env_log_c", env_log, 1, NT * P * FRAMES, + NT * P * FRAMES); + octave_save_float(fout, "rx_timing_log_c", rx_timing_log, 1, FRAMES, FRAMES); + octave_save_complex(fout, "rx_symbols_log_c", (COMP *)rx_symbols_log, + (FDMDV_NC + 1), FRAMES, FRAMES); + octave_save_complex(fout, "phase_difference_log_c", + (COMP *)phase_difference_log, (FDMDV_NC + 1), FRAMES, + FRAMES); + octave_save_float(fout, "sig_est_log_c", (float *)sig_est_log, (FDMDV_NC + 1), + FRAMES, FRAMES); + octave_save_float(fout, "noise_est_log_c", (float *)noise_est_log, + (FDMDV_NC + 1), FRAMES, FRAMES); + octave_save_int(fout, "rx_bits_log_c", rx_bits_log, 1, + FDMDV_BITS_PER_FRAME * FRAMES); + octave_save_float(fout, "foff_fine_log_c", foff_fine_log, 1, FRAMES, FRAMES); + octave_save_int(fout, "sync_bit_log_c", sync_bit_log, 1, FRAMES); + octave_save_int(fout, "sync_log_c", sync_log, 1, FRAMES); + octave_save_int(fout, "nin_log_c", nin_log, 1, FRAMES); + fclose(fout); + + fdmdv_destroy(fdmdv); + + return 0; +} diff --git a/unittest/tfifo.c b/unittest/tfifo.c new file mode 100644 index 0000000..bedda4f --- /dev/null +++ b/unittest/tfifo.c @@ -0,0 +1,103 @@ +/* + tfifo.c + David Rowe + Nov 19 2012 + + Tests FIFOs, in particular thread safety. +*/ + +#include <assert.h> +#include <pthread.h> +#include <stdio.h> + +#include "codec2_fifo.h" + +#define FIFO_SZ 1024 +#define WRITE_SZ 10 +#define READ_SZ 8 +#define N_MAX 100 +#define LOOPS 1000000 + +int run_thread = 1; +struct FIFO *f; + +void writer(void); +void *writer_thread(void *data); +pthread_mutex_t mutex = PTHREAD_MUTEX_INITIALIZER; + +#define USE_THREADS +//#define USE_MUTEX + +int main() { + pthread_t awriter_thread; + int i, j; + short read_buf[READ_SZ]; + int n_out = 0; + int success; + + f = codec2_fifo_create(FIFO_SZ); +#ifdef USE_THREADS + pthread_create(&awriter_thread, NULL, writer_thread, NULL); +#endif + + for (i = 0; i < LOOPS;) { +#ifndef USE_THREADS + writer(); +#endif + +#ifdef USE_MUTEX + pthread_mutex_lock(&mutex); +#endif + success = (codec2_fifo_read(f, read_buf, READ_SZ) == 0); +#ifdef USE_MUTEX + pthread_mutex_unlock(&mutex); +#endif + + if (success) { + for (j = 0; j < READ_SZ; j++) { + if (read_buf[j] != n_out) { + printf("error: %d %d\n", read_buf[j], n_out); + return (1); + } + n_out++; + if (n_out == N_MAX) n_out = 0; + } + i++; + } + } + +#ifdef USE_THREADS + run_thread = 0; + pthread_join(awriter_thread, NULL); +#endif + + printf("%d loops tested OK\n", LOOPS); + return 0; +} + +int n_in = 0; + +void writer(void) { + short write_buf[WRITE_SZ]; + int i; + + if ((FIFO_SZ - codec2_fifo_used(f)) > WRITE_SZ) { + for (i = 0; i < WRITE_SZ; i++) { + write_buf[i] = n_in++; + if (n_in == N_MAX) n_in = 0; + } +#ifdef USE_MUTEX + pthread_mutex_lock(&mutex); +#endif + codec2_fifo_write(f, write_buf, WRITE_SZ); + pthread_mutex_unlock(&mutex); + } +} + +void *writer_thread(void *data) { + while (run_thread) { + writer(); + } + + return NULL; +} diff --git a/unittest/tfmfsk.c b/unittest/tfmfsk.c new file mode 100644 index 0000000..d6818e1 --- /dev/null +++ b/unittest/tfmfsk.c @@ -0,0 +1,200 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfmfsk.c + AUTHOR......: Brady O'Brien + DATE CREATED: 8 February 2016 + + C test driver for fmfsk_mod and fmfsk_demod in fmfsk.c. Reads a file with +input bits/rf and spits out modulated/demoduladed samples and a dump of internal + state. To run unit test, see octave/tfmfsk.m + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2016 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/>. +*/ + +//#define MODEMPROBE_ENABLE + +#include <stdio.h> + +#include "modem_probe.h" + +/* Note: This is a dirty hack to force fsk.c to compile with modem probing + * enabled */ +#include "fmfsk.c" + +#define ST_BITS 10000 +#define ST_FS 48000 +#define ST_RS 2400 +#define ST_EBNO 8 + +#define TEST_SELF_FULL 1 /* No-arg self test */ +#define TEST_MOD 2 /* Test modulator using in and out file */ +#define TEST_DEMOD 3 /* Test demodulator using in and out file */ + +int main(int argc, char *argv[]) { + struct FMFSK *fmfsk; + int Fs, Rs; + FILE *fin, *fout; + + uint8_t *bitbuf = NULL; + float *modbuf = NULL; + uint8_t *bitbufp; + float *modbufp; + + size_t bitbufsize = 0; + size_t modbufsize = 0; + + int test_type; + + int i; + + fin = NULL; + fout = NULL; + + /* Set up full self-test */ + if (argc == 1) { + test_type = TEST_SELF_FULL; + modem_probe_init("fmfsk", "fmfsk_tfmfsk_log.txt"); + Fs = ST_FS; + Rs = ST_RS; + } else if (argc < 7) { + /* Not running any test */ + printf( + "Usage: %s [(M|D) SampleRate BitRate InputFile OutputFile " + "OctaveLogFile]\n", + argv[0]); + exit(1); + } else { + /* Running stim-drivin test */ + /* Mod test */ + if (strcmp(argv[1], "M") == 0 || strcmp(argv[1], "m") == 0) { + test_type = TEST_MOD; + /* Demod test */ + } else if (strcmp(argv[1], "D") == 0 || strcmp(argv[1], "d") == 0) { + test_type = TEST_DEMOD; + } else { + printf("Must specify mod or demod test with M or D\n"); + exit(1); + } + /* Extract parameters */ + Fs = atoi(argv[2]); + Rs = atoi(argv[3]); + + /* Open files */ + fin = fopen(argv[4], "r"); + fout = fopen(argv[5], "w"); + + if (fin == NULL || fout == NULL) { + printf("Couldn't open test vector files\n"); + exit(1); + } + /* Init modem probing */ + modem_probe_init("fmfsk", argv[6]); + } + + srand(1); + + /* set up FSK */ + fmfsk = fmfsk_create(Fs, Rs); + /* Modulate! */ + if (test_type == TEST_MOD || test_type == TEST_SELF_FULL) { + /* Generate random bits for self test */ + if (test_type == TEST_SELF_FULL) { + bitbufsize = ST_BITS; + bitbuf = (uint8_t *)malloc(sizeof(uint8_t) * ST_BITS); + for (i = 0; i < ST_BITS; i++) { + /* Generate a randomish bit */ + bitbuf[i] = (uint8_t)(rand() & 0x01); + } + } else { /* Load bits from a file */ + /* Figure out how many bits are in the input file */ + fseek(fin, 0L, SEEK_END); + bitbufsize = ftell(fin); + fseek(fin, 0L, SEEK_SET); + bitbuf = malloc(sizeof(uint8_t) * bitbufsize); + i = 0; + /* Read in some bits */ + bitbufp = bitbuf; + while (fread(bitbufp, sizeof(uint8_t), fmfsk->nbit, fin) == fmfsk->nbit) { + i++; + bitbufp += fmfsk->nbit; + /* Make sure we don't break the buffer */ + if (i * fmfsk->nbit > bitbufsize) { + bitbuf = + realloc(bitbuf, sizeof(uint8_t) * (bitbufsize + fmfsk->nbit)); + bitbufsize += fmfsk->nbit; + } + } + } + /* Allocate modulation buffer */ + modbuf = (float *)malloc(sizeof(float) * (bitbufsize / fmfsk->nbit) * + fmfsk->N * 4); + modbufsize = (bitbufsize / fmfsk->nbit) * fmfsk->N; + /* Do the modulation */ + modbufp = modbuf; + bitbufp = bitbuf; + while (bitbufp < bitbuf + bitbufsize) { + fmfsk_mod(fmfsk, modbufp, bitbufp); + modbufp += fmfsk->N; + bitbufp += fmfsk->nbit; + } + /* For a mod-only test, write out the result */ + if (test_type == TEST_MOD) { + fwrite(modbuf, sizeof(float), modbufsize, fout); + free(modbuf); + } + /* Free bit buffer */ + free(bitbuf); + } + + /* Add channel imp here */ + + /* Now test the demod */ + if (test_type == TEST_DEMOD || test_type == TEST_SELF_FULL) { + free(modbuf); + modbuf = malloc(sizeof(float) * (fmfsk->N + fmfsk->Ts * 2)); + bitbuf = malloc(sizeof(uint8_t) * fmfsk->nbit); + /* Demod-only test */ + if (test_type == TEST_DEMOD) { + // fprintf(stderr,"%d\n",(fmfsk->N+fmfsk->Ts*2)); + while (fread(modbuf, sizeof(float), fmfsk_nin(fmfsk), fin) == + fmfsk_nin(fmfsk)) { + fmfsk_demod(fmfsk, bitbuf, modbuf); + fwrite(bitbuf, sizeof(uint8_t), fmfsk->nbit, fout); + } + } + /* Demod after channel imp. and mod */ + else { + bitbufp = bitbuf; + modbufp = modbuf; + while (modbufp < modbuf + modbufsize) { + fmfsk_demod(fmfsk, bitbuf, modbuf); + modbufp += fmfsk_nin(fmfsk); + } + } + free(bitbuf); + } + + modem_probe_close(); + if (test_type == TEST_DEMOD || test_type == TEST_MOD) { + fclose(fin); + fclose(fout); + } + fmfsk_destroy(fmfsk); + exit(0); +} diff --git a/unittest/tfreedv_2400A_rawdata.c b/unittest/tfreedv_2400A_rawdata.c new file mode 100644 index 0000000..b81aa47 --- /dev/null +++ b/unittest/tfreedv_2400A_rawdata.c @@ -0,0 +1,113 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfreedv_2400A_rawdata.c + AUTHOR......: Jeroen Vreeken + DATE CREATED: 24 May 2020 + + FreeDV 2400A rawdata test. + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2020 Jeroen Vreeken <[email protected]> + + This library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + This library 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 + Lesser 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 <stdio.h> + +#include "assert.h" +#include "freedv_api.h" + +int main(int argc, char **argv) { + struct freedv *f; + int i; + + printf("freedv_api tests for mode 2400A\n"); + + printf("freedv_open(FREEDV_MODE_2400A) "); + f = freedv_open(FREEDV_MODE_2400A); + assert(f != NULL); + printf("Passed\n"); + + printf("freedv_get_mode() "); + int mode = freedv_get_mode(f); + assert(mode == FREEDV_MODE_2400A); + printf("Passed\n"); + + printf("freedv_get_n_max_modem_samples() "); + int max_samples = freedv_get_n_max_modem_samples(f); + assert(max_samples == 2040); + printf("%d Passed\n", max_samples); + + printf("freedv_get_n_nom_modem_samples() "); + int nom_samples = freedv_get_n_nom_modem_samples(f); + assert(nom_samples == 2000); + printf("%d Passed\n", nom_samples); + + printf("freedv_get_n_speech_samples() "); + int speech_samples = freedv_get_n_speech_samples(f); + assert(speech_samples == 320); + printf("%d Passed\n", speech_samples); + + printf("freedv_get_n_bits_per_codec_frame() "); + int codec_bits = freedv_get_bits_per_codec_frame(f); + assert(codec_bits == 52); + printf("%d Passed\n", codec_bits); + + printf("freedv_get_n_bits_per_modem_frame() "); + int frame_bits = freedv_get_bits_per_modem_frame(f); + assert(frame_bits == 52); + printf("%d Passed\n", frame_bits); + + printf("freedv_rawdatatx()/freedv_rawdatarx() "); + int frames = 0; + { + short mod[nom_samples * 10]; + /* Note: A codec frame is only 6.5 bytes! + so the seventh byte will be half empty! + */ + unsigned char payload[7] = {0x11, 0x22, 0x33, 0x44, 0x55, 0x66, 0x70}; + for (i = 0; i < 10; i++) { + freedv_rawdatatx(f, mod + i * nom_samples, payload); + } + int nin = 0; + for (i = 0; i < nom_samples * 9; i += nin) { + nin = freedv_nin(f); + unsigned char payload_rx[7] = {0}; + int r = freedv_rawdatarx(f, payload_rx, mod + i); + if (r) { + int b; + for (b = 0; b < 7; b++) { + if (payload[b] != payload_rx[b]) { + printf("Received codec bits 0x%02x do not match expected 0x%02x\n", + payload_rx[b], payload[b]); + } + } + frames++; + } + } + } + if (!frames) { + printf("Did not decode any frames successfully\n"); + goto fail; + } + + printf("Tests passed\n"); + return 0; +fail: + printf("Test failed\n"); + return 1; +} diff --git a/unittest/tfreedv_2400B_rawdata.c b/unittest/tfreedv_2400B_rawdata.c new file mode 100644 index 0000000..f515bb5 --- /dev/null +++ b/unittest/tfreedv_2400B_rawdata.c @@ -0,0 +1,113 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfreedv_2400B_rawdata.c + AUTHOR......: Jeroen Vreeken + DATE CREATED: 24 May 2020 + + FreeDV 2400B rawdata test. + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2020 Jeroen Vreeken <[email protected]> + + This library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + This library 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 + Lesser 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 <stdio.h> + +#include "assert.h" +#include "freedv_api.h" + +int main(int argc, char **argv) { + struct freedv *f; + int i; + + printf("freedv_api tests for mode 2400B\n"); + + printf("freedv_open(FREEDV_MODE_2400B) "); + f = freedv_open(FREEDV_MODE_2400B); + assert(f != NULL); + printf("Passed\n"); + + printf("freedv_get_mode() "); + int mode = freedv_get_mode(f); + assert(mode == FREEDV_MODE_2400B); + printf("Passed\n"); + + printf("freedv_get_n_max_modem_samples() "); + int max_samples = freedv_get_n_max_modem_samples(f); + assert(max_samples == 1930); + printf("%d Passed\n", max_samples); + + printf("freedv_get_n_nom_modem_samples() "); + int nom_samples = freedv_get_n_nom_modem_samples(f); + assert(nom_samples == 1920); + printf("%d Passed\n", nom_samples); + + printf("freedv_get_n_speech_samples() "); + int speech_samples = freedv_get_n_speech_samples(f); + assert(speech_samples == 320); + printf("%d Passed\n", speech_samples); + + printf("freedv_get_n_bits_per_codec_frame() "); + int codec_bits = freedv_get_bits_per_codec_frame(f); + assert(codec_bits == 52); + printf("%d Passed\n", codec_bits); + + printf("freedv_get_n_bits_per_modem_frame() "); + int frame_bits = freedv_get_bits_per_modem_frame(f); + assert(frame_bits == 52); + printf("%d Passed\n", frame_bits); + + printf("freedv_rawdatatx()/freedv_rawdatarx() "); + int frames = 0; + { + short mod[nom_samples * 10]; + /* Note: A codec frame is only 6.5 bytes! + so the seventh byte will be half empty! + */ + unsigned char payload[7] = {0x11, 0x22, 0x33, 0x44, 0x55, 0x66, 0x70}; + for (i = 0; i < 10; i++) { + freedv_rawdatatx(f, mod + i * nom_samples, payload); + } + int nin = 0; + for (i = 0; i < nom_samples * 9; i += nin) { + nin = freedv_nin(f); + unsigned char payload_rx[7] = {0}; + int r = freedv_rawdatarx(f, payload_rx, mod + i); + if (r) { + int b; + for (b = 0; b < 7; b++) { + if (payload[b] != payload_rx[b]) { + printf("Received codec bits 0x%02x do not match expected 0x%02x\n", + payload_rx[b], payload[b]); + } + } + frames++; + } + } + } + if (!frames) { + printf("Did not decode any frames successfully\n"); + goto fail; + } + + printf("Tests passed\n"); + return 0; +fail: + printf("Test failed\n"); + return 1; +} diff --git a/unittest/tfreedv_800XA_rawdata.c b/unittest/tfreedv_800XA_rawdata.c new file mode 100644 index 0000000..2b34c41 --- /dev/null +++ b/unittest/tfreedv_800XA_rawdata.c @@ -0,0 +1,147 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfreedv_800XA_rawdata.c + AUTHOR......: Jeroen Vreeken + DATE CREATED: 24 May 2020 + + FreeDV 800XA rawdata test. + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2020 Jeroen Vreeken <[email protected]> + + This library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + This library 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 + Lesser 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 <stdio.h> + +#include "assert.h" +#include "freedv_api.h" + +int main(int argc, char **argv) { + struct freedv *f; + int i; + + printf("freedv_api tests for mode 800XA\n"); + + printf("freedv_open(FREEDV_MODE_800XA) "); + f = freedv_open(FREEDV_MODE_800XA); + assert(f != NULL); + printf("Passed\n"); + + printf("freedv_get_mode() "); + int mode = freedv_get_mode(f); + assert(mode == FREEDV_MODE_800XA); + printf("Passed\n"); + + printf("freedv_get_n_max_modem_samples() "); + int max_samples = freedv_get_n_max_modem_samples(f); + assert(max_samples == 660); + printf("%d Passed\n", max_samples); + + printf("freedv_get_n_nom_modem_samples() "); + int nom_samples = freedv_get_n_nom_modem_samples(f); + assert(nom_samples == 640); + printf("%d Passed\n", nom_samples); + + printf("freedv_get_n_speech_samples() "); + int speech_samples = freedv_get_n_speech_samples(f); + assert(speech_samples == 640); + printf("%d Passed\n", speech_samples); + + printf("freedv_get_n_bits_per_codec_frame() "); + int codec_bits = freedv_get_bits_per_codec_frame(f); + assert(codec_bits == 28); + printf("%d Passed\n", codec_bits); + + printf("freedv_get_n_bits_per_modem_frame() "); + int frame_bits = freedv_get_bits_per_modem_frame(f); + assert(frame_bits == 56); + printf("%d Passed\n", frame_bits); + + /* Note: A codec frame is only 3.5 bytes! + so the fourth and eight bytes will be half empty! + */ + unsigned char payload[8] = {0x12, 0x34, 0x56, 0x70, 0x89, 0xab, 0xcd, 0xe0}; + unsigned char payload_tx[7] = {0x12, 0x34, 0x56, 0x78, 0x9a, 0xbc, 0xde}; + + printf("freedv_codec_frames_from_rawdata() "); + unsigned char codec_frames[8] = {0}; + freedv_codec_frames_from_rawdata(f, codec_frames, payload_tx); + int fails = 0; + for (i = 0; i < 8; i++) { + if (codec_frames[i] != payload[i]) { + printf("byte %d: 0x%02x does not match expected 0x%02x\n", i, + codec_frames[i], payload[i]); + fails++; + } + } + if (fails) goto fail; + printf("Passed\n"); + + printf("freedv_rawdata_from_codec_frames() "); + unsigned char rawdata[7] = {0}; + freedv_rawdata_from_codec_frames(f, rawdata, payload); + fails = 0; + for (i = 0; i < 7; i++) { + if (rawdata[i] != payload_tx[i]) { + printf("byte %d: 0x%02x does not match expected 0x%02x\n", i, rawdata[i], + payload_tx[i]); + fails++; + } + } + if (fails) goto fail; + printf("Passed\n"); + + printf("freedv_rawdatatx()/freedv_rawdatarx() "); + int frames = 0; + fails = 0; + { + short mod[nom_samples * 10]; + for (i = 0; i < 10; i++) { + freedv_rawdatatx(f, mod + i * nom_samples, payload_tx); + } + int nin = 0; + for (i = 0; i < nom_samples * 9; i += nin) { + nin = freedv_nin(f); + unsigned char payload_rx[8] = {0}; + int r = freedv_rawdatarx(f, payload_rx, mod + i); + if (r == 7) { + int b; + for (b = 0; b < 7; b++) { + if (payload_tx[b] != payload_rx[b]) { + printf("Received codec bits 0x%02x do not match expected 0x%02x\n", + payload_rx[b], payload_tx[b]); + fails++; + } + } + frames++; + } + } + } + if (!frames) { + printf("Did not decode any frames successfully\n"); + goto fail; + } + if (fails) goto fail; + printf("Passed\n"); + + printf("Tests passed\n"); + return 0; +fail: + printf("Test failed\n"); + return 1; +} diff --git a/unittest/tfreedv_data_channel.c b/unittest/tfreedv_data_channel.c new file mode 100644 index 0000000..e94cadc --- /dev/null +++ b/unittest/tfreedv_data_channel.c @@ -0,0 +1,282 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfreedv_data_channel + AUTHOR......: Jeroen Vreeken + DATE CREATED: May 3 2016 + + Tests for the data channel code. + Data channel frame behaviour is tested with test vectors. + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2016 Jeroen Vreeken + + 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 <stdio.h> +#include <string.h> + +#include "freedv_data_channel.h" + +unsigned char test_header[] = {0x11, 0x22, 0x33, 0x44, 0x55, 0x66}; +unsigned char bcast_header[] = {0xff, 0xff, 0xff, 0xff, 0xff, 0xff}; + +struct testvec { + char *testname; + + unsigned char *data; + size_t data_size; + + size_t frame_size; + + unsigned char *frame_data; + size_t frame_data_size; + + unsigned char *flags; +} testvec[] = { + { + "Regular packet, does not match header and no broadcast", + (unsigned char[]){0x01, 0x02, 0x03, 0x04, 0x05, 0x06, 0x07, 0x08, 0x09, + 0x0a, 0x0b, 0x0c, 0x0d, 0x0e, 0x0f, 0x10, 0x11, 0x12}, + 0x12, + 8, + (unsigned char[]){0x07, 0x08, 0x09, 0x0a, 0x0b, 0x0c, 0x01, + 0x02, 0x03, 0x04, 0x05, 0x06, 0x0d, 0x0e, + 0x0f, 0x10, 0x11, 0x12, 0x47, 0x6e}, + 0x14, + (unsigned char[]){0x00, 0x00, 0x04}, + }, + { + "Header", + NULL, + 0, + 8, + (unsigned char[]){0x11, 0x22, 0x33, 0x44, 0x55, 0x66, 0x5a, 0x60}, + 0x08, + (unsigned char[]){0x08}, + }, + { + "Broadcast packet", + (unsigned char[]){0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0x11, 0x22, 0x33, + 0x44, 0x55, 0x66, 0x05, 0x06, 0x07, 0x08, 0x09, 0x0a, + 0x0b, 0x0c, 0x0d, 0x0e, 0x0f, 0x10, 0x11}, + 0x19, + 8, + (unsigned char[]){0x05, 0x06, 0x07, 0x08, 0x09, 0x0a, 0x0b, 0x0c, 0x0d, + 0x0e, 0x0f, 0x10, 0x11, 0x3c, 0xbe}, + 0x0f, + (unsigned char[]){0xc0, 0x07}, + }, + { + "Broadcast packet, header does not match", + (unsigned char[]){0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xaa, 0x22, 0xbb, + 0xcc, 0xdd, 0xee, 0x05, 0x06, 0x07, 0x08, 0x09, 0x0a, + 0x0b, 0x0c, 0x0d, 0x0e, 0x0f, 0x10, 0x11}, + 0x19, + 8, + (unsigned char[]){0xaa, 0x22, 0xbb, 0xcc, 0xdd, 0xee, 0x05, + 0x06, 0x07, 0x08, 0x09, 0x0a, 0x0b, 0x0c, + 0x0d, 0x0e, 0x0f, 0x10, 0x11, 0x1a, 0x68}, + 0x15, + (unsigned char[]){0x40, 0x00, 0x05}, + }, + { + "Header 6 bytes", + NULL, + 0, + 6, + (unsigned char[]){0x11, 0x22, 0x33, 0x44, 0x55, 0x66}, + 0x06, + (unsigned char[]){0x2f}, + }, + { + "Broadcast packet (6 byte frames)", + (unsigned char[]){0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0x11, 0x22, 0x33, + 0x44, 0x55, 0x66, 0x05, 0x06, 0x07, 0x08, 0x09, 0x0a, + 0x0b, 0x0c, 0x0d, 0x0e, 0x0f, 0x10, 0x11}, + 0x19, + 6, + (unsigned char[]){0x05, 0x06, 0x07, 0x08, 0x09, 0x0a, 0x0b, 0x0c, 0x0d, + 0x0e, 0x0f, 0x10, 0x11, 0x3c, 0xbe}, + 0x0f, + (unsigned char[]){0xc0, 0x00, 0x03}, + }, + { + "Broadcast packet, header does not match (6 byte frames)", + (unsigned char[]){0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xaa, 0x22, 0xbb, + 0xcc, 0xdd, 0xee, 0x05, 0x06, 0x07, 0x08, 0x09, 0x0a, + 0x0b, 0x0c, 0x0d, 0x0e, 0x0f, 0x10, 0x11}, + 0x19, + 6, + (unsigned char[]){0xaa, 0x22, 0xbb, 0xcc, 0xdd, 0xee, 0x05, + 0x06, 0x07, 0x08, 0x09, 0x0a, 0x0b, 0x0c, + 0x0d, 0x0e, 0x0f, 0x10, 0x11, 0x1a, 0x68}, + 0x15, + (unsigned char[]){0x40, 0x00, 0x00, 0x03}, + }, +}; + +static int ret = 0; +static int vector = 0; +static size_t frame_data_pos = 0; +static int rx_done = 0; + +void *tx_cb_arg = (void *)0xaa55; +void *rx_cb_arg = (void *)0xbb44; + +void tfreedv_data_callback_tx(void *arg, unsigned char *packet, size_t *size) { + if (tx_cb_arg != arg) { + ret++; + printf("FAIL: %s called with wrong argument value\n", __func__); + } + printf("--------------------------------------\n"); + printf("TX callback called for test %zd bytes data for test %d:\n'%s'\n", + testvec[vector].data_size, vector, testvec[vector].testname); + + memcpy(packet, testvec[vector].data, testvec[vector].data_size); + *size = testvec[vector].data_size; + + return; +} + +void tfreedv_data_callback_rx(void *arg, unsigned char *packet, size_t size) { + if (rx_cb_arg != arg) { + ret++; + printf("FAIL: %s called with wrong argument value\n", __func__); + } + printf("RX callback called with %zd bytes\n", size); + + if (testvec[vector].data_size) { + size_t data_size = testvec[vector].data_size; + if (data_size != size) { + printf("FAIL: Received size does not match test vector: %zd != %zd\n", + size, data_size); + ret++; + } else { + size_t i; + for (i = 0; i < data_size; i++) { + if (packet[i] != testvec[vector].data[i]) { + printf("FAIL: byte %zd does not match 0x%02x != 0x%02x\n", i, + packet[i], testvec[vector].data[i]); + ret++; + } + } + } + } else { + if (size != 12) { + printf("FAIL: Received header is not 12 bytes: %zd\n", size); + ret++; + } else { + if (memcmp(packet, bcast_header, 6)) { + printf("FAIL: Header is not a broadcast!\n"); + ret++; + } + if (memcmp(packet + 6, test_header, 6)) { + printf("FAIL: Header does not match!\n"); + ret++; + } + } + } + + rx_done = 1; +} + +int main(int argc, char **argv) { + struct freedv_data_channel *fdc; + + fdc = freedv_data_channel_create(); + + freedv_data_set_header(fdc, test_header); + freedv_data_set_cb_tx(fdc, tfreedv_data_callback_tx, tx_cb_arg); + freedv_data_set_cb_rx(fdc, tfreedv_data_callback_rx, rx_cb_arg); + + while (vector < sizeof(testvec) / sizeof(struct testvec)) { + size_t frame_size = testvec[vector].frame_size; + unsigned char frame[frame_size]; + int from, bcast, crc, end; + int i; + size_t check_size; + unsigned char flags; + int nr_frames; + + freedv_data_channel_tx_frame(fdc, frame, frame_size, &from, &bcast, &crc, + &end); + + check_size = frame_size; + if (frame_data_pos + check_size > testvec[vector].frame_data_size) + check_size = testvec[vector].frame_data_size - frame_data_pos; + + flags = from * 0x80 + bcast * 0x40 + crc * 0x20 + end; + printf("0x%02x:", flags); + for (i = 0; i < check_size; i++) { + if (frame[i] != testvec[vector].frame_data[frame_data_pos + i]) { + printf(" [0x%02x!=0x%02x]", frame[i], + testvec[vector].frame_data[frame_data_pos + i]); + ret++; + } else { + printf(" 0x%02x", frame[i]); + } + } + printf("\n"); + + if (flags != testvec[vector].flags[frame_data_pos / frame_size]) { + printf("FAIL: Flags byte does not match 0x%02x != 0x%02x\n", flags, + testvec[vector].flags[frame_data_pos / frame_size]); + ret++; + } + + freedv_data_channel_rx_frame(fdc, frame, frame_size, from, bcast, crc, end); + + frame_data_pos += frame_size; + + nr_frames = freedv_data_get_n_tx_frames(fdc, frame_size); + + if (frame_data_pos >= testvec[vector].frame_data_size) { + if (nr_frames) { + printf("FAIL: nr_frames is not zero: %d\n", nr_frames); + ret++; + } + vector++; + frame_data_pos = 0; + if (!rx_done) { + printf("FAIL: RX callback not executed\n"); + ret++; + } + rx_done = 0; + } else { + int vec_frames = (testvec[vector].frame_data_size - frame_data_pos); + vec_frames /= frame_size; + vec_frames++; + if (nr_frames != vec_frames) { + printf("FAIL: nr_frames != vec_frames: %d != %d\n", nr_frames, + vec_frames); + ret++; + } + } + } + + freedv_data_channel_destroy(fdc); + + printf("--------------------------------------\n"); + printf("tfreedv_data_channel test result: "); + if (ret) { + printf("Failed %d\n", ret); + } else { + printf("Passed\n"); + } + + return ret; +} diff --git a/unittest/tfsk.c b/unittest/tfsk.c new file mode 100644 index 0000000..337a350 --- /dev/null +++ b/unittest/tfsk.c @@ -0,0 +1,234 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfsk.c + AUTHOR......: Brady O'Brien + DATE CREATED: 20 January 2016 + + C test driver for fsk_mod and fsk_demod in fsk.c. Reads a file with input + bits/rf and spits out modulated/demoduladed samples and a dump of internal + state. To run unit test, see octave/tfsk.m + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2016 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/>. +*/ + +//#define MODEMPROBE_ENABLE + +#include <stdio.h> + +#include "modem_probe.h" + +/* Note: This is a dirty hack to force fsk.c to compile with modem probing + * enabled */ +#include "fsk.c" + +#define ST_BITS 10000 +#define ST_FS 8000 +#define ST_RS 100 +#define ST_F1 1200 +#define ST_Fs 400 +#define ST_EBNO 8 +#define ST_M 2 + +#define TEST_SELF_FULL 1 /* No-arg self test */ +#define TEST_MOD 2 /* Test modulator using in and out file */ +#define TEST_MOD_H 3 /* Test modulator using in and out file */ +#define TEST_DEMOD 4 /* Test demodulator using in and out file */ +#define TEST_DEMOD_H 5 /* Test demodulator using in and out file */ + +int main(int argc, char *argv[]) { + struct FSK *fsk; + int Fs, Rs, f1, fs, M, lock_nin; + FILE *fin, *fout; + + uint8_t *bitbuf = NULL; + float *modbuf = NULL; + uint8_t *bitbufp; + float *modbufp; + + size_t bitbufsize = 0; + size_t modbufsize = 0; + + int test_type; + + int i; + + fin = NULL; + fout = NULL; + + /* Set up full self-test */ + if (argc == 1) { + test_type = TEST_SELF_FULL; + modem_probe_init("fsk2", "fsk2_tfsk_log.txt"); + Fs = ST_FS; + Rs = ST_RS; + f1 = ST_F1; + fs = ST_Fs; + M = ST_M; + lock_nin = 0; + } else if (argc < 10) { + /* Not running any test */ + printf( + "Usage: %s [(M|D|DX) Mode TXFreq1 TXFreqSpace SampleRate SymbolRate " + "lock_nin InputFile OutputFile OctaveLogFile]\n", + argv[0]); + exit(1); + } else { + /* Running stim-drivin test */ + /* Mod test */ + if (strcmp(argv[1], "MX") == 0) { + test_type = TEST_MOD_H; + } else if (strcmp(argv[1], "M") == 0 || strcmp(argv[1], "m") == 0) { + test_type = TEST_MOD; + /* Demod test */ + } else if (strcmp(argv[1], "DX") == 0) { + test_type = TEST_DEMOD_H; + } else if (strcmp(argv[1], "D") == 0 || strcmp(argv[1], "d") == 0) { + test_type = TEST_DEMOD; + } else { + printf("Must specify mod or demod test with M or D\n"); + exit(1); + } + /* Extract parameters */ + M = atoi(argv[2]); + f1 = atoi(argv[3]); + fs = atoi(argv[4]); + Fs = atoi(argv[5]); + Rs = atoi(argv[6]); + lock_nin = atoi(argv[7]); + + /* Open files */ + fin = fopen(argv[8], "r"); + fout = fopen(argv[9], "w"); + + if (fin == NULL || fout == NULL) { + printf("Couldn't open test vector files\n"); + exit(1); + } + /* Init modem probing */ + modem_probe_init("fsk", argv[10]); + } + + srand(1); + /* set up FSK */ + if (test_type == TEST_DEMOD_H || test_type == TEST_MOD_H) { + fsk = fsk_create_hbr(Fs, Rs, M, 10, FSK_DEFAULT_NSYM, f1, fs); + if (test_type == TEST_DEMOD_H) + test_type = TEST_DEMOD; + else + test_type = TEST_MOD; + } else { + fsk = fsk_create(Fs, Rs, M, f1, fs); + } + fsk_set_freq_est_limits(fsk, 300, 2800); + fsk->lock_nin = lock_nin; + + /* Modulate! */ + if (test_type == TEST_MOD || test_type == TEST_SELF_FULL) { + /* Generate random bits for self test */ + if (test_type == TEST_SELF_FULL) { + bitbufsize = ST_BITS; + bitbuf = (uint8_t *)malloc(sizeof(uint8_t) * ST_BITS); + for (i = 0; i < ST_BITS; i++) { + /* Generate a randomish bit */ + bitbuf[i] = (uint8_t)(rand() & 0x01); + } + } else { /* Load bits from a file */ + /* Figure out how many bits are in the input file */ + fseek(fin, 0L, SEEK_END); + bitbufsize = ftell(fin); + fseek(fin, 0L, SEEK_SET); + bitbuf = malloc(sizeof(uint8_t) * bitbufsize); + i = 0; + /* Read in some bits */ + bitbufp = bitbuf; + while (fread(bitbufp, sizeof(uint8_t), fsk->Nbits, fin) == fsk->Nbits) { + i++; + bitbufp += fsk->Nbits; + /* Make sure we don't break the buffer */ + if (i * fsk->Nbits > bitbufsize) { + bitbuf = realloc(bitbuf, sizeof(uint8_t) * (bitbufsize + fsk->Nbits)); + bitbufsize += fsk->Nbits; + } + } + } + /* Allocate modulation buffer */ + modbuf = + (float *)malloc(sizeof(float) * (bitbufsize / fsk->Nbits) * fsk->N * 4); + modbufsize = (bitbufsize / fsk->Nbits) * fsk->N; + /* Do the modulation */ + modbufp = modbuf; + bitbufp = bitbuf; + while (bitbufp < bitbuf + bitbufsize) { + fsk_mod(fsk, modbufp, bitbufp, fsk->Nbits); + modbufp += fsk->N; + bitbufp += fsk->Nbits; + } + /* For a mod-only test, write out the result */ + if (test_type == TEST_MOD) { + fwrite(modbuf, sizeof(float), modbufsize, fout); + free(modbuf); + } + /* Free bit buffer */ + free(bitbuf); + } + + /* Now test the demod */ + if (test_type == TEST_DEMOD || test_type == TEST_SELF_FULL) { + free(modbuf); + modbuf = malloc(sizeof(float) * (fsk->N + fsk->Ts * 2)); + bitbuf = malloc(sizeof(uint8_t) * fsk->Nbits); + /* Demod-only test */ + if (test_type == TEST_DEMOD) { + while (fread(modbuf, sizeof(float), fsk_nin(fsk), fin) == fsk_nin(fsk)) { + int n = fsk_nin(fsk); + COMP modbuf_comp[n]; + for (i = 0; i < n; i++) { + modbuf_comp[i].real = modbuf[i]; + modbuf_comp[i].imag = 0.0; + } + fsk_demod(fsk, bitbuf, modbuf_comp); + fwrite(bitbuf, sizeof(uint8_t), fsk->Nbits, fout); + } + } + /* Demod after channel imp. and mod */ + else { + bitbufp = bitbuf; + modbufp = modbuf; + while (modbufp < modbuf + modbufsize) { + int n = fsk_nin(fsk); + COMP modbuf_comp[n]; + for (i = 0; i < n; i++) { + modbuf_comp[i].real = modbuf[i]; + modbuf_comp[i].imag = 0.0; + } + fsk_demod(fsk, bitbuf, modbuf_comp); + modbufp += fsk_nin(fsk); + } + } + free(bitbuf); + } + + modem_probe_close(); + if (test_type == TEST_DEMOD || test_type == TEST_MOD) { + fclose(fin); + fclose(fout); + } + fsk_destroy(fsk); + exit(0); +} diff --git a/unittest/tfsk_llr.c b/unittest/tfsk_llr.c new file mode 100644 index 0000000..5929452 --- /dev/null +++ b/unittest/tfsk_llr.c @@ -0,0 +1,60 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tfsk_llr.c + AUTHOR......: David Rowe + DATE CREATED: July 2020 + + Simple test program for 4FSK LLR routines. + +\*---------------------------------------------------------------------------*/ + +#include <math.h> +#include <stdio.h> + +#include "mpdecode_core.h" + +#define M 4 +#define BPS 2 +#define NSYM 5 +#define V_EST 2 +#define SNR_EST 10 + +/* Generated test vectors with: + + octave:100> init_cml('~/cml/'); + octave:101> rx_filt=[1 0 0 0; 0 1 0 0; 0 0 1 0; 0 0 0 1] + octave:102> symL = DemodFSK(rx_filt,10,1); -Somap(symL) +*/ + +/* one col per symbol: + 0 1 2 3 4 */ +float rx_filt[] = { + 1.0, 0.0, 0.0, 0.0, 1.0, /* filter 0 */ + 0.0, 1.0, 0.0, 0.0, 0.0, /* filter 1 */ + 0.0, 0.0, 1.0, 0.0, 0.0, /* filter 2 */ + 0.0, 0.0, 0.0, 1.0, 0.0 /* filter 3 */ +}; + +float llrs_target[] = {7.3252, 7.3252, /* bit 0, 1 */ + 7.3252, -7.3252, /* 2, 3, ... */ + -7.3252, 7.3252, -7.3252, -7.3252, 7.3252, 7.3252}; + +int main(void) { + float llrs[BPS * NSYM] = {0}; + + fsk_rx_filt_to_llrs(llrs, rx_filt, V_EST, SNR_EST, M, NSYM); + + float error = 0.0; + for (int i = 0; i < NSYM * BPS; i++) { + fprintf(stderr, "% f\n", llrs[i]); + error += pow(llrs[i] - llrs_target[i], 2.0); + } + + if (error < 1E-3) { + fprintf(stderr, "PASS\n"); + return 0; + } else { + fprintf(stderr, "FAIL\n"); + return 1; + } +} diff --git a/unittest/thash.c b/unittest/thash.c new file mode 100644 index 0000000..d936728 --- /dev/null +++ b/unittest/thash.c @@ -0,0 +1,18 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: thash.c + AUTHOR......: David Rowe + DATE CREATED: July 2020 + + Simple test program for freeDV API get hash function + +\*---------------------------------------------------------------------------*/ + +#include <stdio.h> + +#include "freedv_api.h" + +int main(void) { + printf("%s\n", freedv_get_hash()); + return 0; +} diff --git a/unittest/tnewamp1.c b/unittest/tnewamp1.c new file mode 100644 index 0000000..297f12d --- /dev/null +++ b/unittest/tnewamp1.c @@ -0,0 +1,295 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tnewamp1.c + AUTHOR......: David Rowe + DATE CREATED: Jan 2017 + + Tests for the C version of the newamp1 amplitude modelling used for + 700c. This program outputs a file of Octave vectors that are loaded + and automatically tested against the Octave version of the modem by + the Octave script tnewamp1.m + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2017 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 "codec2_fft.h" +#include "defines.h" +#include "dump.h" +#include "newamp1.h" +#include "nlp.h" +#include "octave.h" +#include "quantise.h" +#include "sine.h" + +#define FRAMES 300 + +int main(int argc, char *argv[]) { + int Fs = 8000; + C2CONST c2const = c2const_create(Fs, N_S); + int n_samp = c2const.n_samp; + int m_pitch = c2const.m_pitch; + short buf[n_samp]; /* input/output buffer */ + float Sn[m_pitch]; /* float input speech samples */ + COMP Sw[FFT_ENC]; /* DFT of Sn[] */ + codec2_fft_cfg fft_fwd_cfg; /* fwd FFT states */ + float w[m_pitch]; /* time domain hamming window */ + float W[FFT_ENC]; /* DFT of w[] */ + MODEL model; + void *nlp_states; + codec2_fft_cfg phase_fft_fwd_cfg, phase_fft_inv_cfg; + float pitch, prev_f0; + int i, m, f, k; + + if (argc != 2) { + printf("usage: ./tnewamp1 RawFile\n"); + exit(1); + } + nlp_states = nlp_create(&c2const); + prev_f0 = 1.0 / P_MAX_S; + fft_fwd_cfg = codec2_fft_alloc(FFT_ENC, 0, NULL, NULL); + make_analysis_window(&c2const, fft_fwd_cfg, w, W); + + phase_fft_fwd_cfg = codec2_fft_alloc(NEWAMP1_PHASE_NFFT, 0, NULL, NULL); + phase_fft_inv_cfg = codec2_fft_alloc(NEWAMP1_PHASE_NFFT, 1, NULL, NULL); + + for (i = 0; i < m_pitch; i++) { + Sn[i] = 1.0; + } + + int K = 20; + float rate_K_sample_freqs_kHz[K]; + float model_octave[FRAMES][MAX_AMP + 2]; // model params in matrix format, + // useful for C <-> Octave + float rate_K_surface[FRAMES][K]; // rate K vecs for each frame, form a + // surface that makes pretty graphs + float rate_K_surface_no_mean[FRAMES][K]; // mean removed surface + float rate_K_surface_no_mean_[FRAMES][K]; // quantised mean removed surface + float mean[FRAMES]; + float mean_[FRAMES]; + float rate_K_surface_[FRAMES][K]; // quantised rate K vecs for each frame + float interpolated_surface_[FRAMES][K]; // dec/interpolated surface + // int voicing[FRAMES]; + int voicing_[FRAMES]; + float model_octave_[FRAMES][MAX_AMP + 2]; + COMP H[FRAMES][MAX_AMP]; + int indexes[FRAMES][NEWAMP1_N_INDEXES]; + float se = 0.0; + float eq[K]; + + for (k = 0; k < K; k++) eq[k] = 0.0; + + for (f = 0; f < FRAMES; f++) { + for (m = 0; m < MAX_AMP + 2; m++) { + model_octave[f][m] = 0.0; + model_octave_[f][m] = 0.0; + } + for (m = 0; m < MAX_AMP; m++) { + H[f][m].real = 0.0; + H[f][m].imag = 0.0; + } + for (k = 0; k < K; k++) interpolated_surface_[f][k] = 0.0; + voicing_[f] = 0; + } + + mel_sample_freqs_kHz(rate_K_sample_freqs_kHz, K, ftomel(200.0), + ftomel(3700.0)); + + // for(int k=0; k<K; k++) + // printf("k: %d sf: %f\n", k, rate_K_sample_freqs_kHz[k]); + + FILE *fin = fopen(argv[1], "rb"); + if (fin == NULL) { + fprintf(stderr, "Problem opening hts1.raw\n"); + exit(1); + } + + int M = 4; + + for (f = 0; f < FRAMES; f++) { + assert(fread(buf, sizeof(short), n_samp, fin) == n_samp); + + /* shift buffer of input samples, and insert new samples */ + + for (i = 0; i < m_pitch - n_samp; i++) { + Sn[i] = Sn[i + n_samp]; + } + for (i = 0; i < n_samp; i++) { + Sn[i + m_pitch - n_samp] = buf[i]; + } + + /* Estimate Sinusoidal Model Parameters ----------------------*/ + + nlp(nlp_states, Sn, n_samp, &pitch, Sw, W, &prev_f0); + model.Wo = TWO_PI / pitch; + + dft_speech(&c2const, fft_fwd_cfg, Sw, Sn, w); + two_stage_pitch_refinement(&c2const, &model, Sw); + estimate_amplitudes(&model, Sw, W, 1); + est_voicing_mbe(&c2const, &model, Sw, W); + // voicing[f] = model.voiced; + + /* newamp1 processing ----------------------------------------*/ + + newamp1_model_to_indexes(&c2const, &indexes[f][0], &model, + &rate_K_surface[f][0], rate_K_sample_freqs_kHz, K, + &mean[f], &rate_K_surface_no_mean[f][0], + &rate_K_surface_no_mean_[f][0], &se, eq, 0); + + newamp1_indexes_to_rate_K_vec( + &rate_K_surface_[f][0], &rate_K_surface_no_mean_[f][0], + rate_K_sample_freqs_kHz, K, &mean_[f], &indexes[f][0], NULL, 1); + +#ifdef VERBOSE + fprintf(stderr, "f: %d Wo: %4.3f L: %d v: %d\n", f, model.Wo, model.L, + model.voiced); + if ((f % M) == 0) { + for (i = 0; i < 5; i++) { + fprintf(stderr, " %5.3f", rate_K_surface_[f][i]); + } + fprintf(stderr, "\n"); + fprintf(stderr, " %d %d %d %d\n", indexes[f][0], indexes[f][1], + indexes[f][2], indexes[f][3]); + } +#endif + /* log vectors */ + + model_octave[f][0] = model.Wo; + model_octave[f][1] = model.L; + for (m = 1; m <= model.L; m++) { + model_octave[f][m + 1] = model.A[m]; + } + } + + /* Decoder */ + + MODEL model__[M]; + float prev_rate_K_vec_[K]; + COMP HH[M][MAX_AMP + 1]; + float Wo_left; + int voicing_left; + + /* initial conditions */ + + for (k = 0; k < K; k++) prev_rate_K_vec_[k] = 0.0; + voicing_left = 0; + Wo_left = 2.0 * M_PI / 100.0; + + /* decoder runs on every M-th frame, 25Hz frame rate, offset at + start is to minimise processing delay (thanks Jeroen!) */ + + fprintf(stderr, "\n"); + for (f = M - 1; f < FRAMES; f += M) { + float a_interpolated_surface_[M][K]; + newamp1_indexes_to_model( + &c2const, model__, (COMP *)HH, (float *)a_interpolated_surface_, + prev_rate_K_vec_, &Wo_left, &voicing_left, rate_K_sample_freqs_kHz, K, + phase_fft_fwd_cfg, phase_fft_inv_cfg, &indexes[f][0], NULL, 1); + +#ifdef VERBOSE + fprintf(stderr, "f: %d\n", f); + fprintf(stderr, " %d %d %d %d\n", indexes[f][0], indexes[f][1], + indexes[f][2], indexes[f][3]); + for (i = 0; i < M; i++) { + fprintf(stderr, " Wo: %4.3f L: %d v: %d\n", model__[i].Wo, model__[i].L, + model__[i].voiced); + } + fprintf(stderr, " rate_K_vec: "); + for (i = 0; i < 5; i++) { + fprintf(stderr, "%5.3f ", prev_rate_K_vec_[i]); + } + fprintf(stderr, "\n"); + fprintf(stderr, " Am H:\n"); + + for (m = 0; m < M; m++) { + fprintf(stderr, " "); + for (i = 1; i <= 5; i++) { + fprintf(stderr, "%5.1f (%5.3f %5.3f) ", model__[m].A[i], HH[m][i].real, + HH[m][i].imag); + } + fprintf(stderr, "\n"); + } + + fprintf(stderr, "\n\n"); +#endif + + // if (f == 7) + // exit(0); + + /* with f == 0, we don't store output, but memories are updated, helps to + match what happens in Codec 2 mode */ + + if (f >= M) { + for (i = 0; i < M; i++) { + for (k = 0; k < K; k++) { + interpolated_surface_[f - M + i][k] = a_interpolated_surface_[i][k]; + } + } + + /* store test vectors */ + + for (i = f - M, m = 0; i < f; i++, m++) { + model_octave_[i][0] = model__[m].Wo; + model_octave_[i][1] = model__[m].L; + voicing_[i] = model__[m].voiced; + } + + int j; + for (i = f - M, j = 0; i < f; i++, j++) { + for (m = 1; m <= model__[j].L; m++) { + model_octave_[i][m + 1] = model__[j].A[m]; + H[i][m - 1] = HH[j][m]; // aH[m]; + } + } + } + } + + fclose(fin); + + /* save vectors in Octave format */ + + FILE *fout = fopen("tnewamp1_out.txt", "wt"); + assert(fout != NULL); + fprintf(fout, "# Created by tnewamp1.c\n"); + octave_save_float(fout, "rate_K_surface_c", (float *)rate_K_surface, FRAMES, + K, K); + octave_save_float(fout, "mean_c", (float *)mean, 1, FRAMES, 1); + octave_save_float(fout, "eq_c", eq, 1, K, K); + octave_save_float(fout, "rate_K_surface_no_mean_c", + (float *)rate_K_surface_no_mean, FRAMES, K, K); + octave_save_float(fout, "rate_K_surface_no_mean__c", + (float *)rate_K_surface_no_mean_, FRAMES, K, K); + octave_save_float(fout, "mean__c", (float *)mean_, FRAMES, 1, 1); + octave_save_float(fout, "rate_K_surface__c", (float *)rate_K_surface_, FRAMES, + K, K); + octave_save_float(fout, "interpolated_surface__c", + (float *)interpolated_surface_, FRAMES, K, K); + octave_save_float(fout, "model_c", (float *)model_octave, FRAMES, MAX_AMP + 2, + MAX_AMP + 2); + octave_save_float(fout, "model__c", (float *)model_octave_, FRAMES, + MAX_AMP + 2, MAX_AMP + 2); + octave_save_int(fout, "voicing__c", (int *)voicing_, 1, FRAMES); + octave_save_complex(fout, "H_c", (COMP *)H, FRAMES, MAX_AMP, MAX_AMP); + fclose(fout); + + printf( + "Done! Now run\n octave:1> " + "tnewamp1(\"../path/to/build_linux/src/hts1a\", " + "\"../path/to/build_linux/unittest\")\n"); + return 0; +} diff --git a/unittest/tofdm.c b/unittest/tofdm.c new file mode 100644 index 0000000..b4d60b7 --- /dev/null +++ b/unittest/tofdm.c @@ -0,0 +1,626 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tofdm.c + AUTHORS.....: David Rowe & Steve Sampson + DATE CREATED: June 2017 + + Tests for the C version of the OFDM modem. This program outputs a + file of Octave vectors that are loaded and automatically tested + against the Octave version of the modem by the Octave script tofdm.m + +\*---------------------------------------------------------------------------*/ + +/* + Copyright (C) 2017 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 <complex.h> +#include <getopt.h> +#include <math.h> +#include <stdint.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "HRA_112_112.h" /* generated by ldpc_fsk_lib.m:ldpc_decode() */ +#include "codec2_ofdm.h" +#include "comp_prim.h" +#include "mpdecode_core.h" +#include "octave.h" +#include "ofdm_internal.h" +#include "test_bits_ofdm.h" + +#define NFRAMES 10 +#define SAMPLE_CLOCK_OFFSET_PPM 100 +#define FOFF_HZ 0.5f + +#define ASCALE (2E5f * 1.1491f / 2.0f) /* scale from shorts back to floats */ + +#define CODED_BITSPERFRAME 224 /* number of LDPC codeword bits/frame */ + +/* QPSK constellation for symbol likelihood calculations */ + +static COMP S_matrix[] = { + {1.0f, 0.0f}, {0.0f, 1.0f}, {0.0f, -1.0f}, {-1.0f, 0.0f}}; + +/* constants we use a lot and don't want to have to deference all the time */ + +static float ofdm_tx_centre; /* TX Center frequency */ +static float ofdm_rx_centre; /* RX Center frequency */ +static float ofdm_fs; /* Sample rate */ +static float ofdm_ts; /* Symbol cycle time */ +static float ofdm_rs; /* Symbol rate */ +static float ofdm_tcp; /* Cyclic prefix duration */ +static float + ofdm_timing_mx_thresh; /* See 700D Part 4 Acquisition blog post and + ofdm_dev.m routines for how this was set */ + +static int ofdm_nc; /* NS-1 data symbols between pilots */ +static int ofdm_ns; +static int ofdm_bps; /* Bits per symbol */ +static int ofdm_m; /* duration of each symbol in samples */ +static int ofdm_ncp; /* duration of CP in samples */ + +static int ofdm_ftwindowwidth; +static int ofdm_bitsperframe; +static int ofdm_rowsperframe; +static int ofdm_samplesperframe; +static int ofdm_samplespersymbol; +static int ofdm_max_samplesperframe; +static int ofdm_nrxbuf; +static int + ofdm_ntxtbits; /* reserve bits/frame for auxiliary text information */ +static int ofdm_nuwbits; /* Unique word, used for positive indication of lock */ + +/*---------------------------------------------------------------------------*\ + + FUNCTION....: fs_offset() + AUTHOR......: David Rowe + DATE CREATED: May 2015 + + Simulates small Fs offset between mod and demod. + (Note: Won't work with float, works OK with double) + +\*---------------------------------------------------------------------------*/ + +static int fs_offset(COMP out[], COMP in[], int n, float sample_rate_ppm) { + double f; + double tin = 0.0; + double step = 1.0 + sample_rate_ppm / 1E6; + int t1, t2; + int tout = 0; + + while (tin < (double)(n - 1)) { + t1 = (int)floor(tin); + t2 = (int)ceil(tin); + assert(t2 < n); + f = tin - (double)t1; + + out[tout].real = + ((double)1.0 - f) * (double)in[t1].real + f * (double)in[t2].real; + out[tout].imag = + ((double)1.0 - f) * (double)in[t1].imag + f * (double)in[t2].imag; + + tin += step; + tout++; + // printf("tin: %f tout: %d f: %f\n", tin, tout, f); + } + + return tout; +} + +/*---------------------------------------------------------------------------*\ + + FUNCTION....: freq_shift() + AUTHOR......: David Rowe + DATE CREATED: 26/4/2012 + + Frequency shift modem signal. The use of complex input and output allows + single sided frequency shifting (no images). + +\*---------------------------------------------------------------------------*/ + +static void freq_shift(COMP rx_fdm_fcorr[], COMP rx_fdm[], float foff, + COMP *foff_phase_rect, int nin) { + float temp = (TAU * foff / ofdm_fs); + COMP foff_rect = {cosf(temp), sinf(temp)}; + int i; + + for (i = 0; i < nin; i++) { + *foff_phase_rect = cmult(*foff_phase_rect, foff_rect); + rx_fdm_fcorr[i] = cmult(rx_fdm[i], *foff_phase_rect); + } + + /* normalise digital oscillator as the magnitude can drift over time */ + + float mag = cabsolute(*foff_phase_rect); + foff_phase_rect->real /= mag; + foff_phase_rect->imag /= mag; +} + +int main(int argc, char *argv[]) { + int opt_Nc = 0; + int ldpc_enable = 1; + struct OFDM *ofdm; + struct OFDM_CONFIG *ofdm_config; + + static struct option long_options[] = {{"nc", required_argument, 0, 'n'}, + {"noldpc", no_argument, 0, 'l'}, + {0, 0, 0, 0}}; + + int opt_index = 0; + char c; + + while ((c = getopt_long(argc, argv, "n:l", long_options, &opt_index)) != -1) { + switch (c) { + case 'n': + opt_Nc = atoi(optarg); + fprintf(stderr, "Nc = %d\n", opt_Nc); + break; + case 'l': + ldpc_enable = 0; + fprintf(stderr, "LDPC disabled\n"); + break; + default: + fprintf(stderr, + "usage: %s [Options]:\n [-l --noldpc]\n [-n --nc " + "NumberoFCarriers]\n", + argv[0]); + exit(1); + } + } + + // init once to get a copy of default config params + + ofdm = ofdm_create(NULL); + assert(ofdm != NULL); + struct OFDM_CONFIG ofdm_config_default; + memcpy(&ofdm_config_default, ofdm_get_config_param(ofdm), + sizeof(struct OFDM_CONFIG)); + ofdm_destroy(ofdm); + + // now do a little customisation on default config, and re-create modem + // instance + + if (opt_Nc) ofdm_config_default.nc = opt_Nc; + // printf("ofdm_create() 2\n"); + ofdm = ofdm_create(&ofdm_config_default); + assert(ofdm != NULL); + ofdm_config = ofdm_get_config_param(ofdm); + ofdm_set_tx_bpf(ofdm, false); + + // same levels as Octave sim + ofdm->amp_scale = 1.0; + + // make local copies for convenience + ofdm_tx_centre = ofdm_config->tx_centre; + ofdm_rx_centre = ofdm_config->rx_centre; + ofdm_fs = ofdm_config->fs; + ofdm_ts = ofdm_config->ts; + ofdm_rs = ofdm_config->rs; + ofdm_tcp = ofdm_config->tcp; + ofdm_timing_mx_thresh = ofdm_config->timing_mx_thresh; + ofdm_nc = ofdm_config->nc; + ofdm_ns = ofdm_config->ns; + ofdm_bps = ofdm_config->bps; + ofdm_m = (int)(ofdm_config->fs / ofdm_config->rs); + ofdm_ncp = (int)(ofdm_config->tcp * ofdm_config->fs); + ofdm_ftwindowwidth = ofdm_config->ftwindowwidth; + ofdm_bitsperframe = ofdm_get_bits_per_frame(ofdm); + ofdm_rowsperframe = ofdm_bitsperframe / (ofdm_config->nc * ofdm_config->bps); + ofdm_samplesperframe = ofdm_get_samples_per_frame(ofdm); + ofdm_samplespersymbol = (ofdm->m + ofdm->ncp); + ofdm_max_samplesperframe = ofdm_get_max_samples_per_frame(ofdm); + ofdm_nrxbuf = ofdm->nrxbuf; + ofdm_ntxtbits = ofdm_config->txtbits; + ofdm_nuwbits = ofdm_config->nuwbits; + + int tx_bits[ofdm_samplesperframe]; + COMP tx[ofdm_samplesperframe]; /* one frame of tx samples */ + + int rx_bits[ofdm_bitsperframe]; /* one frame of rx bits */ + printf("Nc = %d ofdm_bitsperframe: %d\n", ofdm_nc, ofdm_bitsperframe); + + /* log arrays */ + + int tx_bits_log[ofdm_bitsperframe * NFRAMES]; + COMP tx_log[ofdm_samplesperframe * NFRAMES]; + COMP rx_log[ofdm_samplesperframe * NFRAMES]; + COMP rxbuf_in_log[ofdm_max_samplesperframe * NFRAMES]; + COMP rxbuf_log[ofdm_nrxbuf * NFRAMES]; + COMP rx_sym_log[(ofdm_ns + 3) * NFRAMES][ofdm_nc + 2]; + float phase_est_pilot_log[ofdm_rowsperframe * NFRAMES][ofdm_nc]; + COMP rx_np_log[ofdm_rowsperframe * ofdm_nc * NFRAMES]; + float rx_amp_log[ofdm_rowsperframe * ofdm_nc * NFRAMES]; + float foff_hz_log[NFRAMES]; + int rx_bits_log[ofdm_bitsperframe * NFRAMES]; + int timing_est_log[NFRAMES]; + int timing_valid_log[NFRAMES]; + float timing_mx_log[NFRAMES]; + float coarse_foff_est_hz_log[NFRAMES]; + int sample_point_log[NFRAMES]; + float symbol_likelihood_log[(CODED_BITSPERFRAME / ofdm_bps) * + (1 << ofdm_bps) * NFRAMES]; + float bit_likelihood_log[CODED_BITSPERFRAME * NFRAMES]; + int detected_data_log[CODED_BITSPERFRAME * NFRAMES]; + float mean_amp_log[NFRAMES]; + float snr_log[NFRAMES]; + + FILE *fout; + int f, i, j; + + /* set up LDPC code */ + + struct LDPC ldpc; + + ldpc.max_iter = HRA_112_112_MAX_ITER; + ldpc.dec_type = 0; + ldpc.q_scale_factor = 1; + ldpc.r_scale_factor = 1; + ldpc.CodeLength = HRA_112_112_CODELENGTH; + ldpc.NumberParityBits = HRA_112_112_NUMBERPARITYBITS; + ldpc.NumberRowsHcols = HRA_112_112_NUMBERROWSHCOLS; + ldpc.max_row_weight = HRA_112_112_MAX_ROW_WEIGHT; + ldpc.max_col_weight = HRA_112_112_MAX_COL_WEIGHT; + ldpc.H_rows = (uint16_t *)HRA_112_112_H_rows; + ldpc.H_cols = (uint16_t *)HRA_112_112_H_cols; + + /* Main Loop + * ---------------------------------------------------------------------*/ + + for (f = 0; f < NFRAMES; f++) { + /* --------------------------------------------------------*\ + Mod + \*---------------------------------------------------------*/ + + /* See CML startup code in tofdm.m */ + + for (i = 0; i < ofdm_nuwbits; i++) { + tx_bits[i] = ofdm->tx_uw[i]; + } + for (i = ofdm_nuwbits; i < ofdm_nuwbits + ofdm_ntxtbits; i++) { + tx_bits[i] = 0; + } + + if (ldpc_enable) { + unsigned char ibits[HRA_112_112_NUMBERROWSHCOLS]; + unsigned char pbits[HRA_112_112_NUMBERPARITYBITS]; + + assert(HRA_112_112_NUMBERROWSHCOLS == ldpc.CodeLength / 2); + for (i = 0; i < ldpc.CodeLength / 2; i++) { + ibits[i] = payload_data_bits[i]; + } + encode(&ldpc, ibits, pbits); + for (j = 0, i = ofdm_nuwbits + ofdm_ntxtbits; j < ldpc.CodeLength / 2; + i++, j++) { + tx_bits[i] = ibits[j]; + } + for (j = 0; j < ldpc.CodeLength / 2; i++, j++) { + tx_bits[i] = pbits[j]; + } + assert(i == ofdm_bitsperframe); + } else { + int Npayloadbits = ofdm_bitsperframe - (ofdm_nuwbits + ofdm_ntxtbits); + uint16_t r[Npayloadbits]; + uint8_t payload_bits[Npayloadbits]; + + ofdm_rand(r, Npayloadbits); + for (i = 0; i < Npayloadbits; i++) payload_bits[i] = r[i] > 16384; + uint8_t txt_bits[ofdm_ntxtbits]; + for (i = 0; i < ofdm_ntxtbits; i++) txt_bits[i] = 0; + + uint8_t tx_bits_char[ofdm_bitsperframe]; + ofdm_assemble_qpsk_modem_packet(ofdm, tx_bits_char, payload_bits, + txt_bits); + for (i = 0; i < ofdm_bitsperframe; i++) tx_bits[i] = tx_bits_char[i]; + } + + ofdm_mod(ofdm, (COMP *)tx, tx_bits); + + /* tx vector logging */ + + memcpy(&tx_bits_log[ofdm_bitsperframe * f], tx_bits, + sizeof(int) * ofdm_bitsperframe); + memcpy(&tx_log[ofdm_samplesperframe * f], tx, + sizeof(COMP) * ofdm_samplesperframe); + } + + /* --------------------------------------------------------*\ + Channel + \*---------------------------------------------------------*/ + + fs_offset(rx_log, tx_log, ofdm_samplesperframe * NFRAMES, + SAMPLE_CLOCK_OFFSET_PPM); + + COMP foff_phase_rect = {1.0f, 0.0f}; + + freq_shift(rx_log, rx_log, FOFF_HZ, &foff_phase_rect, + ofdm_samplesperframe * NFRAMES); + + /* --------------------------------------------------------*\ + Demod + \*---------------------------------------------------------*/ + + /* Init/pre-load rx with ideal timing so we can test with timing estimation + * disabled */ + + int Nsam = ofdm_samplesperframe * NFRAMES; + int prx = 0; + int nin = ofdm_samplesperframe + 2 * ofdm_samplespersymbol; + + int lnew; + COMP rxbuf_in[ofdm_max_samplesperframe]; + +#define FRONT_LOAD +#ifdef FRONT_LOAD + for (i = 0; i < nin; i++, prx++) { + ofdm->rxbuf[ofdm_nrxbuf - nin + i] = + rx_log[prx].real + rx_log[prx].imag * I; + } +#endif + + int nin_tot = 0; + + /* disable estimators for initial testing */ + + ofdm_set_verbose(ofdm, false); + ofdm_set_timing_enable(ofdm, true); + ofdm_set_foff_est_enable(ofdm, true); + ofdm_set_phase_est_enable(ofdm, true); + +//#define TESTING_FILE +#ifdef TESTING_FILE + FILE *fin = fopen("~/codec2-dev/octave/ofdm_test.raw", "rb"); + assert(fin != NULL); + int Nbitsperframe = ofdm_bitsperframe; + int Nmaxsamperframe = ofdm_max_samplesperframe; + short rx_scaled[Nmaxsamperframe]; +#endif + + /* start this with something sensible otherwise LDPC decode fails in tofdm.m + */ + + ofdm->mean_amp = 1.0; + + for (f = 0; f < NFRAMES; f++) { + /* For initial testing, timing est is off, so nin is always + fixed. TODO: we need a constant for rxbuf_in[] size that + is the maximum possible nin */ + + nin = ofdm_get_nin(ofdm); + assert(nin <= ofdm_max_samplesperframe); + + /* 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. */ + + if ((Nsam - prx) < nin) { + lnew = Nsam - prx; + } else { + lnew = nin; + } + // printf("nin: %d prx: %d lnew: %d\n", nin, prx, lnew); + for (i = 0; i < nin; i++) { + rxbuf_in[i].real = 0.0; + rxbuf_in[i].imag = 0.0; + } + + if (lnew) { + for (i = 0; i < lnew; i++, prx++) { + rxbuf_in[i] = rx_log[prx]; + } + } + assert(prx <= ofdm_max_samplesperframe * NFRAMES); + +#ifdef TESTING_FILE + fread(rx_scaled, sizeof(short), nin, fin); + + for (i = 0; i < nin; i++) { + rxbuf_in[i].real = (float)rx_scaled[i] / ASCALE; + rxbuf_in[i].imag = 0.0; + } +#endif + + /* uncoded OFDM modem ---------------------------------------*/ + + ofdm_demod(ofdm, rx_bits, rxbuf_in); + +#ifdef TESTING_FILE + int Nerrs = 0; + for (i = 0; i < Nbitsperframe; i++) { + if (test_bits_ofdm[i] != rx_bits[i]) { + Nerrs++; + } + } + printf("f: %d Nerr: %d\n", f, Nerrs); +#endif + + float symbol_likelihood[(CODED_BITSPERFRAME / ofdm_bps) * (1 << ofdm_bps)]; + float bit_likelihood[CODED_BITSPERFRAME]; + uint8_t out_char[CODED_BITSPERFRAME]; + + if (ldpc_enable) { + /* LDPC functions --------------------------------------*/ + + float EsNo = 10; + + /* first few symbols are used for UW and txt bits, find start of (224,112) + * LDPC codeword */ + + assert((ofdm_nuwbits + ofdm_ntxtbits + CODED_BITSPERFRAME) == + ofdm_bitsperframe); + + COMP ldpc_codeword_symbols[(CODED_BITSPERFRAME / ofdm_bps)]; + + for (i = 0, j = (ofdm_nuwbits + ofdm_ntxtbits) / ofdm_bps; + i < (CODED_BITSPERFRAME / ofdm_bps); i++, j++) { + ldpc_codeword_symbols[i].real = crealf(ofdm->rx_np[j]); + ldpc_codeword_symbols[i].imag = cimagf(ofdm->rx_np[j]); + } + + float *ldpc_codeword_symbol_amps = + &ofdm->rx_amp[(ofdm_nuwbits + ofdm_ntxtbits) / ofdm_bps]; + + Demod2D(symbol_likelihood, ldpc_codeword_symbols, S_matrix, EsNo, + ldpc_codeword_symbol_amps, ofdm->mean_amp, + CODED_BITSPERFRAME / ofdm_bps); + Somap(bit_likelihood, symbol_likelihood, 1 << ofdm_bps, ofdm_bps, + CODED_BITSPERFRAME / ofdm_bps); + + float llr[CODED_BITSPERFRAME]; + int parityCheckCount; + + // fprintf(stderr, "\n"); + for (i = 0; i < CODED_BITSPERFRAME; i++) { + llr[i] = -bit_likelihood[i]; + // fprintf(stderr, "%f ", llr[i]); + } + + // fprintf(stderr, "\n"); + + run_ldpc_decoder(&ldpc, out_char, llr, &parityCheckCount); + /* + fprintf(stderr, "iter: %d parityCheckCount: %d\n", iter, + parityCheckCount); for(i=0; i<CODED_BITSPERFRAME; i++) { fprintf(stderr, + "%d ", out_char[i]); + } + */ + } + + /* rx vector logging -----------------------------------*/ + + assert(nin_tot < ofdm_samplesperframe * NFRAMES); + memcpy(&rxbuf_in_log[nin_tot], rxbuf_in, sizeof(COMP) * nin); + nin_tot += nin; + + for (i = 0; i < ofdm_nrxbuf; i++) { + rxbuf_log[ofdm_nrxbuf * f + i].real = crealf(ofdm->rxbuf[i]); + rxbuf_log[ofdm_nrxbuf * f + i].imag = cimagf(ofdm->rxbuf[i]); + } + + for (i = 0; i < (ofdm_ns + 3); i++) { + for (j = 0; j < (ofdm_nc + 2); j++) { + rx_sym_log[(ofdm_ns + 3) * f + i][j].real = crealf(ofdm->rx_sym[i][j]); + rx_sym_log[(ofdm_ns + 3) * f + i][j].imag = cimagf(ofdm->rx_sym[i][j]); + } + } + + /* note corrected phase (rx no phase) is one big linear array for frame */ + + for (i = 0; i < ofdm_rowsperframe * ofdm_nc; i++) { + rx_np_log[ofdm_rowsperframe * ofdm_nc * f + i].real = + crealf(ofdm->rx_np[i]); + rx_np_log[ofdm_rowsperframe * ofdm_nc * f + i].imag = + cimagf(ofdm->rx_np[i]); + } + + /* note phase/amp ests the same for each col, but check them all anyway */ + + for (i = 0; i < ofdm_rowsperframe; i++) { + for (j = 0; j < ofdm_nc; j++) { + phase_est_pilot_log[ofdm_rowsperframe * f + i][j] = + ofdm->aphase_est_pilot_log[ofdm_nc * i + j]; + rx_amp_log[ofdm_rowsperframe * ofdm_nc * f + ofdm_nc * i + j] = + ofdm->rx_amp[ofdm_nc * i + j]; + } + } + + foff_hz_log[f] = ofdm->foff_est_hz; + timing_est_log[f] = ofdm->timing_est + 1; /* offset by 1 to match Octave */ + timing_valid_log[f] = ofdm->timing_valid; + timing_mx_log[f] = ofdm->timing_mx; + coarse_foff_est_hz_log[f] = ofdm->coarse_foff_est_hz; + sample_point_log[f] = + ofdm->sample_point + 1; /* offset by 1 to match Octave */ + float EsNodB = ofdm_esno_est_calc(ofdm->rx_np, ofdm_rowsperframe * ofdm_nc); + snr_log[f] = ofdm_snr_from_esno(ofdm, EsNodB); + mean_amp_log[f] = ofdm->mean_amp; + + memcpy(&rx_bits_log[ofdm_bitsperframe * f], rx_bits, sizeof(rx_bits)); + + if (ldpc_enable) { + for (i = 0; i < (CODED_BITSPERFRAME / ofdm_bps) * (1 << ofdm_bps); i++) { + symbol_likelihood_log[(CODED_BITSPERFRAME / ofdm_bps) * + (1 << ofdm_bps) * f + + i] = symbol_likelihood[i]; + } + for (i = 0; i < CODED_BITSPERFRAME; i++) { + bit_likelihood_log[CODED_BITSPERFRAME * f + i] = bit_likelihood[i]; + detected_data_log[CODED_BITSPERFRAME * f + i] = out_char[i]; + } + } + } + + /*---------------------------------------------------------*\ + Dump logs to Octave file for evaluation + by tofdm.m Octave script + \*---------------------------------------------------------*/ + + fout = fopen("tofdm_out.txt", "wt"); + assert(fout != NULL); + fprintf(fout, "# Created by tofdm.c\n"); + octave_save_complex(fout, "pilot_samples_c", (COMP *)ofdm->pilot_samples, 1, + ofdm_samplespersymbol, ofdm_samplespersymbol); + octave_save_int(fout, "tx_bits_log_c", tx_bits_log, 1, + ofdm_bitsperframe * NFRAMES); + octave_save_complex(fout, "tx_log_c", (COMP *)tx_log, 1, + ofdm_samplesperframe * NFRAMES, + ofdm_samplesperframe * NFRAMES); + octave_save_complex(fout, "rx_log_c", (COMP *)rx_log, 1, + ofdm_samplesperframe * NFRAMES, + ofdm_samplesperframe * NFRAMES); + octave_save_complex(fout, "rxbuf_in_log_c", (COMP *)rxbuf_in_log, 1, nin_tot, + nin_tot); + octave_save_complex(fout, "rxbuf_log_c", (COMP *)rxbuf_log, 1, + ofdm_nrxbuf * NFRAMES, ofdm_nrxbuf * NFRAMES); + octave_save_complex(fout, "rx_sym_log_c", (COMP *)rx_sym_log, + (ofdm_ns + 3) * NFRAMES, ofdm_nc + 2, ofdm_nc + 2); + octave_save_float(fout, "phase_est_pilot_log_c", (float *)phase_est_pilot_log, + ofdm_rowsperframe * NFRAMES, ofdm_nc, ofdm_nc); + octave_save_float(fout, "rx_amp_log_c", (float *)rx_amp_log, 1, + ofdm_rowsperframe * ofdm_nc * NFRAMES, + ofdm_rowsperframe * ofdm_nc * NFRAMES); + octave_save_float(fout, "foff_hz_log_c", foff_hz_log, NFRAMES, 1, 1); + octave_save_int(fout, "timing_est_log_c", timing_est_log, NFRAMES, 1); + octave_save_int(fout, "timing_valid_log_c", timing_valid_log, NFRAMES, 1); + octave_save_float(fout, "timing_mx_log_c", timing_mx_log, NFRAMES, 1, 1); + octave_save_float(fout, "coarse_foff_est_hz_log_c", coarse_foff_est_hz_log, + NFRAMES, 1, 1); + octave_save_int(fout, "sample_point_log_c", sample_point_log, NFRAMES, 1); + octave_save_complex(fout, "rx_np_log_c", (COMP *)rx_np_log, 1, + ofdm_rowsperframe * ofdm_nc * NFRAMES, + ofdm_rowsperframe * ofdm_nc * NFRAMES); + octave_save_int(fout, "rx_bits_log_c", rx_bits_log, 1, + ofdm_bitsperframe * NFRAMES); + octave_save_float(fout, "symbol_likelihood_log_c", symbol_likelihood_log, + (CODED_BITSPERFRAME / ofdm_bps) * (1 << ofdm_bps) * NFRAMES, + 1, 1); + octave_save_float(fout, "bit_likelihood_log_c", bit_likelihood_log, + CODED_BITSPERFRAME * NFRAMES, 1, 1); + octave_save_int(fout, "detected_data_log_c", detected_data_log, 1, + CODED_BITSPERFRAME * NFRAMES); + octave_save_float(fout, "snr_log_c", snr_log, NFRAMES, 1, 1); + octave_save_float(fout, "mean_amp_log_c", mean_amp_log, NFRAMES, 1, 1); + fclose(fout); +#ifdef TESTING_FILE + fclose(fin); +#endif + + ofdm_destroy(ofdm); + + return 0; +} diff --git a/unittest/tofdm_acq.c b/unittest/tofdm_acq.c new file mode 100644 index 0000000..cfb0f99 --- /dev/null +++ b/unittest/tofdm_acq.c @@ -0,0 +1,94 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tofdm_acq.c + AUTHORS.....: David Rowe + DATE CREATED: Mar 2021 + + Tests for the acquistion (sync) parts of the C version of the OFDM modem. + This program outputs a file of Octave vectors that are loaded and + automatically tested against the Octave version of the modem by the Octave + script tofdm_acq.m + +\*---------------------------------------------------------------------------*/ + +#include <assert.h> +#include <math.h> +#include <stdint.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "codec2_ofdm.h" +#include "octave.h" +#include "ofdm_internal.h" + +#define MAX_FRAMES 500 + +int main(int argc, char *argv[]) { + struct OFDM *ofdm; + struct OFDM_CONFIG ofdm_config; + + ofdm_init_mode("datac0", &ofdm_config); + ofdm = ofdm_create(&ofdm_config); + ofdm->data_mode = "burst"; + ofdm->verbose = 2; + ofdm->timing_mx_thresh = 0.15; + ofdm->postambledetectoren = 1; + assert(ofdm != NULL); + + int nin = ofdm_get_nin(ofdm); + int rxbufst = ofdm->rxbufst; + + FILE *fin = fopen(argv[1], "rb"); + assert(fin != NULL); + short rx_scaled[ofdm_get_max_samples_per_frame(ofdm)]; + int f = 0; + + float timing_mx_log[MAX_FRAMES]; + int ct_est_log[MAX_FRAMES]; + float foff_est_log[MAX_FRAMES]; + int timing_valid_log[MAX_FRAMES]; + int nin_log[MAX_FRAMES]; + + while (fread(rx_scaled, sizeof(short), nin, fin) == nin) { + fprintf(stderr, "%3d ", f); + ofdm_sync_search_shorts(ofdm, rx_scaled, ofdm->amp_scale / 2.0f); + + if (f < MAX_FRAMES) { + timing_mx_log[f] = ofdm->timing_mx; + ct_est_log[f] = ofdm->ct_est; + foff_est_log[f] = ofdm->foff_est_hz; + timing_valid_log[f] = ofdm->timing_valid; + nin_log[f] = ofdm->nin; + } + f++; + + // reset these to defaults, as they get modified when timing_valid asserted + ofdm->nin = nin; + ofdm->rxbufst = rxbufst; + } + fclose(fin); + + /*---------------------------------------------------------*\ + Dump logs to Octave file for evaluation + by tofdm_acq.m Octave script + \*---------------------------------------------------------*/ + + FILE *fout = fopen("tofdm_acq_out.txt", "wt"); + assert(fout != NULL); + fprintf(fout, "# Created by tofdm_acq.c\n"); + octave_save_complex(fout, "tx_preamble_c", (COMP *)ofdm->tx_preamble, 1, + ofdm->samplesperframe, ofdm->samplesperframe); + octave_save_complex(fout, "tx_postamble_c", (COMP *)ofdm->tx_postamble, 1, + ofdm->samplesperframe, ofdm->samplesperframe); + octave_save_float(fout, "timing_mx_log_c", timing_mx_log, 1, f, f); + octave_save_float(fout, "foff_est_log_c", foff_est_log, 1, f, f); + octave_save_int(fout, "ct_est_log_c", ct_est_log, 1, f); + octave_save_int(fout, "timing_valid_log_c", timing_valid_log, 1, f); + octave_save_int(fout, "nin_log_c", nin_log, 1, f); + fclose(fout); + + ofdm_destroy(ofdm); + + return 0; +} diff --git a/unittest/tqam16.c b/unittest/tqam16.c new file mode 100644 index 0000000..9e4c8d0 --- /dev/null +++ b/unittest/tqam16.c @@ -0,0 +1,35 @@ +/*---------------------------------------------------------------------------*\ + + FILE........: tqam16.c + AUTHOR......: David Rowe + DATE CREATED: August 2020 + + Simple sanity check for QAM16 symbol mapping. + +\*---------------------------------------------------------------------------*/ + +#include <stdio.h> +#include <string.h> + +#include "ofdm_internal.h" + +int main(void) { + int c; + for (c = 0; c < 16; c++) { + int tx_bits[4], rx_bits[4]; + for (int i = 0; i < 4; i++) tx_bits[i] = (c >> (3 - i)) & 0x1; + complex float symbol = qam16_mod(tx_bits); + qam16_demod(symbol, rx_bits); + if (memcmp(tx_bits, rx_bits, 4)) { + fprintf(stderr, "FAIL on %d!\ntx_bits: ", c); + for (int i = 0; i < 4; i++) fprintf(stderr, "%d ", tx_bits[i]); + fprintf(stderr, "%f %f\nrx_bits: ", creal(symbol), cimag(symbol)); + for (int i = 0; i < 4; i++) fprintf(stderr, "%d ", rx_bits[i]); + fprintf(stderr, "%f %f\n", creal(symbol), cimag(symbol)); + return 1; + } + } + + fprintf(stderr, "%d tested OK...\nPASS!\n", c); + return 0; +} diff --git a/unittest/tquisk_filter.c b/unittest/tquisk_filter.c new file mode 100644 index 0000000..1879dff --- /dev/null +++ b/unittest/tquisk_filter.c @@ -0,0 +1,46 @@ +/* + tquisk_filter.c + + Unit test for complex band pass filters in src/filter.c + + cd codec2/build_linux + ./misc/mksine - 1500 2 | unittest/tquisk_filter | aplay + + By adjusting the frequency you can audibly test filter response. +*/ + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> + +#include "filter.h" +#include "filter_coef.h" + +#define N 159 /* processing buffer size (odd number deliberate) */ +#define CENTRE 1500.0 +#define FS 8000.0 + +int main() { + short buf_short[N]; + complex float buf[N]; + struct quisk_cfFilter *bpf; + int i; + + bpf = malloc(sizeof(struct quisk_cfFilter)); + assert(bpf != NULL); + quisk_filt_cfInit(bpf, filtP200S400, sizeof(filtP200S400) / sizeof(float)); + quisk_cfTune(bpf, CENTRE / FS); + + while (fread(buf_short, sizeof(short), N, stdin) == N) { + for (i = 0; i < N; i++) buf[i] = buf_short[i]; + quisk_ccfFilter(buf, buf, N, bpf); + /* we only output the real part in this test */ + for (i = 0; i < N; i++) buf_short[i] = creal(buf[i]); + fwrite(buf_short, sizeof(short), N, stdout); + } + + quisk_filt_destroy(bpf); + free(bpf); + return 0; +} diff --git a/unittest/tvq_mbest.c b/unittest/tvq_mbest.c new file mode 100644 index 0000000..1ea52a9 --- /dev/null +++ b/unittest/tvq_mbest.c @@ -0,0 +1,32 @@ +/* + tvq_mbest.c + David Rowe Dec 2019 + + Generate some test vectors to exercise misc/vq_mbest.c +*/ + +#include <assert.h> +#include <stdio.h> +#include <stdlib.h> + +void write_float_file(char fn[], float *values, int n) { + FILE *f = fopen(fn, "wb"); + assert(f != NULL); + assert(fwrite(values, sizeof(float), n, f) == n); + fclose(f); +} + +int main(void) { + /* we're only interested in searching the inner 2 values, outer elements + should be ignored */ + float target[] = {0.0, 1.0, 1.0, 0.0}; + write_float_file("target.f32", target, 4); + float vq1[] = { + 1.0, 0.9, 0.9, 1.0, /* this will be a better match on first stage */ + 2.0, 0.8, 0.8, 2.0}; /* but after second stage should choose this */ + write_float_file("vq1.f32", vq1, 8); + float vq2[] = {10.0, 0.3, 0.3, 10.0, 20.0, + 0.2, 0.2, 20.0}; /* 0.8+0.2 == 1.0 so best 2nd stage entry */ + write_float_file("vq2.f32", vq2, 8); + return 0; +} diff --git a/unittest/vq_mbest.c b/unittest/vq_mbest.c new file mode 100644 index 0000000..62b7e95 --- /dev/null +++ b/unittest/vq_mbest.c @@ -0,0 +1,300 @@ +/* + vq_mbest.c + David Rowe Dec 2019 + + Utility to perform a mbest VQ search on vectors from stdin, sending + quantised vectors to stdout. +*/ + +#include <assert.h> +#include <getopt.h> +#include <limits.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "mbest.h" + +#define MAX_K 20 +#define MAX_ENTRIES 4096 +#define MAX_STAGES 5 + +void quant_mbest(float vec_out[], int indexes[], float vec_in[], int num_stages, + float vqw[], float vq[], int m[], int k, int mbest_survivors); + +int verbose = 0; + +int main(int argc, char *argv[]) { + float vq[MAX_STAGES * MAX_K * MAX_ENTRIES]; + float vqw[MAX_STAGES * MAX_K * MAX_ENTRIES]; + int m[MAX_STAGES]; + int k = 0, mbest_survivors = 1, num_stages = 0; + char fnames[256], fn[256], *comma, *p; + FILE *fq; + float lower = -1E32; + int st = -1; + int en = -1; + int num = INT_MAX; + int output_vec_usage = 0; + + int o = 0; + int opt_idx = 0; + while (o != -1) { + static struct option long_opts[] = {{"k", required_argument, 0, 'k'}, + {"quant", required_argument, 0, 'q'}, + {"mbest", required_argument, 0, 'm'}, + {"lower", required_argument, 0, 'l'}, + {"verbose", required_argument, 0, 'v'}, + {"st", required_argument, 0, 't'}, + {"en", required_argument, 0, 'e'}, + {"num", required_argument, 0, 'n'}, + {"vec_usage", no_argument, 0, 'u'}, + {0, 0, 0, 0}}; + + o = getopt_long(argc, argv, "hk:q:m:vt:e:n:u", long_opts, &opt_idx); + switch (o) { + case 'k': + k = atoi(optarg); + assert(k <= MAX_K); + break; + case 'q': + /* load up list of comma delimited file names */ + strcpy(fnames, optarg); + p = fnames; + num_stages = 0; + do { + assert(num_stages < MAX_STAGES); + strcpy(fn, p); + comma = strchr(fn, ','); + if (comma) { + *comma = 0; + p = comma + 1; + } + /* load quantiser file */ + fprintf(stderr, "stage: %d loading %s ... ", num_stages, fn); + fq = fopen(fn, "rb"); + if (fq == NULL) { + fprintf(stderr, "Couldn't open: %s\n", fn); + exit(1); + } + /* count how many entries m of dimension k are in this VQ file */ + m[num_stages] = 0; + float dummy[k]; + while (fread(dummy, sizeof(float), k, fq) == (size_t)k) + m[num_stages]++; + assert(m[num_stages] <= MAX_ENTRIES); + fprintf(stderr, "%d entries of vectors width %d\n", m[num_stages], k); + /* now load VQ into memory */ + rewind(fq); + int rd = fread(&vq[num_stages * k * MAX_ENTRIES], sizeof(float), + m[num_stages] * k, fq); + assert(rd == m[num_stages] * k); + num_stages++; + fclose(fq); + } while (comma); + break; + case 'm': + mbest_survivors = atoi(optarg); + fprintf(stderr, "mbest_survivors = %d\n", mbest_survivors); + break; + case 'n': + num = atoi(optarg); + break; + case 'l': + lower = atof(optarg); + break; + case 't': + st = atoi(optarg); + break; + case 'e': + en = atoi(optarg); + break; + case 'u': + output_vec_usage = 1; + break; + case 'v': + verbose = 1; + break; + help: + fprintf(stderr, "\n"); + fprintf(stderr, + "usage: %s -k dimension -q vq1.f32,vq2.f32,.... [Options]\n", + argv[0]); + fprintf(stderr, "\n"); + fprintf(stderr, + "input vectors on stdin, output quantised vectors on stdout\n"); + fprintf(stderr, "\n"); + fprintf(stderr, + "--lower lowermeanLimit Only count vectors with average " + "above this level in distortion calculations\n"); + fprintf(stderr, + "--mbest N number of survivors at each stage, " + "set to 0 for standard VQ search\n"); + fprintf(stderr, + "--st Kst start vector element for error " + "calculation (default 0)\n"); + fprintf(stderr, + "--en Ken end vector element for error " + "calculation (default K-1)\n"); + fprintf(stderr, + "--num numToProcess number of vectors to quantise " + "(default to EOF)\n"); + fprintf(stderr, + "--vec_usage Output a record of how many times " + "each vector is used\n"); + fprintf(stderr, "-v Verbose\n"); + exit(1); + } + } + + if ((num_stages == 0) || (k == 0)) goto help; + + /* default to measuring error on entire vector */ + if (st == -1) st = 0; + if (en == -1) en = k - 1; + + float w[k]; + for (int i = 0; i < st; i++) w[i] = 0.0; + for (int i = st; i <= en; i++) w[i] = 1.0; + for (int i = en + 1; i < k; i++) w[i] = 0.0; + + /* apply weighting to codebook (rather than in search) */ + memcpy(vqw, vq, sizeof(vq)); + for (int s = 0; s < num_stages; s++) { + mbest_precompute_weight(&vqw[s * k * MAX_ENTRIES], w, k, m[s]); + } + + int indexes[num_stages], nvecs = 0; + int vec_usage[m[0]]; + for (int i = 0; i < m[0]; i++) vec_usage[i] = 0; + float target[k], quantised[k]; + float sqe = 0.0; + while (fread(&target, sizeof(float), k, stdin) && (nvecs < num)) { + for (int i = 0; i < k; i++) target[i] *= w[i]; + int dont_count = 0; + /* optional clamping to lower limit or mean */ + float mean = 0.0; + for (int i = 0; i < k; i++) mean += target[i]; + mean /= k; + float difference = mean - lower; + if (difference < 0.0) { + /* bring target up to lower clamping limit */ + for (int i = 0; i < k; i++) target[i] += -difference; + dont_count = 1; + } + quant_mbest(quantised, indexes, target, num_stages, vqw, vq, m, k, + mbest_survivors); + if (dont_count == 0) { + for (int i = st; i <= en; i++) sqe += pow(target[i] - quantised[i], 2.0); + } + fwrite(&quantised, sizeof(float), k, stdout); + nvecs++; + // count number f time each vector is used (just for first stage) + vec_usage[indexes[0]]++; + } + + fprintf(stderr, "MSE: %4.2f\n", sqe / (nvecs * (en - st + 1))); + + if (output_vec_usage) { + for (int i = 0; i < m[0]; i++) fprintf(stderr, "%d\n", vec_usage[i]); + } + + return 0; +} + +// print vector debug function + +void pv(char s[], float v[], int k) { + int i; + if (verbose) { + fprintf(stderr, "%s", s); + for (i = 0; i < k; i++) fprintf(stderr, "%4.2f ", v[i]); + fprintf(stderr, "\n"); + } +} + +// mbest algorithm version, backported from LPCNet/src + +void quant_mbest(float vec_out[], int indexes[], float vec_in[], int num_stages, + float vqw[], float vq[], int m[], int k, int mbest_survivors) { + float err[k], se1; + int i, j, s, s1, ind; + + struct MBEST *mbest_stage[num_stages]; + int index[num_stages]; + float target[k]; + + for (i = 0; i < num_stages; i++) { + mbest_stage[i] = mbest_create(mbest_survivors); + index[i] = 0; + } + + se1 = 0.0; + for (i = 0; i < k; i++) { + err[i] = vec_in[i]; + se1 += err[i] * err[i]; + } + se1 /= k; + + /* now quantise err[] using multi-stage mbest search, preserving + mbest_survivors at each stage */ + + mbest_search(vqw, err, k, m[0], mbest_stage[0], index); + if (verbose) mbest_print("Stage 1:", mbest_stage[0]); + + for (s = 1; s < num_stages; s++) { + /* for each candidate in previous stage, try to find best vector in next + * stage */ + for (j = 0; j < mbest_survivors; j++) { + /* indexes that lead us this far */ + for (s1 = 0; s1 < s; s1++) { + index[s1 + 1] = mbest_stage[s - 1]->list[j].index[s1]; + } + /* target is residual err[] vector given path to this candidate */ + for (i = 0; i < k; i++) target[i] = err[i]; + for (s1 = 0; s1 < s; s1++) { + ind = index[s - s1]; + if (verbose) + fprintf(stderr, " s: %d s1: %d s-s1: %d ind: %d\n", s, s1, s - s1, + ind); + for (i = 0; i < k; i++) { + target[i] -= vqw[s1 * k * MAX_ENTRIES + ind * k + i]; + } + } + pv(" target: ", target, k); + mbest_search(&vqw[s * k * MAX_ENTRIES], target, k, m[s], mbest_stage[s], + index); + } + char str[80]; + sprintf(str, "Stage %d:", s + 1); + if (verbose) mbest_print(str, mbest_stage[s]); + } + + for (s = 0; s < num_stages; s++) { + indexes[s] = mbest_stage[num_stages - 1]->list[0].index[num_stages - 1 - s]; + } + + /* OK put it all back together using best survivor */ + for (i = 0; i < k; i++) vec_out[i] = 0.0; + for (s = 0; s < num_stages; s++) { + int ind = indexes[s]; + float se2 = 0.0; + for (i = 0; i < k; i++) { + err[i] -= vqw[s * k * MAX_ENTRIES + ind * k + i]; + vec_out[i] += vq[s * k * MAX_ENTRIES + ind * k + i]; + se2 += err[i] * err[i]; + } + se2 /= k; + pv(" err: ", err, k); + if (verbose) fprintf(stderr, " se2: %f\n", se2); + } + pv(" vec_out: ", vec_out, k); + + pv("\n vec_in: ", vec_in, k); + pv(" vec_out: ", vec_out, k); + pv(" err: ", err, k); + if (verbose) fprintf(stderr, " se1: %f\n", se1); + + for (i = 0; i < num_stages; i++) mbest_destroy(mbest_stage[i]); +} |