Converting allocation table to 1/32 bit/sample resolution.
[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 3
139
140 static inline void 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    SAVE_STACK;
148
149    logM = LM<<BITRES;
150    lo = 0;
151    hi = 1<<ALLOC_STEPS;
152    for (i=0;i<ALLOC_STEPS;i++)
153    {
154       int mid = (lo+hi)>>1;
155       psum = 0;
156       for (j=start;j<end;j++)
157          psum += (((1<<ALLOC_STEPS)-mid)*bits1[j] + mid*bits2[j])>>ALLOC_STEPS;
158       if (psum > (total<<BITRES))
159          hi = mid;
160       else
161          lo = mid;
162    }
163    psum = 0;
164    /*printf ("interp bisection gave %d\n", lo);*/
165    for (j=start;j<end;j++)
166    {
167       bits[j] = (((1<<ALLOC_STEPS)-lo)*bits1[j] + lo*bits2[j])>>ALLOC_STEPS;
168       psum += bits[j];
169    }
170    /* Allocate the remaining bits */
171    {
172       int left, perband;
173       left = (total<<BITRES)-psum;
174       perband = left/(end-start);
175       for (j=start;j<end;j++)
176          bits[j] += perband;
177       left = left-end*perband;
178       for (j=start;j<start+left;j++)
179          bits[j]++;
180    }
181    for (j=start;j<end;j++)
182    {
183       int N0, N, den;
184       int offset;
185       int NClogN;
186
187       N0 = m->eBands[j+1]-m->eBands[j];
188       N=N0<<LM;
189       NClogN = N*C*(m->logN[j] + logM);
190
191       /* Compensate for the extra DoF in stereo */
192       den=(C*N+ ((C==2 && N>2) ? 1 : 0));
193
194       /* Offset for the number of fine bits by log2(N)/2 + FINE_OFFSET
195          compared to their "fair share" of total/N */
196       offset = (NClogN>>1)-N*C*FINE_OFFSET;
197
198       /* N=2 is the only point that doesn't match the curve */
199       if (N==2)
200          offset += N*C<<BITRES>>2;
201
202       /* Changing the offset for allocating the second and third fine energy bit */
203       if (bits[j] + offset < den*2<<BITRES)
204          offset += NClogN>>2;
205       else if (bits[j] + offset < den*3<<BITRES)
206          offset += NClogN>>3;
207
208       /* Divide with rounding */
209       ebits[j] = (bits[j] + offset + (den<<(BITRES-1))) / (den<<BITRES);
210
211       /* If we rounded down, make it a candidate for final fine energy pass */
212       fine_priority[j] = ebits[j]*(den<<BITRES) >= bits[j]+offset;
213
214       /* For N=1, all bits go to fine energy except for a single sign bit */
215       if (N==1)
216          ebits[j] = (bits[j]/C >> BITRES)-1;
217       /* Make sure the first bit is spent on fine energy */
218       if (ebits[j] < 1)
219          ebits[j] = 1;
220
221       /* Make sure not to bust */
222       if (C*ebits[j] > (bits[j]>>BITRES))
223          ebits[j] = bits[j]/C >> BITRES;
224
225       /* More than that is useless because that's about as far as PVQ can go */
226       if (ebits[j]>7)
227          ebits[j]=7;
228
229       /* The other bits are assigned to PVQ */
230       bits[j] -= C*ebits[j]<<BITRES;
231       if (bits[j] < 0)
232          bits[j] = 0;
233    }
234    RESTORE_STACK;
235 }
236
237 void compute_allocation(const CELTMode *m, int start, int end, int *offsets, int total, int *pulses, int *ebits, int *fine_priority, int _C, int LM)
238 {
239    int lo, hi, len, j;
240    const int C = CHANNELS(_C);
241    VARDECL(int, bits1);
242    VARDECL(int, bits2);
243    SAVE_STACK;
244    
245    len = m->nbEBands;
246    ALLOC(bits1, len, int);
247    ALLOC(bits2, len, int);
248
249    lo = 0;
250    hi = m->nbAllocVectors - 1;
251    while (hi-lo != 1)
252    {
253       int psum = 0;
254       int mid = (lo+hi) >> 1;
255       for (j=start;j<end;j++)
256       {
257          int N = m->eBands[j+1]-m->eBands[j];
258          bits1[j] = ((C*N*m->allocVectors[mid*len+j]<<LM>>2) + offsets[j]);
259          psum += bits1[j];
260          /*printf ("%d ", bits[j]);*/
261       }
262       /*printf ("\n");*/
263       if (psum > (total<<BITRES))
264          hi = mid;
265       else
266          lo = mid;
267       /*printf ("lo = %d, hi = %d\n", lo, hi);*/
268    }
269    /*printf ("interp between %d and %d\n", lo, hi);*/
270    for (j=start;j<end;j++)
271    {
272       int N = m->eBands[j+1]-m->eBands[j];
273       bits1[j] = (C*N*m->allocVectors[lo*len+j]<<LM>>2) + offsets[j];
274       bits2[j] = (C*N*m->allocVectors[hi*len+j]<<LM>>2) + offsets[j];
275    }
276    interp_bits2pulses(m, start, end, bits1, bits2, total, pulses, ebits, fine_priority, len, C, LM);
277    RESTORE_STACK;
278 }
279