32-bit fixes
[opus.git] / libcelt / bands.c
index d348b30..70b126b 100644 (file)
@@ -225,17 +225,17 @@ static void stereo_split(celt_norm *X, celt_norm *Y, int N)
    for (j=0;j<N;j++)
    {
       celt_norm r, l;
-      l = MULT16_16_Q15(QCONST16(.70711f,15), X[j]);
-      r = MULT16_16_Q15(QCONST16(.70711f,15), Y[j]);
+      l = MULT16_16_Q15(QCONST16(.70710678f,15), X[j]);
+      r = MULT16_16_Q15(QCONST16(.70710678f,15), Y[j]);
       X[j] = l+r;
       Y[j] = r-l;
    }
 }
 
-static void stereo_merge(celt_norm *X, celt_norm *Y, celt_word16 mid, celt_word16 side, int N)
+static void stereo_merge(celt_norm *X, celt_norm *Y, celt_word16 mid, int N)
 {
    int j;
-   celt_word32 xp=0;
+   celt_word32 xp=0, side=0;
    celt_word32 El, Er;
 #ifdef FIXED_POINT
    int kl, kr;
@@ -244,12 +244,15 @@ static void stereo_merge(celt_norm *X, celt_norm *Y, celt_word16 mid, celt_word1
 
    /* Compute the norm of X+Y and X-Y as |X|^2 + |Y|^2 +/- sum(xy) */
    for (j=0;j<N;j++)
+   {
       xp = MAC16_16(xp, X[j], Y[j]);
+      side = MAC16_16(side, Y[j], Y[j]);
+   }
    /* mid and side are in Q15, not Q14 like X and Y */
    mid = SHR32(mid, 1);
-   side = SHR32(side, 1);
-   El = MULT16_16(mid, mid) + MULT16_16(side, side) - 2*xp;
-   Er = MULT16_16(mid, mid) + MULT16_16(side, side) + 2*xp;
+   //side = SHR32(side, 1);
+   El = MULT16_16(mid, mid) + side - 2*xp;
+   Er = MULT16_16(mid, mid) + side + 2*xp;
    if (Er < EPSILON)
       Er = EPSILON;
    if (El < EPSILON)
@@ -387,7 +390,18 @@ void measure_norm_mse(const CELTMode *m, float *X, float *X0, float *bandE, floa
 
 #endif
 
-static void interleave_vector(celt_norm *X, int N0, int stride)
+/* Indexing table for converting from natural Hadamard to ordery Hadamard
+   This is essentially a bit-reversed Gray, on top of which we've added
+   an inversion of the order because we want the DC at the end rather than
+   the beginning. The lines are for N=2, 4, 8, 16 */
+static const int ordery_table[] = {
+       1,  0,
+       3,  0,  2,  1,
+       7,  0,  4,  3,  6,  1,  5,  2,
+      15,  0,  8,  7, 12,  3, 11,  4, 14,  1,  9,  6, 13,  2, 10,  5,
+};
+
+static void deinterleave_hadamard(celt_norm *X, int N0, int stride, int hadamard)
 {
    int i,j;
    VARDECL(celt_norm, tmp);
@@ -395,15 +409,25 @@ static void interleave_vector(celt_norm *X, int N0, int stride)
    SAVE_STACK;
    N = N0*stride;
    ALLOC(tmp, N, celt_norm);
-   for (i=0;i<stride;i++)
-      for (j=0;j<N0;j++)
-         tmp[j*stride+i] = X[i*N0+j];
+   if (hadamard)
+   {
+      const int *ordery = ordery_table+stride-2;
+      for (i=0;i<stride;i++)
+      {
+         for (j=0;j<N0;j++)
+            tmp[ordery[i]*N0+j] = X[j*stride+i];
+      }
+   } else {
+      for (i=0;i<stride;i++)
+         for (j=0;j<N0;j++)
+            tmp[i*N0+j] = X[j*stride+i];
+   }
    for (j=0;j<N;j++)
       X[j] = tmp[j];
    RESTORE_STACK;
 }
 
-static void deinterleave_vector(celt_norm *X, int N0, int stride)
+static void interleave_hadamard(celt_norm *X, int N0, int stride, int hadamard)
 {
    int i,j;
    VARDECL(celt_norm, tmp);
@@ -411,9 +435,17 @@ static void deinterleave_vector(celt_norm *X, int N0, int stride)
    SAVE_STACK;
    N = N0*stride;
    ALLOC(tmp, N, celt_norm);
-   for (i=0;i<stride;i++)
-      for (j=0;j<N0;j++)
-         tmp[i*N0+j] = X[j*stride+i];
+   if (hadamard)
+   {
+      const int *ordery = ordery_table+stride-2;
+      for (i=0;i<stride;i++)
+         for (j=0;j<N0;j++)
+            tmp[j*stride+i] = X[ordery[i]*N0+j];
+   } else {
+      for (i=0;i<stride;i++)
+         for (j=0;j<N0;j++)
+            tmp[j*stride+i] = X[i*N0+j];
+   }
    for (j=0;j<N;j++)
       X[j] = tmp[j];
    RESTORE_STACK;
@@ -427,8 +459,8 @@ void haar1(celt_norm *X, int N0, int stride)
       for (j=0;j<N0;j++)
       {
          celt_norm tmp1, tmp2;
-         tmp1 = MULT16_16_Q15(QCONST16(.7070678f,15), X[stride*2*j+i]);
-         tmp2 = MULT16_16_Q15(QCONST16(.7070678f,15), X[stride*(2*j+1)+i]);
+         tmp1 = MULT16_16_Q15(QCONST16(.70710678f,15), X[stride*2*j+i]);
+         tmp2 = MULT16_16_Q15(QCONST16(.70710678f,15), X[stride*(2*j+1)+i]);
          X[stride*2*j+i] = tmp1 + tmp2;
          X[stride*(2*j+1)+i] = tmp1 - tmp2;
       }
@@ -487,6 +519,9 @@ static void quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_
    int recombine=0;
    int inv = 0;
    celt_word16 mid=0, side=0;
+   int longBlocks;
+
+   longBlocks = B0==1;
 
    N_B /= B;
    N_B0 = N_B;
@@ -565,9 +600,9 @@ static void quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_
       if (B0>1)
       {
          if (encode)
-            deinterleave_vector(X, N_B, B0);
+            deinterleave_hadamard(X, N_B, B0, longBlocks);
          if (lowband)
-            deinterleave_vector(lowband, N_B, B0);
+            deinterleave_hadamard(lowband, N_B, B0, longBlocks);
       }
    }
 
@@ -613,7 +648,7 @@ static void quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_
 
          /* Entropy coding of the angle. We use a uniform pdf for the
             first stereo split but a triangular one for the rest. */
-         if (stereo || B>1)
+         if (stereo || B0>1)
          {
             if (encode)
                ec_enc_uint((ec_enc*)ec, itheta, qn+1);
@@ -765,9 +800,15 @@ static void quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_
          int next_level=0;
 
          /* Give more bits to low-energy MDCTs than they would otherwise deserve */
-         if (B>1 && !stereo && itheta > 8192)
-            delta -= delta>>(1+level);
-
+         if (B0>1 && !stereo)
+         {
+            if (itheta > 8192)
+               /* Rough approximation for pre-echo masking */
+               delta -= delta>>(4-LM);
+            else
+               /* Corresponds to a forward-masking slope of 1.5 dB per 10 ms */
+               delta = IMIN(0, delta + (N<<BITRES>>(5-LM)));
+         }
          mbits = (b-qalloc-delta)/2;
          if (mbits > b-qalloc)
             mbits = b-qalloc;
@@ -834,7 +875,7 @@ static void quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_
                   for (j=0;j<N;j++)
                   {
                      *seed = lcg_rand(*seed);
-                     X[j] = (int)(*seed)>>20;
+                     X[j] = (celt_int32)(*seed)>>20;
                   }
                } else {
                   /* Folded spectrum */
@@ -853,7 +894,7 @@ static void quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_
       if (stereo)
       {
          if (N!=2)
-            stereo_merge(X, Y, mid, side, N);
+            stereo_merge(X, Y, mid, N);
          if (inv)
          {
             int j;
@@ -866,7 +907,7 @@ static void quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_
 
          /* Undo the sample reorganization going from time order to frequency order */
          if (B0>1)
-            interleave_vector(X, N_B, B0);
+            interleave_hadamard(X, N_B, B0, longBlocks);
 
          /* Undo time-freq changes that we did earlier */
          N_B = N_B0;
@@ -903,7 +944,8 @@ void quant_all_bands(int encode, const CELTMode *m, int start, int end,
       int shortBlocks, int spread, int dual_stereo, int intensity, int *tf_res, int resynth,
       int total_bits, void *ec, int LM, int codedBands)
 {
-   int i, balance;
+   int i;
+   celt_int32 balance;
    celt_int32 remaining_bits;
    const celt_int16 * restrict eBands = m->eBands;
    celt_norm * restrict norm, * restrict norm2;
@@ -957,10 +999,10 @@ void quant_all_bands(int encode, const CELTMode *m, int start, int end,
    lowband_offset = -1;
    for (i=start;i<end;i++)
    {
-      int tell;
+      celt_int32 tell;
       int b;
       int N;
-      int curr_balance;
+      celt_int32 curr_balance;
       int effective_lowband=-1;
       celt_norm * restrict X, * restrict Y;
       int tf_change=0;
@@ -979,22 +1021,19 @@ void quant_all_bands(int encode, const CELTMode *m, int start, int end,
       /* Compute how many bits we want to allocate to this band */
       if (i != start)
          balance -= tell;
-      remaining_bits = (total_bits<<BITRES)-tell-1;
+      remaining_bits = ((celt_int32)total_bits<<BITRES)-tell-1;
       if (i <= codedBands-1)
       {
          curr_balance = (codedBands-i);
          if (curr_balance > 3)
             curr_balance = 3;
          curr_balance = balance / curr_balance;
-         b = IMIN(remaining_bits+1,pulses[i]+curr_balance);
+         b = IMIN(16384, IMIN(remaining_bits+1,pulses[i]+curr_balance));
          if (b<0)
             b = 0;
       } else {
          b = 0;
       }
-      /* Prevents ridiculous bit depths */
-      if (b > C*16*N<<BITRES)
-         b = C*16*N<<BITRES;
 
       if (M*eBands[i]-N >= M*eBands[start] && (update_lowband || lowband_offset==-1))
             lowband_offset = M*eBands[i];