Bands are now (again) a multiple of the number of short MDCTs
[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_FRAME_SIZE:
65          *value = mode->mdctSize;
66          break;
67       case CELT_GET_LOOKAHEAD:
68          *value = mode->overlap;
69          break;
70       case CELT_GET_BITSTREAM_VERSION:
71          *value = CELT_BITSTREAM_VERSION;
72          break;
73       case CELT_GET_SAMPLE_RATE:
74          *value = mode->Fs;
75          break;
76       default:
77          return CELT_UNIMPLEMENTED;
78    }
79    return CELT_OK;
80 }
81
82 #ifndef STATIC_MODES
83
84 /* Defining 25 critical bands for the full 0-20 kHz audio bandwidth
85    Taken from http://ccrma.stanford.edu/~jos/bbt/Bark_Frequency_Scale.html */
86 #define BARK_BANDS 25
87 static const celt_int16 bark_freq[BARK_BANDS+1] = {
88       0,   100,   200,   300,   400,
89     510,   630,   770,   920,  1080,
90    1270,  1480,  1720,  2000,  2320,
91    2700,  3150,  3700,  4400,  5300,
92    6400,  7700,  9500, 12000, 15500,
93   20000};
94
95 /* This allocation table is per critical band. When creating a mode, the bits get added together 
96    into the codec bands, which are sometimes larger than one critical band at low frequency */
97
98 #ifdef STDIN_TUNING
99 int BITALLOC_SIZE;
100 int *band_allocation;
101 #else
102 #define BITALLOC_SIZE 12
103 static const int band_allocation[BARK_BANDS*BITALLOC_SIZE] = 
104    /* 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  */
105    {  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*/
106       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*/
107       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*/
108       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*/
109       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*/
110       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*/
111       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*/
112       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*/
113       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*/
114       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*/
115       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*/
116       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*/
117    };
118 #endif
119
120 static celt_int16 *compute_ebands(celt_int32 Fs, int frame_size, int nbShortMdcts, int *nbEBands)
121 {
122    int min_bins = 2;
123    celt_int16 *eBands;
124    int i, res, min_width, lin, low, high, nBark, offset=0;
125
126    if (min_bins < nbShortMdcts)
127       min_bins = nbShortMdcts;
128    res = (Fs+frame_size)/(2*frame_size);
129    min_width = min_bins*res;
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] >= min_width)
139          break;
140
141    low = (bark_freq[lin]+res*min_bins/2)/(res*min_bins);
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] = min_bins*i;
152    /* Spacing follows critical bands */
153    for (i=0;i<high;i++)
154    {
155       int target = bark_freq[lin+i];
156       eBands[i+low] = (target+(offset+res*nbShortMdcts)/2)/(res*nbShortMdcts)*nbShortMdcts;
157       offset = eBands[i+low]*res - target;
158    }
159    /* Enforce the minimum spacing at the boundary */
160    for (i=0;i<*nbEBands;i++)
161       if (eBands[i] < min_bins*i)
162          eBands[i] = min_bins*i;
163    eBands[*nbEBands] = (bark_freq[nBark]+res*nbShortMdcts/2)/(res*nbShortMdcts)*nbShortMdcts;
164    eBands[*nbEBands+1] = frame_size;
165    if (eBands[*nbEBands] > eBands[*nbEBands+1])
166       eBands[*nbEBands] = eBands[*nbEBands+1];
167    for (i=1;i<*nbEBands-1;i++)
168    {
169       if (eBands[i+1]-eBands[i] < eBands[i]-eBands[i-1])
170       {
171          eBands[i] -= (2*eBands[i]-eBands[i-1]-eBands[i+1]+nbShortMdcts)/(2*nbShortMdcts)*nbShortMdcts;
172       }
173    }
174    /*for (i=0;i<*nbEBands+1;i++)
175       printf ("%d ", eBands[i]);
176    printf ("\n");
177    exit(1);*/
178    /* FIXME: Remove last band if too small */
179    return eBands;
180 }
181
182 static void compute_allocation_table(CELTMode *mode, int res)
183 {
184    int i, j, nBark;
185    celt_int16 *allocVectors;
186
187    /* Find the number of critical bands supported by our sampling rate */
188    for (nBark=1;nBark<BARK_BANDS;nBark++)
189     if (bark_freq[nBark+1]*2 >= mode->Fs)
190        break;
191
192    mode->nbAllocVectors = BITALLOC_SIZE;
193    allocVectors = celt_alloc(sizeof(celt_int16)*(BITALLOC_SIZE*mode->nbEBands));
194    if (allocVectors==NULL)
195       return;
196    /* Compute per-codec-band allocation from per-critical-band matrix */
197    for (i=0;i<BITALLOC_SIZE;i++)
198    {
199       celt_int32 current = 0;
200       int eband = 0;
201       for (j=0;j<nBark;j++)
202       {
203          int edge, low, high;
204          celt_int32 alloc;
205
206          alloc = mode->mdctSize*band_allocation[i*BARK_BANDS+j];
207          low = bark_freq[j];
208          high = bark_freq[j+1];
209
210          edge = mode->eBands[eband+1]*res;
211          while (edge <= high && eband < mode->nbEBands)
212          {
213             celt_int32 num;
214             int den, bits;
215             num = alloc * (edge-low);
216             den = high-low;
217             /* Divide with rounding */
218             bits = (2*num+den)/(2*den);
219             allocVectors[i*mode->nbEBands+eband] = (current+bits+128)>>8;
220
221             /* Remove the part of the band we just allocated */
222             low = edge;
223             alloc -= bits;
224
225             /* Move to next eband */
226             current = 0;
227             eband++;
228             edge = mode->eBands[eband+1]*res;
229          }
230          current += alloc;
231       }
232       if (eband < mode->nbEBands)
233          allocVectors[i*mode->nbEBands+eband] = (current+128)>>8;
234    }
235    mode->allocVectors = allocVectors;
236 }
237
238 #endif /* STATIC_MODES */
239
240 CELTMode *celt_mode_create(celt_int32 Fs, int frame_size, int *error)
241 {
242    int i;
243 #ifdef STDIN_TUNING
244    scanf("%d ", &MIN_BINS);
245    scanf("%d ", &BITALLOC_SIZE);
246    band_allocation = celt_alloc(sizeof(int)*BARK_BANDS*BITALLOC_SIZE);
247    for (i=0;i<BARK_BANDS*BITALLOC_SIZE;i++)
248    {
249       scanf("%d ", band_allocation+i);
250    }
251 #endif
252 #ifdef STATIC_MODES
253    const CELTMode *m = NULL;
254    CELTMode *mode=NULL;
255    ALLOC_STACK;
256 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
257    if (global_stack==NULL)
258    {
259       celt_free(global_stack);
260       goto failure;
261    }
262 #endif 
263    for (i=0;i<TOTAL_MODES;i++)
264    {
265       if (Fs == static_mode_list[i]->Fs &&
266           frame_size == static_mode_list[i]->mdctSize)
267       {
268          m = static_mode_list[i];
269          break;
270       }
271    }
272    if (m == NULL)
273    {
274       celt_warning("Mode not included as part of the static modes");
275       if (error)
276          *error = CELT_BAD_ARG;
277       return NULL;
278    }
279    mode = (CELTMode*)celt_alloc(sizeof(CELTMode));
280    if (mode==NULL)
281       goto failure;
282    CELT_COPY(mode, m, 1);
283    mode->marker_start = MODEPARTIAL;
284 #else
285    int res;
286    CELTMode *mode=NULL;
287    celt_word16 *window;
288    celt_int16 *logN;
289    ALLOC_STACK;
290 #if !defined(VAR_ARRAYS) && !defined(USE_ALLOCA)
291    if (global_stack==NULL)
292    {
293       celt_free(global_stack);
294       goto failure;
295    }
296 #endif 
297
298    /* The good thing here is that permutation of the arguments will automatically be invalid */
299    
300    if (Fs < 32000 || Fs > 96000)
301    {
302       celt_warning("Sampling rate must be between 32 kHz and 96 kHz");
303       if (error)
304          *error = CELT_BAD_ARG;
305       return NULL;
306    }
307    if (frame_size < 64 || frame_size > 1024 || frame_size%2!=0)
308    {
309       celt_warning("Only even frame sizes from 64 to 1024 are supported");
310       if (error)
311          *error = CELT_BAD_ARG;
312       return NULL;
313    }
314    res = (Fs+frame_size)/(2*frame_size);
315    
316    mode = celt_alloc(sizeof(CELTMode));
317    if (mode==NULL)
318       goto failure;
319    mode->marker_start = MODEPARTIAL;
320    mode->Fs = Fs;
321    mode->mdctSize = frame_size;
322    mode->ePredCoef = QCONST16(.8f,15);
323
324    if (frame_size > 640 && (frame_size%16)==0)
325    {
326      mode->nbShortMdcts = 8;
327    } else if (frame_size > 384 && (frame_size%8)==0)
328    {
329      mode->nbShortMdcts = 4;
330    } else if (frame_size > 384 && (frame_size%10)==0)
331    {
332      mode->nbShortMdcts = 5;
333    } else if (frame_size > 256 && (frame_size%6)==0)
334    {
335      mode->nbShortMdcts = 3;
336    } else if (frame_size > 256 && (frame_size%8)==0)
337    {
338      mode->nbShortMdcts = 4;
339    } else if (frame_size > 64 && (frame_size%4)==0)
340    {
341      mode->nbShortMdcts = 2;
342    } else if (frame_size > 128 && (frame_size%6)==0)
343    {
344      mode->nbShortMdcts = 3;
345    } else
346    {
347      mode->nbShortMdcts = 1;
348    }
349
350    mode->eBands = compute_ebands(Fs, frame_size, mode->nbShortMdcts, &mode->nbEBands);
351    if (mode->eBands==NULL)
352       goto failure;
353
354    mode->pitchEnd = 4000*(celt_int32)frame_size/Fs;
355    
356    /* Overlap must be divisible by 4 */
357    if (mode->nbShortMdcts > 1)
358       mode->overlap = ((frame_size/mode->nbShortMdcts)>>2)<<2; 
359    else
360       mode->overlap = (frame_size>>3)<<2;
361
362    compute_allocation_table(mode, res);
363    if (mode->allocVectors==NULL)
364       goto failure;
365    
366    window = (celt_word16*)celt_alloc(mode->overlap*sizeof(celt_word16));
367    if (window==NULL)
368       goto failure;
369
370 #ifndef FIXED_POINT
371    for (i=0;i<mode->overlap;i++)
372       window[i] = Q15ONE*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap));
373 #else
374    for (i=0;i<mode->overlap;i++)
375       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))));
376 #endif
377    mode->window = window;
378
379    mode->bits = (const celt_int16 **)compute_alloc_cache(mode, 1);
380    if (mode->bits==NULL)
381       goto failure;
382
383    logN = (celt_int16*)celt_alloc(mode->nbEBands*sizeof(celt_int16));
384    if (logN==NULL)
385       goto failure;
386
387    for (i=0;i<mode->nbEBands;i++)
388       logN[i] = log2_frac(mode->eBands[i+1]-mode->eBands[i], BITRES);
389    mode->logN = logN;
390 #endif /* !STATIC_MODES */
391
392    clt_mdct_init(&mode->mdct, 2*mode->mdctSize);
393
394    mode->shortMdctSize = mode->mdctSize/mode->nbShortMdcts;
395    clt_mdct_init(&mode->shortMdct, 2*mode->shortMdctSize);
396    mode->shortWindow = mode->window;
397    mode->prob = quant_prob_alloc(mode);
398    if ((mode->mdct.trig==NULL) || (mode->shortMdct.trig==NULL)
399 #ifndef ENABLE_TI_DSPLIB55
400         || (mode->mdct.kfft==NULL) || (mode->shortMdct.kfft==NULL)
401 #endif
402         || (mode->prob==NULL))
403      goto failure;
404
405    mode->marker_start = MODEVALID;
406    mode->marker_end   = MODEVALID;
407    if (error)
408       *error = CELT_OK;
409    return mode;
410 failure: 
411    if (error)
412       *error = CELT_INVALID_MODE;
413    if (mode!=NULL)
414       celt_mode_destroy(mode);
415    return NULL;
416 }
417
418 void celt_mode_destroy(CELTMode *mode)
419 {
420    int i;
421    const celt_int16 *prevPtr = NULL;
422    if (mode == NULL)
423    {
424       celt_warning("NULL passed to celt_mode_destroy");
425       return;
426    }
427
428    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
429    {
430       celt_warning("Freeing a mode which has already been freed"); 
431       return;
432    }
433
434    if (mode->marker_start != MODEVALID && mode->marker_start != MODEPARTIAL)
435    {
436       celt_warning("This is not a valid CELT mode structure");
437       return;  
438    }
439    mode->marker_start = MODEFREED;
440 #ifndef STATIC_MODES
441    if (mode->bits!=NULL)
442    {
443       for (i=0;i<mode->nbEBands;i++)
444       {
445          if (mode->bits[i] != prevPtr)
446          {
447             prevPtr = mode->bits[i];
448             celt_free((int*)mode->bits[i]);
449           }
450       }
451    }   
452    celt_free((celt_int16**)mode->bits);
453    celt_free((celt_int16*)mode->eBands);
454    celt_free((celt_int16*)mode->allocVectors);
455    
456    celt_free((celt_word16*)mode->window);
457    celt_free((celt_int16*)mode->logN);
458
459 #endif
460    clt_mdct_clear(&mode->mdct);
461    clt_mdct_clear(&mode->shortMdct);
462    quant_prob_free(mode->prob);
463    mode->marker_end = MODEFREED;
464    celt_free((CELTMode *)mode);
465 }
466
467 int check_mode(const CELTMode *mode)
468 {
469    if (mode==NULL)
470       return CELT_INVALID_MODE;
471    if (mode->marker_start == MODEVALID && mode->marker_end == MODEVALID)
472       return CELT_OK;
473    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
474       celt_warning("Using a mode that has already been freed");
475    else
476       celt_warning("This is not a valid CELT mode");
477    return CELT_INVALID_MODE;
478 }