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