Minor tweaks to the max allocation
[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, int skip_start,
144       const int *bits1, const int *bits2, const int *thresh, int total, int skip_rsv,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 stereo;
153    int codedBands=-1;
154    int alloc_floor;
155    int left, percoeff;
156    int done;
157    int balance;
158    SAVE_STACK;
159
160    alloc_floor = C<<BITRES;
161    stereo = C>1;
162
163    logM = LM<<BITRES;
164    lo = 0;
165    hi = 1<<ALLOC_STEPS;
166    for (i=0;i<ALLOC_STEPS;i++)
167    {
168       int mid = (lo+hi)>>1;
169       psum = 0;
170       done = 0;
171       for (j=end;j-->start;)
172       {
173          int tmp = bits1[j] + (mid*bits2[j]>>ALLOC_STEPS);
174          if (tmp >= thresh[j] || done)
175          {
176             done = 1;
177             /* Don't allocate more than we can actually use */
178             psum += IMIN(tmp, 64*C<<BITRES<<LM);
179          } else {
180             if (tmp >= alloc_floor)
181                psum += alloc_floor;
182          }
183       }
184       if (psum > total)
185          hi = mid;
186       else
187          lo = mid;
188    }
189    psum = 0;
190    /*printf ("interp bisection gave %d\n", lo);*/
191    done = 0;
192    for (j=end;j-->start;)
193    {
194       int tmp = bits1[j] + (lo*bits2[j]>>ALLOC_STEPS);
195       if (tmp < thresh[j] && !done)
196       {
197          if (tmp >= alloc_floor)
198             tmp = alloc_floor;
199          else
200             tmp = 0;
201       } else
202          done = 1;
203       /* Don't allocate more than we can actually use */
204       tmp = IMIN(tmp, 64*C<<BITRES<<LM);
205       bits[j] = tmp;
206       psum += tmp;
207    }
208
209    /* Decide which bands to skip, working backwards from the end. */
210    for (codedBands=end;;codedBands--)
211    {
212       int band_width;
213       int band_bits;
214       int rem;
215       j = codedBands-1;
216       /*Figure out how many left-over bits we would be adding to this band.
217         This can include bits we've stolen back from higher, skipped bands.*/
218       left = total-psum;
219       percoeff = left/(m->eBands[codedBands]-m->eBands[start]);
220       left -= (m->eBands[codedBands]-m->eBands[start])*percoeff;
221       /* Never skip the first band, nor a band that has been boosted by
222           dynalloc.
223          In the first case, we'd be coding a bit to signal we're going to waste
224           all the other bits.
225          In the second case, we'd be coding a bit to redistribute all the bits
226           we just signaled should be cocentrated in this band. */
227       if (j<=skip_start)
228       {
229          /* Give the bit we reserved to end skipping back to this band. */
230          bits[j] += skip_rsv;
231          break;
232       }
233       rem = IMAX(left-(m->eBands[j]-m->eBands[start]),0);
234       band_width = m->eBands[codedBands]-m->eBands[j];
235       band_bits = bits[j] + percoeff*band_width + rem;
236       /*Only code a skip decision if we're above the threshold for this band.
237         Otherwise it is force-skipped.
238         This ensures that we have enough bits to code the skip flag.*/
239       if (band_bits >= IMAX(thresh[j], alloc_floor+(1<<BITRES)))
240       {
241          if (encode)
242          {
243             /*This if() block is the only part of the allocation function that
244                is not a mandatory part of the bitstream: any bands we choose to
245                skip here must be explicitly signaled.*/
246             /*Choose a threshold with some hysteresis to keep bands from
247                fluctuating in and out.*/
248             if (band_bits > ((j<prev?7:9)*band_width<<LM<<BITRES)>>4)
249             {
250                ec_enc_bit_logp((ec_enc *)ec, 1, 1);
251                break;
252             }
253             ec_enc_bit_logp((ec_enc *)ec, 0, 1);
254          } else if (ec_dec_bit_logp((ec_dec *)ec, 1)) {
255             break;
256          }
257          /*We used a bit to skip this band.*/
258          psum += 1<<BITRES;
259          band_bits -= 1<<BITRES;
260       }
261       /*Reclaim the bits originally allocated to this band.*/
262       psum -= bits[j];
263       if (band_bits >= alloc_floor)
264       {
265          /*If we have enough for a fine energy bit per channel, use it.*/
266          psum += alloc_floor;
267          bits[j] = alloc_floor;
268       } else {
269          /*Otherwise this band gets nothing at all.*/
270          bits[j] = 0;
271       }
272    }
273
274    /* Allocate the remaining bits */
275    if (codedBands>start) {
276       for (j=start;j<codedBands;j++)
277          bits[j] += percoeff*(m->eBands[j+1]-m->eBands[j]);
278       for (j=start;j<codedBands;j++)
279       {
280          int tmp = IMIN(left, m->eBands[j+1]-m->eBands[j]);
281          bits[j] += tmp;
282          left -= tmp;
283       }
284    }
285    /*for (j=0;j<end;j++)printf("%d ", bits[j]);printf("\n");*/
286
287    balance = 0;
288    for (j=start;j<codedBands;j++)
289    {
290       int N0, N, den;
291       int offset;
292       int NClogN;
293
294       celt_assert(bits[j] >= 0);
295       N0 = m->eBands[j+1]-m->eBands[j];
296       N=N0<<LM;
297
298       if (N>1)
299       {
300          NClogN = N*C*(m->logN[j] + logM);
301
302          /* Compensate for the extra DoF in stereo */
303          den=(C*N+ ((C==2 && N>2) ? 1 : 0));
304
305          /* Offset for the number of fine bits by log2(N)/2 + FINE_OFFSET
306             compared to their "fair share" of total/N */
307          offset = (NClogN>>1)-N*C*FINE_OFFSET;
308
309          /* N=2 is the only point that doesn't match the curve */
310          if (N==2)
311             offset += N*C<<BITRES>>2;
312
313          /* Changing the offset for allocating the second and third
314              fine energy bit */
315          if (bits[j] + offset < den*2<<BITRES)
316             offset += NClogN>>2;
317          else if (bits[j] + offset < den*3<<BITRES)
318             offset += NClogN>>3;
319
320          /* Divide with rounding */
321          ebits[j] = IMAX(0, (bits[j] + offset + (den<<(BITRES-1))) / (den<<BITRES));
322
323          /* If we rounded down, make it a candidate for final
324              fine energy pass */
325          fine_priority[j] = ebits[j]*(den<<BITRES) >= bits[j]+offset;
326
327          /* Make sure not to bust */
328          if (C*ebits[j] > (bits[j]>>BITRES))
329             ebits[j] = bits[j] >> stereo >> BITRES;
330
331          /* More than 8 is useless because that's about as far as PVQ can go */
332          if (ebits[j]>8)
333             ebits[j]=8;
334
335       } else {
336          /* For N=1, all bits go to fine energy except for a single sign bit */
337          ebits[j] = IMIN(IMAX(0,(bits[j] >> stereo >> BITRES)-1),7);
338          fine_priority[j] = (ebits[j]+1)*C<<BITRES >= (bits[j]-balance);
339          /* N=1 bands can't take advantage of the re-balancing in
340              quant_all_bands() because they don't have shape, only fine energy.
341             Instead, do the re-balancing here.*/
342          balance = IMAX(0,bits[j] - ((ebits[j]+1)*C<<BITRES));
343          if (j+1<codedBands)
344          {
345             bits[j] -= balance;
346             bits[j+1] += balance;
347          }
348       }
349
350       /* The other bits are assigned to PVQ */
351       bits[j] -= C*ebits[j]<<BITRES;
352       celt_assert(bits[j] >= 0);
353       celt_assert(ebits[j] >= 0);
354    }
355    /* The skipped bands use all their bits for fine energy. */
356    for (;j<end;j++)
357    {
358       ebits[j] = bits[j] >> stereo >> BITRES;
359       celt_assert(C*ebits[j]<<BITRES == bits[j]);
360       bits[j] = 0;
361       fine_priority[j] = ebits[j]<1;
362    }
363    RESTORE_STACK;
364    return codedBands;
365 }
366
367 int compute_allocation(const CELTMode *m, int start, int end, const int *offsets, int alloc_trim,
368       int total, int *pulses, int *ebits, int *fine_priority, int _C, int LM, void *ec, int encode, int prev)
369 {
370    int lo, hi, len, j;
371    const int C = CHANNELS(_C);
372    int codedBands;
373    int skip_start;
374    int skip_rsv;
375    VARDECL(int, bits1);
376    VARDECL(int, bits2);
377    VARDECL(int, thresh);
378    VARDECL(int, trim_offset);
379    SAVE_STACK;
380    
381    total = IMAX(total, 0);
382    len = m->nbEBands;
383    skip_start = start;
384    /* Reserve a bit to signal the end of manually skipped bands. */
385    skip_rsv = total >= 1<<BITRES ? 1<<BITRES : 0;
386    total -= skip_rsv;
387    ALLOC(bits1, len, int);
388    ALLOC(bits2, len, int);
389    ALLOC(thresh, len, int);
390    ALLOC(trim_offset, len, int);
391
392    for (j=start;j<end;j++)
393    {
394       /* Below this threshold, we're sure not to allocate any PVQ bits */
395       thresh[j] = IMAX((C)<<BITRES, (3*(m->eBands[j+1]-m->eBands[j])<<LM<<BITRES)>>4);
396       /* Tilt of the allocation curve */
397       trim_offset[j] = C*(m->eBands[j+1]-m->eBands[j])*(alloc_trim-5-LM)*(m->nbEBands-j-1)
398             <<(LM+BITRES)>>6;
399       /* Giving less resolution to single-coefficient bands because they get
400          more benefit from having one coarse value per coefficient*/
401       if ((m->eBands[j+1]-m->eBands[j])<<LM==1)
402          trim_offset[j] -= C<<BITRES;
403    }
404    lo = 1;
405    hi = m->nbAllocVectors - 2;
406    do
407    {
408       int done = 0;
409       int psum = 0;
410       int mid = (lo+hi) >> 1;
411       for (j=end;j-->start;)
412       {
413          int N = m->eBands[j+1]-m->eBands[j];
414          bits1[j] = C*N*m->allocVectors[mid*len+j]<<LM>>2;
415          if (bits1[j] > 0)
416             bits1[j] = IMAX(0, bits1[j] + trim_offset[j]);
417          bits1[j] += offsets[j];
418          if (bits1[j] >= thresh[j] || done)
419          {
420             done = 1;
421             /* Don't allocate more than we can actually use */
422             psum += IMIN(bits1[j], 64*C<<BITRES<<LM);
423          } else {
424             if (bits1[j] >= C<<BITRES)
425                psum += C<<BITRES;
426          }
427       }
428       if (psum > total)
429          hi = mid - 1;
430       else
431          lo = mid + 1;
432       /*printf ("lo = %d, hi = %d\n", lo, hi);*/
433    }
434    while (lo <= hi);
435    hi = lo--;
436    /*printf ("interp between %d and %d\n", lo, hi);*/
437    for (j=start;j<end;j++)
438    {
439       int N = m->eBands[j+1]-m->eBands[j];
440       bits1[j] = C*N*m->allocVectors[lo*len+j]<<LM>>2;
441       bits2[j] = C*N*m->allocVectors[hi*len+j]<<LM>>2;
442       if (bits1[j] > 0)
443          bits1[j] = IMAX(0, bits1[j] + trim_offset[j]);
444       if (bits2[j] > 0)
445          bits2[j] = IMAX(0, bits2[j] + trim_offset[j]);
446       if (lo > 0)
447          bits1[j] += offsets[j];
448       bits2[j] += offsets[j];
449       if (offsets[j]>0)
450          skip_start = j;
451       bits2[j] -= bits1[j];
452    }
453    codedBands = interp_bits2pulses(m, start, end, skip_start, bits1, bits2, thresh,
454          total, skip_rsv, pulses, ebits, fine_priority, len, C, LM, ec, encode, prev);
455    RESTORE_STACK;
456    return codedBands;
457 }
458