nothing to see here
[opus.git] / libcelt / modes.c
1 /* Copyright (c) 2007-2008 CSIRO
2    Copyright (c) 2007-2009 Xiph.Org Foundation
3    Copyright (c) 2008 Gregory Maxwell 
4    Written by Jean-Marc Valin and Gregory Maxwell */
5 /*
6    Redistribution and use in source and binary forms, with or without
7    modification, are permitted provided that the following conditions
8    are met:
9    
10    - Redistributions of source code must retain the above copyright
11    notice, this list of conditions and the following disclaimer.
12    
13    - Redistributions in binary form must reproduce the above copyright
14    notice, this list of conditions and the following disclaimer in the
15    documentation and/or other materials provided with the distribution.
16    
17    - Neither the name of the Xiph.org Foundation nor the names of its
18    contributors may be used to endorse or promote products derived from
19    this software without specific prior written permission.
20    
21    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
22    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
23    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
24    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
25    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
26    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
27    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
28    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
29    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
30    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
31    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32 */
33
34 #ifdef HAVE_CONFIG_H
35 #include "config.h"
36 #endif
37
38 #include "celt.h"
39 #include "modes.h"
40 #include "rate.h"
41 #include "os_support.h"
42 #include "stack_alloc.h"
43 #include "quant_bands.h"
44
45 #ifdef STATIC_MODES
46 #include "static_modes.c"
47 #endif
48
49 #define MODEVALID   0xa110ca7e
50 #define MODEPARTIAL 0x7eca10a1
51 #define MODEFREED   0xb10cf8ee
52
53 #ifndef M_PI
54 #define M_PI 3.141592653
55 #endif
56
57
58 int celt_mode_info(const CELTMode *mode, int request, celt_int32 *value)
59 {
60    if (check_mode(mode) != CELT_OK)
61       return CELT_INVALID_MODE;
62    switch (request)
63    {
64       case CELT_GET_LOOKAHEAD:
65          *value = mode->overlap;
66          break;
67       case CELT_GET_BITSTREAM_VERSION:
68          *value = CELT_BITSTREAM_VERSION;
69          break;
70       case CELT_GET_SAMPLE_RATE:
71          *value = mode->Fs;
72          break;
73       default:
74          return CELT_UNIMPLEMENTED;
75    }
76    return CELT_OK;
77 }
78
79 #ifndef STATIC_MODES
80
81 /* Defining 25 critical bands for the full 0-20 kHz audio bandwidth
82    Taken from http://ccrma.stanford.edu/~jos/bbt/Bark_Frequency_Scale.html */
83 #define BARK_BANDS 25
84 static const celt_int16 bark_freq[BARK_BANDS+1] = {
85       0,   100,   200,   300,   400,
86     510,   630,   770,   920,  1080,
87    1270,  1480,  1720,  2000,  2320,
88    2700,  3150,  3700,  4400,  5300,
89    6400,  7700,  9500, 12000, 15500,
90   20000};
91
92 /* This allocation table is per critical band. When creating a mode, the bits get added together 
93    into the codec bands, which are sometimes larger than one critical band at low frequency */
94
95 #define BITALLOC_SIZE 12
96
97 static const celt_int16 eband5ms[] = {
98        0,  1,  2,  3,  4,  5,  6,  7,  8, 10, 12, 14, 16, 20, 24, 28, 34, 40, 48, 60, 78, 100
99 };
100
101 static const unsigned char band_allocation[] = {
102     /* 0 200 400 600 800  1k 1.2 1.4 1.6  2k 2.4 2.8 3.2  4k 4.8 5.6 6.8  8k 9.6 12k 15.6 */
103       10,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,
104       10,  3,  8,  2,  0,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,
105       10,  6,  8,  6,  5,  4,  3,  2,  7, 10, 11,  9,  7,  3,  1,  0,  0,  0,  0,  0,  0,
106       10, 10, 14, 11, 10,  8,  6,  5, 10, 12, 13, 11,  8,  4,  2,  1,  0,  0,  0,  0,  0,
107       13, 10, 17, 16, 14, 12, 10,  8, 12, 14, 14, 12,  9,  5,  3,  2,  2,  1,  0,  0,  0,
108       17, 21, 23, 26, 24, 20, 17, 16, 17, 18, 16, 14, 11,  6,  3,  2,  2,  1,  1,  0,  0,
109       21, 21, 36, 32, 28, 24, 23, 23, 22, 18, 18, 14, 11,  7,  5,  5,  5,  3,  3,  0,  0,
110       31, 35, 40, 32, 30, 28, 26, 26, 25, 24, 19, 15, 15, 13,  9,  9,  8,  7,  5,  2,  0,
111       42, 46, 46, 37, 35, 34, 33, 32, 34, 35, 32, 31, 27, 24, 23, 23, 18, 14, 11,  7,  0,
112       46, 49, 46, 46, 42, 43, 44, 47, 50, 52, 51, 48, 39, 32, 27, 24, 22, 19, 17, 11,  5,
113       53, 53, 49, 48, 55, 66, 71, 71, 71, 65, 64, 64, 56, 47, 41, 37, 31, 24, 20, 16, 10,
114       60, 64, 74, 74, 87,103,106,102,101,100,101, 95, 80, 69, 63, 55, 47, 36, 26, 21, 15,
115 };
116
117 static celt_int16 *compute_ebands(celt_int32 Fs, int frame_size, int res, int *nbEBands)
118 {
119    celt_int16 *eBands;
120    int i, lin, low, high, nBark, offset=0;
121
122    if (Fs == 400*(celt_int32)frame_size && Fs >= 40000)
123    {
124       *nbEBands = sizeof(eband5ms)/sizeof(eband5ms[0])-1;
125       eBands = celt_alloc(sizeof(celt_int16)*(*nbEBands+2));
126       for (i=0;i<*nbEBands+2;i++)
127          eBands[i] = eband5ms[i];
128       eBands[*nbEBands+1] = frame_size;
129       return eBands;
130    }
131    /* Find the number of critical bands supported by our sampling rate */
132    for (nBark=1;nBark<BARK_BANDS;nBark++)
133     if (bark_freq[nBark+1]*2 >= Fs)
134        break;
135
136    /* Find where the linear part ends (i.e. where the spacing is more than min_width */
137    for (lin=0;lin<nBark;lin++)
138       if (bark_freq[lin+1]-bark_freq[lin] >= res)
139          break;
140
141    low = (bark_freq[lin]+res/2)/res;
142    high = nBark-lin;
143    *nbEBands = low+high;
144    eBands = celt_alloc(sizeof(celt_int16)*(*nbEBands+2));
145    
146    if (eBands==NULL)
147       return NULL;
148    
149    /* Linear spacing (min_width) */
150    for (i=0;i<low;i++)
151       eBands[i] = i;
152    if (low>0)
153       offset = eBands[low-1]*res - bark_freq[lin-1];
154    /* Spacing follows critical bands */
155    for (i=0;i<high;i++)
156    {
157       int target = bark_freq[lin+i];
158       eBands[i+low] = (target+(offset+res)/2)/res;
159       offset = eBands[i+low]*res - target;
160    }
161    /* Enforce the minimum spacing at the boundary */
162    for (i=0;i<*nbEBands;i++)
163       if (eBands[i] < i)
164          eBands[i] = i;
165    eBands[*nbEBands] = (bark_freq[nBark]+res/2)/res;
166    eBands[*nbEBands+1] = frame_size;
167    if (eBands[*nbEBands] > eBands[*nbEBands+1])
168       eBands[*nbEBands] = eBands[*nbEBands+1];
169    for (i=1;i<*nbEBands-1;i++)
170    {
171       if (eBands[i+1]-eBands[i] < eBands[i]-eBands[i-1])
172       {
173          eBands[i] -= (2*eBands[i]-eBands[i-1]-eBands[i+1])/2;
174       }
175    }
176    /*for (i=0;i<=*nbEBands+1;i++)
177       printf ("%d ", eBands[i]);
178    printf ("\n");
179    exit(1);*/
180    /* FIXME: Remove last band if too small */
181    return eBands;
182 }
183
184 static void compute_allocation_table(CELTMode *mode, int res)
185 {
186    int i, j, nBark;
187    unsigned char *allocVectors;
188    int maxBands = sizeof(eband5ms)/sizeof(eband5ms[0])-1;
189
190    mode->nbAllocVectors = BITALLOC_SIZE;
191    allocVectors = celt_alloc(sizeof(unsigned char)*(BITALLOC_SIZE*mode->nbEBands));
192    if (allocVectors==NULL)
193       return;
194
195    /* Check for standard mode */
196    if (mode->Fs == 400*(celt_int32)mode->shortMdctSize && mode->Fs >= 40000)
197    {
198       for (i=0;i<BITALLOC_SIZE*mode->nbEBands;i++)
199          allocVectors[i] = band_allocation[i];
200       mode->allocVectors = allocVectors;
201       return;
202    }
203
204    /* If not the standard mode, interpolate */
205
206    /* Find the number of critical bands supported by our sampling rate */
207    for (nBark=1;nBark<maxBands;nBark++)
208     if (eband5ms[j+1]*400 >= mode->Fs)
209        break;
210
211    /* Compute per-codec-band allocation from per-critical-band matrix */
212    for (i=0;i<BITALLOC_SIZE;i++)
213    {
214       celt_int32 current = 0;
215       int eband = 0;
216       for (j=0;j<nBark;j++)
217       {
218          int edge, low, high;
219          celt_int32 alloc;
220          alloc = band_allocation[i*maxBands + j]*(mode->eBands[eband+1]-mode->eBands[eband])<<4;
221          low = eband5ms[j]*200;
222          high = eband5ms[j+1]*200;
223          edge = mode->eBands[eband+1]*res;
224          while (edge <= high && eband < mode->nbEBands)
225          {
226             celt_int32 num;
227             int den, bits;
228             int N = (mode->eBands[eband+1]-mode->eBands[eband]);
229             num = alloc * (edge-low);
230             den = high-low;
231             /* Divide with rounding */
232             bits = (2*num+den)/(2*den);
233             allocVectors[i*mode->nbEBands+eband] = (2*(current+bits)+(N<<4))/(2*N<<4);
234             /* Remove the part of the band we just allocated */
235             low = edge;
236             alloc -= bits;
237
238             /* Move to next eband */
239             current = 0;
240             eband++;
241             edge = mode->eBands[eband+1]*res;
242          }
243          current += alloc;
244       }
245       if (eband < mode->nbEBands)
246       {
247          int N = (mode->eBands[eband+1]-mode->eBands[eband]);
248          allocVectors[i*mode->nbEBands+eband] = (2*current+(N<<4))/(2*N<<4);
249       }
250    }
251    /*printf ("\n");
252    for (i=0;i<BITALLOC_SIZE;i++)
253    {
254       for (j=0;j<mode->nbEBands;j++)
255          printf ("%d ", allocVectors[i*mode->nbEBands+j]);
256       printf ("\n");
257    }
258    exit(0);*/
259
260    mode->allocVectors = allocVectors;
261 }
262
263 #endif /* STATIC_MODES */
264
265 CELTMode *celt_mode_create(celt_int32 Fs, int frame_size, int *error)
266 {
267    int i;
268 #ifdef STDIN_TUNING
269    scanf("%d ", &MIN_BINS);
270    scanf("%d ", &BITALLOC_SIZE);
271    band_allocation = celt_alloc(sizeof(int)*BARK_BANDS*BITALLOC_SIZE);
272    for (i=0;i<BARK_BANDS*BITALLOC_SIZE;i++)
273    {
274       scanf("%d ", band_allocation+i);
275    }
276 #endif
277 #ifdef STATIC_MODES
278    const CELTMode *m = NULL;
279    CELTMode *mode=NULL;
280    ALLOC_STACK;
281 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
282    if (global_stack==NULL)
283    {
284       celt_free(global_stack);
285       goto failure;
286    }
287 #endif 
288    for (i=0;i<TOTAL_MODES;i++)
289    {
290       if (Fs == static_mode_list[i]->Fs &&
291           frame_size == static_mode_list[i]->mdctSize)
292       {
293          m = static_mode_list[i];
294          break;
295       }
296    }
297    if (m == NULL)
298    {
299       celt_warning("Mode not included as part of the static modes");
300       if (error)
301          *error = CELT_BAD_ARG;
302       return NULL;
303    }
304    mode = (CELTMode*)celt_alloc(sizeof(CELTMode));
305    if (mode==NULL)
306       goto failure;
307    CELT_COPY(mode, m, 1);
308    mode->marker_start = MODEPARTIAL;
309 #else
310    int res;
311    CELTMode *mode=NULL;
312    celt_word16 *window;
313    celt_int16 *logN;
314    ALLOC_STACK;
315 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
316    if (global_stack==NULL)
317    {
318       celt_free(global_stack);
319       goto failure;
320    }
321 #endif 
322
323    /* The good thing here is that permutation of the arguments will automatically be invalid */
324    
325    if (Fs < 32000 || Fs > 96000)
326    {
327       celt_warning("Sampling rate must be between 32 kHz and 96 kHz");
328       if (error)
329          *error = CELT_BAD_ARG;
330       return NULL;
331    }
332    if (frame_size < 64 || frame_size > 1024 || frame_size%2!=0)
333    {
334       celt_warning("Only even frame sizes from 64 to 1024 are supported");
335       if (error)
336          *error = CELT_BAD_ARG;
337       return NULL;
338    }
339    
340    mode = celt_alloc(sizeof(CELTMode));
341    if (mode==NULL)
342       goto failure;
343    mode->marker_start = MODEPARTIAL;
344    mode->Fs = Fs;
345    mode->mdctSize = frame_size;
346    mode->ePredCoef = QCONST16(.8f,15);
347
348    if (frame_size >= 640 && (frame_size%16)==0)
349    {
350      mode->nbShortMdcts = 8;
351    } else if (frame_size >= 320 && (frame_size%8)==0)
352    {
353      mode->nbShortMdcts = 4;
354    } else if (frame_size >= 160 && (frame_size%4)==0)
355    {
356      mode->nbShortMdcts = 2;
357    } else
358    {
359      mode->nbShortMdcts = 1;
360    }
361
362    mode->shortMdctSize = mode->mdctSize/mode->nbShortMdcts;
363    res = (mode->Fs+mode->shortMdctSize)/(2*mode->shortMdctSize);
364
365    mode->eBands = compute_ebands(Fs, mode->shortMdctSize, res, &mode->nbEBands);
366    if (mode->eBands==NULL)
367       goto failure;
368
369    mode->pitchEnd = 4000*(celt_int32)mode->shortMdctSize/Fs;
370    
371    /* Overlap must be divisible by 4 */
372    if (mode->nbShortMdcts > 1)
373       mode->overlap = (mode->shortMdctSize>>2)<<2;
374    else
375       mode->overlap = (frame_size>>3)<<2;
376
377
378    compute_allocation_table(mode, res);
379    if (mode->allocVectors==NULL)
380       goto failure;
381    
382    window = (celt_word16*)celt_alloc(mode->overlap*sizeof(celt_word16));
383    if (window==NULL)
384       goto failure;
385
386 #ifndef FIXED_POINT
387    for (i=0;i<mode->overlap;i++)
388       window[i] = Q15ONE*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap));
389 #else
390    for (i=0;i<mode->overlap;i++)
391       window[i] = MIN32(32767,floor(.5+32768.*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap))));
392 #endif
393    mode->window = window;
394
395    mode->bits = mode->_bits+1;
396    for (i=0;(1<<i)<=mode->nbShortMdcts;i++)
397    {
398       mode->bits[i] = (const celt_int16 **)compute_alloc_cache(mode, 1<<i);
399       if (mode->bits[i]==NULL)
400          goto failure;
401    }
402    mode->bits[-1] = (const celt_int16 **)compute_alloc_cache(mode, 0);
403    if (mode->bits[-1]==NULL)
404       goto failure;
405
406    logN = (celt_int16*)celt_alloc(mode->nbEBands*sizeof(celt_int16));
407    if (logN==NULL)
408       goto failure;
409
410    for (i=0;i<mode->nbEBands;i++)
411       logN[i] = log2_frac(mode->eBands[i+1]-mode->eBands[i], BITRES);
412    mode->logN = logN;
413 #endif /* !STATIC_MODES */
414
415    for (i=0;(1<<i)<=mode->nbShortMdcts;i++)
416    {
417       clt_mdct_init(&mode->mdct[i], 2*mode->shortMdctSize<<i);
418       if ((mode->mdct[i].trig==NULL)
419 #ifndef ENABLE_TI_DSPLIB55
420            || (mode->mdct[i].kfft==NULL)
421 #endif
422       )
423         goto failure;
424    }
425    mode->prob = quant_prob_alloc(mode);
426    if (mode->prob==NULL)
427      goto failure;
428
429    mode->marker_start = MODEVALID;
430    mode->marker_end   = MODEVALID;
431    if (error)
432       *error = CELT_OK;
433    return mode;
434 failure: 
435    if (error)
436       *error = CELT_INVALID_MODE;
437    if (mode!=NULL)
438       celt_mode_destroy(mode);
439    return NULL;
440 }
441
442 void celt_mode_destroy(CELTMode *mode)
443 {
444    int i, m;
445    const celt_int16 *prevPtr = NULL;
446    if (mode == NULL)
447    {
448       celt_warning("NULL passed to celt_mode_destroy");
449       return;
450    }
451
452    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
453    {
454       celt_warning("Freeing a mode which has already been freed"); 
455       return;
456    }
457
458    if (mode->marker_start != MODEVALID && mode->marker_start != MODEPARTIAL)
459    {
460       celt_warning("This is not a valid CELT mode structure");
461       return;  
462    }
463    mode->marker_start = MODEFREED;
464 #ifndef STATIC_MODES
465    for (m=0;(1<<m)<=mode->nbShortMdcts;m++)
466    {
467       if (mode->bits[m]!=NULL)
468       {
469          for (i=0;i<mode->nbEBands;i++)
470          {
471             if (mode->bits[m][i] != prevPtr)
472             {
473                prevPtr = mode->bits[m][i];
474                celt_free((int*)mode->bits[m][i]);
475             }
476          }
477       }
478       celt_free((celt_int16**)mode->bits[m]);
479    }
480    if (mode->bits[-1]!=NULL)
481    {
482       for (i=0;i<mode->nbEBands;i++)
483       {
484          if (mode->bits[-1][i] != prevPtr)
485          {
486             prevPtr = mode->bits[-1][i];
487             celt_free((int*)mode->bits[-1][i]);
488          }
489       }
490    }
491    celt_free((celt_int16**)mode->bits[-1]);
492
493    celt_free((celt_int16*)mode->eBands);
494    celt_free((celt_int16*)mode->allocVectors);
495    
496    celt_free((celt_word16*)mode->window);
497    celt_free((celt_int16*)mode->logN);
498
499 #endif
500    for (i=0;(1<<i)<=mode->nbShortMdcts;i++)
501       clt_mdct_clear(&mode->mdct[i]);
502
503    quant_prob_free(mode->prob);
504    mode->marker_end = MODEFREED;
505    celt_free((CELTMode *)mode);
506 }
507
508 int check_mode(const CELTMode *mode)
509 {
510    if (mode==NULL)
511       return CELT_INVALID_MODE;
512    if (mode->marker_start == MODEVALID && mode->marker_end == MODEVALID)
513       return CELT_OK;
514    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
515       celt_warning("Using a mode that has already been freed");
516    else
517       celt_warning("This is not a valid CELT mode");
518    return CELT_INVALID_MODE;
519 }