Merged compute_allocation_table() and compute_energy_allocation_table()
[opus.git] / libcelt / modes.c
1 /* (C) 2007-2008 Jean-Marc Valin, CSIRO
2 */
3 /*
4    Redistribution and use in source and binary forms, with or without
5    modification, are permitted provided that the following conditions
6    are met:
7    
8    - Redistributions of source code must retain the above copyright
9    notice, this list of conditions and the following disclaimer.
10    
11    - Redistributions in binary form must reproduce the above copyright
12    notice, this list of conditions and the following disclaimer in the
13    documentation and/or other materials provided with the distribution.
14    
15    - Neither the name of the Xiph.org Foundation nor the names of its
16    contributors may be used to endorse or promote products derived from
17    this software without specific prior written permission.
18    
19    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
20    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
21    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
22    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
23    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
24    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
25    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
26    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
27    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
28    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
29    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30 */
31
32 #ifdef HAVE_CONFIG_H
33 #include "config.h"
34 #endif
35
36 #include "celt.h"
37 #include "modes.h"
38 #include "rate.h"
39 #include "os_support.h"
40 #include "stack_alloc.h"
41 #include "quant_bands.h"
42
43 #ifdef STATIC_MODES
44 #include "static_modes.c"
45 #endif
46
47 #define MODEVALID 0xa110ca7e
48 #define MODEFREED 0xb10cf8ee
49
50 #ifndef M_PI
51 #define M_PI 3.141592653
52 #endif
53
54
55 int celt_mode_info(const CELTMode *mode, int request, celt_int32_t *value)
56 {
57    switch (request)
58    {
59       case CELT_GET_FRAME_SIZE:
60          *value = mode->mdctSize;
61          break;
62       case CELT_GET_LOOKAHEAD:
63          *value = mode->overlap;
64          break;
65       case CELT_GET_NB_CHANNELS:
66          *value = mode->nbChannels;
67          break;
68       case CELT_GET_BITSTREAM_VERSION:
69          *value = CELT_BITSTREAM_VERSION;
70          break;
71       default:
72          return CELT_BAD_ARG;
73    }
74    return CELT_OK;
75 }
76
77 #ifndef STATIC_MODES
78
79 #define PBANDS 8
80
81 #ifdef STDIN_TUNING
82 int MIN_BINS;
83 #else
84 #define MIN_BINS 3
85 #endif
86
87 /* Defining 25 critical bands for the full 0-20 kHz audio bandwidth
88    Taken from http://ccrma.stanford.edu/~jos/bbt/Bark_Frequency_Scale.html */
89 #define BARK_BANDS 25
90 static const celt_int16_t bark_freq[BARK_BANDS+1] = {
91       0,   100,   200,   300,   400,
92     510,   630,   770,   920,  1080,
93    1270,  1480,  1720,  2000,  2320,
94    2700,  3150,  3700,  4400,  5300,
95    6400,  7700,  9500, 12000, 15500,
96   20000};
97
98 static const celt_int16_t pitch_freq[PBANDS+1] ={0, 345, 689, 1034, 1378, 2067, 3273, 5340, 6374};
99
100 /* This allocation table is per critical band. When creating a mode, the bits get added together 
101    into the codec bands, which are sometimes larger than one critical band at low frequency */
102
103 #ifdef STDIN_TUNING
104 int BITALLOC_SIZE;
105 int *band_allocation;
106 #else
107 #define BITALLOC_SIZE 11
108 static const int band_allocation[BARK_BANDS*BITALLOC_SIZE] = 
109    {  2,  2,  1,  1,  2,  2,  1,  1,  1,  1,  1,  1,  1,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,
110       2,  2,  2,  1,  2,  2,  2,  2,  2,  1,  2,  2,  4,  5,  7,  7,  7,  5,  4,  0,  0,  0,  0,  0,  0,
111       2,  2,  2,  2,  3,  2,  2,  2,  2,  2,  3,  3,  5,  6,  8,  8,  8,  6,  5,  4,  0,  0,  0,  0,  0,
112       3,  2,  2,  2,  3,  3,  2,  3,  2,  3,  4,  4,  6,  7,  9,  9,  9,  7,  6,  5,  5,  5,  0,  0,  0,
113       3,  3,  2,  2,  3,  3,  3,  3,  3,  4,  4,  5,  7,  9, 10, 10, 10,  9,  6,  5,  5,  5,  5,  1,  0,
114       4,  3,  3,  3,  3,  3,  3,  3,  4,  4,  6,  7,  7,  9, 11, 10, 10,  9,  9,  8, 11, 10, 10,  1,  1,
115       5,  5,  4,  4,  5,  5,  5,  5,  6,  6,  8,  8, 10, 12, 15, 15, 13, 12, 12, 12, 18, 18, 16, 10,  1,
116       6,  6,  6,  6,  6,  6,  7,  7,  9,  9, 11, 12, 13, 18, 22, 23, 24, 25, 28, 30, 35, 35, 35, 35, 15,
117       7,  7,  7,  7,  7,  7, 10, 10, 10, 13, 14, 18, 20, 24, 28, 32, 32, 35, 38, 38, 42, 50, 59, 54, 31,
118       8,  8,  8,  8,  8,  9, 10, 12, 14, 20, 22, 25, 28, 30, 35, 42, 46, 50, 55, 60, 62, 62, 62, 62, 62,
119      12, 12, 12, 12, 12, 13, 15, 18, 22, 30, 32, 35, 40, 45, 55, 62, 66, 70, 85, 90, 92, 92, 92, 92, 92,
120    };
121 #endif
122
123 static celt_int16_t *compute_ebands(celt_int32_t Fs, int frame_size, int *nbEBands)
124 {
125    celt_int16_t *eBands;
126    int i, res, min_width, lin, low, high;
127    res = (Fs+frame_size)/(2*frame_size);
128    min_width = MIN_BINS*res;
129    /*printf ("min_width = %d\n", min_width);*/
130
131    /* Find where the linear part ends (i.e. where the spacing is more than min_width */
132    for (lin=0;lin<BARK_BANDS;lin++)
133       if (bark_freq[lin+1]-bark_freq[lin] >= min_width)
134          break;
135    
136    /*printf ("lin = %d (%d Hz)\n", lin, bark_freq[lin]);*/
137    low = ((bark_freq[lin]/res)+(MIN_BINS-1))/MIN_BINS;
138    high = BARK_BANDS-lin;
139    *nbEBands = low+high;
140    eBands = celt_alloc(sizeof(celt_int16_t)*(*nbEBands+2));
141    
142    /* Linear spacing (min_width) */
143    for (i=0;i<low;i++)
144       eBands[i] = MIN_BINS*i;
145    /* Spacing follows critical bands */
146    for (i=0;i<high;i++)
147       eBands[i+low] = (bark_freq[lin+i]+res/2)/res;
148    /* Enforce the minimum spacing at the boundary */
149    for (i=0;i<*nbEBands;i++)
150       if (eBands[i] < MIN_BINS*i)
151          eBands[i] = MIN_BINS*i;
152    eBands[*nbEBands] = (bark_freq[BARK_BANDS]+res/2)/res;
153    eBands[*nbEBands+1] = frame_size;
154    if (eBands[*nbEBands] > eBands[*nbEBands+1])
155       eBands[*nbEBands] = eBands[*nbEBands+1];
156    
157    /* FIXME: Remove last band if too small */
158    /*for (i=0;i<*nbEBands+2;i++)
159       printf("%d ", eBands[i]);
160    printf ("\n");*/
161    return eBands;
162 }
163
164 static void compute_pbands(CELTMode *mode, int res)
165 {
166    int i;
167    celt_int16_t *pBands;
168    pBands=celt_alloc(sizeof(celt_int16_t)*(PBANDS+2));
169    mode->nbPBands = PBANDS;
170    for (i=0;i<PBANDS+1;i++)
171    {
172       pBands[i] = (pitch_freq[i]+res/2)/res;
173       if (pBands[i] < mode->eBands[i])
174          pBands[i] = mode->eBands[i];
175    }
176    pBands[PBANDS+1] = mode->eBands[mode->nbEBands+1];
177    for (i=1;i<mode->nbPBands+1;i++)
178    {
179       int j;
180       for (j=0;j<mode->nbEBands;j++)
181          if (mode->eBands[j] <= pBands[i] && mode->eBands[j+1] > pBands[i])
182             break;
183       /*printf ("%d %d\n", i, j);*/
184       if (mode->eBands[j] != pBands[i])
185       {
186          if (pBands[i]-mode->eBands[j] < mode->eBands[j+1]-pBands[i] && 
187              mode->eBands[j] != pBands[i-1])
188             pBands[i] = mode->eBands[j];
189          else
190             pBands[i] = mode->eBands[j+1];
191       }
192    }
193    /*for (i=0;i<mode->nbPBands+2;i++)
194       printf("%d ", pBands[i]);
195    printf ("\n");*/
196    mode->pBands = pBands;
197    mode->pitchEnd = pBands[PBANDS];
198 }
199
200 static void compute_allocation_table(CELTMode *mode, int res)
201 {
202    int i, j, eband;
203    celt_int16_t *allocVectors, *allocEnergy;
204    const int C = CHANNELS(mode);
205
206    mode->nbAllocVectors = BITALLOC_SIZE;
207    allocVectors = celt_alloc(sizeof(celt_int16_t)*(BITALLOC_SIZE*mode->nbEBands));
208    allocEnergy = celt_alloc(sizeof(celt_int16_t)*(mode->nbAllocVectors*(mode->nbEBands+1)));
209    /* Compute per-codec-band allocation from per-critical-band matrix */
210    for (i=0;i<BITALLOC_SIZE;i++)
211    {
212       eband = 0;
213       for (j=0;j<BARK_BANDS;j++)
214       {
215          int edge, low;
216          celt_int32_t alloc;
217          edge = mode->eBands[eband+1]*res;
218          alloc = band_allocation[i*BARK_BANDS+j];
219          alloc = alloc*C*mode->mdctSize/256;
220          if (edge < bark_freq[j+1])
221          {
222             int num, den;
223             num = alloc * (edge-bark_freq[j]);
224             den = bark_freq[j+1]-bark_freq[j];
225             low = (num+den/2)/den;
226             allocVectors[i*mode->nbEBands+eband] += low;
227             eband++;
228             allocVectors[i*mode->nbEBands+eband] += alloc-low;
229          } else {
230             allocVectors[i*mode->nbEBands+eband] += alloc;
231          }
232       }
233    }
234    /* Compute fine energy resolution and update the pulse allocation table to subtract that */
235    for (i=0;i<mode->nbAllocVectors;i++)
236    {
237       int sum = 0;
238       int min_bits = 1;
239       if (allocVectors[i*mode->nbEBands]>12)
240          min_bits = 2;
241       if (allocVectors[i*mode->nbEBands]>24)
242          min_bits = 3;
243       for (j=0;j<mode->nbEBands;j++)
244       {
245          allocEnergy[i*(mode->nbEBands+1)+j] = allocVectors[i*mode->nbEBands+j]
246                                          / (C*(mode->eBands[j+1]-mode->eBands[j]));
247          if (allocEnergy[i*(mode->nbEBands+1)+j]<min_bits)
248             allocEnergy[i*(mode->nbEBands+1)+j] = min_bits;
249          if (allocEnergy[i*(mode->nbEBands+1)+j]>7)
250             allocEnergy[i*(mode->nbEBands+1)+j] = 7;
251          /* The bits used for fine allocation can't be used for pulses */
252          allocVectors[i*mode->nbEBands+j] -= C*allocEnergy[i*(mode->nbEBands+1)+j];
253          sum += allocEnergy[i*(mode->nbEBands+1)+j];
254       }
255       allocEnergy[i*(mode->nbEBands+1)+mode->nbEBands] = sum;
256    }
257    mode->energy_alloc = allocEnergy;
258    mode->allocVectors = allocVectors;
259 }
260
261 #endif /* STATIC_MODES */
262
263 CELTMode *celt_mode_create(celt_int32_t Fs, int channels, int frame_size, int *error)
264 {
265    int i;
266 #ifdef STDIN_TUNING
267    scanf("%d ", &MIN_BINS);
268    scanf("%d ", &BITALLOC_SIZE);
269    band_allocation = celt_alloc(sizeof(int)*BARK_BANDS*BITALLOC_SIZE);
270    for (i=0;i<BARK_BANDS*BITALLOC_SIZE;i++)
271    {
272       scanf("%d ", band_allocation+i);
273    }
274 #endif
275 #ifdef STATIC_MODES
276    const CELTMode *m = NULL;
277    CELTMode *mode=NULL;
278    ALLOC_STACK;
279    for (i=0;i<TOTAL_MODES;i++)
280    {
281       if (Fs == static_mode_list[i]->Fs &&
282           channels == static_mode_list[i]->nbChannels &&
283           frame_size == static_mode_list[i]->mdctSize &&
284           lookahead == static_mode_list[i]->overlap)
285       {
286          m = static_mode_list[i];
287          break;
288       }
289    }
290    if (m == NULL)
291    {
292       celt_warning("Mode not included as part of the static modes");
293       if (error)
294          *error = CELT_BAD_ARG;
295       return NULL;
296    }
297    mode = (CELTMode*)celt_alloc(sizeof(CELTMode));
298    CELT_COPY(mode, m, 1);
299 #else
300    int res;
301    CELTMode *mode;
302    celt_word16_t *window;
303    ALLOC_STACK;
304
305    /* The good thing here is that permutation of the arguments will automatically be invalid */
306    
307    if (Fs < 32000 || Fs > 64000)
308    {
309       celt_warning("Sampling rate must be between 32 kHz and 64 kHz");
310       if (error)
311          *error = CELT_BAD_ARG;
312       return NULL;
313    }
314    if (channels < 0 || channels > 2)
315    {
316       celt_warning("Only mono and stereo supported");
317       if (error)
318          *error = CELT_BAD_ARG;
319       return NULL;
320    }
321    if (frame_size < 64 || frame_size > 512 || frame_size%2!=0)
322    {
323       celt_warning("Only even frame sizes between 64 and 512 are supported");
324       if (error)
325          *error = CELT_BAD_ARG;
326       return NULL;
327    }
328    res = (Fs+frame_size)/(2*frame_size);
329    
330    mode = celt_alloc(sizeof(CELTMode));
331    mode->Fs = Fs;
332    mode->mdctSize = frame_size;
333    mode->nbChannels = channels;
334    mode->eBands = compute_ebands(Fs, frame_size, &mode->nbEBands);
335    compute_pbands(mode, res);
336    mode->ePredCoef = QCONST16(.8f,15);
337    
338    if (frame_size <= 64)
339    {
340       mode->nbShortMdcts = 1;
341    } else if (frame_size <= 256)
342    {
343       mode->nbShortMdcts = 2;
344    } else if (frame_size <= 384)
345    {
346       mode->nbShortMdcts = 3;
347    } else {
348       mode->nbShortMdcts = 4;
349    }
350    if (mode->nbShortMdcts > 1)
351       mode->overlap = frame_size/mode->nbShortMdcts;
352    else
353       mode->overlap = frame_size/2;
354    
355    compute_allocation_table(mode, res);
356    /*printf ("%d bands\n", mode->nbEBands);*/
357    
358    window = (celt_word16_t*)celt_alloc(mode->overlap*sizeof(celt_word16_t));
359
360 #ifndef FIXED_POINT
361    for (i=0;i<mode->overlap;i++)
362       window[i] = Q15ONE*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap));
363 #else
364    for (i=0;i<mode->overlap;i++)
365       window[i] = MIN32(32767,32768.*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap)));
366 #endif
367    mode->window = window;
368
369    mode->bits = (const celt_int16_t **)compute_alloc_cache(mode, 1);
370
371    mode->bits_stereo = NULL;
372 #ifndef SHORTCUTS
373    psydecay_init(&mode->psy, MAX_PERIOD/2, mode->Fs);
374 #endif
375    
376    mode->marker_start = MODEVALID;
377    mode->marker_end = MODEVALID;
378 #endif /* !STATIC_MODES */
379    mdct_init(&mode->mdct, 2*mode->mdctSize);
380    mode->fft = pitch_state_alloc(MAX_PERIOD);
381
382    mode->shortMdctSize = mode->mdctSize/mode->nbShortMdcts;
383    mdct_init(&mode->shortMdct, 2*mode->shortMdctSize);
384    mode->shortWindow = mode->window;
385
386    mode->prob = quant_prob_alloc(mode);
387    
388    if (mode->nbChannels>=2)
389       mode->bits_stereo = (const celt_int16_t **)compute_alloc_cache(mode, mode->nbChannels);
390
391    if (error)
392       *error = CELT_OK;
393    return mode;
394 }
395
396 void celt_mode_destroy(CELTMode *mode)
397 {
398 #ifndef STATIC_MODES
399    int i;
400    const celt_int16_t *prevPtr = NULL;
401    for (i=0;i<mode->nbEBands;i++)
402    {
403       if (mode->bits[i] != prevPtr)
404       {
405          prevPtr = mode->bits[i];
406          celt_free((int*)mode->bits[i]);
407       }
408    }
409    celt_free((int**)mode->bits);
410    if (check_mode(mode) != CELT_OK)
411       return;
412    celt_free((int*)mode->eBands);
413    celt_free((int*)mode->pBands);
414    celt_free((int*)mode->allocVectors);
415    
416    celt_free((celt_word16_t*)mode->window);
417
418    mode->marker_start = MODEFREED;
419    mode->marker_end = MODEFREED;
420 #ifndef SHORTCUTS
421    psydecay_clear(&mode->psy);
422 #endif
423 #endif
424    mdct_clear(&mode->mdct);
425    pitch_state_free(mode->fft);
426    quant_prob_free(mode->prob);
427    celt_free((celt_int16_t *)mode->energy_alloc);
428    celt_free((CELTMode *)mode);
429 }
430
431 int check_mode(const CELTMode *mode)
432 {
433    if (mode->marker_start == MODEVALID && mode->marker_end == MODEVALID)
434       return CELT_OK;
435    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
436       celt_warning("Using a mode that has already been freed");
437    else
438       celt_warning("This is not a valid CELT mode");
439    return CELT_INVALID_MODE;
440 }