Setting fine_priority for skipped bands
[opus.git] / libcelt / rate.c
1 /* Copyright (c) 2007-2008 CSIRO
2    Copyright (c) 2007-2009 Xiph.Org Foundation
3    Written by Jean-Marc Valin */
4 /*
5    Redistribution and use in source and binary forms, with or without
6    modification, are permitted provided that the following conditions
7    are met:
8    
9    - Redistributions of source code must retain the above copyright
10    notice, this list of conditions and the following disclaimer.
11    
12    - Redistributions in binary form must reproduce the above copyright
13    notice, this list of conditions and the following disclaimer in the
14    documentation and/or other materials provided with the distribution.
15    
16    - Neither the name of the Xiph.org Foundation nor the names of its
17    contributors may be used to endorse or promote products derived from
18    this software without specific prior written permission.
19    
20    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
23    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
24    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
25    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
26    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
27    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
28    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
29    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
30    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 */
32
33 #ifdef HAVE_CONFIG_H
34 #include "config.h"
35 #endif
36
37 #include <math.h>
38 #include "modes.h"
39 #include "cwrs.h"
40 #include "arch.h"
41 #include "os_support.h"
42
43 #include "entcode.h"
44 #include "rate.h"
45
46
47 #ifndef STATIC_MODES
48
49 /*Determines if V(N,K) fits in a 32-bit unsigned integer.
50   N and K are themselves limited to 15 bits.*/
51 static int fits_in32(int _n, int _k)
52 {
53    static const celt_int16 maxN[15] = {
54       32767, 32767, 32767, 1476, 283, 109,  60,  40,
55        29,  24,  20,  18,  16,  14,  13};
56    static const celt_int16 maxK[15] = {
57       32767, 32767, 32767, 32767, 1172, 238,  95,  53,
58        36,  27,  22,  18,  16,  15,  13};
59    if (_n>=14)
60    {
61       if (_k>=14)
62          return 0;
63       else
64          return _n <= maxN[_k];
65    } else {
66       return _k <= maxK[_n];
67    }
68 }
69
70 void compute_pulse_cache(CELTMode *m, int LM)
71 {
72    int i;
73    int curr=0;
74    int nbEntries=0;
75    int entryN[100], entryK[100], entryI[100];
76    const celt_int16 *eBands = m->eBands;
77    PulseCache *cache = &m->cache;
78    celt_int16 *cindex;
79    unsigned char *bits;
80
81    cindex = celt_alloc(sizeof(cache->index[0])*m->nbEBands*(LM+2));
82    cache->index = cindex;
83
84    /* Scan for all unique band sizes */
85    for (i=0;i<=LM+1;i++)
86    {
87       int j;
88       for (j=0;j<m->nbEBands;j++)
89       {
90          int k;
91          int N = (eBands[j+1]-eBands[j])<<i>>1;
92          cindex[i*m->nbEBands+j] = -1;
93          /* Find other bands that have the same size */
94          for (k=0;k<=i;k++)
95          {
96             int n;
97             for (n=0;n<m->nbEBands && (k!=i || n<j);n++)
98             {
99                if (N == (eBands[n+1]-eBands[n])<<k>>1)
100                {
101                   cindex[i*m->nbEBands+j] = cindex[k*m->nbEBands+n];
102                   break;
103                }
104             }
105          }
106          if (cache->index[i*m->nbEBands+j] == -1 && N!=0)
107          {
108             int K;
109             entryN[nbEntries] = N;
110             K = 0;
111             while (fits_in32(N,get_pulses(K+1)) && K<MAX_PSEUDO)
112                K++;
113             entryK[nbEntries] = K;
114             cindex[i*m->nbEBands+j] = curr;
115             entryI[nbEntries] = curr;
116
117             curr += K+1;
118             nbEntries++;
119          }
120       }
121    }
122    bits = celt_alloc(sizeof(unsigned char)*curr);
123    cache->bits = bits;
124    cache->size = curr;
125    /* Compute the cache for all unique sizes */
126    for (i=0;i<nbEntries;i++)
127    {
128       int j;
129       unsigned char *ptr = bits+entryI[i];
130       celt_int16 tmp[MAX_PULSES+1];
131       get_required_bits(tmp, entryN[i], get_pulses(entryK[i]), BITRES);
132       for (j=1;j<=entryK[i];j++)
133          ptr[j] = tmp[get_pulses(j)]-1;
134       ptr[0] = entryK[i];
135    }
136 }
137
138 #endif /* !STATIC_MODES */
139
140
141 #define ALLOC_STEPS 6
142
143 static inline int interp_bits2pulses(const CELTMode *m, int start, int end,
144       int *bits1, int *bits2, const int *thresh, int total, int *bits,
145       int *ebits, int *fine_priority, int len, int _C, int LM, void *ec, int encode, int prev)
146 {
147    int psum;
148    int lo, hi;
149    int i, j;
150    int logM;
151    const int C = CHANNELS(_C);
152    int codedBands=-1;
153    int alloc_floor;
154    int left, percoeff;
155    int done;
156    int balance;
157    SAVE_STACK;
158
159    alloc_floor = C<<BITRES;
160
161    logM = LM<<BITRES;
162    lo = 0;
163    hi = 1<<ALLOC_STEPS;
164    for (i=0;i<ALLOC_STEPS;i++)
165    {
166       int mid = (lo+hi)>>1;
167       psum = 0;
168       done = 0;
169       for (j=end;j-->start;)
170       {
171          int tmp = bits1[j] + (mid*bits2[j]>>ALLOC_STEPS);
172          if (tmp >= thresh[j] || done)
173          {
174             done = 1;
175             /* Don't allocate more than we can actually use */
176             psum += IMIN(tmp, 64*C<<BITRES<<LM);
177          } else {
178             if (tmp >= alloc_floor)
179                psum += alloc_floor;
180          }
181       }
182       if (psum > total)
183          hi = mid;
184       else
185          lo = mid;
186    }
187    psum = 0;
188    /*printf ("interp bisection gave %d\n", lo);*/
189    done = 0;
190    for (j=end;j-->start;)
191    {
192       int tmp = bits1[j] + (lo*bits2[j]>>ALLOC_STEPS);
193       if (tmp < thresh[j] && !done)
194       {
195          if (tmp >= alloc_floor)
196             tmp = alloc_floor;
197          else
198             tmp = 0;
199       } else
200          done = 1;
201       /* Don't allocate more than we can actually use */
202       tmp = IMIN(tmp, 64*C<<BITRES<<LM);
203       bits[j] = tmp;
204       psum += tmp;
205    }
206
207    /* Decide which bands to skip, working backwards from the end. */
208    for (codedBands=end;;codedBands--)
209    {
210       int band_width;
211       int band_bits;
212       int rem;
213       j = codedBands-1;
214       /*Figure out how many left-over bits we would be adding to this band.
215         This can include bits we've stolen back from higher, skipped bands.*/
216       left = total-psum;
217       percoeff = left/(m->eBands[codedBands]-m->eBands[start]);
218       left -= (m->eBands[codedBands]-m->eBands[start])*percoeff;
219       /* Never skip the first band: we'd be coding a bit to signal that we're
220           going to waste all the other bits.
221          This means we won't be using the extra bit we reserved to signal the
222           end of manual skipping, but that will get added back in by
223           quant_all_bands().*/
224       if (j<=start)
225          break;
226       rem = IMAX(left-(m->eBands[j]-m->eBands[start]),0);
227       band_width = m->eBands[codedBands]-m->eBands[j];
228       band_bits = bits[j] + percoeff*band_width + rem;
229       /*Only code a skip decision if we're above the threshold for this band.
230         Otherwise it is force-skipped.
231         This ensures that we have enough bits to code the skip flag.*/
232       if (band_bits >= IMAX(thresh[j], alloc_floor+(1<<BITRES)))
233       {
234          if (encode)
235          {
236             /*This if() block is the only part of the allocation function that
237                is not a mandatory part of the bitstream: any bands we choose to
238                skip here must be explicitly signaled.*/
239             /*Choose a threshold with some hysteresis to keep bands from
240                fluctuating in and out.*/
241             if (band_bits > ((j<prev?7:9)*band_width<<LM<<BITRES)>>4)
242             {
243                ec_enc_bit_prob((ec_enc *)ec, 1, 32768);
244                break;
245             }
246             ec_enc_bit_prob((ec_enc *)ec, 0, 32768);
247          } else if (ec_dec_bit_prob((ec_dec *)ec, 32768)) {
248             break;
249          }
250          /*We used a bit to skip this band.*/
251          psum += 1<<BITRES;
252          band_bits -= 1<<BITRES;
253       }
254       /*Reclaim the bits originally allocated to this band.*/
255       psum -= bits[j];
256       if (band_bits >= alloc_floor)
257       {
258          /*If we have enough for a fine energy bit per channel, use it.*/
259          psum += alloc_floor;
260          bits[j] = alloc_floor;
261       } else {
262          /*Otherwise this band gets nothing at all.*/
263          bits[j] = 0;
264       }
265    }
266
267    /* Allocate the remaining bits */
268    if (codedBands>start) {
269       for (j=start;j<codedBands;j++)
270          bits[j] += percoeff*(m->eBands[j+1]-m->eBands[j]);
271       for (j=start;j<codedBands;j++)
272       {
273          int tmp = IMIN(left, m->eBands[j+1]-m->eBands[j]);
274          bits[j] += tmp;
275          left -= tmp;
276       }
277    }
278    /*for (j=0;j<end;j++)printf("%d ", bits[j]);printf("\n");*/
279
280    balance = 0;
281    for (j=start;j<codedBands;j++)
282    {
283       int N0, N, den;
284       int offset;
285       int NClogN;
286
287       celt_assert(bits[j] >= 0);
288       N0 = m->eBands[j+1]-m->eBands[j];
289       N=N0<<LM;
290
291       if (N>1)
292       {
293          NClogN = N*C*(m->logN[j] + logM);
294
295          /* Compensate for the extra DoF in stereo */
296          den=(C*N+ ((C==2 && N>2) ? 1 : 0));
297
298          /* Offset for the number of fine bits by log2(N)/2 + FINE_OFFSET
299             compared to their "fair share" of total/N */
300          offset = (NClogN>>1)-N*C*FINE_OFFSET;
301
302          /* N=2 is the only point that doesn't match the curve */
303          if (N==2)
304             offset += N*C<<BITRES>>2;
305
306          /* Changing the offset for allocating the second and third
307              fine energy bit */
308          if (bits[j] + offset < den*2<<BITRES)
309             offset += NClogN>>2;
310          else if (bits[j] + offset < den*3<<BITRES)
311             offset += NClogN>>3;
312
313          /* Divide with rounding */
314          ebits[j] = IMAX(0, (bits[j] + offset + (den<<(BITRES-1))) / (den<<BITRES));
315
316          /* If we rounded down, make it a candidate for final
317              fine energy pass */
318          fine_priority[j] = ebits[j]*(den<<BITRES) >= bits[j]+offset;
319
320          /* Make sure not to bust */
321          if (C*ebits[j] > (bits[j]>>BITRES))
322             ebits[j] = bits[j]/C >> BITRES;
323
324          /* More than 7 is useless because that's about as far as PVQ can go */
325          if (ebits[j]>7)
326             ebits[j]=7;
327
328       } else {
329          /* For N=1, all bits go to fine energy except for a single sign bit */
330          ebits[j] = IMIN(IMAX(0,(bits[j]/C >> BITRES)-1),7);
331          fine_priority[j] = (ebits[j]+1)*C<<BITRES >= (bits[j]-balance);
332          /* N=1 bands can't take advantage of the re-balancing in
333              quant_all_bands() because they don't have shape, only fine energy.
334             Instead, do the re-balancing here.*/
335          balance = IMAX(0,bits[j] - ((ebits[j]+1)*C<<BITRES));
336          if (j+1<codedBands)
337          {
338             bits[j] -= balance;
339             bits[j+1] += balance;
340          }
341       }
342
343       /* The other bits are assigned to PVQ */
344       bits[j] -= C*ebits[j]<<BITRES;
345       celt_assert(bits[j] >= 0);
346       celt_assert(ebits[j] >= 0);
347    }
348    /* The skipped bands use all their bits for fine energy. */
349    for (;j<end;j++)
350    {
351       ebits[j] = bits[j]/C >> BITRES;
352       celt_assert(C*ebits[j]<<BITRES == bits[j]);
353       bits[j] = 0;
354       fine_priority[j] = ebits[j]<1;
355    }
356    RESTORE_STACK;
357    return codedBands;
358 }
359
360 int compute_allocation(const CELTMode *m, int start, int end, int *offsets, int alloc_trim,
361       int total, int *pulses, int *ebits, int *fine_priority, int _C, int LM, void *ec, int encode, int prev)
362 {
363    int lo, hi, len, j;
364    const int C = CHANNELS(_C);
365    int codedBands;
366    VARDECL(int, bits1);
367    VARDECL(int, bits2);
368    VARDECL(int, thresh);
369    VARDECL(int, trim_offset);
370    SAVE_STACK;
371    
372    total = IMAX(total, 0);
373    len = m->nbEBands;
374    ALLOC(bits1, len, int);
375    ALLOC(bits2, len, int);
376    ALLOC(thresh, len, int);
377    ALLOC(trim_offset, len, int);
378
379    /* Below this threshold, we're sure not to allocate any PVQ bits */
380    for (j=start;j<end;j++)
381       thresh[j] = IMAX((C)<<BITRES, (3*(m->eBands[j+1]-m->eBands[j])<<LM<<BITRES)>>4);
382    /* Tilt of the allocation curve */
383    for (j=start;j<end;j++)
384       trim_offset[j] = C*(m->eBands[j+1]-m->eBands[j])*(alloc_trim-5-LM)*(m->nbEBands-j-1)
385             <<(LM+BITRES)>>6;
386
387    lo = 0;
388    hi = m->nbAllocVectors - 1;
389    while (hi-lo != 1)
390    {
391       int psum = 0;
392       int mid = (lo+hi) >> 1;
393       for (j=start;j<end;j++)
394       {
395          int N = m->eBands[j+1]-m->eBands[j];
396          bits1[j] = C*N*m->allocVectors[mid*len+j]<<LM>>2;
397          if (bits1[j] > 0)
398             bits1[j] += trim_offset[j];
399          if (bits1[j] < 0)
400             bits1[j] = 0;
401          bits1[j] += offsets[j];
402          if (bits1[j] >= thresh[j])
403             psum += bits1[j];
404          else if (bits1[j] >= C<<BITRES)
405             psum += C<<BITRES;
406
407          /*printf ("%d ", bits[j]);*/
408       }
409       /*printf ("\n");*/
410       if (psum > total)
411          hi = mid;
412       else
413          lo = mid;
414       /*printf ("lo = %d, hi = %d\n", lo, hi);*/
415    }
416    /*printf ("interp between %d and %d\n", lo, hi);*/
417    for (j=start;j<end;j++)
418    {
419       int N = m->eBands[j+1]-m->eBands[j];
420       bits1[j] = (C*N*m->allocVectors[lo*len+j]<<LM>>2);
421       bits2[j] = (C*N*m->allocVectors[hi*len+j]<<LM>>2) - bits1[j];
422       if (bits1[j] > 0)
423          bits1[j] += trim_offset[j];
424       if (bits1[j] < 0)
425          bits1[j] = 0;
426       bits1[j] += offsets[j];
427    }
428    codedBands = interp_bits2pulses(m, start, end, bits1, bits2, thresh,
429          total, pulses, ebits, fine_priority, len, C, LM, ec, encode, prev);
430    RESTORE_STACK;
431    return codedBands;
432 }
433