Preventing bands from being coded at a rate below (for now) 3/8 bit/sample
[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    for (i=0;i<=LM+1;i++)
85    {
86       int j;
87       for (j=0;j<m->nbEBands;j++)
88       {
89          int k;
90          int N = (eBands[j+1]-eBands[j])<<i>>1;
91          cindex[i*m->nbEBands+j] = -1;
92          for (k=0;k<=i;k++)
93          {
94             int n;
95             for (n=0;n<m->nbEBands && (k!=i || n<j);n++)
96             {
97                if (N == (eBands[n+1]-eBands[n])<<k>>1)
98                {
99                   cindex[i*m->nbEBands+j] = cindex[k*m->nbEBands+n];
100                   break;
101                }
102             }
103          }
104          if (cache->index[i*m->nbEBands+j] == -1 && N!=0)
105          {
106             int K;
107             entryN[nbEntries] = N;
108             K = 0;
109             while (fits_in32(N,get_pulses(K+1)) && K<MAX_PSEUDO-1)
110                K++;
111             entryK[nbEntries] = K;
112             cindex[i*m->nbEBands+j] = curr;
113             entryI[nbEntries] = curr;
114
115             curr += K+1;
116             nbEntries++;
117          }
118       }
119    }
120    bits = celt_alloc(sizeof(unsigned char)*curr);
121    cache->bits = bits;
122    cache->size = curr;
123    for (i=0;i<nbEntries;i++)
124    {
125       int j;
126       unsigned char *ptr = bits+entryI[i];
127       celt_int16 tmp[MAX_PULSES];
128       get_required_bits(tmp, entryN[i], get_pulses(entryK[i]), BITRES);
129       for (j=1;j<=entryK[i];j++)
130          ptr[j] = tmp[get_pulses(j)]-1;
131       ptr[0] = entryK[i];
132    }
133 }
134
135 #endif /* !STATIC_MODES */
136
137
138 #define ALLOC_STEPS 6
139
140 static inline int interp_bits2pulses(const CELTMode *m, int start, int end, int *bits1, int *bits2, int total, int *bits, int *ebits, int *fine_priority, int len, int _C, int LM)
141 {
142    int psum;
143    int lo, hi;
144    int i, j;
145    int logM;
146    const int C = CHANNELS(_C);
147    int codedBands=-1;
148    VARDECL(int, thresh);
149    SAVE_STACK;
150
151    ALLOC(thresh, len, int);
152
153    /* Threshold: don't allow any band to go below 3/8 bit/sample */
154    for (j=start;j<end;j++)
155       thresh[j] = 3*(C*(m->eBands[j+1]-m->eBands[j])<<LM<<BITRES)>>3;
156    logM = LM<<BITRES;
157    lo = 0;
158    hi = 1<<ALLOC_STEPS;
159    for (i=0;i<ALLOC_STEPS;i++)
160    {
161       int mid = (lo+hi)>>1;
162       psum = 0;
163       for (j=start;j<end;j++)
164       {
165          int tmp = bits1[j] + (mid*bits2[j]>>ALLOC_STEPS);
166          if (tmp >= thresh[j])
167             psum += tmp;
168          else if (tmp >= 1<<BITRES)
169             psum += 1<<BITRES;
170       }
171       if (psum > (total<<BITRES))
172          hi = mid;
173       else
174          lo = mid;
175    }
176    psum = 0;
177    /*printf ("interp bisection gave %d\n", lo);*/
178    for (j=start;j<end;j++)
179    {
180       int tmp = bits1[j] + (lo*bits2[j]>>ALLOC_STEPS);
181       if (tmp >= thresh[j])
182       {
183          bits[j] = tmp;
184          codedBands = j;
185       } else if (tmp >= 1<<BITRES)
186          bits[j] = 1<<BITRES;
187       else
188          bits[j] = 0;
189       psum += bits[j];
190    }
191    codedBands++;
192    /* Allocate the remaining bits */
193    {
194       int left, perband;
195       left = (total<<BITRES)-psum;
196       perband = left/(codedBands-start);
197       for (j=start;j<codedBands;j++)
198          bits[j] += perband;
199       left = left-codedBands*perband;
200       for (j=start;j<start+left;j++)
201          bits[j]++;
202    }
203    for (j=start;j<end;j++)
204    {
205       int N0, N, den;
206       int offset;
207       int NClogN;
208
209       N0 = m->eBands[j+1]-m->eBands[j];
210       N=N0<<LM;
211       NClogN = N*C*(m->logN[j] + logM);
212
213       /* Compensate for the extra DoF in stereo */
214       den=(C*N+ ((C==2 && N>2) ? 1 : 0));
215
216       /* Offset for the number of fine bits by log2(N)/2 + FINE_OFFSET
217          compared to their "fair share" of total/N */
218       offset = (NClogN>>1)-N*C*FINE_OFFSET;
219
220       /* N=2 is the only point that doesn't match the curve */
221       if (N==2)
222          offset += N*C<<BITRES>>2;
223
224       /* Changing the offset for allocating the second and third fine energy bit */
225       if (bits[j] + offset < den*2<<BITRES)
226          offset += NClogN>>2;
227       else if (bits[j] + offset < den*3<<BITRES)
228          offset += NClogN>>3;
229
230       /* Divide with rounding */
231       ebits[j] = (bits[j] + offset + (den<<(BITRES-1))) / (den<<BITRES);
232
233       /* If we rounded down, make it a candidate for final fine energy pass */
234       fine_priority[j] = ebits[j]*(den<<BITRES) >= bits[j]+offset;
235
236       /* For N=1, all bits go to fine energy except for a single sign bit */
237       if (N==1)
238          ebits[j] = (bits[j]/C >> BITRES)-1;
239       /* Make sure the first bit is spent on fine energy */
240       if (ebits[j] < 1)
241          ebits[j] = 1;
242
243       /* Make sure not to bust */
244       if (C*ebits[j] > (bits[j]>>BITRES))
245          ebits[j] = bits[j]/C >> BITRES;
246
247       /* More than that is useless because that's about as far as PVQ can go */
248       if (ebits[j]>7)
249          ebits[j]=7;
250
251       /* The other bits are assigned to PVQ */
252       bits[j] -= C*ebits[j]<<BITRES;
253       if (bits[j] < 0)
254          bits[j] = 0;
255    }
256    RESTORE_STACK;
257    return codedBands;
258 }
259
260 int compute_allocation(const CELTMode *m, int start, int end, int *offsets, int total, int *pulses, int *ebits, int *fine_priority, int _C, int LM)
261 {
262    int lo, hi, len, j;
263    const int C = CHANNELS(_C);
264    int codedBands;
265    VARDECL(int, bits1);
266    VARDECL(int, bits2);
267    SAVE_STACK;
268    
269    len = m->nbEBands;
270    ALLOC(bits1, len, int);
271    ALLOC(bits2, len, int);
272
273    lo = 0;
274    hi = m->nbAllocVectors - 1;
275    while (hi-lo != 1)
276    {
277       int psum = 0;
278       int mid = (lo+hi) >> 1;
279       for (j=start;j<end;j++)
280       {
281          int N = m->eBands[j+1]-m->eBands[j];
282          bits1[j] = ((C*N*m->allocVectors[mid*len+j]<<LM>>2) + offsets[j]);
283          psum += bits1[j];
284          /*printf ("%d ", bits[j]);*/
285       }
286       /*printf ("\n");*/
287       if (psum > (total<<BITRES))
288          hi = mid;
289       else
290          lo = mid;
291       /*printf ("lo = %d, hi = %d\n", lo, hi);*/
292    }
293    /*printf ("interp between %d and %d\n", lo, hi);*/
294    for (j=start;j<end;j++)
295    {
296       int N = m->eBands[j+1]-m->eBands[j];
297       bits1[j] = (C*N*m->allocVectors[lo*len+j]<<LM>>2);
298       bits2[j] = (C*N*m->allocVectors[hi*len+j]<<LM>>2) - bits1[j];
299       bits1[j] += offsets[j];
300    }
301    codedBands = interp_bits2pulses(m, start, end, bits1, bits2, total, pulses, ebits, fine_priority, len, C, LM);
302    RESTORE_STACK;
303    return codedBands;
304 }
305