Bunch of fixes for frames of 2.5 ms.
[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 #ifdef STDIN_TUNING
96 int BITALLOC_SIZE;
97 int *band_allocation;
98 #else
99 #define BITALLOC_SIZE 12
100 static const unsigned char band_allocation[BARK_BANDS*BITALLOC_SIZE] =
101    /* 0 100 200 300 400 510 630 770 920 1k  1.2 1.5 1.7 2k  2.3 2.7 3.1 3.7 4.4 5.3 6.4 7.7 9.5 12k 15k  */
102    {  4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, /*0*/
103       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, /*1*/
104       2,  2,  2,  1,  2,  2,  2,  2,  2,  2,  2,  2,  4,  5,  7,  7,  7,  5,  4,  0,  0,  0,  0,  0,  0, /*2*/
105       2,  2,  2,  2,  3,  3,  3,  3,  3,  3,  3,  3,  5,  6,  8,  8,  8,  6,  5,  4,  0,  0,  0,  0,  0, /*3*/
106       3,  2,  2,  2,  3,  4,  4,  4,  4,  4,  4,  4,  6,  7,  9,  9,  9,  7,  6,  5,  5,  5,  0,  0,  0, /*4*/
107       3,  3,  3,  4,  4,  5,  6,  6,  6,  6,  6,  7,  7,  9, 10, 10, 10,  9,  6,  5,  5,  5,  5,  1,  0, /*5*/
108       4,  3,  3,  4,  6,  7,  7,  7,  7,  7,  8,  9,  9,  9, 11, 10, 10,  9,  9,  8, 11, 10, 10,  1,  0, /*6*/
109       5,  5,  5,  6,  7,  7,  7,  7,  8,  8,  9, 10, 10, 12, 12, 11, 11, 17, 12, 15, 15, 20, 18, 10,  1, /*7*/
110       6,  7,  7,  7,  8,  8,  8,  8,  9, 10, 11, 12, 14, 17, 18, 21, 22, 27, 29, 39, 37, 38, 40, 35,  1, /*8*/
111       7,  7,  7,  8,  8,  8, 10, 10, 10, 13, 14, 18, 20, 24, 28, 32, 32, 35, 38, 38, 42, 50, 59, 54, 31, /*9*/
112       8,  8,  8,  8,  8,  9, 10, 12, 14, 20, 22, 25, 28, 30, 35, 42, 46, 50, 55, 60, 62, 62, 72, 82, 62, /*10*/
113       9,  9,  9, 10, 12, 13, 15, 18, 22, 30, 32, 35, 40, 45, 55, 62, 66, 70, 85, 90, 92, 92, 92,102, 92, /*11*/
114    };
115 #endif
116
117 static const celt_int16 eband5ms[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 10, 12, 14, 16, 20, 24, 28, 34, 40, 48, 60, 78, 100};
118
119 static celt_int16 *compute_ebands(celt_int32 Fs, int frame_size, int res, int *nbEBands)
120 {
121    celt_int16 *eBands;
122    int i, lin, low, high, nBark, offset=0;
123
124    if (Fs == 400*(celt_int32)frame_size && Fs >= 40000)
125    {
126       *nbEBands = sizeof(eband5ms)/sizeof(eband5ms[0])-1;
127       eBands = celt_alloc(sizeof(celt_int16)*(*nbEBands+2));
128       for (i=0;i<*nbEBands+2;i++)
129          eBands[i] = eband5ms[i];
130       eBands[*nbEBands+1] = frame_size;
131       return eBands;
132    }
133    /* Find the number of critical bands supported by our sampling rate */
134    for (nBark=1;nBark<BARK_BANDS;nBark++)
135     if (bark_freq[nBark+1]*2 >= Fs)
136        break;
137
138    /* Find where the linear part ends (i.e. where the spacing is more than min_width */
139    for (lin=0;lin<nBark;lin++)
140       if (bark_freq[lin+1]-bark_freq[lin] >= res)
141          break;
142
143    low = (bark_freq[lin]+res/2)/res;
144    high = nBark-lin;
145    *nbEBands = low+high;
146    eBands = celt_alloc(sizeof(celt_int16)*(*nbEBands+2));
147    
148    if (eBands==NULL)
149       return NULL;
150    
151    /* Linear spacing (min_width) */
152    for (i=0;i<low;i++)
153       eBands[i] = i;
154    if (low>0)
155       offset = eBands[low-1]*res - bark_freq[lin-1];
156    /* Spacing follows critical bands */
157    for (i=0;i<high;i++)
158    {
159       int target = bark_freq[lin+i];
160       eBands[i+low] = (target+(offset+res)/2)/res;
161       offset = eBands[i+low]*res - target;
162    }
163    /* Enforce the minimum spacing at the boundary */
164    for (i=0;i<*nbEBands;i++)
165       if (eBands[i] < i)
166          eBands[i] = i;
167    eBands[*nbEBands] = (bark_freq[nBark]+res/2)/res;
168    eBands[*nbEBands+1] = frame_size;
169    if (eBands[*nbEBands] > eBands[*nbEBands+1])
170       eBands[*nbEBands] = eBands[*nbEBands+1];
171    for (i=1;i<*nbEBands-1;i++)
172    {
173       if (eBands[i+1]-eBands[i] < eBands[i]-eBands[i-1])
174       {
175          eBands[i] -= (2*eBands[i]-eBands[i-1]-eBands[i+1])/2;
176       }
177    }
178    /*for (i=0;i<=*nbEBands+1;i++)
179       printf ("%d ", eBands[i]);
180    printf ("\n");
181    exit(1);*/
182    /* FIXME: Remove last band if too small */
183    return eBands;
184 }
185
186 static void compute_allocation_table(CELTMode *mode, int res)
187 {
188    int i, j, nBark;
189    celt_int16 *allocVectors;
190
191    /* Find the number of critical bands supported by our sampling rate */
192    for (nBark=1;nBark<BARK_BANDS;nBark++)
193     if (bark_freq[nBark+1]*2 >= mode->Fs)
194        break;
195
196    mode->nbAllocVectors = BITALLOC_SIZE;
197    allocVectors = celt_alloc(sizeof(celt_int16)*(BITALLOC_SIZE*mode->nbEBands));
198    if (allocVectors==NULL)
199       return;
200    /* Compute per-codec-band allocation from per-critical-band matrix */
201    for (i=0;i<BITALLOC_SIZE;i++)
202    {
203       celt_int32 current = 0;
204       int eband = 0;
205       for (j=0;j<nBark;j++)
206       {
207          int edge, low, high;
208          celt_int32 alloc;
209
210          alloc = mode->shortMdctSize*band_allocation[i*BARK_BANDS+j];
211          low = bark_freq[j];
212          high = bark_freq[j+1];
213
214          edge = mode->eBands[eband+1]*res;
215          while (edge <= high && eband < mode->nbEBands)
216          {
217             celt_int32 num;
218             int den, bits;
219             num = alloc * (edge-low);
220             den = high-low;
221             /* Divide with rounding */
222             bits = (2*num+den)/(2*den);
223             allocVectors[i*mode->nbEBands+eband] = (current+bits+128)>>(8-BITRES);
224             /* Remove one bit from every band -- FIXME: this is just a temporary hack*/
225             allocVectors[i*mode->nbEBands+eband] -= 1<<BITRES;
226             if (allocVectors[i*mode->nbEBands+eband]<0)
227                allocVectors[i*mode->nbEBands+eband]=0;
228             /* Remove the part of the band we just allocated */
229             low = edge;
230             alloc -= bits;
231
232             /* Move to next eband */
233             current = 0;
234             eband++;
235             edge = mode->eBands[eband+1]*res;
236          }
237          current += alloc;
238       }
239       if (eband < mode->nbEBands)
240       {
241          allocVectors[i*mode->nbEBands+eband] = (current+128)>>(8-BITRES);
242          /* Same hack as above FIXME: again */
243          allocVectors[i*mode->nbEBands+eband] -= 1<<BITRES;
244          if (allocVectors[i*mode->nbEBands+eband]<0)
245             allocVectors[i*mode->nbEBands+eband]=0;
246       }
247    }
248    mode->allocVectors = allocVectors;
249 }
250
251 #endif /* STATIC_MODES */
252
253 CELTMode *celt_mode_create(celt_int32 Fs, int frame_size, int *error)
254 {
255    int i;
256 #ifdef STDIN_TUNING
257    scanf("%d ", &MIN_BINS);
258    scanf("%d ", &BITALLOC_SIZE);
259    band_allocation = celt_alloc(sizeof(int)*BARK_BANDS*BITALLOC_SIZE);
260    for (i=0;i<BARK_BANDS*BITALLOC_SIZE;i++)
261    {
262       scanf("%d ", band_allocation+i);
263    }
264 #endif
265 #ifdef STATIC_MODES
266    const CELTMode *m = NULL;
267    CELTMode *mode=NULL;
268    ALLOC_STACK;
269 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
270    if (global_stack==NULL)
271    {
272       celt_free(global_stack);
273       goto failure;
274    }
275 #endif 
276    for (i=0;i<TOTAL_MODES;i++)
277    {
278       if (Fs == static_mode_list[i]->Fs &&
279           frame_size == static_mode_list[i]->mdctSize)
280       {
281          m = static_mode_list[i];
282          break;
283       }
284    }
285    if (m == NULL)
286    {
287       celt_warning("Mode not included as part of the static modes");
288       if (error)
289          *error = CELT_BAD_ARG;
290       return NULL;
291    }
292    mode = (CELTMode*)celt_alloc(sizeof(CELTMode));
293    if (mode==NULL)
294       goto failure;
295    CELT_COPY(mode, m, 1);
296    mode->marker_start = MODEPARTIAL;
297 #else
298    int res;
299    CELTMode *mode=NULL;
300    celt_word16 *window;
301    celt_int16 *logN;
302    ALLOC_STACK;
303 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
304    if (global_stack==NULL)
305    {
306       celt_free(global_stack);
307       goto failure;
308    }
309 #endif 
310
311    /* The good thing here is that permutation of the arguments will automatically be invalid */
312    
313    if (Fs < 32000 || Fs > 96000)
314    {
315       celt_warning("Sampling rate must be between 32 kHz and 96 kHz");
316       if (error)
317          *error = CELT_BAD_ARG;
318       return NULL;
319    }
320    if (frame_size < 64 || frame_size > 1024 || frame_size%2!=0)
321    {
322       celt_warning("Only even frame sizes from 64 to 1024 are supported");
323       if (error)
324          *error = CELT_BAD_ARG;
325       return NULL;
326    }
327    
328    mode = celt_alloc(sizeof(CELTMode));
329    if (mode==NULL)
330       goto failure;
331    mode->marker_start = MODEPARTIAL;
332    mode->Fs = Fs;
333    mode->mdctSize = frame_size;
334    mode->ePredCoef = QCONST16(.8f,15);
335
336    if (frame_size >= 640 && (frame_size%16)==0)
337    {
338      mode->nbShortMdcts = 8;
339    } else if (frame_size >= 320 && (frame_size%8)==0)
340    {
341      mode->nbShortMdcts = 4;
342    } else if (frame_size >= 160 && (frame_size%4)==0)
343    {
344      mode->nbShortMdcts = 2;
345    } else
346    {
347      mode->nbShortMdcts = 1;
348    }
349
350    mode->shortMdctSize = mode->mdctSize/mode->nbShortMdcts;
351    res = (mode->Fs+mode->shortMdctSize)/(2*mode->shortMdctSize);
352
353    mode->eBands = compute_ebands(Fs, mode->shortMdctSize, res, &mode->nbEBands);
354    if (mode->eBands==NULL)
355       goto failure;
356
357    mode->pitchEnd = 4000*(celt_int32)mode->shortMdctSize/Fs;
358    
359    /* Overlap must be divisible by 4 */
360    if (mode->nbShortMdcts > 1)
361       mode->overlap = (mode->shortMdctSize>>2)<<2;
362    else
363       mode->overlap = (frame_size>>3)<<2;
364
365
366    compute_allocation_table(mode, res);
367    if (mode->allocVectors==NULL)
368       goto failure;
369    
370    window = (celt_word16*)celt_alloc(mode->overlap*sizeof(celt_word16));
371    if (window==NULL)
372       goto failure;
373
374 #ifndef FIXED_POINT
375    for (i=0;i<mode->overlap;i++)
376       window[i] = Q15ONE*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap));
377 #else
378    for (i=0;i<mode->overlap;i++)
379       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))));
380 #endif
381    mode->window = window;
382
383    mode->bits = mode->_bits+1;
384    for (i=0;(1<<i)<=mode->nbShortMdcts;i++)
385    {
386       mode->bits[i] = (const celt_int16 **)compute_alloc_cache(mode, 1<<i);
387       if (mode->bits[i]==NULL)
388          goto failure;
389    }
390    mode->bits[-1] = (const celt_int16 **)compute_alloc_cache(mode, 0);
391    if (mode->bits[-1]==NULL)
392       goto failure;
393
394    logN = (celt_int16*)celt_alloc(mode->nbEBands*sizeof(celt_int16));
395    if (logN==NULL)
396       goto failure;
397
398    for (i=0;i<mode->nbEBands;i++)
399       logN[i] = log2_frac(mode->eBands[i+1]-mode->eBands[i], BITRES);
400    mode->logN = logN;
401 #endif /* !STATIC_MODES */
402
403    for (i=0;(1<<i)<=mode->nbShortMdcts;i++)
404    {
405       clt_mdct_init(&mode->mdct[i], 2*mode->shortMdctSize<<i);
406       if ((mode->mdct[i].trig==NULL)
407 #ifndef ENABLE_TI_DSPLIB55
408            || (mode->mdct[i].kfft==NULL)
409 #endif
410       )
411         goto failure;
412    }
413    mode->prob = quant_prob_alloc(mode);
414    if (mode->prob==NULL)
415      goto failure;
416
417    mode->marker_start = MODEVALID;
418    mode->marker_end   = MODEVALID;
419    if (error)
420       *error = CELT_OK;
421    return mode;
422 failure: 
423    if (error)
424       *error = CELT_INVALID_MODE;
425    if (mode!=NULL)
426       celt_mode_destroy(mode);
427    return NULL;
428 }
429
430 void celt_mode_destroy(CELTMode *mode)
431 {
432    int i, m;
433    const celt_int16 *prevPtr = NULL;
434    if (mode == NULL)
435    {
436       celt_warning("NULL passed to celt_mode_destroy");
437       return;
438    }
439
440    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
441    {
442       celt_warning("Freeing a mode which has already been freed"); 
443       return;
444    }
445
446    if (mode->marker_start != MODEVALID && mode->marker_start != MODEPARTIAL)
447    {
448       celt_warning("This is not a valid CELT mode structure");
449       return;  
450    }
451    mode->marker_start = MODEFREED;
452 #ifndef STATIC_MODES
453    for (m=0;(1<<m)<=mode->nbShortMdcts;m++)
454    {
455       if (mode->bits[m]!=NULL)
456       {
457          for (i=0;i<mode->nbEBands;i++)
458          {
459             if (mode->bits[m][i] != prevPtr)
460             {
461                prevPtr = mode->bits[m][i];
462                celt_free((int*)mode->bits[m][i]);
463             }
464          }
465       }
466       celt_free((celt_int16**)mode->bits[m]);
467    }
468    if (mode->bits[-1]!=NULL)
469    {
470       for (i=0;i<mode->nbEBands;i++)
471       {
472          if (mode->bits[-1][i] != prevPtr)
473          {
474             prevPtr = mode->bits[-1][i];
475             celt_free((int*)mode->bits[-1][i]);
476          }
477       }
478    }
479    celt_free((celt_int16**)mode->bits[-1]);
480
481    celt_free((celt_int16*)mode->eBands);
482    celt_free((celt_int16*)mode->allocVectors);
483    
484    celt_free((celt_word16*)mode->window);
485    celt_free((celt_int16*)mode->logN);
486
487 #endif
488    for (i=0;(1<<i)<=mode->nbShortMdcts;i++)
489       clt_mdct_clear(&mode->mdct[i]);
490
491    quant_prob_free(mode->prob);
492    mode->marker_end = MODEFREED;
493    celt_free((CELTMode *)mode);
494 }
495
496 int check_mode(const CELTMode *mode)
497 {
498    if (mode==NULL)
499       return CELT_INVALID_MODE;
500    if (mode->marker_start == MODEVALID && mode->marker_end == MODEVALID)
501       return CELT_OK;
502    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
503       celt_warning("Using a mode that has already been freed");
504    else
505       celt_warning("This is not a valid CELT mode");
506    return CELT_INVALID_MODE;
507 }