More negative bits prevention
[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, int *skip, 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    SAVE_STACK;
154
155    logM = LM<<BITRES;
156    lo = 0;
157    hi = 1<<ALLOC_STEPS;
158    for (i=0;i<ALLOC_STEPS;i++)
159    {
160       int mid = (lo+hi)>>1;
161       psum = 0;
162       for (j=start;j<end;j++)
163       {
164          int tmp = bits1[j] + (mid*bits2[j]>>ALLOC_STEPS);
165          /* Don't allocate more than we can actually use */
166          if (tmp >= thresh[j])
167             psum += tmp;
168          else if (tmp >= C<<BITRES)
169             psum += C<<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 >= C<<BITRES)
186          bits[j] = C<<BITRES;
187       else
188          bits[j] = 0;
189       /* Don't allocate more than we can actually use */
190       if (bits[j] > 64*C<<BITRES<<LM)
191          bits[j] = 64*C<<BITRES<<LM;
192       psum += bits[j];
193    }
194    codedBands++;
195    if (*skip==-1)
196    {
197       *skip=0;
198       for (j=codedBands-1;j>=0;j--)
199       {
200          if ((bits[j] > (7*(m->eBands[j+1]-m->eBands[j])<<LM<<BITRES)>>4 && j<prev)
201                || (bits[j] > (9*(m->eBands[j+1]-m->eBands[j])<<LM<<BITRES)>>4))
202             break;
203          else
204             (*skip)++;
205       }
206       *skip = IMIN(*skip, codedBands-start-1);
207    }
208    for (i=0;i<*skip;i++)
209    {
210       int alloc_floor;
211       if (C==2)
212          alloc_floor = 4<<BITRES; /* Fine energy plus intensity stereo */
213       else
214          alloc_floor = 1<<BITRES;
215       if (bits[codedBands-1] >= alloc_floor)
216       {
217          psum = psum - bits[codedBands-1] + alloc_floor;
218          bits[codedBands-1] = alloc_floor;
219       } else if (bits[codedBands-1] < C<<BITRES)
220       {
221          psum = psum - bits[codedBands-1];
222          bits[codedBands-1] = 0;
223       }
224       codedBands--;
225    }
226    /* Allocate the remaining bits */
227    if (codedBands) {
228       int left, perband;
229       left = (total<<BITRES)-psum;
230       perband = left/(m->eBands[codedBands]-m->eBands[start]);
231       for (j=start;j<codedBands;j++)
232          bits[j] += perband*(m->eBands[j+1]-m->eBands[j]);
233       left = left-(m->eBands[codedBands]-m->eBands[start])*perband;
234       for (j=start;j<codedBands;j++)
235       {
236          int tmp = IMIN(left, m->eBands[j+1]-m->eBands[j]);
237          bits[j] += tmp;
238          left -= tmp;
239       }
240    }
241    /*for (j=0;j<end;j++)printf("%d ", bits[j]);printf("\n");*/
242    for (j=start;j<end;j++)
243    {
244       int N0, N, den;
245       int offset;
246       int NClogN;
247
248       N0 = m->eBands[j+1]-m->eBands[j];
249       N=N0<<LM;
250       NClogN = N*C*(m->logN[j] + logM);
251
252       /* Compensate for the extra DoF in stereo */
253       den=(C*N+ ((C==2 && N>2) ? 1 : 0));
254
255       /* Offset for the number of fine bits by log2(N)/2 + FINE_OFFSET
256          compared to their "fair share" of total/N */
257       offset = (NClogN>>1)-N*C*FINE_OFFSET;
258
259       /* N=2 is the only point that doesn't match the curve */
260       if (N==2)
261          offset += N*C<<BITRES>>2;
262
263       /* Changing the offset for allocating the second and third fine energy bit */
264       if (bits[j] + offset < den*2<<BITRES)
265          offset += NClogN>>2;
266       else if (bits[j] + offset < den*3<<BITRES)
267          offset += NClogN>>3;
268
269       /* Divide with rounding */
270       ebits[j] = (bits[j] + offset + (den<<(BITRES-1))) / (den<<BITRES);
271
272       /* If we rounded down, make it a candidate for final fine energy pass */
273       fine_priority[j] = ebits[j]*(den<<BITRES) >= bits[j]+offset;
274
275       /* For N=1, all bits go to fine energy except for a single sign bit */
276       if (N==1)
277       {
278          ebits[j] = IMAX(0,(bits[j]/C >> BITRES)-1);
279          fine_priority[j] = (ebits[j]+1)*C<<BITRES >= bits[j];
280       }
281       /* Make sure not to bust */
282       if (C*ebits[j] > (bits[j]>>BITRES))
283          ebits[j] = bits[j]/C >> BITRES;
284
285       /* More than that is useless because that's about as far as PVQ can go */
286       if (ebits[j]>7)
287          ebits[j]=7;
288
289       /* The other bits are assigned to PVQ */
290       bits[j] -= C*ebits[j]<<BITRES;
291       celt_assert(bits[j] >= 0);
292    }
293    RESTORE_STACK;
294    return codedBands;
295 }
296
297 int compute_allocation(const CELTMode *m, int start, int end, int *offsets, int alloc_trim,
298       int total, int *pulses, int *ebits, int *fine_priority, int _C, int LM, int *skip, int prev)
299 {
300    int lo, hi, len, j;
301    const int C = CHANNELS(_C);
302    int codedBands;
303    VARDECL(int, bits1);
304    VARDECL(int, bits2);
305    VARDECL(int, thresh);
306    VARDECL(int, trim_offset);
307    SAVE_STACK;
308    
309    len = m->nbEBands;
310    ALLOC(bits1, len, int);
311    ALLOC(bits2, len, int);
312    ALLOC(thresh, len, int);
313    ALLOC(trim_offset, len, int);
314
315    /* Below this threshold, we don't allocate any PVQ bits */
316    for (j=start;j<end;j++)
317       thresh[j] = ((C+1)<<BITRES) + (3*(m->eBands[j+1]-m->eBands[j])<<LM<<BITRES)>>4;
318    /* Tilt of the allocation curve */
319    for (j=start;j<end;j++)
320       trim_offset[j] = C*(m->eBands[j+1]-m->eBands[j])*(alloc_trim-5-LM)*(m->nbEBands-j-1)
321             <<(LM+BITRES)>>6;
322
323    lo = 0;
324    hi = m->nbAllocVectors - 1;
325    while (hi-lo != 1)
326    {
327       int psum = 0;
328       int mid = (lo+hi) >> 1;
329       for (j=start;j<end;j++)
330       {
331          int N = m->eBands[j+1]-m->eBands[j];
332          bits1[j] = C*N*m->allocVectors[mid*len+j]<<LM>>2;
333          if (bits1[j] > 0)
334             bits1[j] += trim_offset[j];
335          if (bits1[j] < 0)
336             bits1[j] = 0;
337          bits1[j] += offsets[j];
338          if (bits1[j] >= thresh[j])
339             psum += bits1[j];
340          else if (bits1[j] >= C<<BITRES)
341             psum += C<<BITRES;
342
343          /*printf ("%d ", bits[j]);*/
344       }
345       /*printf ("\n");*/
346       if (psum > (total<<BITRES))
347          hi = mid;
348       else
349          lo = mid;
350       /*printf ("lo = %d, hi = %d\n", lo, hi);*/
351    }
352    /*printf ("interp between %d and %d\n", lo, hi);*/
353    for (j=start;j<end;j++)
354    {
355       int N = m->eBands[j+1]-m->eBands[j];
356       bits1[j] = (C*N*m->allocVectors[lo*len+j]<<LM>>2);
357       bits2[j] = (C*N*m->allocVectors[hi*len+j]<<LM>>2) - bits1[j];
358       if (bits1[j] > 0)
359          bits1[j] += trim_offset[j];
360       if (bits1[j] < 0)
361          bits1[j] = 0;
362       bits1[j] += offsets[j];
363    }
364    codedBands = interp_bits2pulses(m, start, end, bits1, bits2, thresh,
365          total, pulses, ebits, fine_priority, len, C, LM, skip, prev);
366    RESTORE_STACK;
367    return codedBands;
368 }
369