feature calculation for BWE added

TheGrumpySnail01 · TheGrumpySnail01 · commit 702e5b26e907 · 2025-03-30T16:34:39.000-07:00
diff --git a/dnn/bwe_demo.c b/dnn/bwe_demo.c
@@ -36,6 +36,7 @@
 #include "lpcnet.h"
 #include "os_support.h"
 #include "cpu_support.h"
+#include "osce_features.h"
 
 
 void usage(void) {
@@ -60,7 +61,40 @@ int main(int argc, char **argv) {
         exit(1);
     }
 
-    printf("BWE coming soon!\n");
+    printf("Feature calculation with signal (100 * (n % 90)) - 8900\n");
+    int n = 0, i;
+    opus_int16 frame[160];
+    int frame_counter = 0;
+    float features[32 + 2 * 41];
+
+    for (frame_counter = 0; frame_counter < 10; frame_counter ++)
+    {
+        for (i = 0; i < 160; i ++ )
+        {
+            frame[i] = 100 * n++ - 8900;
+            n = n % 90;
+        }
+
+        osce_bwe_calculate_features(features, frame, 160);
+
+        printf("frame[%d]\n", frame_counter);
+        printf("lmspec: ");
+        for (i = 0; i < 32; i ++)
+        {
+            printf(" %f ", features[i]);
+        }
+        printf("\nphasediff: ");
+        for (;i < 32 + 2 * 41; i ++)
+        {
+            printf(" %f ", features[i]);
+        }
+        printf("\n\n");
+
+        fwrite(frame, sizeof(frame[0]), 160, fout);
+
+    }
+
+
 
     fclose(fin);
     fclose(fout);
diff --git a/dnn/osce_features.c b/dnn/osce_features.c
@@ -68,6 +68,13 @@ static const int center_bins_noisy[18] = {
     136,    160
 };
 
+static const int center_bins_bwe[32] = {
+      0,   5,  10,  15,  20,  25,  30,  35,
+     40,  45,  50,  55,  60,  65,  70,  75,
+     80,  85,  90,  95, 100, 105, 110, 115,
+    120, 125, 130, 135, 140, 145, 150, 160
+};
+
 static const float band_weights_clean[64] = {
      0.666666666667f,     0.400000000000f,     0.333333333333f,     0.400000000000f,
      0.500000000000f,     0.400000000000f,     0.333333333333f,     0.400000000000f,
@@ -95,6 +102,17 @@ static const float band_weights_noisy[18] = {
      0.041666666667f,     0.080000000000f
 };
 
+static const float band_weights_bwe[32] = {
+    0.333333333, 0.200000000, 0.200000000, 0.200000000,
+    0.200000000, 0.200000000, 0.200000000, 0.200000000,
+    0.200000000, 0.200000000, 0.200000000, 0.200000000,
+    0.200000000, 0.200000000, 0.200000000, 0.200000000,
+    0.200000000, 0.200000000, 0.200000000, 0.200000000,
+    0.200000000, 0.200000000, 0.200000000, 0.200000000,
+    0.200000000, 0.200000000, 0.200000000, 0.200000000,
+    0.200000000, 0.200000000, 0.133333333, 0.181818182
+};
+
 static float osce_window[OSCE_SPEC_WINDOW_SIZE] = {
      0.004908718808f,     0.014725683311f,     0.024541228523f,     0.034354408400f,     0.044164277127f,
      0.053969889210f,     0.063770299562f,     0.073564563600f,     0.083351737332f,     0.093130877450f,
@@ -440,6 +458,98 @@ void osce_calculate_features(
 }
 
 
+#define OSCE_BWE_MAX_INSTAFREQ_BIN 40
+#define OSCE_BWE_HALF_WINDOW_SIZE 160
+#define OSCE_BWE_WINDOW_SIZE (2 * (OSCE_BWE_HALF_WINDOW_SIZE))
+#define OSCE_BWE_NUM_BANDS 32
+
+void osce_bwe_calculate_features(
+    float                       *features,                      /* O    input features                              */
+    const opus_int16            xq[],                           /* I    Decoded speech                              */
+    int                         num_samples                     /* I    number of input samples                     */
+ )
+ {
+    static int init;
+    static float signal_history[OSCE_BWE_HALF_WINDOW_SIZE];
+    static float last_spec[2 * OSCE_BWE_MAX_INSTAFREQ_BIN + 2] = {0};
+
+    int n, k;
+    kiss_fft_cpx fft_buffer[OSCE_BWE_WINDOW_SIZE];
+    float spec[2 * OSCE_BWE_MAX_INSTAFREQ_BIN + 2];
+    float buffer[OSCE_BWE_WINDOW_SIZE];
+    float mag_spec[OSCE_SPEC_NUM_FREQS];
+    float *lmspec, *instafreq;
+
+    /* OSCE_BWE_WINDOW_SIZE == 320 is a hard requirement */
+    celt_assert(num_samples == OSCE_BWE_HALF_WINDOW_SIZE && OSCE_BWE_WINDOW_SIZE == 320);
+
+    if (init == 0)
+    {
+        /* ToDo: fix python and retrain */
+        for (k = 0; k <= OSCE_BWE_MAX_INSTAFREQ_BIN; k ++)
+        {
+            last_spec[2*k] = 1e-9;
+        }
+        init = 1;
+    }
+
+    lmspec = features;
+    instafreq = features + OSCE_BWE_NUM_BANDS;
+
+    OPUS_COPY(buffer, signal_history, OSCE_BWE_HALF_WINDOW_SIZE);
+    for (n = 0; n < num_samples; n++)
+    {
+        buffer[n + OSCE_BWE_HALF_WINDOW_SIZE] = (float) xq[n] / (1U<<15);
+    }
+
+    /* update signal history buffer */
+    OPUS_COPY(signal_history, buffer + OSCE_BWE_HALF_WINDOW_SIZE, OSCE_BWE_HALF_WINDOW_SIZE);
+
+    /* apply window */
+    for (n = 0; n < OSCE_BWE_WINDOW_SIZE; n ++)
+    {
+        buffer[n] *= osce_window[n];
+    }
+
+    /* DFT */
+    forward_transform(fft_buffer, buffer);
+
+    /* instafreq */
+    for (k = 0; k <= OSCE_BWE_MAX_INSTAFREQ_BIN; k++)
+    {
+
+        float aux_r, aux_i, aux_abs;
+        float re1, re2, im1, im2;
+        spec[2*k]   = OSCE_BWE_WINDOW_SIZE * fft_buffer[k].r + 1e-9; /* ToDo: remove 1e-9 from python code*/
+        spec[2*k+1] = OSCE_BWE_WINDOW_SIZE * fft_buffer[k].i;
+        re1 = spec[2*k];
+        im1 = spec[2*k+1];
+        re2 = last_spec[2*k];
+        im2 = last_spec[2*k+1];
+        aux_r = re1 * re2 + im1 * im2;
+        aux_i = im1 * re2 - re1 * im2;
+        aux_abs = sqrt(aux_r * aux_r + aux_i * aux_i);
+        instafreq[k] = aux_r / (aux_abs + 1e-9);
+        instafreq[k + OSCE_BWE_MAX_INSTAFREQ_BIN + 1] = aux_i / (aux_abs + 1e-9);
+    }
+
+    /* erb-scale magnitude spectrogram */
+    for (k = 0; k < OSCE_SPEC_NUM_FREQS; k ++)
+    {
+        mag_spec[k] = OSCE_BWE_WINDOW_SIZE * sqrt(fft_buffer[k].r * fft_buffer[k].r + fft_buffer[k].i * fft_buffer[k].i);
+    }
+
+    apply_filterbank(lmspec, mag_spec, center_bins_bwe, band_weights_bwe, OSCE_BWE_NUM_BANDS);
+
+    for (k = 0; k < OSCE_BWE_NUM_BANDS; k++)
+    {
+        lmspec[k] = log(lmspec[k] + 1e-9);
+    }
+
+    /* update instafreq buffer */
+    OPUS_COPY(last_spec, spec, 2 * OSCE_BWE_MAX_INSTAFREQ_BIN + 2);
+}
+
 void osce_cross_fade_10ms(float *x_enhanced, float *x_in, int length)
 {
     int i;
diff --git a/dnn/osce_features.h b/dnn/osce_features.h
@@ -44,6 +44,11 @@ void osce_calculate_features(
     opus_int32                  num_bits                        /* I    Size of SILK payload in bits                */
 );
 
+void osce_bwe_calculate_features(
+   float                       *features,                      /* O    input features                              */
+   const opus_int16            xq[],                           /* I    Decoded speech                              */
+   int                         num_samples                     /* I    number of input samples                     */
+);
 
 void osce_cross_fade_10ms(float *x_enhanced, float *x_in, int length);
 
diff --git a/dnn/torch/osce/scripts/bwe_extract_filterbank.py b/dnn/torch/osce/scripts/bwe_extract_filterbank.py
@@ -0,0 +1,29 @@
+import argparse
+import sys
+sys.path.append('./')
+
+import torch
+from utils.spec import create_filter_bank
+import numpy as np
+
+
+parser = argparse.ArgumentParser()
+parser.add_argument('checkpoint', type=str)
+
+
+
+if __name__ == "__main__":
+    args = parser.parse_args()
+
+    c = torch.load(args.checkpoint, map_location='cpu')
+
+    num_bands = c['setup']['data']['spec_num_bands']
+    fb, center_bins = create_filter_bank(num_bands, n_fft=320, fs=16000, scale='erb', round_center_bins=True, normalize=False, return_center_bins=True)
+    weights = 1/fb.sum(axis=-1)
+
+    print(f"center_bins:")
+
+    print("".join([f"{int(cb):4d}," for cb in center_bins]))
+
+    print(f"band_weights:")
+    print("".join([f" {w:1.9f}," for w in weights]))
diff --git a/dnn/torch/osce/scripts/print_features.py b/dnn/torch/osce/scripts/print_features.py
@@ -0,0 +1,22 @@
+import argparse
+import sys
+sys.path.append('./')
+
+import torch
+from utils.bwe_features import load_inference_data
+import numpy as np
+
+parser = argparse.ArgumentParser()
+parser.add_argument('testsignal', type=str)
+
+if __name__ == "__main__":
+    args = parser.parse_args()
+
+    _, features = load_inference_data(args.testsignal)
+
+    N = features.shape[0]
+
+    for n in range(N):
+        print(f"frame[{n}]")
+        print(f"lmspec: {features[n, :32]}")
+        print(f"freqs: {features[n,32:]}")
diff --git a/dnn/torch/osce/utils/spec.py b/dnn/torch/osce/utils/spec.py
@@ -65,7 +65,7 @@ def gen_filterbank(N, Fs=16000, keep_size=False):
     RE = RE/norm[:, np.newaxis]
     return torch.from_numpy(RE)
 
-def create_filter_bank(num_bands, n_fft=320, fs=16000, scale='bark', round_center_bins=False, return_upper=False, normalize=False):
+def create_filter_bank(num_bands, n_fft=320, fs=16000, scale='bark', round_center_bins=False, return_upper=False, normalize=False, return_center_bins=False):
 
     f0 = 0
     num_bins = n_fft // 2 + 1
@@ -111,6 +111,9 @@ def create_filter_bank(num_bands, n_fft=320, fs=16000, scale='bark', round_cente
     if normalize:
         filter_bank = filter_bank / np.sum(filter_bank, axis=1).reshape(-1, 1)
 
+    if return_center_bins:
+        return filter_bank, center_bins
+
     return filter_bank
 
 
@@ -232,7 +235,7 @@ def instafreq(x, frame_size, max_bin, window=None):
     X = np.fft.fft(x_unfold, n=frame_size, axis=-1)
 
     # instantaneus frequency
-    X_trunc = X[..., :max_bin + 1] + 1e-9
+    X_trunc = X[..., :max_bin + 1]  + 1e-9
     Y = X_trunc[1:] * np.conj(X_trunc[:-1])
     Y = Y / (np.abs(Y) + 1e-9)
 
