From ff74e396e4616be054e596101f60132aab828d04 Mon Sep 17 00:00:00 2001
From: Jean-Marc Valin <Jean-Marc.Valin@csiro.au>
Date: Wed, 27 Feb 2008 15:35:43 +1100
Subject: [PATCH] fixed-point: converted compute_pitch_gain() and removed the
 energy-based weighting that didn't seem to help anyway.

---
 libcelt/arch.h  |  5 ++++-
 libcelt/bands.c | 37 +++++++++++++++----------------------
 2 files changed, 19 insertions(+), 23 deletions(-)

diff --git a/libcelt/arch.h b/libcelt/arch.h
index 7bc1e790..73ebef5e 100644
--- a/libcelt/arch.h
+++ b/libcelt/arch.h
@@ -67,12 +67,15 @@ typedef float celt_mask_t;
 #define NORM_SCALING_1 (1.f/16384.f)
 #define ENER_SCALING 16384.f
 #define ENER_SCALING_1 (1.f/16384.f)
+
 #define PGAIN_SCALING 32768.f
 #define PGAIN_SCALING_1 (1.f/32768.f)
+#define PGAIN_SHIFT 15
 
 #define DB_SCALING 256.f
 #define DB_SCALING_1 (1.f/256.f)
 
+#define EPSILON 1
 #define VERY_SMALL 0
 #define VERY_LARGE32 ((celt_word32_t)2147483647)
 #define VERY_LARGE16 ((celt_word16_t)32767)
@@ -121,7 +124,7 @@ typedef float celt_mask_t;
 #define DB_SCALING 1.f
 #define DB_SCALING_1 1.f
 
-
+#define EPSILON 1e-15
 #define VERY_SMALL 1e-15f
 #define VERY_LARGE32 1e15f
 #define VERY_LARGE16 1e15f
diff --git a/libcelt/bands.c b/libcelt/bands.c
index 200af28e..7573c3a2 100644
--- a/libcelt/bands.c
+++ b/libcelt/bands.c
@@ -163,36 +163,29 @@ void compute_pitch_gain(const CELTMode *m, celt_norm_t *X, celt_norm_t *P, celt_
    int i, B;
    const int *eBands = m->eBands;
    const int *pBands = m->pBands;
-   VARDECL(float *w);
    B = m->nbMdctBlocks*m->nbChannels;
-   ALLOC(w, B*eBands[m->nbEBands], float);
-   for (i=0;i<m->nbEBands;i++)
-   {
-      int j;
-      for (j=B*eBands[i];j<B*eBands[i+1];j++)
-         w[j] = bank[i]*ENER_SCALING_1;
-   }
-
    
    for (i=0;i<m->nbPBands;i++)
    {
-      float Sxy=0;
-      float Sxx = 0;
+      celt_word32_t Sxy=0, Sxx=0;
       int j;
-      float gain;
+      /* We know we're not going to overflow because Sxx can't be more than 1 (Q28) */
       for (j=B*pBands[i];j<B*pBands[i+1];j++)
       {
-         Sxy += 1.f*X[j]*P[j]*w[j];
-         Sxx += 1.f*X[j]*X[j]*w[j];
+         Sxy = MAC16_16(Sxy, X[j], P[j]);
+         Sxx = MAC16_16(Sxx, X[j], X[j]);
       }
-      gain = Sxy/(1e-10*NORM_SCALING*NORM_SCALING+Sxx);
-      if (gain > 1.f)
-         gain = 1.f;
-      if (gain < 0.0f)
-         gain = 0.0f;
-      /* We need to be a bit conservative, otherwise residual doesn't quantise well */
-      gain *= .9f;
-      gains[i] = PGAIN_SCALING*gain;
+      /* No negative gain allowed */
+      if (Sxy < 0)
+         Sxy = 0;
+      /* Not sure how that would happen, just making sure */
+      if (Sxy > Sxx)
+         Sxy = Sxx;
+      /* We need to be a bit conservative (multiply gain by 0.9), otherwise the
+         residual doesn't quantise well */
+      Sxy = MULT16_32_Q15(QCONST16(.9f, 15), Sxy);
+      /* gain = Sxy/Sxx */
+      gains[i] = DIV32_16(Sxy,ADD32(SHR32(Sxx, PGAIN_SHIFT),EPSILON));
       /*printf ("%f ", 1-sqrt(1-gain*gain));*/
    }
    /*if(rand()%10==0)