Some tuning improvements, more to come
[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 EXPORT 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       default:
69          return CELT_BAD_ARG;
70    }
71    return CELT_OK;
72 }
73
74 #ifndef STATIC_MODES
75
76 #define PBANDS 8
77
78 #ifdef STDIN_TUNING
79 int MIN_BINS;
80 #else
81 #define MIN_BINS 3
82 #endif
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_t 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 static const celt_int16_t pitch_freq[PBANDS+1] ={0, 345, 689, 1034, 1378, 2067, 3273, 5340, 6374};
96
97 /* This allocation table is per critical band. When creating a mode, the bits get added together 
98    into the codec bands, which are sometimes larger than one critical band at low frequency */
99
100 #ifdef STDIN_TUNING
101 int BITALLOC_SIZE;
102 int *band_allocation;
103 #else
104 #define BITALLOC_SIZE 10
105 static const int band_allocation[BARK_BANDS*BITALLOC_SIZE] = 
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,
107       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,
108       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,
109       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,
110       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,
111       4,  4,  3,  3,  3,  3,  3,  3,  3,  4,  5,  6,  7,  9, 10, 10, 10,  9,  9,  8, 11, 10, 10,  1,  1,
112       6,  6,  5,  5,  5,  5,  5,  5,  5,  6,  7,  8, 10, 12, 13, 13, 13, 12, 12, 12, 18, 18, 16, 10,  1,
113       7,  7,  6,  6,  7,  7,  7,  7,  8,  9, 10, 11, 13, 18, 20, 22, 25, 22, 25, 30, 30, 35, 35, 35, 15,
114       9,  9,  9,  9, 10, 10, 10, 10, 11, 12, 18, 22, 25, 24, 23, 27, 29, 27, 35, 40, 42, 50, 59, 54, 31,
115      16, 16, 16, 16, 16, 18, 20, 22, 24, 28, 30, 32, 35, 34, 33, 37, 39, 40, 45, 50, 52, 60, 60, 60, 60,
116 };
117 #endif
118
119 static celt_int16_t *compute_ebands(celt_int32_t Fs, int frame_size, int *nbEBands)
120 {
121    celt_int16_t *eBands;
122    int i, res, min_width, lin, low, high;
123    res = (Fs+frame_size)/(2*frame_size);
124    min_width = MIN_BINS*res;
125    /*printf ("min_width = %d\n", min_width);*/
126
127    /* Find where the linear part ends (i.e. where the spacing is more than min_width */
128    for (lin=0;lin<BARK_BANDS;lin++)
129       if (bark_freq[lin+1]-bark_freq[lin] >= min_width)
130          break;
131    
132    /*printf ("lin = %d (%d Hz)\n", lin, bark_freq[lin]);*/
133    low = ((bark_freq[lin]/res)+(MIN_BINS-1))/MIN_BINS;
134    high = BARK_BANDS-lin;
135    *nbEBands = low+high;
136    eBands = celt_alloc(sizeof(celt_int16_t)*(*nbEBands+2));
137    
138    /* Linear spacing (min_width) */
139    for (i=0;i<low;i++)
140       eBands[i] = MIN_BINS*i;
141    /* Spacing follows critical bands */
142    for (i=0;i<high;i++)
143       eBands[i+low] = (bark_freq[lin+i]+res/2)/res;
144    /* Enforce the minimum spacing at the boundary */
145    for (i=0;i<*nbEBands;i++)
146       if (eBands[i] < MIN_BINS*i)
147          eBands[i] = MIN_BINS*i;
148    eBands[*nbEBands] = (bark_freq[BARK_BANDS]+res/2)/res;
149    eBands[*nbEBands+1] = frame_size;
150    if (eBands[*nbEBands] > eBands[*nbEBands+1])
151       eBands[*nbEBands] = eBands[*nbEBands+1];
152    
153    /* FIXME: Remove last band if too small */
154    /*for (i=0;i<*nbEBands+2;i++)
155       printf("%d ", eBands[i]);
156    printf ("\n");*/
157    return eBands;
158 }
159
160 static void compute_pbands(CELTMode *mode, int res)
161 {
162    int i;
163    celt_int16_t *pBands;
164    pBands=celt_alloc(sizeof(celt_int16_t)*(PBANDS+2));
165    mode->nbPBands = PBANDS;
166    for (i=0;i<PBANDS+1;i++)
167    {
168       pBands[i] = (pitch_freq[i]+res/2)/res;
169       if (pBands[i] < mode->eBands[i])
170          pBands[i] = mode->eBands[i];
171    }
172    pBands[PBANDS+1] = mode->eBands[mode->nbEBands+1];
173    for (i=1;i<mode->nbPBands+1;i++)
174    {
175       int j;
176       for (j=0;j<mode->nbEBands;j++)
177          if (mode->eBands[j] <= pBands[i] && mode->eBands[j+1] > pBands[i])
178             break;
179       /*printf ("%d %d\n", i, j);*/
180       if (mode->eBands[j] != pBands[i])
181       {
182          if (pBands[i]-mode->eBands[j] < mode->eBands[j+1]-pBands[i] && 
183              mode->eBands[j] != pBands[i-1])
184             pBands[i] = mode->eBands[j];
185          else
186             pBands[i] = mode->eBands[j+1];
187       }
188    }
189    /*for (i=0;i<mode->nbPBands+2;i++)
190       printf("%d ", pBands[i]);
191    printf ("\n");*/
192    mode->pBands = pBands;
193    mode->pitchEnd = pBands[PBANDS];
194 }
195
196 static void compute_allocation_table(CELTMode *mode, int res)
197 {
198    int i, j, eband;
199    celt_int16_t *allocVectors;
200    
201    mode->nbAllocVectors = BITALLOC_SIZE;
202    allocVectors = celt_alloc(sizeof(celt_int16_t)*(BITALLOC_SIZE*mode->nbEBands));
203    for (i=0;i<BITALLOC_SIZE;i++)
204    {
205       eband = 0;
206       for (j=0;j<BARK_BANDS;j++)
207       {
208          int edge, low;
209          edge = mode->eBands[eband+1]*res;
210          if (edge < bark_freq[j+1])
211          {
212             int num, den;
213             num = band_allocation[i*BARK_BANDS+j] * (edge-bark_freq[j]);
214             den = bark_freq[j+1]-bark_freq[j];
215             low = (num+den/2)/den;
216             allocVectors[i*mode->nbEBands+eband] += low;
217             eband++;
218             allocVectors[i*mode->nbEBands+eband] += band_allocation[i*BARK_BANDS+j]-low;
219          } else {
220             allocVectors[i*mode->nbEBands+eband] += band_allocation[i*BARK_BANDS+j];
221          }
222       }
223    }
224    /*for (i=0;i<BITALLOC_SIZE;i++)
225    {
226       for (j=0;j<mode->nbEBands;j++)
227          printf ("%2d ", allocVectors[i*mode->nbEBands+j]);
228       printf ("\n");
229    }*/
230    mode->allocVectors = allocVectors;
231 }
232
233 #endif /* STATIC_MODES */
234
235 CELTMode EXPORT *celt_mode_create(celt_int32_t Fs, int channels, int frame_size, int lookahead, int *error)
236 {
237    int i;
238 #ifdef STDIN_TUNING
239    scanf("%d ", &MIN_BINS);
240    scanf("%d ", &BITALLOC_SIZE);
241    band_allocation = celt_alloc(sizeof(int)*BARK_BANDS*BITALLOC_SIZE);
242    for (i=0;i<BARK_BANDS*BITALLOC_SIZE;i++)
243    {
244       scanf("%d ", band_allocation+i);
245    }
246 #endif
247 #ifdef STATIC_MODES
248    const CELTMode *m = NULL;
249    CELTMode *mode=NULL;
250    ALLOC_STACK;
251    for (i=0;i<TOTAL_MODES;i++)
252    {
253       if (Fs == static_mode_list[i]->Fs &&
254           channels == static_mode_list[i]->nbChannels &&
255           frame_size == static_mode_list[i]->mdctSize &&
256           lookahead == static_mode_list[i]->overlap)
257       {
258          m = static_mode_list[i];
259          break;
260       }
261    }
262    if (m == NULL)
263    {
264       celt_warning("Mode not included as part of the static modes");
265       if (error)
266          *error = CELT_BAD_ARG;
267       return NULL;
268    }
269    mode = (CELTMode*)celt_alloc(sizeof(CELTMode));
270    CELT_COPY(mode, m, 1);
271 #else
272    int res;
273    CELTMode *mode;
274    celt_word16_t *window;
275    ALLOC_STACK;
276
277    /* The good thing here is that permutation of the arguments will automatically be invalid */
278    
279    if (Fs < 32000 || Fs > 64000)
280    {
281       celt_warning("Sampling rate must be between 32 kHz and 64 kHz");
282       if (error)
283          *error = CELT_BAD_ARG;
284       return NULL;
285    }
286    if (channels < 0 || channels > 2)
287    {
288       celt_warning("Only mono and stereo supported");
289       if (error)
290          *error = CELT_BAD_ARG;
291       return NULL;
292    }
293    if (frame_size < 64 || frame_size > 256 || frame_size%2!=0)
294    {
295       celt_warning("Only even frame sizes between 64 and 256 are supported");
296       if (error)
297          *error = CELT_BAD_ARG;
298       return NULL;
299    }
300    if (lookahead < 32 || lookahead > frame_size)
301    {
302       celt_warning("The overlap must be between 32 and the frame size");
303       if (error)
304          *error = CELT_BAD_ARG;
305       return NULL;
306    }
307    res = (Fs+frame_size)/(2*frame_size);
308    
309    mode = celt_alloc(sizeof(CELTMode));
310    mode->Fs = Fs;
311    mode->overlap = lookahead;
312    mode->mdctSize = frame_size;
313    mode->nbChannels = channels;
314    mode->eBands = compute_ebands(Fs, frame_size, &mode->nbEBands);
315    compute_pbands(mode, res);
316    mode->ePredCoef = QCONST16(.8f,15);
317    
318    compute_allocation_table(mode, res);
319    /*printf ("%d bands\n", mode->nbEBands);*/
320    
321    window = (celt_word16_t*)celt_alloc(mode->overlap*sizeof(celt_word16_t));
322
323 #ifndef FIXED_POINT
324    for (i=0;i<mode->overlap;i++)
325       window[i] = Q15ONE*sin(.5*M_PI* sin(.5*M_PI*(i+.5)/mode->overlap) * sin(.5*M_PI*(i+.5)/mode->overlap));
326 #else
327    for (i=0;i<mode->overlap;i++)
328       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)));
329 #endif
330    mode->window = window;
331
332    compute_alloc_cache(mode);
333
334 #ifndef SHORTCUTS
335    psydecay_init(&mode->psy, MAX_PERIOD/2, mode->Fs);
336 #endif
337    
338    mode->marker_start = MODEVALID;
339    mode->marker_end = MODEVALID;
340 #endif /* !STATIC_MODES */
341    mdct_init(&mode->mdct, 2*mode->mdctSize);
342    mode->fft = pitch_state_alloc(MAX_PERIOD);
343
344    mode->prob = quant_prob_alloc(mode);
345    if (error)
346       *error = CELT_OK;
347    return mode;
348 }
349
350 void EXPORT celt_mode_destroy(CELTMode *mode)
351 {
352 #ifndef STATIC_MODES
353    int i;
354    const celt_int16_t *prevPtr = NULL;
355    for (i=0;i<mode->nbEBands;i++)
356    {
357       if (mode->bits[i] != prevPtr)
358       {
359          prevPtr = mode->bits[i];
360          celt_free((int*)mode->bits[i]);
361       }
362    }
363    celt_free((int**)mode->bits);
364    if (check_mode(mode) != CELT_OK)
365       return;
366    celt_free((int*)mode->eBands);
367    celt_free((int*)mode->pBands);
368    celt_free((int*)mode->allocVectors);
369    
370    celt_free((celt_word16_t*)mode->window);
371
372    mode->marker_start = MODEFREED;
373    mode->marker_end = MODEFREED;
374 #ifndef SHORTCUTS
375    psydecay_clear(&mode->psy);
376 #endif
377 #endif
378    mdct_clear(&mode->mdct);
379    pitch_state_free(mode->fft);
380    quant_prob_free(mode->prob);
381    celt_free((CELTMode *)mode);
382 }
383
384 int check_mode(const CELTMode *mode)
385 {
386    if (mode->marker_start == MODEVALID && mode->marker_end == MODEVALID)
387       return CELT_OK;
388    if (mode->marker_start == MODEFREED || mode->marker_end == MODEFREED)
389       celt_warning("Using a mode that has already been freed");
390    else
391       celt_warning("This is not a valid CELT mode");
392    return CELT_INVALID_MODE;
393 }