86543936a43ebd6dad5ace4cd575604831002164
[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;
204    
205    mode->nbAllocVectors = BITALLOC_SIZE;
206    allocVectors = celt_alloc(sizeof(celt_int16_t)*(BITALLOC_SIZE*mode->nbEBands));
207    for (i=0;i<BITALLOC_SIZE;i++)
208    {
209       eband = 0;
210       for (j=0;j<BARK_BANDS;j++)
211       {
212          int edge, low;
213          celt_int32_t alloc;
214          edge = mode->eBands[eband+1]*res;
215          alloc = band_allocation[i*BARK_BANDS+j];
216          if (mode->nbChannels == 2)
217             alloc = alloc*3*mode->mdctSize/512;
218          else
219             alloc = alloc*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    /*for (i=0;i<BITALLOC_SIZE;i++)
235    {
236       for (j=0;j<mode->nbEBands;j++)
237          printf ("%2d ", allocVectors[i*mode->nbEBands+j]);
238       printf ("\n");
239    }*/
240    mode->allocVectors = allocVectors;
241 }
242
243 #endif /* STATIC_MODES */
244
245 static void compute_energy_allocation_table(CELTMode *mode)
246 {
247    int i, j;
248    celt_int16_t *alloc;
249    const int C = CHANNELS(mode);
250    
251    alloc = celt_alloc(sizeof(celt_int16_t)*(mode->nbAllocVectors*(mode->nbEBands+1)));
252    for (i=0;i<mode->nbAllocVectors;i++)
253    {
254       int sum = 0;
255       int min_bits = 1;
256       if (mode->allocVectors[i*mode->nbEBands]>12)
257          min_bits = 2;
258       if (mode->allocVectors[i*mode->nbEBands]>24)
259          min_bits = 3;
260       for (j=0;j<mode->nbEBands;j++)
261       {
262          alloc[i*(mode->nbEBands+1)+j] = mode->allocVectors[i*mode->nbEBands+j]
263                                          / (C*(mode->eBands[j+1]-mode->eBands[j]));
264          if (alloc[i*(mode->nbEBands+1)+j]<min_bits)
265             alloc[i*(mode->nbEBands+1)+j] = min_bits;
266          if (alloc[i*(mode->nbEBands+1)+j]>7)
267             alloc[i*(mode->nbEBands+1)+j] = 7;
268          sum += alloc[i*(mode->nbEBands+1)+j];
269          /*printf ("%d ", alloc[i*(mode->nbEBands+1)+j]);*/
270          /*printf ("%f ", mode->allocVectors[i*mode->nbEBands+j]*1.f/(mode->eBands[j+1]-mode->eBands[j]-1));*/
271       }
272       alloc[i*(mode->nbEBands+1)+mode->nbEBands] = sum;
273       /*printf ("\n");*/
274    }
275    mode->energy_alloc = alloc;
276 }
277
278 CELTMode *celt_mode_create(celt_int32_t Fs, int channels, int frame_size, int *error)
279 {
280    int i;
281 #ifdef STDIN_TUNING
282    scanf("%d ", &MIN_BINS);
283    scanf("%d ", &BITALLOC_SIZE);
284    band_allocation = celt_alloc(sizeof(int)*BARK_BANDS*BITALLOC_SIZE);
285    for (i=0;i<BARK_BANDS*BITALLOC_SIZE;i++)
286    {
287       scanf("%d ", band_allocation+i);
288    }
289 #endif
290 #ifdef STATIC_MODES
291    const CELTMode *m = NULL;
292    CELTMode *mode=NULL;
293    ALLOC_STACK;
294    for (i=0;i<TOTAL_MODES;i++)
295    {
296       if (Fs == static_mode_list[i]->Fs &&
297           channels == static_mode_list[i]->nbChannels &&
298           frame_size == static_mode_list[i]->mdctSize &&
299           lookahead == static_mode_list[i]->overlap)
300       {
301          m = static_mode_list[i];
302          break;
303       }
304    }
305    if (m == NULL)
306    {
307       celt_warning("Mode not included as part of the static modes");
308       if (error)
309          *error = CELT_BAD_ARG;
310       return NULL;
311    }
312    mode = (CELTMode*)celt_alloc(sizeof(CELTMode));
313    CELT_COPY(mode, m, 1);
314 #else
315    int res;
316    CELTMode *mode;
317    celt_word16_t *window;
318    ALLOC_STACK;
319
320    /* The good thing here is that permutation of the arguments will automatically be invalid */
321    
322    if (Fs < 32000 || Fs > 64000)
323    {
324       celt_warning("Sampling rate must be between 32 kHz and 64 kHz");
325       if (error)
326          *error = CELT_BAD_ARG;
327       return NULL;
328    }
329    if (channels < 0 || channels > 2)
330    {
331       celt_warning("Only mono and stereo supported");
332       if (error)
333          *error = CELT_BAD_ARG;
334       return NULL;
335    }
336    if (frame_size < 64 || frame_size > 512 || frame_size%2!=0)
337    {
338       celt_warning("Only even frame sizes between 64 and 512 are supported");
339       if (error)
340          *error = CELT_BAD_ARG;
341       return NULL;
342    }
343    res = (Fs+frame_size)/(2*frame_size);
344    
345    mode = celt_alloc(sizeof(CELTMode));
346    mode->Fs = Fs;
347    mode->mdctSize = frame_size;
348    mode->nbChannels = channels;
349    mode->eBands = compute_ebands(Fs, frame_size, &mode->nbEBands);
350    compute_pbands(mode, res);
351    mode->ePredCoef = QCONST16(.8f,15);
352    
353    if (frame_size <= 64)
354    {
355       mode->nbShortMdcts = 1;
356    } else if (frame_size <= 256)
357    {
358       mode->nbShortMdcts = 2;
359    } else if (frame_size <= 384)
360    {
361       mode->nbShortMdcts = 3;
362    } else {
363       mode->nbShortMdcts = 4;
364    }
365    if (mode->nbShortMdcts > 1)
366       mode->overlap = frame_size/mode->nbShortMdcts;
367    else
368       mode->overlap = frame_size/2;
369    
370    compute_allocation_table(mode, res);
371    /*printf ("%d bands\n", mode->nbEBands);*/
372    
373    window = (celt_word16_t*)celt_alloc(mode->overlap*sizeof(celt_word16_t));
374
375 #ifndef FIXED_POINT
376    for (i=0;i<mode->overlap;i++)
377       window[i] = Q15ONE*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap));
378 #else
379    for (i=0;i<mode->overlap;i++)
380       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)));
381 #endif
382    mode->window = window;
383
384    mode->bits = (const celt_int16_t **)compute_alloc_cache(mode, 1);
385
386    mode->bits_stereo = NULL;
387 #ifndef SHORTCUTS
388    psydecay_init(&mode->psy, MAX_PERIOD/2, mode->Fs);
389 #endif
390    
391    mode->marker_start = MODEVALID;
392    mode->marker_end = MODEVALID;
393 #endif /* !STATIC_MODES */
394    mdct_init(&mode->mdct, 2*mode->mdctSize);
395    mode->fft = pitch_state_alloc(MAX_PERIOD);
396
397    mode->shortMdctSize = mode->mdctSize/mode->nbShortMdcts;
398    mdct_init(&mode->shortMdct, 2*mode->shortMdctSize);
399    mode->shortWindow = mode->window;
400
401    mode->prob = quant_prob_alloc(mode);
402    compute_energy_allocation_table(mode);
403    
404    if (mode->nbChannels>=2)
405       mode->bits_stereo = (const celt_int16_t **)compute_alloc_cache(mode, mode->nbChannels);
406
407    if (error)
408       *error = CELT_OK;
409    return mode;
410 }
411
412 void celt_mode_destroy(CELTMode *mode)
413 {
414 #ifndef STATIC_MODES
415    int i;
416    const celt_int16_t *prevPtr = NULL;
417    for (i=0;i<mode->nbEBands;i++)
418    {
419       if (mode->bits[i] != prevPtr)
420       {
421          prevPtr = mode->bits[i];
422          celt_free((int*)mode->bits[i]);
423       }
424    }
425    celt_free((int**)mode->bits);
426    if (check_mode(mode) != CELT_OK)
427       return;
428    celt_free((int*)mode->eBands);
429    celt_free((int*)mode->pBands);
430    celt_free((int*)mode->allocVectors);
431    
432    celt_free((celt_word16_t*)mode->window);
433
434    mode->marker_start = MODEFREED;
435    mode->marker_end = MODEFREED;
436 #ifndef SHORTCUTS
437    psydecay_clear(&mode->psy);
438 #endif
439 #endif
440    mdct_clear(&mode->mdct);
441    pitch_state_free(mode->fft);
442    quant_prob_free(mode->prob);
443    celt_free((celt_int16_t *)mode->energy_alloc);
444    celt_free((CELTMode *)mode);
445 }
446
447 int check_mode(const CELTMode *mode)
448 {
449    if (mode->marker_start == MODEVALID && mode->marker_end == MODEVALID)
450       return CELT_OK;
451    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
452       celt_warning("Using a mode that has already been freed");
453    else
454       celt_warning("This is not a valid CELT mode");
455    return CELT_INVALID_MODE;
456 }