xref: /plan9/sys/src/games/mp3enc/util.c (revision 8f5875f3e9b20916b4c52ad4336922bc8653eb7b)
1 /*
2  *	lame utility library source file
3  *
4  *	Copyright (c) 1999 Albert L Faber
5  *
6  * This library is free software; you can redistribute it and/or
7  * modify it under the terms of the GNU Library General Public
8  * License as published by the Free Software Foundation; either
9  * version 2 of the License, or (at your option) any later version.
10  *
11  * This library is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.	 See the GNU
14  * Library General Public License for more details.
15  *
16  * You should have received a copy of the GNU Library General Public
17  * License along with this library; if not, write to the
18  * Free Software Foundation, Inc., 59 Temple Place - Suite 330,
19  * Boston, MA 02111-1307, USA.
20  */
21 
22 /* $Id: util.c,v 1.67 2001/03/20 00:42:56 markt Exp $ */
23 
24 #ifdef HAVE_CONFIG_H
25 # include <config.h>
26 #endif
27 
28 #define PRECOMPUTE
29 
30 #include "util.h"
31 #include <ctype.h>
32 #include <assert.h>
33 #include <stdarg.h>
34 
35 #if defined(__FreeBSD__) && !defined(__alpha__)
36 # include <machine/floatingpoint.h>
37 #endif
38 
39 #ifdef WITH_DMALLOC
40 #include <dmalloc.h>
41 #endif
42 
43 /***********************************************************************
44 *
45 *  Global Function Definitions
46 *
47 ***********************************************************************/
48 /*empty and close mallocs in gfc */
49 
freegfc(lame_internal_flags * const gfc)50 void  freegfc ( lame_internal_flags* const gfc )   /* bit stream structure */
51 {
52     int  i;
53 
54 #ifdef KLEMM_44
55     if (gfc->resample_in != NULL) {
56         resample_close(gfc->resample_in);
57         gfc->resample_in = NULL;
58     }
59     free(gfc->mfbuf[0]);
60     free(gfc->mfbuf[1]);
61 #endif
62 
63     for ( i = 0 ; i <= 2*BPC; i++ )
64         if ( gfc->blackfilt[i] != NULL ) {
65             free ( gfc->blackfilt[i] );
66 	    gfc->blackfilt[i] = NULL;
67 	}
68     if ( gfc->inbuf_old[0] ) {
69         free ( gfc->inbuf_old[0] );
70 	gfc->inbuf_old[0] = NULL;
71     }
72     if ( gfc->inbuf_old[1] ) {
73         free ( gfc->inbuf_old[1] );
74 	gfc->inbuf_old[1] = NULL;
75     }
76 
77     if ( gfc->bs.buf != NULL ) {
78         free ( gfc->bs.buf );
79         gfc->bs.buf = NULL;
80     }
81 
82     if ( gfc->VBR_seek_table.bag ) {
83         free ( gfc->VBR_seek_table.bag );
84     }
85     if ( gfc->ATH ) {
86         free ( gfc->ATH );
87     }
88     free ( gfc );
89 }
90 
ATHformula_old(FLOAT8 f)91 FLOAT8 ATHformula_old(FLOAT8 f)
92 {
93   FLOAT8 ath;
94   f /= 1000;  // convert to khz
95   f  = Max(0.01, f);
96   f  = Min(18.0, f);
97 
98   /* from Painter & Spanias, 1997 */
99   /* minimum: (i=77) 3.3kHz = -5db */
100   ath =    3.640 * pow(f,-0.8)
101          - 6.500 * exp(-0.6*pow(f-3.3,2.0))
102          + 0.001 * pow(f,4.0);
103   return ath;
104 }
105 
ATHformula_GB(FLOAT8 f)106 FLOAT8 ATHformula_GB(FLOAT8 f)
107 {
108   FLOAT8 ath;
109   f /= 1000;  // convert to khz
110   f  = Max(0.01, f);
111   f  = Min(18.0, f);
112 
113   /* from Painter & Spanias, 1997 */
114   /* modified by Gabriel Bouvigne to better fit to the reality */
115   ath =    3.640 * pow(f,-0.8)
116          - 6.800 * exp(-0.6*pow(f-3.4,2.0))
117          + 6.000 * exp(-0.15*pow(f-8.7,2.0))
118          + 0.6* 0.001 * pow(f,4.0);
119   return ath;
120 }
121 
ATHformula_GBtweak(FLOAT8 f)122 FLOAT8 ATHformula_GBtweak(FLOAT8 f)
123 {
124   FLOAT8 ath;
125   f /= 1000;  // convert to khz
126   f  = Max(0.01, f);
127   f  = Min(18.0, f);
128 
129   /* from Painter & Spanias, 1997 */
130   /* modified by Gabriel Bouvigne to better fit to the reality */
131   ath =    3.640 * pow(f,-0.8)
132          - 6.800 * exp(-0.6*pow(f-3.4,2.0))
133          + 6.000 * exp(-0.15*pow(f-8.7,2.0))
134          + 0.57* 0.001 * pow(f,4.0) //0.57 to maximize HF importance
135          + 6; //std --athlower -6 for
136   return ath;
137 }
138 
139 
140 /*
141  *  Klemm 1994 and 1997. Experimental data. Sorry, data looks a little bit
142  *  dodderly. Data below 30 Hz is extrapolated from other material, above 18
143  *  kHz the ATH is limited due to the original purpose (too much noise at
144  *  ATH is not good even if it's theoretically inaudible).
145  */
146 
ATHformula_Frank(FLOAT8 freq)147 FLOAT8  ATHformula_Frank( FLOAT8 freq )
148 {
149     /*
150      * one value per 100 cent = 1
151      * semitone = 1/4
152      * third = 1/12
153      * octave = 1/40 decade
154      * rest is linear interpolated, values are currently in decibel rel. 20 �Pa
155      */
156     static FLOAT tab [] = {
157         /*    10.0 */  96.69, 96.69, 96.26, 95.12,
158         /*    12.6 */  93.53, 91.13, 88.82, 86.76,
159         /*    15.8 */  84.69, 82.43, 79.97, 77.48,
160         /*    20.0 */  74.92, 72.39, 70.00, 67.62,
161         /*    25.1 */  65.29, 63.02, 60.84, 59.00,
162         /*    31.6 */  57.17, 55.34, 53.51, 51.67,
163         /*    39.8 */  50.04, 48.12, 46.38, 44.66,
164         /*    50.1 */  43.10, 41.73, 40.50, 39.22,
165         /*    63.1 */  37.23, 35.77, 34.51, 32.81,
166         /*    79.4 */  31.32, 30.36, 29.02, 27.60,
167         /*   100.0 */  26.58, 25.91, 24.41, 23.01,
168         /*   125.9 */  22.12, 21.25, 20.18, 19.00,
169         /*   158.5 */  17.70, 16.82, 15.94, 15.12,
170         /*   199.5 */  14.30, 13.41, 12.60, 11.98,
171         /*   251.2 */  11.36, 10.57,  9.98,  9.43,
172         /*   316.2 */   8.87,  8.46,  7.44,  7.12,
173         /*   398.1 */   6.93,  6.68,  6.37,  6.06,
174         /*   501.2 */   5.80,  5.55,  5.29,  5.02,
175         /*   631.0 */   4.75,  4.48,  4.22,  3.98,
176         /*   794.3 */   3.75,  3.51,  3.27,  3.22,
177         /*  1000.0 */   3.12,  3.01,  2.91,  2.68,
178         /*  1258.9 */   2.46,  2.15,  1.82,  1.46,
179         /*  1584.9 */   1.07,  0.61,  0.13, -0.35,
180         /*  1995.3 */  -0.96, -1.56, -1.79, -2.35,
181         /*  2511.9 */  -2.95, -3.50, -4.01, -4.21,
182         /*  3162.3 */  -4.46, -4.99, -5.32, -5.35,
183         /*  3981.1 */  -5.13, -4.76, -4.31, -3.13,
184         /*  5011.9 */  -1.79,  0.08,  2.03,  4.03,
185         /*  6309.6 */   5.80,  7.36,  8.81, 10.22,
186         /*  7943.3 */  11.54, 12.51, 13.48, 14.21,
187         /* 10000.0 */  14.79, 13.99, 12.85, 11.93,
188         /* 12589.3 */  12.87, 15.19, 19.14, 23.69,
189         /* 15848.9 */  33.52, 48.65, 59.42, 61.77,
190         /* 19952.6 */  63.85, 66.04, 68.33, 70.09,
191         /* 25118.9 */  70.66, 71.27, 71.91, 72.60,
192     };
193     FLOAT8    freq_log;
194     unsigned  index;
195 
196     if ( freq <    10. ) freq =    10.;
197     if ( freq > 29853. ) freq = 29853.;
198 
199     freq_log = 40. * log10 (0.1 * freq);   /* 4 steps per third, starting at 10 Hz */
200     index    = (unsigned) freq_log;
201     assert ( index < sizeof(tab)/sizeof(*tab) );
202     return tab [index] * (1 + index - freq_log) + tab [index+1] * (freq_log - index);
203 }
204 
ATHformula(FLOAT8 f,lame_global_flags * gfp)205 FLOAT8 ATHformula(FLOAT8 f,lame_global_flags *gfp)
206 {
207   switch(gfp->ATHtype)
208     {
209     case 0:
210       return ATHformula_old(f);
211     case 1:
212       return ATHformula_Frank(f);
213     case 2:
214       return ATHformula_GB(f);
215     case 3:
216       return ATHformula_GBtweak(f);
217     }
218 
219   return ATHformula_Frank(f);
220 }
221 
222 /* see for example "Zwicker: Psychoakustik, 1982; ISBN 3-540-11401-7 */
freq2bark(FLOAT8 freq)223 FLOAT8 freq2bark(FLOAT8 freq)
224 {
225   /* input: freq in hz  output: barks */
226     if (freq<0) freq=0;
227     freq = freq * 0.001;
228     return 13.0*atan(.76*freq) + 3.5*atan(freq*freq/(7.5*7.5));
229 }
230 
231 /* see for example "Zwicker: Psychoakustik, 1982; ISBN 3-540-11401-7 */
freq2cbw(FLOAT8 freq)232 FLOAT8 freq2cbw(FLOAT8 freq)
233 {
234   /* input: freq in hz  output: critical band width */
235     freq = freq * 0.001;
236     return 25+75*pow(1+1.4*(freq*freq),0.69);
237 }
238 
239 
240 
241 
242 
243 
244 /***********************************************************************
245  * compute bitsperframe and mean_bits for a layer III frame
246  **********************************************************************/
getframebits(lame_global_flags * gfp,int * bitsPerFrame,int * mean_bits)247 void getframebits(lame_global_flags *gfp, int *bitsPerFrame, int *mean_bits)
248 {
249   lame_internal_flags *gfc=gfp->internal_flags;
250   int  whole_SpF;  /* integral number of Slots per Frame without padding */
251   int  bit_rate;
252 
253   /* get bitrate in kbps [?] */
254   if (gfc->bitrate_index)
255     bit_rate = bitrate_table[gfp->version][gfc->bitrate_index];
256   else
257     bit_rate = gfp->brate;
258   assert ( bit_rate <= 550 );
259 
260   // bytes_per_frame = bitrate * 1000 / ( gfp->out_samplerate / (gfp->version == 1  ?  1152  :  576 )) / 8;
261   // bytes_per_frame = bitrate * 1000 / gfp->out_samplerate * (gfp->version == 1  ?  1152  :  576 ) / 8;
262   // bytes_per_frame = bitrate * ( gfp->version == 1  ?  1152/8*1000  :  576/8*1000 ) / gfp->out_samplerate;
263 
264   whole_SpF = (gfp->version+1)*72000*bit_rate / gfp->out_samplerate;
265 
266   // There must be somewhere code toggling gfc->padding on and off
267 
268   /* main encoding routine toggles padding on and off */
269   /* one Layer3 Slot consists of 8 bits */
270   *bitsPerFrame = 8 * (whole_SpF + gfc->padding);
271 
272   // sideinfo_len
273   *mean_bits = (*bitsPerFrame - 8*gfc->sideinfo_len) / gfc->mode_gr;
274 }
275 
276 
277 
278 
279 #define ABS(A) (((A)>0) ? (A) : -(A))
280 
FindNearestBitrate(int bRate,int version,int samplerate)281 int FindNearestBitrate(
282 int bRate,        /* legal rates from 32 to 448 */
283 int version,      /* MPEG-1 or MPEG-2 LSF */
284 int samplerate)   /* convert bitrate in kbps to index */
285 {
286     int  bitrate = 0;
287     int  i;
288 
289     for ( i = 1; i <= 14; i++ )
290         if ( ABS (bitrate_table[version][i] - bRate) < ABS (bitrate - bRate) )
291             bitrate = bitrate_table [version] [i];
292 
293     return bitrate;
294 }
295 
296 
297 /* map frequency to a valid MP3 sample frequency
298  *
299  * Robert.Hegemann@gmx.de 2000-07-01
300  */
map2MP3Frequency(int freq)301 int map2MP3Frequency(int freq)
302 {
303     if (freq <=  8000) return  8000;
304     if (freq <= 11025) return 11025;
305     if (freq <= 12000) return 12000;
306     if (freq <= 16000) return 16000;
307     if (freq <= 22050) return 22050;
308     if (freq <= 24000) return 24000;
309     if (freq <= 32000) return 32000;
310     if (freq <= 44100) return 44100;
311 
312     return 48000;
313 }
314 
BitrateIndex(int bRate,int version,int samplerate)315 int BitrateIndex(
316 int bRate,        /* legal rates from 32 to 448 kbps */
317 int version,      /* MPEG-1 or MPEG-2/2.5 LSF */
318 int samplerate)   /* convert bitrate in kbps to index */
319 {
320     int  i;
321 
322     for ( i = 0; i <= 14; i++)
323         if ( bitrate_table [version] [i] == bRate )
324             return i;
325 
326     return -1;
327 }
328 
329 /* convert samp freq in Hz to index */
330 
SmpFrqIndex(int sample_freq,int * const version)331 int SmpFrqIndex ( int sample_freq, int* const version )
332 {
333     switch ( sample_freq ) {
334     case 44100: *version = 1; return  0;
335     case 48000: *version = 1; return  1;
336     case 32000: *version = 1; return  2;
337     case 22050: *version = 0; return  0;
338     case 24000: *version = 0; return  1;
339     case 16000: *version = 0; return  2;
340     case 11025: *version = 0; return  0;
341     case 12000: *version = 0; return  1;
342     case  8000: *version = 0; return  2;
343     default:    *version = 0; return -1;
344     }
345 }
346 
347 
348 /*****************************************************************************
349 *
350 *  End of bit_stream.c package
351 *
352 *****************************************************************************/
353 
354 /* reorder the three short blocks By Takehiro TOMINAGA */
355 /*
356   Within each scalefactor band, data is given for successive
357   time windows, beginning with window 0 and ending with window 2.
358   Within each window, the quantized values are then arranged in
359   order of increasing frequency...
360 */
freorder(int scalefac_band[],FLOAT8 ix_orig[576])361 void freorder(int scalefac_band[],FLOAT8 ix_orig[576]) {
362   int i,sfb, window, j=0;
363   FLOAT8 ix[576];
364   for (sfb = 0; sfb < SBMAX_s; sfb++) {
365     int start = scalefac_band[sfb];
366     int end   = scalefac_band[sfb + 1];
367     for (window = 0; window < 3; window++) {
368       for (i = start; i < end; ++i) {
369 	ix[j++] = ix_orig[3*i+window];
370       }
371     }
372   }
373   memcpy(ix_orig,ix,576*sizeof(FLOAT8));
374 }
375 
376 
377 
378 
379 
380 
381 
382 #ifndef KLEMM_44
383 
384 
385 /* resampling via FIR filter, blackman window */
blackman(FLOAT8 x,FLOAT8 fcn,int l)386 inline static FLOAT8 blackman(FLOAT8 x,FLOAT8 fcn,int l)
387 {
388   /* This algorithm from:
389 SIGNAL PROCESSING ALGORITHMS IN FORTRAN AND C
390 S.D. Stearns and R.A. David, Prentice-Hall, 1992
391   */
392   FLOAT8 bkwn,x2;
393   FLOAT8 wcn = (PI * fcn);
394 
395   x /= l;
396   if (x<0) x=0;
397   if (x>1) x=1;
398   x2 = x - .5;
399 
400   bkwn = 0.42 - 0.5*cos(2*x*PI)  + 0.08*cos(4*x*PI);
401   if (fabs(x2)<1e-9) return wcn/PI;
402   else
403     return  (  bkwn*sin(l*wcn*x2)  / (PI*l*x2)  );
404 
405 
406 }
407 
408 /* gcd - greatest common divisor */
409 /* Joint work of Euclid and M. Hendry */
410 
gcd(int i,int j)411 int gcd ( int i, int j )
412 {
413 //    assert ( i > 0  &&  j > 0 );
414     return j ? gcd(j, i % j) : i;
415 }
416 
417 
418 
419 /* copy in new samples from in_buffer into mfbuf, with resampling & scaling
420    if necessary.  n_in = number of samples from the input buffer that
421    were used.  n_out = number of samples copied into mfbuf  */
422 
fill_buffer(lame_global_flags * gfp,sample_t * mfbuf[2],sample_t * in_buffer[2],int nsamples,int * n_in,int * n_out)423 void fill_buffer(lame_global_flags *gfp,
424 		 sample_t *mfbuf[2],
425 		 sample_t *in_buffer[2],
426 		 int nsamples, int *n_in, int *n_out)
427 {
428     lame_internal_flags *gfc = gfp->internal_flags;
429     int ch,i;
430 
431     /* copy in new samples into mfbuf, with resampling if necessary */
432     if (gfc->resample_ratio != 1.0) {
433 	for (ch = 0; ch < gfc->channels_out; ch++) {
434 	    *n_out =
435 		fill_buffer_resample(gfp, &mfbuf[ch][gfc->mf_size],
436 				     gfp->framesize, in_buffer[ch],
437 				     nsamples, n_in, ch);
438 	}
439     }
440     else {
441 	*n_out = Min(gfp->framesize, nsamples);
442 	*n_in = *n_out;
443 	for (i = 0; i < *n_out; ++i) {
444 	    mfbuf[0][gfc->mf_size + i] = in_buffer[0][i];
445 	    if (gfc->channels_out == 2)
446 		mfbuf[1][gfc->mf_size + i] = in_buffer[1][i];
447 	}
448     }
449 
450     /* user selected scaling of the samples */
451     if (gfp->scale != 0) {
452 	for (i=0 ; i<*n_out; ++i) {
453 	    mfbuf[0][gfc->mf_size+i] *= gfp->scale;
454 	    if (gfc->channels_out == 2)
455 		mfbuf[1][gfc->mf_size + i] *= gfp->scale;
456 	}
457     }
458 
459 }
460 
461 
462 
463 
fill_buffer_resample(lame_global_flags * gfp,sample_t * outbuf,int desired_len,sample_t * inbuf,int len,int * num_used,int ch)464 int fill_buffer_resample(
465        lame_global_flags *gfp,
466        sample_t *outbuf,
467        int desired_len,
468        sample_t *inbuf,
469        int len,
470        int *num_used,
471        int ch)
472 {
473 
474 
475   lame_internal_flags *gfc=gfp->internal_flags;
476   int BLACKSIZE;
477   FLOAT8 offset,xvalue;
478   int i,j=0,k;
479   int filter_l;
480   FLOAT8 fcn,intratio;
481   FLOAT *inbuf_old;
482   int bpc;   /* number of convolution functions to pre-compute */
483   bpc = gfp->out_samplerate/gcd(gfp->out_samplerate,gfp->in_samplerate);
484   if (bpc>BPC) bpc = BPC;
485 
486   intratio=( fabs(gfc->resample_ratio - floor(.5+gfc->resample_ratio)) < .0001 );
487   fcn = 1.00/gfc->resample_ratio;
488   if (fcn>1.00) fcn=1.00;
489   filter_l = gfp->quality < 7 ? 31 : 7;
490   filter_l = 31;
491   if (0==filter_l % 2 ) --filter_l;/* must be odd */
492   filter_l += intratio;            /* unless resample_ratio=int, it must be even */
493 
494 
495   BLACKSIZE = filter_l+1;  /* size of data needed for FIR */
496 
497   if ( gfc->fill_buffer_resample_init == 0 ) {
498     gfc->inbuf_old[0]=calloc(BLACKSIZE,sizeof(gfc->inbuf_old[0][0]));
499     gfc->inbuf_old[1]=calloc(BLACKSIZE,sizeof(gfc->inbuf_old[0][0]));
500     for (i=0; i<=2*bpc; ++i)
501       gfc->blackfilt[i]=calloc(BLACKSIZE,sizeof(gfc->blackfilt[0][0]));
502 
503     gfc->itime[0]=0;
504     gfc->itime[1]=0;
505 
506     /* precompute blackman filter coefficients */
507     for ( j = 0; j <= 2*bpc; j++ ) {
508         FLOAT8 sum = 0.;
509         offset = (j-bpc) / (2.*bpc);
510         for ( i = 0; i <= filter_l; i++ )
511             sum +=
512 	    gfc->blackfilt[j][i]  = blackman(i-offset,fcn,filter_l);
513 	for ( i = 0; i <= filter_l; i++ )
514 	  gfc->blackfilt[j][i] /= sum;
515     }
516     gfc->fill_buffer_resample_init = 1;
517   }
518 
519   inbuf_old=gfc->inbuf_old[ch];
520 
521   /* time of j'th element in inbuf = itime + j/ifreq; */
522   /* time of k'th element in outbuf   =  j/ofreq */
523   for (k=0;k<desired_len;k++) {
524     FLOAT time0;
525     int joff;
526 
527     time0 = k*gfc->resample_ratio;       /* time of k'th output sample */
528     j = floor( time0 -gfc->itime[ch]  );
529 
530     /* check if we need more input data */
531     if ((filter_l + j - filter_l/2) >= len) break;
532 
533     /* blackman filter.  by default, window centered at j+.5(filter_l%2) */
534     /* but we want a window centered at time0.   */
535     offset = ( time0 -gfc->itime[ch] - (j + .5*(filter_l%2)));
536     assert(fabs(offset)<=.500001);
537 
538     /* find the closest precomputed window for this offset: */
539     joff = floor((offset*2*bpc) + bpc +.5);
540 
541     xvalue = 0.;
542     for (i=0 ; i<=filter_l ; ++i) {
543       int j2 = i+j-filter_l/2;
544       int y;
545       assert(j2<len);
546       assert(j2+BLACKSIZE >= 0);
547       y = (j2<0) ? inbuf_old[BLACKSIZE+j2] : inbuf[j2];
548 #define PRECOMPUTE
549 #ifdef PRECOMPUTE
550       xvalue += y*gfc->blackfilt[joff][i];
551 #else
552       xvalue += y*blackman(i-offset,fcn,filter_l);  /* very slow! */
553 #endif
554     }
555     outbuf[k]=xvalue;
556   }
557 
558 
559   /* k = number of samples added to outbuf */
560   /* last k sample used data from [j-filter_l/2,j+filter_l-filter_l/2]  */
561 
562   /* how many samples of input data were used:  */
563   *num_used = Min(len,filter_l+j-filter_l/2);
564 
565   /* adjust our input time counter.  Incriment by the number of samples used,
566    * then normalize so that next output sample is at time 0, next
567    * input buffer is at time itime[ch] */
568   gfc->itime[ch] += *num_used - k*gfc->resample_ratio;
569 
570   /* save the last BLACKSIZE samples into the inbuf_old buffer */
571   if (*num_used >= BLACKSIZE) {
572       for (i=0;i<BLACKSIZE;i++)
573 	  inbuf_old[i]=inbuf[*num_used + i -BLACKSIZE];
574   }else{
575       /* shift in *num_used samples into inbuf_old  */
576        int n_shift = BLACKSIZE-*num_used;  /* number of samples to shift */
577 
578        /* shift n_shift samples by *num_used, to make room for the
579 	* num_used new samples */
580        for (i=0; i<n_shift; ++i )
581 	   inbuf_old[i] = inbuf_old[i+ *num_used];
582 
583        /* shift in the *num_used samples */
584        for (j=0; i<BLACKSIZE; ++i, ++j )
585 	   inbuf_old[i] = inbuf[j];
586 
587        assert(j==*num_used);
588   }
589   return k;  /* return the number samples created at the new samplerate */
590 }
591 
592 
593 #endif /* ndef KLEMM_44 */
594 
595 
596 
597 /***********************************************************************
598 *
599 *  Message Output
600 *
601 ***********************************************************************/
lame_debugf(const lame_internal_flags * gfc,const char * format,...)602 void  lame_debugf (const lame_internal_flags *gfc, const char* format, ... )
603 {
604     va_list  args;
605 
606     va_start ( args, format );
607 
608     if ( gfc->report.debugf != NULL ) {
609         gfc->report.debugf( format, args );
610     } else {
611         (void) vfprintf ( stderr, format, args );
612         fflush ( stderr );      /* an debug function should flush immediately */
613     }
614 
615     va_end   ( args );
616 }
617 
618 
lame_msgf(const lame_internal_flags * gfc,const char * format,...)619 void  lame_msgf (const lame_internal_flags *gfc, const char* format, ... )
620 {
621     va_list  args;
622 
623     va_start ( args, format );
624 
625     if ( gfc->report.msgf != NULL ) {
626         gfc->report.msgf( format, args );
627     } else {
628         (void) vfprintf ( stderr, format, args );
629         fflush ( stderr );     /* we print to stderr, so me may want to flush */
630     }
631 
632     va_end   ( args );
633 }
634 
635 
lame_errorf(const lame_internal_flags * gfc,const char * format,...)636 void  lame_errorf (const lame_internal_flags *gfc, const char* format, ... )
637 {
638     va_list  args;
639 
640     va_start ( args, format );
641 
642     if ( gfc->report.errorf != NULL ) {
643         gfc->report.errorf( format, args );
644     } else {
645         (void) vfprintf ( stderr, format, args );
646         fflush   ( stderr );    /* an error function should flush immediately */
647     }
648 
649     va_end   ( args );
650 }
651 
652 
653 
654 /***********************************************************************
655  *
656  *      routines to detect CPU specific features like 3DNow, MMX, SIMD
657  *
658  *  donated by Frank Klemm
659  *  added Robert Hegemann 2000-10-10
660  *
661  ***********************************************************************/
662 
has_i387(void)663 int  has_i387 ( void )
664 {
665 #ifdef HAVE_NASM
666     return 1;
667 #else
668     return 0;   /* don't know, assume not */
669 #endif
670 }
671 
has_MMX(void)672 int  has_MMX ( void )
673 {
674 #ifdef HAVE_NASM
675     extern int has_MMX_nasm ( void );
676     return has_MMX_nasm ();
677 #else
678     return 0;   /* don't know, assume not */
679 #endif
680 }
681 
has_3DNow(void)682 int  has_3DNow ( void )
683 {
684 #ifdef HAVE_NASM
685     extern int has_3DNow_nasm ( void );
686     return has_3DNow_nasm ();
687 #else
688     return 0;   /* don't know, assume not */
689 #endif
690 }
691 
has_SIMD(void)692 int  has_SIMD ( void )
693 {
694 #ifdef HAVE_NASM
695     extern int has_SIMD_nasm ( void );
696     return has_SIMD_nasm ();
697 #else
698     return 0;   /* don't know, assume not */
699 #endif
700 }
701 
has_SIMD2(void)702 int  has_SIMD2 ( void )
703 {
704 #ifdef HAVE_NASM
705     extern int has_SIMD2_nasm ( void );
706     return has_SIMD2_nasm ();
707 #else
708     return 0;   /* don't know, assume not */
709 #endif
710 }
711 
712 /***********************************************************************
713  *
714  *  some simple statistics
715  *
716  *  bitrate index 0: free bitrate -> not allowed in VBR mode
717  *  : bitrates, kbps depending on MPEG version
718  *  bitrate index 15: forbidden
719  *
720  *  mode_ext:
721  *  0:  LR
722  *  1:  LR-i
723  *  2:  MS
724  *  3:  MS-i
725  *
726  ***********************************************************************/
727 
updateStats(lame_internal_flags * const gfc)728 void updateStats( lame_internal_flags * const gfc )
729 {
730     assert ( gfc->bitrate_index < 16u );
731     assert ( gfc->mode_ext      <  4u );
732 
733     /* count bitrate indices */
734     gfc->bitrate_stereoMode_Hist [gfc->bitrate_index] [4] ++;
735 
736     /* count 'em for every mode extension in case of 2 channel encoding */
737     if (gfc->channels_out == 2)
738         gfc->bitrate_stereoMode_Hist [gfc->bitrate_index] [gfc->mode_ext]++;
739 }
740 
741 
742 
743 /*  caution: a[] will be resorted!!
744  */
select_kth_int(int a[],int N,int k)745 int select_kth_int(int a[], int N, int k)
746 {
747     int i, j, l, r, v, w;
748 
749     l = 0;
750     r = N-1;
751     while (r > l) {
752         v = a[r];
753         i = l-1;
754         j = r;
755         for (;;) {
756             while (a[++i] < v) /*empty*/;
757             while (a[--j] > v) /*empty*/;
758             if (i >= j)
759                 break;
760             /* swap i and j */
761             w = a[i];
762             a[i] = a[j];
763             a[j] = w;
764         }
765         /* swap i and r */
766         w = a[i];
767         a[i] = a[r];
768         a[r] = w;
769         if (i >= k)
770             r = i-1;
771         if (i <= k)
772             l = i+1;
773     }
774     return a[k];
775 }
776 
777 
778 
disable_FPE(void)779 void disable_FPE(void) {
780 /* extremly system dependent stuff, move to a lib to make the code readable */
781 /*==========================================================================*/
782 
783     /*
784      *  Disable floating point exceptions
785      */
786 #if defined(__FreeBSD__) && !defined(__alpha__)
787     {
788         /* seet floating point mask to the Linux default */
789         fp_except_t mask;
790         mask = fpgetmask();
791         /* if bit is set, we get SIGFPE on that error! */
792         fpsetmask(mask & ~(FP_X_INV | FP_X_DZ));
793         /*  DEBUGF("FreeBSD mask is 0x%x\n",mask); */
794     }
795 #endif
796 
797 #if defined(__riscos__) && !defined(ABORTFP)
798     /* Disable FPE's under RISC OS */
799     /* if bit is set, we disable trapping that error! */
800     /*   _FPE_IVO : invalid operation */
801     /*   _FPE_DVZ : divide by zero */
802     /*   _FPE_OFL : overflow */
803     /*   _FPE_UFL : underflow */
804     /*   _FPE_INX : inexact */
805     DisableFPETraps(_FPE_IVO | _FPE_DVZ | _FPE_OFL);
806 #endif
807 
808     /*
809      *  Debugging stuff
810      *  The default is to ignore FPE's, unless compiled with -DABORTFP
811      *  so add code below to ENABLE FPE's.
812      */
813 
814 #if defined(ABORTFP)
815 #if defined(_MSC_VER)
816     {
817 #include <float.h>
818         unsigned int mask;
819         mask = _controlfp(0, 0);
820         mask &= ~(_EM_OVERFLOW | _EM_UNDERFLOW | _EM_ZERODIVIDE | _EM_INVALID);
821         mask = _controlfp(mask, _MCW_EM);
822     }
823 #elif defined(__CYGWIN__)
824 #  define _FPU_GETCW(cw) __asm__ ("fnstcw %0" : "=m" (*&cw))
825 #  define _FPU_SETCW(cw) __asm__ ("fldcw %0" : : "m" (*&cw))
826 
827 #  define _EM_INEXACT     0x00000020 /* inexact (precision) */
828 #  define _EM_UNDERFLOW   0x00000010 /* underflow */
829 #  define _EM_OVERFLOW    0x00000008 /* overflow */
830 #  define _EM_ZERODIVIDE  0x00000004 /* zero divide */
831 #  define _EM_INVALID     0x00000001 /* invalid */
832     {
833         unsigned int mask;
834         _FPU_GETCW(mask);
835         /* Set the FPU control word to abort on most FPEs */
836         mask &= ~(_EM_OVERFLOW | _EM_ZERODIVIDE | _EM_INVALID);
837         _FPU_SETCW(mask);
838     }
839 # elif defined(__linux__)
840     {
841 
842 #  include <fpu_control.h>
843 #  ifndef _FPU_GETCW
844 #  define _FPU_GETCW(cw) __asm__ ("fnstcw %0" : "=m" (*&cw))
845 #  endif
846 #  ifndef _FPU_SETCW
847 #  define _FPU_SETCW(cw) __asm__ ("fldcw %0" : : "m" (*&cw))
848 #  endif
849 
850         /*
851          * Set the Linux mask to abort on most FPE's
852          * if bit is set, we _mask_ SIGFPE on that error!
853          *  mask &= ~( _FPU_MASK_IM | _FPU_MASK_ZM | _FPU_MASK_OM | _FPU_MASK_UM );
854          */
855 
856         unsigned int mask;
857         _FPU_GETCW(mask);
858         mask &= ~(_FPU_MASK_IM | _FPU_MASK_ZM | _FPU_MASK_OM);
859         _FPU_SETCW(mask);
860     }
861 #endif
862 #endif /* ABORTFP */
863 }
864 
865 
866 /* end of util.c */
867