First cut of c2wideband.c and its wideband_enc_dec() function.
authorphilayres <philayres@01035d8c-6547-0410-b346-abe4f91aad63>
Wed, 19 Jul 2017 13:08:57 +0000 (13:08 +0000)
committerphilayres <philayres@01035d8c-6547-0410-b346-abe4f91aad63>
Wed, 19 Jul 2017 13:08:57 +0000 (13:08 +0000)
A unittest/tc2wideband.c file has been created to test some of the core calculation functions, plus to run wideband_enc_dec() with random data in order to track down buffer overflows, etc.
At this commit, the unit tests in this new file run. The core functions operate as expected. wideband_enc_dec() runs on a single block of 16 frames without crashing, and has been stepped through to validate that most code is working as expected. Real data has not been used at this point in time.
Created generate_wideband_map.c, which creates c2wideband_map.h from the Octave text file ../Octave/c2wideband_map
codec2.c, c2dec.c, c2enc.c have been adapted to include the new wideband option.
Makefiles have been adapted and appear to work.

git-svn-id: https://svn.code.sf.net/p/freetel/code@3314 01035d8c-6547-0410-b346-abe4f91aad63

15 files changed:
codec2-dev/octave/c2wideband_batch.m
codec2-dev/octave/c2wideband_map [new file with mode: 0644]
codec2-dev/raw/speech_orig_16k.wav [new file with mode: 0644]
codec2-dev/src/CMakeLists.txt
codec2-dev/src/c2dec.c
codec2-dev/src/c2enc.c
codec2-dev/src/c2wideband.c [new file with mode: 0644]
codec2-dev/src/c2wideband.h [new file with mode: 0644]
codec2-dev/src/c2wideband_map.h [new file with mode: 0644]
codec2-dev/src/codebook/c2wideband_map.txt [new file with mode: 0644]
codec2-dev/src/codec2.c
codec2-dev/src/codec2.h
codec2-dev/src/generate_wideband_map.c [new file with mode: 0644]
codec2-dev/unittest/CMakeLists.txt
codec2-dev/unittest/tc2wideband.c [new file with mode: 0644]

index c4ba7147a7d913fe89a0eb9754e843d1804e339f..55dac99d4fe795c2fda2287359dac1bfe42e822a 100644 (file)
@@ -197,7 +197,8 @@ function [model_ rate_K_surface] = experiment_rate_K_dct2(model, plots=1)
   % break into blocks of (Nt time samples) x (K freq samples)
 
   Nblocks = floor(frames/(Nt*dec));
-  %printf("frames: %d Nblocks: %d\n", frames, Nblocks);
+  
+  printf("frames: %d Nblocks: %d\n", frames, Nblocks);
 
   % map that defines order we read out and quantise DCT coeffs
   % TODO: for C port we need an Octave function to write Map to a C
diff --git a/codec2-dev/octave/c2wideband_map b/codec2-dev/octave/c2wideband_map
new file mode 100644 (file)
index 0000000..946aa19
--- /dev/null
@@ -0,0 +1,8 @@
+ 1.00000000e+00 3.00000000e+00 7.00000000e+00 8.00000000e+00 1.10000000e+01 1.50000000e+01 1.90000000e+01 2.60000000e+01 2.10000000e+01 2.40000000e+01 2.00000000e+01 3.00000000e+01 3.80000000e+01 4.80000000e+01 2.90000000e+01 3.20000000e+01 4.20000000e+01 6.40000000e+01 6.30000000e+01 5.40000000e+01 5.60000000e+01 5.00000000e+01 7.20000000e+01 9.10000000e+01 7.80000000e+01 6.70000000e+01 5.70000000e+01 7.00000000e+01 7.60000000e+01 1.04000000e+02
+ 2.00000000e+00 5.00000000e+00 1.20000000e+01 2.20000000e+01 2.30000000e+01 3.30000000e+01 3.10000000e+01 4.30000000e+01 3.40000000e+01 3.90000000e+01 4.40000000e+01 4.10000000e+01 4.90000000e+01 6.00000000e+01 5.10000000e+01 7.70000000e+01 9.70000000e+01 9.00000000e+01 1.14000000e+02 8.10000000e+01 1.21000000e+02 1.13000000e+02 8.00000000e+01 9.90000000e+01 1.18000000e+02 1.09000000e+02 8.90000000e+01 1.02000000e+02 1.07000000e+02 1.19000000e+02
+ 4.00000000e+00 9.00000000e+00 1.60000000e+01 2.70000000e+01 4.00000000e+01 4.50000000e+01 4.70000000e+01 5.50000000e+01 4.60000000e+01 6.10000000e+01 6.50000000e+01 6.60000000e+01 7.50000000e+01 8.30000000e+01 6.80000000e+01 1.06000000e+02 9.40000000e+01 1.27000000e+02 1.36000000e+02 1.38000000e+02 1.37000000e+02 1.41000000e+02 8.80000000e+01 1.50000000e+02 1.59000000e+02 1.47000000e+02 1.33000000e+02 1.28000000e+02 1.90000000e+02 1.31000000e+02
+ 6.00000000e+00 1.80000000e+01 2.80000000e+01 5.20000000e+01 8.40000000e+01 6.20000000e+01 7.30000000e+01 7.90000000e+01 7.10000000e+01 9.80000000e+01 8.50000000e+01 9.30000000e+01 8.70000000e+01 1.12000000e+02 1.15000000e+02 1.10000000e+02 1.29000000e+02 1.99000000e+02 1.42000000e+02 1.70000000e+02 2.18000000e+02 1.62000000e+02 1.74000000e+02 2.05000000e+02 1.43000000e+02 1.64000000e+02 1.72000000e+02 1.75000000e+02 2.02000000e+02 2.16000000e+02
+ 1.00000000e+01 2.50000000e+01 3.50000000e+01 6.90000000e+01 7.40000000e+01 9.20000000e+01 9.60000000e+01 1.35000000e+02 1.11000000e+02 1.56000000e+02 9.50000000e+01 1.05000000e+02 1.26000000e+02 1.34000000e+02 1.63000000e+02 2.15000000e+02 1.49000000e+02 1.48000000e+02 1.53000000e+02 2.28000000e+02 1.80000000e+02 2.30000000e+02 2.11000000e+02 2.13000000e+02 2.40000000e+02 2.24000000e+02 2.09000000e+02 2.22000000e+02 2.26000000e+02 1.73000000e+02
+ 1.30000000e+01 3.70000000e+01 5.30000000e+01 1.08000000e+02 8.60000000e+01 1.17000000e+02 1.03000000e+02 1.00000000e+02 1.24000000e+02 1.22000000e+02 2.01000000e+02 1.92000000e+02 2.03000000e+02 2.19000000e+02 1.91000000e+02 1.69000000e+02 2.20000000e+02 2.14000000e+02 1.39000000e+02 2.04000000e+02 2.36000000e+02 1.79000000e+02 1.81000000e+02 2.00000000e+02 1.93000000e+02 2.21000000e+02 2.35000000e+02 1.87000000e+02 2.08000000e+02 2.17000000e+02
+ 1.40000000e+01 3.60000000e+01 5.90000000e+01 1.16000000e+02 1.30000000e+02 1.78000000e+02 1.40000000e+02 1.32000000e+02 1.57000000e+02 1.52000000e+02 1.68000000e+02 1.60000000e+02 1.23000000e+02 1.97000000e+02 1.83000000e+02 2.34000000e+02 1.44000000e+02 2.23000000e+02 1.65000000e+02 1.67000000e+02 2.31000000e+02 1.96000000e+02 1.71000000e+02 2.33000000e+02 1.82000000e+02 2.25000000e+02 1.88000000e+02 2.37000000e+02 2.27000000e+02 2.10000000e+02
+ 1.70000000e+01 5.80000000e+01 8.20000000e+01 1.01000000e+02 1.20000000e+02 1.86000000e+02 1.46000000e+02 1.25000000e+02 1.51000000e+02 1.94000000e+02 1.61000000e+02 1.84000000e+02 1.58000000e+02 1.89000000e+02 1.77000000e+02 1.95000000e+02 1.45000000e+02 1.66000000e+02 2.06000000e+02 1.54000000e+02 1.98000000e+02 2.12000000e+02 2.38000000e+02 2.29000000e+02 1.85000000e+02 1.76000000e+02 2.39000000e+02 1.55000000e+02 2.07000000e+02 2.32000000e+02
diff --git a/codec2-dev/raw/speech_orig_16k.wav b/codec2-dev/raw/speech_orig_16k.wav
new file mode 100644 (file)
index 0000000..737d9a5
Binary files /dev/null and b/codec2-dev/raw/speech_orig_16k.wav differ
index 7dc9e8c036944f723bce75c1cbaac8135a7ef908..7d758d37233eb8579a4dbeb832b1437e80c6642d 100644 (file)
@@ -66,6 +66,13 @@ set(CODEBOOKSNEWAMP1_ENERGY
     ${D}/newamp1_energy_q.txt
 )
 
+
+# generated wideband map
+set(C2WB_MAP
+    ${D}/c2wideband_map.txt
+)
+
+
 # when crosscompiling we need a native executable
 if(CMAKE_CROSSCOMPILING)
     include(ExternalProject)
@@ -146,6 +153,40 @@ add_custom_command(
     DEPENDS generate_codebook ${CODEBOOKSNEWAMP1_ENERGY}
 )
 
+
+
+# when crosscompiling we need a native executable
+if(CMAKE_CROSSCOMPILING)
+    include(ExternalProject)
+    ExternalProject_Add(codec2_native
+       SOURCE_DIR ${CMAKE_SOURCE_DIR}
+       BUILD_COMMAND $(MAKE) generate_wideband_map
+       INSTALL_COMMAND ${CMAKE_COMMAND} -E copy src/generate_wideband_map ${CMAKE_CURRENT_BINARY_DIR}
+    )
+    add_executable(generate_wideband_map IMPORTED)
+    set_target_properties(generate_wideband_map PROPERTIES
+        IMPORTED_LOCATION ${CMAKE_CURRENT_BINARY_DIR}/generate_wideband_map)
+    add_dependencies(generate_wideband_map codec2_native)
+else(CMAKE_CROSSCOMPILING)
+# Build code generator binaries. These do not get installed.
+    # generate_codebook
+    add_executable(generate_wideband_map generate_wideband_map.c)
+    target_link_libraries(generate_wideband_map ${CMAKE_REQUIRED_LIBRARIES})
+    # Make native builds available for cross-compiling.
+    export(TARGETS generate_wideband_map
+        FILE ${CMAKE_BINARY_DIR}/ImportExecutables.cmake)
+endif(CMAKE_CROSSCOMPILING)
+
+
+# c2wideband_map.h
+add_custom_command(
+    OUTPUT ${CMAKE_CURRENT_BINARY_DIR}/c2wideband_map.h
+    COMMAND generate_wideband_map c2wideband_map ${C2WB_MAP} > ${CMAKE_CURRENT_BINARY_DIR}/c2wideband_map.h
+    DEPENDS generate_wideband_map ${C2WB_MAP}
+)
+
+
+
 #
 # codec2 library sources
 #
@@ -188,6 +229,7 @@ set(CODEC2_SRCS
     freedv_data_channel.c
     varicode.c
     modem_stats.c
+    c2wideband.c
 )
 
 set(CODEC2_PUBLIC_HEADERS
index 57f081543d94ac8bb9008324575601d0fbcd87ff..1a590ef4eae20aa1fd5aafe22180cf8e4e18f2cf 100644 (file)
@@ -98,11 +98,13 @@ int main(int argc, char *argv[])
     else if (strcmp(argv[1],"700") == 0)
        mode = CODEC2_MODE_700;
     else if (strcmp(argv[1],"700B") == 0)
-       mode = CODEC2_MODE_700B;
-     else if (strcmp(argv[1],"700C") == 0)
+       mode = CODEC2_MODE_700B; 
+    else if (strcmp(argv[1],"700C") == 0)
        mode = CODEC2_MODE_700C;
+    else if (strcmp(argv[1],"WB") == 0)
+       mode = CODEC2_MODE_WB;
    else {
-       fprintf(stderr, "Error in mode: %s.  Must be 3200, 2400, 1600, 1400, 1300, 1200, 700, 700B or 700C\n", argv[1]);
+       fprintf(stderr, "Error in mode: %s.  Must be 3200, 2400, 1600, 1400, 1300, 1200, 700, 700B, 700C or WB\n", argv[1]);
        exit(1);
     }
     bit_rate = atoi(argv[1]);
@@ -320,7 +322,7 @@ void print_help(const struct option* long_options, int num_opts, char* argv[])
        int i;
        char *option_parameters;
        fprintf(stderr, "\nc2dec - Codec 2 decoder and bit error simulation program\n"
-               "usage: %s 3200|2400|1600|1400|1300|1200|700|700B InputFile OutputRawFile [OPTIONS]\n\n"
+               "usage: %s 3200|2400|1600|1400|1300|1200|700|700B|WB InputFile OutputRawFile [OPTIONS]\n\n"
                 "Options:\n", argv[0]);
         for(i=0; i<num_opts-1; i++) {
                if(long_options[i].has_arg == no_argument) {
index 713ce2c1ec58f2d620d6f44f128658e32da33759..a1eeb9ae095e85a6950a5f88f596a05ede5c2711 100644 (file)
@@ -46,7 +46,7 @@ int main(int argc, char *argv[])
     int            bit, byte,i;
 
     if (argc < 4) {
-       printf("usage: c2enc 3200|2400|1600|1400|1300|1200|700|700B|700C InputRawspeechFile OutputBitFile [--natural] [--softdec]\n");
+       printf("usage: c2enc 3200|2400|1600|1400|1300|1200|700|700B|700C|WB InputRawspeechFile OutputBitFile [--natural] [--softdec]\n");
        printf("e.g    c2enc 1400 ../raw/hts1a.raw hts1a.c2\n");
        printf("e.g    c2enc 1300 ../raw/hts1a.raw hts1a.c2 --natural\n");
        exit(1);
@@ -70,8 +70,10 @@ int main(int argc, char *argv[])
        mode = CODEC2_MODE_700B;
     else if (strcmp(argv[1],"700C") == 0)
        mode = CODEC2_MODE_700C;
+    else if (strcmp(argv[1],"WB") == 0)
+       mode = CODEC2_MODE_WB;
     else {
-       fprintf(stderr, "Error in mode: %s.  Must be 3200, 2400, 1600, 1400, 1300, 1200, 700, 700B or 700C\n", argv[1]);
+       fprintf(stderr, "Error in mode: %s.  Must be 3200, 2400, 1600, 1400, 1300, 1200, 700, 700B, 700C or WB\n", argv[1]);
        exit(1);
     }
 
diff --git a/codec2-dev/src/c2wideband.c b/codec2-dev/src/c2wideband.c
new file mode 100644 (file)
index 0000000..94cab45
--- /dev/null
@@ -0,0 +1,723 @@
+/*---------------------------------------------------------------------------*\
+
+  FILE........: c2wideband.c
+  AUTHOR......: David Rowe & Phil Ayres
+  DATE CREATED: July 2017
+
+ * TODO - what is this file doing?
+ * 
+\*---------------------------------------------------------------------------*/
+
+/*
+  Copyright David Rowe 2017
+
+  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 <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+#include <unistd.h>
+
+#include "defines.h"
+#include "codec2_fft.h"
+#include "sine.h"
+#include "nlp.h"
+#include "dump.h"
+#include "lpc.h"
+#include "quantise.h"
+#include "phase.h"
+#include "interp.h"
+#include "postfilter.h"
+#include "codec2.h"
+#include "lsp.h"
+#include "codec2_internal.h"
+#include "machdep.h"
+#include "codec2.h"
+#include "newamp1.h"
+
+#include "c2wideband.h"
+#include "c2wideband_map.h"
+
+float mean(float data[], int n);
+int unit_test();
+float std(float data[], int rows);
+void diff_de(int rows, int cols, float D[rows][cols], float E[rows][cols], float diff_de[rows][cols]);
+void array_col_to_row(int rows, int cols, float data[rows][cols], int col, float res[]);
+void std_on_cols(int rows, int cols, float data[rows][cols], float res[]);
+void setup_map(WIDEBAND_MAP * wb_map, int Nt, int K);
+#define C2WB_PLOT
+#ifdef C2WB_PLOT
+static int total_frames;
+static int curr_block_pos;
+#endif
+
+
+//TODO: move all this to the standard codec2.c functions
+#ifdef IGNORE
+
+int main(int argc, char *argv[]) {
+  struct CODEC2 *codec2;
+  FILE *fin;
+  FILE *fout;
+  short *buf;
+  unsigned char *bits;
+  int nsam, nbit, i, r;
+  int errno = -1;
+  for (i = 0; i < 10; i++) {
+    r = codec2_rand();
+    printf("[%d] r = %d\n", i, r);
+  }
+
+  if (argc != 3) {
+    printf("usage: %s InputRawSpeechFile OutputRawSpeechFile\n", argv[0]);
+    exit(1);
+  }
+
+  if ((fin = fopen(argv[1], "rb")) == NULL) {
+    fprintf(stderr, "Error opening input speech file: %s: %s.\n",
+            argv[1], strerror(errno));
+    exit(1);
+  }
+
+  if ((fout = fopen(argv[2], "wb")) == NULL) {
+    fprintf(stderr, "Error opening output speech file: %s: %s.\n",
+            argv[2], strerror(errno));
+    exit(1);
+  }
+
+#ifdef DUMP
+  dump_on("c2demo");
+#endif
+
+  /* Note only one set of Codec 2 states is required for an encoder
+     and decoder pair. */
+
+  codec2 = codec2_create(CODEC2_MODE_WB);
+  nsam = C2WB_SPERF; //codec2_samples_per_frame(codec2);
+  buf = (short*) malloc(nsam * sizeof (short));
+  nbit = C2WB_BPERF; //codec2_bits_per_frame(codec2);
+  bits = (unsigned char*) malloc(nbit * sizeof (char));
+
+  while (fread(buf, sizeof (short), nsam, fin) == (size_t) nsam) {
+    codec2_encode_wb(codec2, bits, buf);
+    codec2_decode_wb(codec2, buf, bits);
+    fwrite(buf, sizeof (short), nsam, fout);
+  }
+
+  free(buf);
+  free(bits);
+  codec2_destroy(codec2);
+
+  fclose(fin);
+  fclose(fout);
+
+  return 0;
+}
+#endif
+
+void codec2_encode_wb(struct CODEC2 *c2, unsigned char * bits, short speech[]) {
+  return;
+}
+
+void codec2_decode_wb(struct CODEC2 *c2, short speech[], const unsigned char * bits) {
+  return;
+}
+
+void calculate_Am_freqs_kHz(float Wo, int L, float Am_freqs_kHz[L]) {
+  //(1:L)*Wo*4/pi;
+
+  for (int i = 0; i < L; i++) {
+    Am_freqs_kHz[i] = i * Wo * 4 / PI;
+  }
+}
+
+
+// Adapted from newamp.m
+
+void resample_const_rate_f_mel(C2CONST *c2const, MODEL * model, float K, float* rate_K_surface, float* rate_K_sample_freqs_kHz) {
+  mel_sample_freqs_kHz(rate_K_sample_freqs_kHz, K, 100, 0.95 * c2const->Fs / 2);
+  resample_const_rate_f(c2const, model, rate_K_surface, rate_K_sample_freqs_kHz, K);
+}
+
+
+// Updates rate_K_vec with the corrected values
+// function [rate_K_vec_corrected orig_error error nasty_error_log nasty_error_m_log] = correct_rate_K_vec(rate_K_vec, rate_K_sample_freqs_kHz, AmdB, AmdB_, K, Wo, L, Fs)
+
+void correct_rate_K_vec(MODEL *model, float rate_K_vec[], float rate_K_sample_freqs_kHz[], float Am_freqs_kHz[], float orig_AmdB[], int K, float Wo, int L, int Fs, float rate_K_vec_corrected[]) {
+
+  /*
+      % aliasing correction --------------------------------------
+
+      % The mel sample rate decreases as frequency increases. Look for
+      % any regions above 1000Hz where we have missed definition of a
+      % spectral peak (formant) due to aliasing.  Adjust the rate K
+      % sample levels to restore peaks.  Theory is that correct
+      % definition of a formant is less important than the frequency of
+      % the formant.  As long as we define a formant in that general
+      % frequency area it will sound OK.
+   */
+
+
+  // this is passed in
+  //Am_freqs_kHz = (1:L)*Wo*4/pi;    
+
+  //% Lets see where we have made an error
+  float error[MAX_AMP + 1];
+  //float orig_error[MAX_AMP+1];
+  float AmdB_[MAX_AMP + 1];
+
+  float nasty_error_freq;
+  float rate_K_prev_sample_kHz;
+  float rate_K_next_sample_kHz;
+
+  int Ncorrections = 3; //% maximum number of rate K samples to correct
+  int error_thresh = 3; //% only worry about errors larger than thresh    
+
+  int closest_k = 0;
+
+  // regenerate the AmdB values from the updated model
+  for (int i = 0; i < MAX_AMP + 1; i++) {
+    AmdB_[i] = 20 * log10(model->A[i]);
+  }
+
+  // calculate error between original AmdB and the new values
+  for (int i = 0; i < MAX_AMP + 1; i++) {
+    int a = orig_AmdB[i] - AmdB_[i];
+    error[i] = a;
+    //  orig_error[i] = a;
+  }
+
+
+
+  //% first 1000Hz is densely sampled so ignore
+  int start_m = floor(L * 1000 / (Fs / 2));
+  for (int i = 0; i < start_m; i++) {
+    error[i] = 0;
+  }
+  //start_m = floor(L*1000/(Fs/2));
+  //error(1:start_m) = 0;  % first 1000Hz is densly sampled so ignore
+
+  // ignore these
+  //float nasty_error_m_log = []; 
+  //float nasty_error_log = [];
+
+
+  // could probably memcpy this, but I need to think about it...
+  // TODO copy??
+  for (int i = 0; i < K; i++) {
+    //rate_K_vec_corrected = rate_K_vec
+    rate_K_vec_corrected[i] = rate_K_vec[i];
+  }
+
+  for (int i = 0; i < Ncorrections; i++) {
+
+    // [mx mx_m] = max(error);
+    int mx_m = 0;
+    int mx = 0;
+    for (int m = 0; m < start_m; m++) {
+      if (error[m] > mx) {
+        mx_m = m;
+        mx = error[m];
+      }
+    }
+
+    if (mx > error_thresh) {
+
+      // ignore these
+      //nasty_error_log = [nasty_error_log mx];
+      //nasty_error_m_log = [nasty_error_m_log mx_m];
+
+      //% find closest rate K sample to nasty error
+
+      nasty_error_freq = (float) mx_m * Wo * Fs / (2 * PI * 1000);
+
+      //[tmp closest_k] = min(abs(rate_K_sample_freqs_kHz - nasty_error_freq));
+      closest_k = -1;
+      float cka = K;
+      for (int k = 0; k < K; k++) {
+        float a = abs(rate_K_sample_freqs_kHz[k] - nasty_error_freq);
+
+        if (closest_k == -1 || a < cka) {
+          closest_k = k;
+          cka = a;
+        }
+      }
+
+      // this sets the value from the original AmdB value set for the model prior to recalculation
+      // TODO - check this is correct
+      rate_K_vec_corrected[closest_k] = orig_AmdB[mx_m];
+
+      //% zero out error in this region and look for another large error region
+
+      // here I'm assuming that the required result is a zero based index
+      // rather than an Octave 1 based index
+      // closest_k is already 0 based in this code
+      // TODO - check this
+      int k = (0 > closest_k - 1 ? 0 : closest_k - 1);
+      //k = max(1, closest_k-1); 
+
+      rate_K_prev_sample_kHz = rate_K_sample_freqs_kHz[k];
+
+
+      // again, assuming the required result is 0 rather than 1 based
+      // TODO - check 
+      k = (K - 1 < closest_k + 1 ? K : closest_k + 1);
+      //k = min(K, closest_k+1); 
+
+      rate_K_next_sample_kHz = rate_K_sample_freqs_kHz[k];
+
+      int st = -1;
+      int st_m = 0;
+      int en = -1;
+      int en_m = 0;
+      for (int m = 0; m < C2WB_K - 1; m++) {
+        //[tmp st_m] = min(abs(Am_freqs_kHz - rate_K_prev_sample_kHz));
+        //[tmp en_m] = min(abs(Am_freqs_kHz - rate_K_next_sample_kHz));
+        int pa = abs(Am_freqs_kHz[m] - rate_K_prev_sample_kHz);
+        int na = abs(abs(Am_freqs_kHz[m] - rate_K_next_sample_kHz));
+        if (st == -1 || pa < st) {
+          st = pa;
+          st_m = m;
+        }
+        if (en == -1 || na < en) {
+          en = na;
+          en_m = m;
+        }
+      }
+      if (closest_k == K)
+        en_m = L;
+
+      for (int m = st_m; m < en_m; m++) {
+        error[m] = 0;
+      }
+      //error(st_m:en_m) = 0;
+    }
+  }
+
+
+
+}
+
+void diff_de(int rows, int cols, float D[rows][cols], float E[rows][cols], float diff_de[rows][cols]) {
+
+
+  for (int col = 0; col < cols; col++) {
+    for (int row = 0; row < rows; row++) {
+      diff_de[row][col] = D[row][col] - E[row][col];
+    }
+  }
+
+}
+
+float mean(float data[], int n) {
+  float sum = 0.0;
+  for (int i = 0; i < n; ++i) {
+    sum += data[i];
+  }
+
+  return sum / n;
+}
+
+void array_col_to_row(int rows, int cols, float data[rows][cols], int col, float res[]) {
+  for (int row = 0; row < rows; row++) {
+    res[row] = data[row][col];
+  }
+}
+
+float std(float data[], int rows) {
+  float standardDeviation = 0.0;
+
+  for (int i = 0; i < rows; ++i) {
+    standardDeviation += pow(data[i] - mean(data, rows), 2);
+  }
+  return sqrt(standardDeviation / rows);
+}
+
+void std_on_cols(int rows, int cols, float data[rows][cols], float res[]) {
+  float row_from_col[cols];
+  for (int col = 0; col < cols; col++) {
+    array_col_to_row(rows, cols, data, col, row_from_col);
+    res[col] = std(row_from_col, rows);
+  }
+}
+
+float mean_std_diff_de(int rows, int cols, float D[rows][cols], float E[rows][cols]) {
+  float matrix_diff_de[rows][cols];
+  float std_array[cols];
+  diff_de(rows, cols, D, E, matrix_diff_de);
+  std_on_cols(rows, cols, matrix_diff_de, std_array);
+  return mean(std_array, cols);
+}
+
+//TODO : everything
+
+void dct2(int rows, int cols, float matrix[], float res[rows][cols]) {
+
+  for (int r = 0; r < rows; r++) {
+    for (int c = 0; c < cols; c++) {
+      res[r][c] = (float) rand() / RAND_MAX;
+    }
+  }
+  /*
+     
+    if m == 1
+      y = dct (x.', n).';
+    elseif n == 1
+      y = dct (x, m);
+    else
+      y = dct (dct (x, m).', n).';
+    endif
+   */
+
+}
+
+//TODO : everything
+
+void idct2(int rows, int cols, float inrks[rows][cols], float rate_K_surface_block[]) {
+  int K = C2WB_K;
+  for (int k = 0; k < K; k++) {
+    rate_K_surface_block[k] = (float) rand() / RAND_MAX;
+  }
+}
+//% Encode/decoder a 160ms block of model parameters
+//% TODO: (i) quantisation of DCT coeffs (ii) break into separate encoder and decoder functions
+//
+//function [model_block_ dct2_sd qn rate_K_surface_block rate_K_surface_block_] = wideband_enc_dec(model_block, rmap, cmap)
+
+void wideband_enc_dec(C2CONST *c2const, int n_block_frames, MODEL model_block[n_block_frames], WIDEBAND_MAP * wb_map,
+        MODEL model_block_[n_block_frames], float dct2_sd[n_block_frames], int * p_qn, float rate_K_surface_block[n_block_frames][C2WB_K], float rate_K_surface_block_[n_block_frames][C2WB_K]) {
+
+  //    c2wideband_const;
+  //
+  //    sim_quant = 1; % used to simulate quantisation, set to 1,2,4, etc
+  //    dist_dB   = 2; % use enough coefficients to get this distortion ond DCT coeffs
+  int sim_quant = 1;
+  float dist_dB = 2;
+
+  //    int Fs = c2const->Fs; 
+  //    int L;   
+  //    float Wo;
+
+  int K = C2WB_K;
+  int Nt = C2WB_NT;
+
+  int rows = Nt;
+  int cols = K;
+  int dec = C2WB_DEC;
+  printf("starting iteration\n");
+  // iterate through the frames in the block
+  for (int f = 0; f < n_block_frames; f++) {
+    printf("processing frame %d of %d\n", f, n_block_frames);
+
+    float rate_K_sample_freqs_kHz[C2WB_K];
+    //MODEL * model = &model_block[f];
+
+
+    //    % Resample variable rate L vectors to fixed length rate K.  We have
+    //    % left high end correction out for now, this is less of an issue
+    //    % with a higher K
+
+
+    //    [rate_K_surface_block rate_K_sample_freqs_kHz] = resample_const_rate_f_mel(model_block, K, Fs);
+    printf("resample_const_rate_f_mel\n");
+    resample_const_rate_f_mel(c2const, &model_block[f], K, rate_K_surface_block[f], rate_K_sample_freqs_kHz);
+
+    //    % decimate down to 20ms time resolution, and DCT
+
+    //    D = dct2(rate_K_surface_block(1:dec:Nt*dec,:));
+
+    float D[rows][cols];
+    float E[rows][cols];
+    printf("dct2\n");
+    dct2(rows, cols, rate_K_surface_block[f], D);
+
+    //    % So D is the 2D block of DCT coeffs at the encoder.  We want to
+    //    % create a quantised version at the "decoder" E.  This loop copies
+    //    % DCTs coeffs from D to E, until we get beneath a distortion
+    //    % threshold.
+    //
+    //    % This is essentially variable rate quantisation, but gives us
+    //    % some idea of the final bit rate.  In practice we will also need
+    //    % to constrain the total number of bits (ie bit rate), and
+    //    % quantise each coefficient.
+    //
+    //    % Turns out than mean SD (across many blocks/frames) is about the
+    //    % same in the DCT domain as the rate K domain.  So we can just
+    //    % measure MSE between D and E to estimate mean SD of the rate K
+    //    % vectors after quantisation.
+
+    //    E = mapped = zeros(Nt,K);                
+    for (int r = 0; r < rows; r++) {
+      memset(E[r], '\0', cols * sizeof (float));
+    }
+
+    //    qn = 0;
+    *p_qn = 0;
+
+    //    adct2_sd = mean(std(D-E));
+
+    float adct2_sd;
+
+
+    adct2_sd = mean_std_diff_de(rows, cols, D, E);
+
+    printf("while adct2_sd > dist_dB (dist_dB=%f)\n", dist_dB);
+    //    while adct2_sd > dist_dB
+    //TODO : validate the addition of  *p_qn limit check is correct
+    while (adct2_sd > dist_dB && *p_qn < C2WB_NT) {
+
+      //      qn++;                        
+      // --- this has moved to the end to cope with 0 base
+      assert(*p_qn < C2WB_NT);
+      int rmapqn = wb_map->rmap[*p_qn];
+      int cmapqn = wb_map->cmap[*p_qn];
+      assert(rmapqn < rows);
+      assert(cmapqn < cols);
+      //      E(rmap(qn), cmap(qn)) = sim_quant*round(D(rmap(qn), cmap(qn))/sim_quant);
+      E[rmapqn][cmapqn] = sim_quant * round(D[rmapqn][cmapqn] / sim_quant);
+      //      adct2_sd = mean(std(D-E));
+      adct2_sd = mean_std_diff_de(rows, cols, D, E);
+
+      //printf("qn %d %f\n", *p_qn, adct2_sd);
+
+      (*p_qn)++;
+      //      %printf("qn %d %f\n", qn, adct2_sd);
+      //    end
+    }
+
+    //    % note neat trick to interpolate to 10ms frames despite dec to 20ms, this means
+    //    % we don't need a separate decode side interpolator.
+
+    //    dct2_sd = mean(std(D-E));
+
+
+    dct2_sd[f] = mean_std_diff_de(rows, cols, D, E);
+
+
+    //    rate_K_surface_block_ = idct2([sqrt(dec)*E; zeros(Nt*(dec-1), K)]);
+    float inrks[rows * dec][K];
+
+    printf("setting inrks\n");
+    for (int r = 0; r < rows; r++) {
+      for (int c = 0; c < cols; c++) {
+        inrks[r][c] = sqrt(dec) * E[r][c];
+      }
+    }
+    for (int r = 0; r < Nt * (dec - 1); r++) {
+      for (int c = 0; c < K; c++) {
+        inrks[r + rows][c] = 0;
+      }
+    }
+
+
+    //TODO ???
+    //= [sqrt(dec)*E; zeros(Nt*(dec-1), K)];
+    printf("idct2\n");
+    idct2(rows, cols, inrks, rate_K_surface_block_[f]);
+
+    //    model_block_ = resample_rate_L(model_block, rate_K_surface_block_, rate_K_sample_freqs_kHz, Fs);        
+    printf("resample_rate_L\n");
+    resample_rate_L(c2const, &model_block[f], rate_K_surface_block_[f], rate_K_sample_freqs_kHz, K);
+    //endfunction 
+  }
+}
+
+void setup_map(WIDEBAND_MAP * wb_map, int Nt, int K) {
+
+  /*
+    % map that defines order we read out and quantise DCT coeffs
+    % TODO: for C port we need an Octave function to write Map to a C
+    % include file
+
+    map = load("c2wideband_map");
+
+    % create arrays to reverse map quantiser_num to r,c Luts
+
+    rmap = cmap = zeros(1,Nt*K);
+    for r=1:Nt
+      for c=1:K
+         quantiser_num = map(r,c);
+         rmap(quantiser_num) = r;
+         cmap(quantiser_num) = c;
+      end
+    end
+
+   */
+
+
+  int quantiser_num;
+
+  memset(wb_map->rmap, '\0', Nt * K * sizeof *wb_map->rmap);
+  memset(wb_map->cmap, '\0', Nt * K * sizeof *wb_map->cmap);
+
+  for (int r = 0; r < Nt; r++) {
+    for (int c = 0; c < K; c++) {
+      quantiser_num = c2wideband_map[r][c];
+      wb_map->rmap[quantiser_num] = r;
+      wb_map->cmap[quantiser_num] = c;
+    }
+  }
+
+}
+
+/*
+% ---------------------------------------------------------------------------------------
+% rate K mel-resampling, high end correction, and DCT experiment workhorse
+
+function [model_ rate_K_surface] = experiment_rate_K_dct2(model, plots=1)
+ */
+
+//TODO - review and produce a unit test for this
+// this is really just a draft at this point
+void experiment_rate_K_dct2(C2CONST *c2const, MODEL model_frames[], int frames, int vq_en, int plots, int* voicing, float *mean_) {
+
+
+  //  newamp;
+  //  c2wideband_const;
+  //  [frames nc] = size(model);
+
+  int n_block_frames = C2WB_NT * C2WB_DEC;
+
+  int K = C2WB_K;
+  int Nt = C2WB_NT;
+  //int dec = C2WB_DEC;    
+  WIDEBAND_MAP wb_map;
+
+
+
+#ifdef C2WB_PLOT
+  total_frames = frames;
+#endif
+
+  //  % break into blocks of (Nt time samples) x (K freq samples)
+
+  // Nblocks = floor(frames/(Nt*dec));
+  // %printf("frames: %d Nblocks: %d\n", frames, Nblocks);
+  //int Nblocks = floor(frames/(Nt*dec));
+
+
+  setup_map(&wb_map, Nt, K);
+
+  // Step through the model in blocks of 8 frames
+  for (int f = 0; f < frames; f += n_block_frames) {
+    MODEL * model_block = &model_frames[f];
+
+#ifdef C2WB_PLOT
+    curr_block_pos = f;
+#endif
+    rate_K_dct2(c2const, n_block_frames, model_block, &wb_map);
+
+
+  }
+
+}
+
+
+
+
+//TODO - review and produce a unit test for this
+// this is really just a draft at this point
+void rate_K_dct2(C2CONST *c2const, int n_block_frames, MODEL model_block[n_block_frames], WIDEBAND_MAP * wb_map) {
+
+  //% per-block processing ----------------------------------------------------
+
+  //  % init a bunch of output variables
+
+  //  rate_K_surface_ = zeros(Nblocks*Nt*dec, K);  
+  // sumnz = zeros(1,Nblocks);
+  // dct2_sd = zeros(1,Nblocks);
+  // model_ = [];
+
+
+  int Nblocks = n_block_frames;
+  int dec = C2WB_DEC;
+  int Nt = C2WB_NT;
+  int K = C2WB_K;
+  int Tf = C2WB_TF;
+  int rate_K_surface_size = Nblocks * Nt * dec;
+
+  float rate_K_surface_block[rate_K_surface_size][K]; // rate K vecs for each frame, form a surface that makes pretty graphs   
+  float rate_K_surface_block_[rate_K_surface_size][K];
+
+#ifdef C2WB_PLOT    
+  float rate_K_surface[total_frames][K];
+  float rate_K_surface_[total_frames][K];
+  MODEL model_[total_frames];
+#endif
+
+  float sumnz[Nblocks];
+  float dct2_sd[Nblocks];
+  int qn;
+
+  MODEL model_block_[n_block_frames];
+
+  //for n=1:Nblocks
+  for (int n = 0; n < Nblocks; n++) {
+    //st = (n-1)*dec*Nt+1; en = st + dec*Nt - 1;
+    int st = (n - 1) * dec * Nt + 1;
+    int en = st + dec * Nt - 1;
+    //%printf("st: %d en: %d\n", st, en);
+
+    //[model_block_ adct2_sd qn rate_K_surface_block rate_K_surface_block_] = wideband_enc_dec(model(st:en,:), rmap, cmap);
+    wideband_enc_dec(c2const, n_block_frames, model_block, wb_map,
+            model_block_, &dct2_sd[n], &qn, &rate_K_surface_block[n], &rate_K_surface_block_[n]
+            );
+
+#ifdef C2WB_PLOT
+    //model_ = [model_; model_block_];
+    model_[curr_block_pos + n] = model_block_[n];
+
+    //% log these for plotting/development
+
+    //rate_K_surface(st:en,:) = rate_K_surface_block;
+    //rate_K_surface_(st:en,:) = rate_K_surface_block_;
+    int pst = 0;
+    for (int p = st; p < en; p++) {
+      for (int k = 0; k < K; k++) {
+        rate_K_surface[p][k] = rate_K_surface_block[n + pst][k];
+        rate_K_surface_[p][k] = rate_K_surface_block_[n + pst][k];
+      }
+      pst++;
+    }
+    //dct2_sd(n) = adct2_sd;            
+    sumnz[n] = (float) qn;
+    // already handled in call
+#endif
+
+    //end
+  }
+
+
+#ifdef C2WB_PLOT
+  printf("average dct spectral distortion: %3.2f dB\n", mean(dct2_sd, Nblocks));
+  printf("mean number of coeffs/DCT: %3.2f/%d\n", mean(sumnz, Nblocks), Nt * K);
+  printf("coeffs/second: %3.2f\n", mean(sumnz, Nblocks) / (Nt * Tf * dec));
+  printf("bits/s: %3.2f\n", 2.9 * mean(sumnz, Nblocks) / (Nt * Tf * dec));
+
+  //TODO
+  //dist = std((rate_K_surface_(1:dec:Nblocks*Nt*dec,:) - rate_K_surface(1:dec:Nblocks*Nt*dec,:))');
+  //float dist = 
+  //if plots
+  //figure(1); clf; plot(dist); title('Rate K SD');
+  //printf("Rate K spectral distortion mean: %3.2f dB var: %3.2f\n", mean(dist), var(dist));
+  //end
+#endif
+
+
+
+}
diff --git a/codec2-dev/src/c2wideband.h b/codec2-dev/src/c2wideband.h
new file mode 100644 (file)
index 0000000..d41379c
--- /dev/null
@@ -0,0 +1,57 @@
+/* 
+ * File:   c2wideband.h
+ * Author: phil
+ *
+ * Created on 02 July 2017, 20:42
+ */
+
+#ifndef C2WIDEBAND_H 
+#define        C2WIDEBAND_H 
+
+#include "defines.h"
+#include "phase.h"
+#include "quantise.h"
+#include "newamp1.h"
+#include "codec2_internal.h"
+
+#define C2WB_K           30  /* rate K vector length */
+#define C2WB_FS       16000 
+
+#define C2WB_NT           8  /* number of blocks in time = 160ms blocks */
+#define C2WB_TF           0.02  /* 20ms frames */
+#define C2WB_DEC          2  /* decimation factor */
+#define C2WB_SPERF       30  /* samples per frame */
+//TODO: decide what this is
+#define C2WB_BPERF      256  /* bits per frame */
+
+
+typedef struct {
+    int rmap[C2WB_K * C2WB_NT];
+    int cmap[C2WB_K * C2WB_NT];
+} WIDEBAND_MAP;
+
+
+
+void codec2_decode_wb(struct CODEC2 *c2, short speech[], const unsigned char * bits);
+
+void calculate_Am_freqs_kHz(float Wo, int L, float p_Am_freqs_kHz[]);
+void resample_const_rate_f_mel(C2CONST *c2const, MODEL * model, float K, float* rate_K_surface, float* rate_K_sample_freqs_kHz);
+void correct_rate_K_vec(MODEL *model, float rate_K_vec[], float rate_K_sample_freqs_kHz[], float Am_freqs_kHz[], float orig_AmdB[],  int K, float Wo, int L, int Fs, float rate_K_vec_corrected[]);
+void batch_rate_K_dct2(C2CONST *c2const, MODEL model_frames[], int frames, int vq_en, int plots, int* voicing, float *mean_);
+void rate_K_dct2(C2CONST *c2const, int n_block_frames, MODEL model_block[n_block_frames], WIDEBAND_MAP * wb_map);
+void wideband_enc_dec(C2CONST *c2const, int n_block_frames, MODEL model_block[], WIDEBAND_MAP * wb_map,
+        MODEL model_block_[], float * p_dct2_sd,  int * p_qn , float rate_K_surface_block[][C2WB_K], float rate_K_surface_block_[][C2WB_K]
+        );
+void codec2_decode_wb(struct CODEC2 *c2, short speech[], const unsigned char * bits);
+void codec2_encode_wb(struct CODEC2 *c2, unsigned char * bits, short speech[]);
+
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* C2WIDEBAND_H */
+
diff --git a/codec2-dev/src/c2wideband_map.h b/codec2-dev/src/c2wideband_map.h
new file mode 100644 (file)
index 0000000..9ea5a3a
--- /dev/null
@@ -0,0 +1,20 @@
+/* THIS IS A GENERATED FILE. Edit generate_wideband_map.c and its input */
+
+/*
+ * This intermediary file and the files that used to create it are under 
+ * The LGPL. See the file COPYING.
+ */
+
+#include "defines.h"
+
+  /* ./codec2-dev/src/codebook/c2wideband_map.txt */
+static const float c2wideband_map[8][30] = {
+{   1,   3,   7,   8,   11,   15,   19,   26,   21,   24,   20,   30,   38,   48,   29,   32,   42,   64,   63,   54,   56,   50,   72,   91,   78,   67,   57,   70,   76,   104 },
+{   2,   5,   12,   22,   23,   33,   31,   43,   34,   39,   44,   41,   49,   60,   51,   77,   97,   90,   114,   81,   121,   113,   80,   99,   118,   109,   89,   102,   107,   119 },
+{   4,   9,   16,   27,   40,   45,   47,   55,   46,   61,   65,   66,   75,   83,   68,   106,   94,   127,   136,   138,   137,   141,   88,   150,   159,   147,   133,   128,   190,   131 },
+{   6,   18,   28,   52,   84,   62,   73,   79,   71,   98,   85,   93,   87,   112,   115,   110,   129,   199,   142,   170,   218,   162,   174,   205,   143,   164,   172,   175,   202,   216 },
+{   10,   25,   35,   69,   74,   92,   96,   135,   111,   156,   95,   105,   126,   134,   163,   215,   149,   148,   153,   228,   180,   230,   211,   213,   240,   224,   209,   222,   226,   173 },
+{   13,   37,   53,   108,   86,   117,   103,   100,   124,   122,   201,   192,   203,   219,   191,   169,   220,   214,   139,   204,   236,   179,   181,   200,   193,   221,   235,   187,   208,   217 },
+{   14,   36,   59,   116,   130,   178,   140,   132,   157,   152,   168,   160,   123,   197,   183,   234,   144,   223,   165,   167,   231,   196,   171,   233,   182,   225,   188,   237,   227,   210 },
+{   17,   58,   82,   101,   120,   186,   146,   125,   151,   194,   161,   184,   158,   189,   177,   195,   145,   166,   206,   154,   198,   212,   238,   229,   185,   176,   239,   155,   207,   232 }
+};
diff --git a/codec2-dev/src/codebook/c2wideband_map.txt b/codec2-dev/src/codebook/c2wideband_map.txt
new file mode 100644 (file)
index 0000000..946aa19
--- /dev/null
@@ -0,0 +1,8 @@
+ 1.00000000e+00 3.00000000e+00 7.00000000e+00 8.00000000e+00 1.10000000e+01 1.50000000e+01 1.90000000e+01 2.60000000e+01 2.10000000e+01 2.40000000e+01 2.00000000e+01 3.00000000e+01 3.80000000e+01 4.80000000e+01 2.90000000e+01 3.20000000e+01 4.20000000e+01 6.40000000e+01 6.30000000e+01 5.40000000e+01 5.60000000e+01 5.00000000e+01 7.20000000e+01 9.10000000e+01 7.80000000e+01 6.70000000e+01 5.70000000e+01 7.00000000e+01 7.60000000e+01 1.04000000e+02
+ 2.00000000e+00 5.00000000e+00 1.20000000e+01 2.20000000e+01 2.30000000e+01 3.30000000e+01 3.10000000e+01 4.30000000e+01 3.40000000e+01 3.90000000e+01 4.40000000e+01 4.10000000e+01 4.90000000e+01 6.00000000e+01 5.10000000e+01 7.70000000e+01 9.70000000e+01 9.00000000e+01 1.14000000e+02 8.10000000e+01 1.21000000e+02 1.13000000e+02 8.00000000e+01 9.90000000e+01 1.18000000e+02 1.09000000e+02 8.90000000e+01 1.02000000e+02 1.07000000e+02 1.19000000e+02
+ 4.00000000e+00 9.00000000e+00 1.60000000e+01 2.70000000e+01 4.00000000e+01 4.50000000e+01 4.70000000e+01 5.50000000e+01 4.60000000e+01 6.10000000e+01 6.50000000e+01 6.60000000e+01 7.50000000e+01 8.30000000e+01 6.80000000e+01 1.06000000e+02 9.40000000e+01 1.27000000e+02 1.36000000e+02 1.38000000e+02 1.37000000e+02 1.41000000e+02 8.80000000e+01 1.50000000e+02 1.59000000e+02 1.47000000e+02 1.33000000e+02 1.28000000e+02 1.90000000e+02 1.31000000e+02
+ 6.00000000e+00 1.80000000e+01 2.80000000e+01 5.20000000e+01 8.40000000e+01 6.20000000e+01 7.30000000e+01 7.90000000e+01 7.10000000e+01 9.80000000e+01 8.50000000e+01 9.30000000e+01 8.70000000e+01 1.12000000e+02 1.15000000e+02 1.10000000e+02 1.29000000e+02 1.99000000e+02 1.42000000e+02 1.70000000e+02 2.18000000e+02 1.62000000e+02 1.74000000e+02 2.05000000e+02 1.43000000e+02 1.64000000e+02 1.72000000e+02 1.75000000e+02 2.02000000e+02 2.16000000e+02
+ 1.00000000e+01 2.50000000e+01 3.50000000e+01 6.90000000e+01 7.40000000e+01 9.20000000e+01 9.60000000e+01 1.35000000e+02 1.11000000e+02 1.56000000e+02 9.50000000e+01 1.05000000e+02 1.26000000e+02 1.34000000e+02 1.63000000e+02 2.15000000e+02 1.49000000e+02 1.48000000e+02 1.53000000e+02 2.28000000e+02 1.80000000e+02 2.30000000e+02 2.11000000e+02 2.13000000e+02 2.40000000e+02 2.24000000e+02 2.09000000e+02 2.22000000e+02 2.26000000e+02 1.73000000e+02
+ 1.30000000e+01 3.70000000e+01 5.30000000e+01 1.08000000e+02 8.60000000e+01 1.17000000e+02 1.03000000e+02 1.00000000e+02 1.24000000e+02 1.22000000e+02 2.01000000e+02 1.92000000e+02 2.03000000e+02 2.19000000e+02 1.91000000e+02 1.69000000e+02 2.20000000e+02 2.14000000e+02 1.39000000e+02 2.04000000e+02 2.36000000e+02 1.79000000e+02 1.81000000e+02 2.00000000e+02 1.93000000e+02 2.21000000e+02 2.35000000e+02 1.87000000e+02 2.08000000e+02 2.17000000e+02
+ 1.40000000e+01 3.60000000e+01 5.90000000e+01 1.16000000e+02 1.30000000e+02 1.78000000e+02 1.40000000e+02 1.32000000e+02 1.57000000e+02 1.52000000e+02 1.68000000e+02 1.60000000e+02 1.23000000e+02 1.97000000e+02 1.83000000e+02 2.34000000e+02 1.44000000e+02 2.23000000e+02 1.65000000e+02 1.67000000e+02 2.31000000e+02 1.96000000e+02 1.71000000e+02 2.33000000e+02 1.82000000e+02 2.25000000e+02 1.88000000e+02 2.37000000e+02 2.27000000e+02 2.10000000e+02
+ 1.70000000e+01 5.80000000e+01 8.20000000e+01 1.01000000e+02 1.20000000e+02 1.86000000e+02 1.46000000e+02 1.25000000e+02 1.51000000e+02 1.94000000e+02 1.61000000e+02 1.84000000e+02 1.58000000e+02 1.89000000e+02 1.77000000e+02 1.95000000e+02 1.45000000e+02 1.66000000e+02 2.06000000e+02 1.54000000e+02 1.98000000e+02 2.12000000e+02 2.38000000e+02 2.29000000e+02 1.85000000e+02 1.76000000e+02 2.39000000e+02 1.55000000e+02 2.07000000e+02 2.32000000e+02
index de79710ca7a6df6d4617e010f553f52bf7f62c37..17c9a3291192bea6c592dd15f5be240da09d441c 100644 (file)
@@ -48,7 +48,7 @@
 #include "machdep.h"
 #include "bpf.h"
 #include "bpfb.h"
-
+#include "c2wideband.h"
 /*---------------------------------------------------------------------------*\
 
                              FUNCTION HEADERS
@@ -103,7 +103,7 @@ struct CODEC2 * codec2_create(int mode)
     struct CODEC2 *c2;
     int            i,l;
 
-    if (!((mode >= 0) && (mode <= CODEC2_MODE_700C))) {
+    if (!((mode >= 0) && (mode <= CODEC2_MODE_WB))) {
         return NULL;
     }  
 
@@ -264,6 +264,9 @@ int codec2_bits_per_frame(struct CODEC2 *c2) {
        return 28;
     if  (c2->mode == CODEC2_MODE_700C)
        return 28;
+    //TODO: verify this
+    if (c2->mode == CODEC2_MODE_WB)
+       return 64;
 
     return 0; /* shouldn't get here */
 }
@@ -298,7 +301,8 @@ int codec2_samples_per_frame(struct CODEC2 *c2) {
        return 320;
     if  (c2->mode == CODEC2_MODE_700C)
        return 320;
-
+    if  (c2->mode == CODEC2_MODE_WB)
+       return 160;
     return 0; /* shouldnt get here */
 }
 
@@ -327,6 +331,9 @@ void codec2_encode(struct CODEC2 *c2, unsigned char *bits, short speech[])
     if (c2->mode == CODEC2_MODE_700C)
        codec2_encode_700c(c2, bits, speech);
 #endif
+    if (c2->mode == CODEC2_MODE_WB)
+       codec2_encode_wb(c2, bits, speech);
+
 }
 
 void codec2_decode(struct CODEC2 *c2, short speech[], const unsigned char *bits)
@@ -359,6 +366,8 @@ void codec2_decode_ber(struct CODEC2 *c2, short speech[], const unsigned char *b
     if (c2->mode == CODEC2_MODE_700C)
        codec2_decode_700c(c2, speech, bits);
 #endif
+    if (c2->mode == CODEC2_MODE_WB)
+       codec2_decode_wb(c2, speech, bits);
 }
 
 
index 0e720f928fa10b67bee95795f377e2f43d292739..98ed32026df302837c92c8700225ca2a9e910481 100644 (file)
@@ -42,6 +42,7 @@
 #define CODEC2_MODE_700  6
 #define CODEC2_MODE_700B 7
 #define CODEC2_MODE_700C 8
+#define CODEC2_MODE_WB 9
 
 struct CODEC2;
 
diff --git a/codec2-dev/src/generate_wideband_map.c b/codec2-dev/src/generate_wideband_map.c
new file mode 100644 (file)
index 0000000..7687ec0
--- /dev/null
@@ -0,0 +1,156 @@
+/*---------------------------------------------------------------------------*\
+
+  FILE........: generate_wideband_map.c
+  AUTHOR......: Phil Ayres
+  DATE CREATED: 17 Jul 2017
+
+  Generate header file containing wideband DCT2 map, runs at compile time.
+  Adapted from generate_codebook.c
+
+\*---------------------------------------------------------------------------*/
+
+/*
+  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 <stdlib.h>
+#include <stdio.h>
+#include <ctype.h>
+#include <math.h>
+#include "c2wideband.h"
+
+static const int Nt = C2WB_NT;
+static const int K = C2WB_K;
+
+static const char usage[] =
+"Usage: %s filename array_name\n"
+"\tCreate C code for wideband DCT2 map.\n";
+
+static const char format[] =
+"The table format must be:\n"
+"\t8 rows by 30 (Nt x K) floating point numbers to fill the specified dimensions.\n";
+
+static const char header[] =
+"/* THIS IS A GENERATED FILE. Edit generate_wideband_map.c and its input */\n\n"
+"/*\n"
+" * This intermediary file and the files that used to create it are under \n"
+" * The LGPL. See the file COPYING.\n"
+" */\n\n"
+"#include \"defines.h\"\n\n";
+
+
+static void
+dump_array(float b[Nt][K])
+{
+  
+  printf("static const float c2wideband_map[%d][%d] = {\n", Nt, K);
+  for (int row = 0; row < Nt; row++ ) {
+      printf("{ ");
+      for (int col = 0; col < K; col++ ) {
+            printf("  %g", b[row][col]);
+            if ( col < K - 1 )
+                printf(", ");
+            else
+                printf(" }");
+      }
+      if ( row < Nt - 1 )
+          printf(",\n");
+      else
+          printf("\n");
+  }
+  printf("};\n");
+}
+
+
+
+float
+get_float(FILE * in, const char * name, char * * cursor, char * buffer,
+ int size)
+{
+  for ( ; ; ) {
+    char *     s = *cursor;
+    char       c;
+
+    while ( (c = *s) != '\0' && !isdigit(c) && c != '-' && c != '.' )
+      s++;
+
+    /* Comments start with "#" and continue to the end of the line. */
+    if ( c != '\0' && c != '#' ) {
+      char *   end = 0;
+      float    f = 0;
+
+      f = strtod(s, &end);
+
+      if ( end != s )
+        *cursor = end;
+      return f;
+    }
+
+    if ( fgets(buffer, size, in) == NULL ) {
+      fprintf(stderr, "%s: Format error. %s\n", name, format);
+      exit(1);
+    }
+    *cursor = buffer;
+  }
+}
+
+static void
+load(FILE * file, const char * name, float b[Nt][K])
+{
+  char                 line[1024];
+  char *               cursor = line;
+
+  *cursor = '\0';
+
+  for (int row = 0; row < Nt; row++ ) {
+      
+      for (int col = 0; col < K; col++ ) {
+        
+            b[row][col]  = get_float(file, name, &cursor, line, sizeof(line));
+      }
+  }
+  
+}
+
+
+int main(int argc, char * * argv)
+{
+  float b[Nt][K];
+  int  i;
+
+  if ( argc < 2 ) {
+    fprintf(stderr, usage, argv[0]);
+    fprintf(stderr, format);
+    exit(1);
+  }
+
+  i=0;
+  FILE *       in = fopen(argv[i + 2], "r");
+
+  if ( in == NULL ) {
+      perror(argv[i + 2]);
+      exit(1);
+  }
+
+  load(in, argv[i + 2], b);
+
+  fclose(in);
+  printf(header);
+  
+  printf("  /* %s */\n", argv[i + 2]);
+  dump_array(b);
+
+  return 0;
+}
index a4b06b9eff2be881d7a62ee13387d04305ce7fbf..2cae7a66c862b37a1935d4d0dc40d286699325ca 100644 (file)
@@ -95,3 +95,6 @@ target_link_libraries(tsrc samplerate)
 
 add_executable(tlininterp tlininterp.c)
 add_executable(tdec tdec.c)
+
+add_executable(tc2wideband tc2wideband.c ../src/c2wideband.c)
+target_link_libraries(tc2wideband codec2)
diff --git a/codec2-dev/unittest/tc2wideband.c b/codec2-dev/unittest/tc2wideband.c
new file mode 100644 (file)
index 0000000..f6c2ab6
--- /dev/null
@@ -0,0 +1,239 @@
+#include <assert.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+#include <unistd.h>
+
+#include "c2wideband.h"
+#include "codec2.h"
+#include "defines.h"
+#include "codec2_fft.h"
+#include "sine.h"
+#include "nlp.h"
+#include "dump.h"
+#include "octave.h"
+#include "newamp1.h"
+#include "quantise.h"
+
+float mean(float data[], int n);
+int unit_test();
+float std(float data[], int rows);
+void diff_de(int rows, int cols, float D[rows][cols], float E[rows][cols], float diff_de[rows][cols]);
+void array_col_to_row(int rows, int cols, float data[rows][cols], int col, float res[]);
+void std_on_cols(int rows, int cols, float data[rows][cols], float res[]);
+float mean_std_diff_de(int rows, int cols, float D[rows][cols], float E[rows][cols]);
+void test_wideband_enc_dec();
+void setup_map(WIDEBAND_MAP * wb_map, int Nt, int K);
+
+char *fn;
+
+void test(char * tfn) {
+  fn = tfn;
+  printf("test function: %s\n", fn);
+}
+
+void test_failed() {
+  printf("Failed to calculate %s.\n", fn);
+  exit(1);
+}
+
+void test_failed_s(char * expected, char * res) {
+
+  printf("Failed to calculate %s.\n", fn);
+
+  printf("expected: %s\ngot: %s\n", expected, res);
+  exit(1);
+}
+
+void test_failed_f(float expected, float res) {
+
+  printf("Failed to calculate %s.\n", fn);
+  printf("expected: %f\ngot: %f\n", expected, res);
+  exit(1);
+}
+
+int main(int argc, char *argv[]) {
+
+  printf("Testing file c2wideband\n");
+
+  test("mean");
+  int n = 5;
+  float in[] = {1.0, 2.0, 4.0, 6.0, 10.0};
+  float res_f;
+  float expect_f;
+
+  res_f = mean(in, n);
+  expect_f = 23.0 / 5;
+  if (res_f != expect_f) {
+    test_failed_f(expect_f, res_f);
+  }
+
+  test("std");
+  res_f = std(in, n);
+  expect_f = 3.2000000000000002;
+  if (res_f != expect_f) {
+    test_failed_f(expect_f, res_f);
+  }
+
+
+  test("diff_de");
+  int rows = 5;
+  int cols = 3;
+  float D[rows][cols];
+  float E[rows][cols];
+  float res_diff_de[rows][cols];
+  float expect_diff_de[rows][cols];
+
+  for (int r = 0; r < rows; r++) {
+    for (int c = 0; c < cols; c++) {
+      float d = rand();
+      float e = rand();
+      D[r][c] = d;
+      E[r][c] = e;
+      expect_diff_de[r][c] = d - e;
+    }
+  }
+
+  diff_de(rows, cols, D, E, res_diff_de);
+
+  if (memcmp(res_diff_de, expect_diff_de, rows * cols * sizeof (float)) != 0) {
+    test_failed();
+  }
+
+
+  test("array_col_to_row");
+  float data[rows][cols];
+  float res_data[rows];
+  float expect_data[cols][rows];
+  for (int r = 0; r < rows; r++) {
+    for (int c = 0; c < cols; c++) {
+      float d = rand();
+      data[r][c] = d;
+      expect_data[c][r] = d;
+    }
+  }
+  for (int c = 0; c < cols; c++) {
+
+    array_col_to_row(rows, cols, data, c, res_data);
+    if (memcmp(res_data, expect_data[c], cols * sizeof (float)) != 0) {
+      test_failed();
+    }
+  }
+
+
+  test("std_on_cols");
+
+  cols = 4;
+
+  float data_std[5][4] = {
+    {1.0, 1.0, 4.0, 6.0},
+    {1.0, 2.0, 5.5, 6.0},
+    {1.0, 4.0, 6.0, 6.0},
+    {1.0, 6.0, 3.3, 6.0},
+    {1.0, 10.0, -0.2, 6.0}
+  };
+
+  float res_data_std[cols];
+  float expect_data_std[] = {0.0, 3.2000000000000002, 2.1903424389807178, 0.0};
+
+  std_on_cols(rows, cols, data_std, res_data_std);
+  if (memcmp(res_data_std, expect_data_std, cols * sizeof (float)) != 0) {
+    test_failed(fn);
+  }
+
+  test("mean_std_diff_de");
+
+  cols = 4;
+
+  float data_std_d[5][4] = {
+    {2.0, 1.0, 5.0, 6.0},
+    {1.0, 2.0, 10.5, 6.0},
+    {4.0, 4.0, 12.0, 6.0},
+    {1.0, 6.0, 4.6, 6.0},
+    {6.5, 10.0, 0.0, 6.0}
+  };
+
+  float data_std_e[5][4] = {
+    {1.0, 0.0, 1.0, 0.0},
+    {0.0, 0.0, 5.0, 0.0},
+    {3.0, 0.0, 6.0, 0.0},
+    {0.0, 0.0, 1.3, 0.0},
+    {5.5, 0.0, 0.2, 0.0}
+  };
+
+  float expected_msd = (3.2000000000000002 + 2.1903424389807178) / cols;
+
+  float res_msd;
+  res_msd = mean_std_diff_de(rows, cols, data_std_d, data_std_e);
+  if (abs(res_msd - expected_msd) > 0.00001) {
+    test_failed_f(res_msd, expected_msd);
+  }
+
+  test_wideband_enc_dec();
+
+
+  return 1;
+}
+
+void test_wideband_enc_dec() {
+
+  test("wideband_enc_dec");
+
+  short *buf;
+  unsigned char *bits;
+  //int            nsam, nbit, i, r;
+  int nbit;
+  int nsam;
+  int n_block_frames = C2WB_NT * C2WB_DEC;
+  int test_me = 123;
+  int K = C2WB_K;
+  int Nt = C2WB_NT;
+  int Fs = C2WB_FS;
+  int dec = C2WB_DEC;
+  int Nblocks = n_block_frames;
+  int rate_K_surface_size = Nblocks * Nt * dec;
+
+  struct CODEC2 * codec2 = codec2_create(CODEC2_MODE_WB);
+
+  C2CONST c2const = c2const_create(Fs);
+  nsam = C2WB_SPERF; //codec2_samples_per_frame(codec2);
+  nbit = C2WB_BPERF; //codec2_bits_per_frame(codec2);
+
+  buf = (short*) malloc(nsam * sizeof (short));
+  bits = (unsigned char*) malloc(nbit * sizeof (char));
+
+  WIDEBAND_MAP wb_map;
+  setup_map(&wb_map, Nt, K);
+
+  MODEL model_block[n_block_frames];
+  MODEL model_block_[n_block_frames];
+  float dct2_sd[Nblocks];
+  int qn;
+  float rate_K_surface_block[rate_K_surface_size][K]; // rate K vecs for each frame, form a surface that makes pretty graphs   
+  float rate_K_surface_block_[rate_K_surface_size][K];
+
+
+  for (int n = 0; n < n_block_frames; n++) {
+    model_block[n].L = 10;
+    model_block[n].Wo = (float) rand() / RAND_MAX;
+    for (int i = 0; i < MAX_AMP + 1; i++) {
+      model_block[n].phi[i] = (float) rand() / RAND_MAX;
+      model_block[n].A[i] = (float) rand() / RAND_MAX;
+    }
+    model_block[n].voiced = round((float) rand() / RAND_MAX);
+  }
+
+  int n = 0;
+
+  printf("setup complete. now calling the function\n");
+
+  wideband_enc_dec(&c2const, n_block_frames, model_block, &wb_map,
+          model_block_, dct2_sd, &qn, &rate_K_surface_block[n], &rate_K_surface_block_[n]);
+
+  codec2_destroy(codec2);
+  free(bits);
+  free(buf);
+  printf("made it to the end\n");
+}
+