All modes based on 2.5 ms short blocks now use the same definition
[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    /* All modes that have 2.5 ms short blocks use the same definition */
123    if (Fs == 400*(celt_int32)frame_size)
124    {
125       *nbEBands = sizeof(eband5ms)/sizeof(eband5ms[0])-1;
126       eBands = celt_alloc(sizeof(celt_int16)*(*nbEBands+2));
127       for (i=0;i<*nbEBands+2;i++)
128          eBands[i] = eband5ms[i];
129       eBands[*nbEBands+1] = frame_size;
130       return eBands;
131    }
132    /* Find the number of critical bands supported by our sampling rate */
133    for (nBark=1;nBark<BARK_BANDS;nBark++)
134     if (bark_freq[nBark+1]*2 >= Fs)
135        break;
136
137    /* Find where the linear part ends (i.e. where the spacing is more than min_width */
138    for (lin=0;lin<nBark;lin++)
139       if (bark_freq[lin+1]-bark_freq[lin] >= res)
140          break;
141
142    low = (bark_freq[lin]+res/2)/res;
143    high = nBark-lin;
144    *nbEBands = low+high;
145    eBands = celt_alloc(sizeof(celt_int16)*(*nbEBands+2));
146    
147    if (eBands==NULL)
148       return NULL;
149    
150    /* Linear spacing (min_width) */
151    for (i=0;i<low;i++)
152       eBands[i] = i;
153    if (low>0)
154       offset = eBands[low-1]*res - bark_freq[lin-1];
155    /* Spacing follows critical bands */
156    for (i=0;i<high;i++)
157    {
158       int target = bark_freq[lin+i];
159       eBands[i+low] = (target+(offset+res)/2)/res;
160       offset = eBands[i+low]*res - target;
161    }
162    /* Enforce the minimum spacing at the boundary */
163    for (i=0;i<*nbEBands;i++)
164       if (eBands[i] < i)
165          eBands[i] = i;
166    eBands[*nbEBands] = (bark_freq[nBark]+res/2)/res;
167    eBands[*nbEBands+1] = frame_size;
168    if (eBands[*nbEBands] > eBands[*nbEBands+1])
169       eBands[*nbEBands] = eBands[*nbEBands+1];
170    for (i=1;i<*nbEBands-1;i++)
171    {
172       if (eBands[i+1]-eBands[i] < eBands[i]-eBands[i-1])
173       {
174          eBands[i] -= (2*eBands[i]-eBands[i-1]-eBands[i+1])/2;
175       }
176    }
177    /*for (i=0;i<=*nbEBands+1;i++)
178       printf ("%d ", eBands[i]);
179    printf ("\n");
180    exit(1);*/
181    /* FIXME: Remove last band if too small */
182    return eBands;
183 }
184
185 static void compute_allocation_table(CELTMode *mode, int res)
186 {
187    int i, j;
188    unsigned char *allocVectors;
189    int maxBands = sizeof(eband5ms)/sizeof(eband5ms[0])-1;
190
191    mode->nbAllocVectors = BITALLOC_SIZE;
192    allocVectors = celt_alloc(sizeof(unsigned char)*(BITALLOC_SIZE*mode->nbEBands));
193    if (allocVectors==NULL)
194       return;
195
196    /* Check for standard mode */
197    if (mode->Fs == 400*(celt_int32)mode->shortMdctSize && mode->Fs >= 40000)
198    {
199       for (i=0;i<BITALLOC_SIZE*mode->nbEBands;i++)
200          allocVectors[i] = band_allocation[i];
201       mode->allocVectors = allocVectors;
202       return;
203    }
204    /* If not the standard mode, interpolate */
205
206    /* Compute per-codec-band allocation from per-critical-band matrix */
207    for (i=0;i<BITALLOC_SIZE;i++)
208    {
209       celt_int32 current = 0;
210       int eband = 0;
211       /* We may be looping over too many bands, but eband will stop being
212          incremented once we reach the last band */
213       for (j=0;j<maxBands;j++)
214       {
215          int edge, low, high;
216          celt_int32 alloc;
217          alloc = band_allocation[i*maxBands + j]*(mode->eBands[eband+1]-mode->eBands[eband])<<4;
218          low = eband5ms[j]*200;
219          high = eband5ms[j+1]*200;
220          edge = mode->eBands[eband+1]*res;
221          while (edge <= high && eband < mode->nbEBands)
222          {
223             celt_int32 num;
224             int den, bits;
225             int N = (mode->eBands[eband+1]-mode->eBands[eband]);
226             num = alloc * (edge-low);
227             den = high-low;
228             /* Divide with rounding */
229             bits = (2*num+den)/(2*den);
230             allocVectors[i*mode->nbEBands+eband] = (2*(current+bits)+(N<<4))/(2*N<<4);
231             /* Remove the part of the band we just allocated */
232             low = edge;
233             alloc -= bits;
234
235             /* Move to next eband */
236             current = 0;
237             eband++;
238             edge = mode->eBands[eband+1]*res;
239          }
240          current += alloc;
241       }
242       if (eband < mode->nbEBands)
243       {
244          int N = (mode->eBands[eband+1]-mode->eBands[eband]);
245          allocVectors[i*mode->nbEBands+eband] = (2*current+(N<<4))/(2*N<<4);
246       }
247    }
248    /*printf ("\n");
249    for (i=0;i<BITALLOC_SIZE;i++)
250    {
251       for (j=0;j<mode->nbEBands;j++)
252          printf ("%d ", allocVectors[i*mode->nbEBands+j]);
253       printf ("\n");
254    }
255    exit(0);*/
256
257    mode->allocVectors = allocVectors;
258 }
259
260 #endif /* STATIC_MODES */
261
262 CELTMode *celt_mode_create(celt_int32 Fs, int frame_size, int *error)
263 {
264    int i;
265    int LM;
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 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
280    if (global_stack==NULL)
281       goto failure;
282 #endif 
283    for (i=0;i<TOTAL_MODES;i++)
284    {
285       if (Fs == static_mode_list[i]->Fs &&
286           frame_size == static_mode_list[i]->shortMdctSize*static_mode_list[i]->nbShortMdcts)
287       {
288          m = static_mode_list[i];
289          break;
290       }
291    }
292    if (m == NULL)
293    {
294       celt_warning("Mode not included as part of the static modes");
295       if (error)
296          *error = CELT_BAD_ARG;
297       return NULL;
298    }
299    mode = (CELTMode*)celt_alloc(sizeof(CELTMode));
300    if (mode==NULL)
301       goto failure;
302    CELT_COPY(mode, m, 1);
303    mode->bits = mode->_bits+1;
304    mode->marker_start = MODEPARTIAL;
305 #else
306    int res;
307    CELTMode *mode=NULL;
308    celt_word16 *window;
309    celt_int16 *logN;
310    ALLOC_STACK;
311 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
312    if (global_stack==NULL)
313       goto failure;
314 #endif 
315
316    /* The good thing here is that permutation of the arguments will automatically be invalid */
317    
318    if (Fs < 8000 || Fs > 96000)
319    {
320       celt_warning("Sampling rate must be between 8 kHz and 96 kHz");
321       if (error)
322          *error = CELT_BAD_ARG;
323       return NULL;
324    }
325    if (frame_size < 40 || frame_size > 1024 || frame_size%2!=0)
326    {
327       celt_warning("Only even frame sizes from 40 to 1024 are supported");
328       if (error)
329          *error = CELT_BAD_ARG;
330       return NULL;
331    }
332    
333    mode = celt_alloc(sizeof(CELTMode));
334    if (mode==NULL)
335       goto failure;
336    mode->marker_start = MODEPARTIAL;
337    mode->Fs = Fs;
338    mode->ePredCoef = QCONST16(.8f,15);
339
340    /* Pre/de-emphasis depends on sampling rate. The "standard" pre-emphasis
341       is defined as A(z) = 1 - 0.85*z^-1 at 48 kHz. Other rates should
342       approximate that. */
343    if(Fs < 12000) /* 8 kHz */
344    {
345       mode->preemph[0] =  QCONST16(.35f, 15);
346       mode->preemph[1] = -QCONST16(.18f, 15);
347       mode->preemph[2] =  QCONST16(.272f, SIG_SHIFT);
348       mode->preemph[3] =  QCONST16(3.6765f, 13);
349    } else if(Fs < 24000) /* 16 kHz */
350    {
351       mode->preemph[0] =  QCONST16(.6f, 15);
352       mode->preemph[1] = -QCONST16(.18f, 15);
353       mode->preemph[2] =  QCONST16(.4425f, SIG_SHIFT);
354       mode->preemph[3] =  QCONST16(2.259887f, 13);
355    } else if(Fs < 40000) /* 32 kHz */
356    {
357       mode->preemph[0] =  QCONST16(.78f, 15);
358       mode->preemph[1] = -QCONST16(.1f, 15);
359       mode->preemph[2] =  QCONST16(.75f, SIG_SHIFT);
360       mode->preemph[3] =  QCONST16(1.33333333f, 13);
361    } else /* 48 kHz */
362    {
363       mode->preemph[0] =  QCONST16(.85f, 15);
364       mode->preemph[1] =  QCONST16(.0f, 15);
365       mode->preemph[2] =  QCONST16(1.f, SIG_SHIFT);
366       mode->preemph[3] =  QCONST16(1.f, 13);
367    }
368
369    if ((celt_int32)frame_size*75 >= Fs && (frame_size%16)==0)
370    {
371      LM = 3;
372    } else if ((celt_int32)frame_size*150 >= Fs && (frame_size%8)==0)
373    {
374      LM = 2;
375    } else if ((celt_int32)frame_size*300 >= Fs && (frame_size%4)==0)
376    {
377      LM = 1;
378    } else
379    {
380      LM = 0;
381    }
382
383    mode->maxLM = LM;
384    mode->nbShortMdcts = 1<<LM;
385    mode->shortMdctSize = frame_size/mode->nbShortMdcts;
386    res = (mode->Fs+mode->shortMdctSize)/(2*mode->shortMdctSize);
387
388    mode->eBands = compute_ebands(Fs, mode->shortMdctSize, res, &mode->nbEBands);
389    if (mode->eBands==NULL)
390       goto failure;
391
392    mode->pitchEnd = 4000*(celt_int32)mode->shortMdctSize/Fs;
393    
394    /* Overlap must be divisible by 4 */
395    if (mode->nbShortMdcts > 1)
396       mode->overlap = (mode->shortMdctSize>>2)<<2;
397    else
398       mode->overlap = (frame_size>>3)<<2;
399
400
401    compute_allocation_table(mode, res);
402    if (mode->allocVectors==NULL)
403       goto failure;
404    
405    window = (celt_word16*)celt_alloc(mode->overlap*sizeof(celt_word16));
406    if (window==NULL)
407       goto failure;
408
409 #ifndef FIXED_POINT
410    for (i=0;i<mode->overlap;i++)
411       window[i] = Q15ONE*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap));
412 #else
413    for (i=0;i<mode->overlap;i++)
414       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))));
415 #endif
416    mode->window = window;
417
418    mode->bits = mode->_bits+1;
419    for (i=0;(1<<i)<=mode->nbShortMdcts;i++)
420    {
421       mode->bits[i] = (const celt_int16 **)compute_alloc_cache(mode, 1<<i);
422       if (mode->bits[i]==NULL)
423          goto failure;
424    }
425    mode->bits[-1] = (const celt_int16 **)compute_alloc_cache(mode, 0);
426    if (mode->bits[-1]==NULL)
427       goto failure;
428
429    logN = (celt_int16*)celt_alloc(mode->nbEBands*sizeof(celt_int16));
430    if (logN==NULL)
431       goto failure;
432
433    for (i=0;i<mode->nbEBands;i++)
434       logN[i] = log2_frac(mode->eBands[i+1]-mode->eBands[i], BITRES);
435    mode->logN = logN;
436 #endif /* !STATIC_MODES */
437
438    clt_mdct_init(&mode->mdct, 2*mode->shortMdctSize*mode->nbShortMdcts, LM);
439    if ((mode->mdct.trig==NULL)
440 #ifndef ENABLE_TI_DSPLIB55
441          || (mode->mdct.kfft==NULL)
442 #endif
443    )
444       goto failure;
445
446    mode->prob = quant_prob_alloc(mode);
447    if (mode->prob==NULL)
448      goto failure;
449
450    mode->marker_start = MODEVALID;
451    mode->marker_end   = MODEVALID;
452    if (error)
453       *error = CELT_OK;
454    return mode;
455 failure: 
456    if (error)
457       *error = CELT_INVALID_MODE;
458    if (mode!=NULL)
459       celt_mode_destroy(mode);
460    return NULL;
461 }
462
463 void celt_mode_destroy(CELTMode *mode)
464 {
465    int i, m;
466    const celt_int16 *prevPtr = NULL;
467    if (mode == NULL)
468    {
469       celt_warning("NULL passed to celt_mode_destroy");
470       return;
471    }
472
473    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
474    {
475       celt_warning("Freeing a mode which has already been freed"); 
476       return;
477    }
478
479    if (mode->marker_start != MODEVALID && mode->marker_start != MODEPARTIAL)
480    {
481       celt_warning("This is not a valid CELT mode structure");
482       return;  
483    }
484    mode->marker_start = MODEFREED;
485 #ifndef STATIC_MODES
486    for (m=0;(1<<m)<=mode->nbShortMdcts;m++)
487    {
488       if (mode->bits[m]!=NULL)
489       {
490          for (i=0;i<mode->nbEBands;i++)
491          {
492             if (mode->bits[m][i] != prevPtr)
493             {
494                prevPtr = mode->bits[m][i];
495                celt_free((int*)mode->bits[m][i]);
496             }
497          }
498       }
499       celt_free((celt_int16**)mode->bits[m]);
500    }
501    if (mode->bits[-1]!=NULL)
502    {
503       for (i=0;i<mode->nbEBands;i++)
504       {
505          if (mode->bits[-1][i] != prevPtr)
506          {
507             prevPtr = mode->bits[-1][i];
508             celt_free((int*)mode->bits[-1][i]);
509          }
510       }
511    }
512    celt_free((celt_int16**)mode->bits[-1]);
513
514    celt_free((celt_int16*)mode->eBands);
515    celt_free((celt_int16*)mode->allocVectors);
516    
517    celt_free((celt_word16*)mode->window);
518    celt_free((celt_int16*)mode->logN);
519
520 #endif
521    clt_mdct_clear(&mode->mdct);
522
523    quant_prob_free(mode->prob);
524    mode->marker_end = MODEFREED;
525    celt_free((CELTMode *)mode);
526 }
527
528 int check_mode(const CELTMode *mode)
529 {
530    if (mode==NULL)
531       return CELT_INVALID_MODE;
532    if (mode->marker_start == MODEVALID && mode->marker_end == MODEVALID)
533       return CELT_OK;
534    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
535       celt_warning("Using a mode that has already been freed");
536    else
537       celt_warning("This is not a valid CELT mode");
538    return CELT_INVALID_MODE;
539 }