libcelt/quant_bands.c

   1 /* (C) 2007 Jean-Marc Valin, CSIRO
   2 */
   3 /*
   4    Redistribution and use in source and binary forms, with or without
   5    modification, are permitted provided that the following conditions
   6    are met:
   7
   8    - Redistributions of source code must retain the above copyright
   9    notice, this list of conditions and the following disclaimer.
  10
  11    - Redistributions in binary form must reproduce the above copyright
  12    notice, this list of conditions and the following disclaimer in the
  13    documentation and/or other materials provided with the distribution.
  14
  15    - Neither the name of the Xiph.org Foundation nor the names of its
  16    contributors may be used to endorse or promote products derived from
  17    this software without specific prior written permission.
  18
  19    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
  20    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
  21    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
  22    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
  23    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
  24    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
  25    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
  26    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
  27    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
  28    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
  29    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  30 */
  31
  32 #ifdef HAVE_CONFIG_H
  33 #include "config.h"
  34 #endif
  35
  36 #include "quant_bands.h"
  37 #include "laplace.h"
  38 #include <math.h>
  39 #include "os_support.h"
  40 #include "arch.h"
  41 #include "mathops.h"
  42
  43 #ifdef FIXED_POINT
  44 const celt_word16_t eMeans[24] = {11520, -2048, -3072, -640, 256, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
  45 #else
  46 const celt_word16_t eMeans[24] = {45.f, -8.f, -12.f, -2.5f, 1.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f};
  47 #endif
  48
  49 /*const int frac[24] = {4, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2};*/
  50 const int frac[24] = {8, 6, 5, 4, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2};
  51
  52 #ifdef FIXED_POINT
  53 static inline celt_ener_t dB2Amp(celt_ener_t dB)
  54 {
  55    celt_ener_t amp;
  56    amp = PSHR32(celt_exp2(MULT16_16_Q14(21771,dB)),2)-QCONST16(.3f, 14);
  57    if (amp < 0)
  58       amp = 0;
  59    return amp;
  60 }
  61
  62 #define DBofTWO 24661
  63 static inline celt_word16_t amp2dB(celt_ener_t amp)
  64 {
  65    /* equivalent to return 6.0207*log2(.3+amp) */
  66    return ROUND(MULT16_16(24661,celt_log2(ADD32(QCONST32(.3f,14),amp))),12);
  67    /* return DB_SCALING*20*log10(.3+ENER_SCALING_1*amp); */
  68 }
  69 #else
  70 static inline celt_ener_t dB2Amp(celt_ener_t dB)
  71 {
  72    celt_ener_t amp;
  73    amp = pow(10, .05*dB)-.3;
  74    if (amp < 0)
  75       amp = 0;
  76    return amp;
  77 }
  78 static inline celt_word16_t amp2dB(celt_ener_t amp)
  79 {
  80    return 20*log10(.3+amp);
  81 }
  82 #endif
  83
  84 static const celt_word16_t base_resolution = QCONST16(6.f,8);
  85
  86 static void quant_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, ec_enc *enc)
  87 {
  88    int i;
  89    int bits;
  90    celt_word16_t prev = 0;
  91    float coef = m->ePredCoef;
  92    VARDECL(celt_word16_t *error);
  93    /* The .7 is a heuristic */
  94    float beta = .7*coef;
  95
  96    ALLOC(error, m->nbEBands, celt_word16_t);
  97    bits = ec_enc_tell(enc, 0);
  98    for (i=0;i<m->nbEBands;i++)
  99    {
 100       int qi;
 101       celt_word16_t q;   /* dB */
 102       celt_word16_t x;   /* dB */
 103       celt_word16_t f;   /* Q8 */
 104       celt_word16_t mean = (1-coef)*eMeans[i];
 105       x = amp2dB(eBands[i]);
 106       f = DIV32_16(SHL32(EXTEND32(x-mean-coef*oldEBands[i]-prev),8),base_resolution);
 107 #ifdef FIXED_POINT
 108       /* Rounding to nearest integer here is really important! */
 109       qi = (f+128)>>8;
 110 #else
 111       qi = (int)floor(.5+f);
 112 #endif
 113       /*ec_laplace_encode(enc, qi, i==0?11192:6192);*/
 114       /*ec_laplace_encode(enc, qi, 8500-i*200);*/
 115       /* If we don't have enough bits to encode all the energy, just assume something safe. */
 116       if (ec_enc_tell(enc, 0) - bits > budget)
 117          qi = -1;
 118       else
 119          ec_laplace_encode(enc, qi, 6000-i*200);
 120       q = qi*base_resolution;
 121       error[i] = f - SHL16(qi,8);
 122
 123       /*printf("%d ", qi);*/
 124       /*printf("%f %f ", pred+prev+q, x);*/
 125       /*printf("%f ", x-pred);*/
 126
 127       oldEBands[i] = mean+coef*oldEBands[i]+prev+q;
 128
 129       prev = mean+prev+(1-beta)*q;
 130    }
 131    /*bits = ec_enc_tell(enc, 0) - bits;*/
 132    /*printf ("%d\n", bits);*/
 133    for (i=0;i<m->nbEBands;i++)
 134    {
 135       int q2;
 136       celt_word16_t offset = (error[i]+QCONST16(.5f,8))*frac[i];
 137       /* FIXME: Instead of giving up without warning, we should degrade everything gracefully */
 138       if (ec_enc_tell(enc, 0) - bits +EC_ILOG(frac[i])> budget)
 139          break;
 140 #ifdef FIXED_POINT
 141       /* Has to be without rounding */
 142       q2 = offset>>8;
 143 #else
 144       q2 = (int)floor(offset);
 145 #endif
 146       if (q2 > frac[i]-1)
 147          q2 = frac[i]-1;
 148       ec_enc_uint(enc, q2, frac[i]);
 149       offset = DIV32_16(SHL16(q2,8)+QCONST16(.5,8),frac[i])-QCONST16(.5f,8);
 150       oldEBands[i] += PSHR32(MULT16_16(DB_SCALING*6,offset),8);
 151       /*printf ("%f ", error[i] - offset);*/
 152    }
 153    for (i=0;i<m->nbEBands;i++)
 154    {
 155       eBands[i] = dB2Amp(oldEBands[i]);
 156    }
 157    /*printf ("%d\n", ec_enc_tell(enc, 0)-9);*/
 158
 159    /*printf ("\n");*/
 160 }
 161
 162 static void unquant_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, ec_dec *dec)
 163 {
 164    int i;
 165    int bits;
 166    celt_word16_t prev = 0;
 167    float coef = m->ePredCoef;
 168    /* The .7 is a heuristic */
 169    float beta = .7*coef;
 170    bits = ec_dec_tell(dec, 0);
 171    for (i=0;i<m->nbEBands;i++)
 172    {
 173       int qi;
 174       celt_word16_t q;
 175       celt_word16_t mean = (1-coef)*eMeans[i];
 176       /* If we didn't have enough bits to encode all the energy, just assume something safe. */
 177       if (ec_dec_tell(dec, 0) - bits > budget)
 178          qi = -1;
 179       else
 180          qi = ec_laplace_decode(dec, 6000-i*200);
 181       q = qi*base_resolution;
 182
 183       /*printf("%d ", qi);*/
 184       /*printf("%f %f ", pred+prev+q, x);*/
 185       /*printf("%f ", x-pred);*/
 186
 187       oldEBands[i] = mean+coef*oldEBands[i]+prev+q;
 188
 189       prev = mean+prev+(1-beta)*q;
 190    }
 191    for (i=0;i<m->nbEBands;i++)
 192    {
 193       int q2;
 194       celt_word16_t offset;
 195       if (ec_dec_tell(dec, 0) - bits +EC_ILOG(frac[i])> budget)
 196          break;
 197       q2 = ec_dec_uint(dec, frac[i]);
 198       offset = DIV32_16(SHL16(q2,8)+QCONST16(.5,8),frac[i])-QCONST16(.5f,8);
 199       oldEBands[i] += PSHR32(MULT16_16(DB_SCALING*6,offset),8);
 200    }
 201    for (i=0;i<m->nbEBands;i++)
 202    {
 203       eBands[i] = dB2Amp(oldEBands[i]);
 204    }
 205    /*printf ("\n");*/
 206 }
 207
 208
 209
 210 void quant_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, ec_enc *enc)
 211 {
 212    int C;
 213
 214    C = m->nbChannels;
 215
 216    if (C==1)
 217       quant_energy_mono(m, eBands, oldEBands, budget, enc);
 218    else
 219 #if 1
 220    {
 221       int c;
 222       VARDECL(celt_ener_t *E);
 223       ALLOC(E, m->nbEBands, celt_ener_t);
 224       for (c=0;c<C;c++)
 225       {
 226          int i;
 227          for (i=0;i<m->nbEBands;i++)
 228             E[i] = eBands[C*i+c];
 229          quant_energy_mono(m, E, oldEBands+c*m->nbEBands, budget/C, enc);
 230          for (i=0;i<m->nbEBands;i++)
 231             eBands[C*i+c] = E[i];
 232       }
 233    }
 234 #else
 235       if (C==2)
 236    {
 237       int i;
 238       int NB = m->nbEBands;
 239       celt_ener_t mid[NB];
 240       celt_ener_t side[NB];
 241       for (i=0;i<NB;i++)
 242       {
 243          //left = eBands[C*i];
 244          //right = eBands[C*i+1];
 245          mid[i] = ENER_SCALING_1*sqrt(eBands[C*i]*eBands[C*i] + eBands[C*i+1]*eBands[C*i+1]);
 246          side[i] = 20*log10((ENER_SCALING_1*eBands[2*i]+.3)/(ENER_SCALING_1*eBands[2*i+1]+.3));
 247          //printf ("%f %f ", mid[i], side[i]);
 248       }
 249       //printf ("\n");
 250       quant_energy_mono(m, mid, oldEBands, enc);
 251       for (i=0;i<NB;i++)
 252          side[i] = pow(10.f,floor(.5f+side[i])/10.f);
 253
 254       //quant_energy_side(m, side, oldEBands+NB, enc);
 255       for (i=0;i<NB;i++)
 256       {
 257          eBands[C*i] = ENER_SCALING*mid[i]*sqrt(side[i]/(1.f+side[i]));
 258          eBands[C*i+1] = ENER_SCALING*mid[i]*sqrt(1.f/(1.f+side[i]));
 259          //printf ("%f %f ", mid[i], side[i]);
 260       }
 261
 262    } else {
 263       celt_fatal("more than 2 channels not supported");
 264    }
 265 #endif
 266 }
 267
 268
 269
 270 void unquant_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, ec_dec *dec)
 271 {
 272    int C;
 273    C = m->nbChannels;
 274
 275    if (C==1)
 276       unquant_energy_mono(m, eBands, oldEBands, budget, dec);
 277    else {
 278       int c;
 279       VARDECL(celt_ener_t *E);
 280       ALLOC(E, m->nbEBands, celt_ener_t);
 281       for (c=0;c<C;c++)
 282       {
 283          int i;
 284          unquant_energy_mono(m, E, oldEBands+c*m->nbEBands, budget/C, dec);
 285          for (i=0;i<m->nbEBands;i++)
 286             eBands[C*i+c] = E[i];
 287       }
 288    }
 289 }