ENCODE.C
上传用户:njqiyou
上传日期:2007-01-08
资源大小:574k
文件大小:53k
源码类别:

mpeg/mp3

开发平台:

C/C++

  1. /**********************************************************************
  2.  * ISO MPEG Audio Subgroup Software Simulation Group (1996)
  3.  * ISO 13818-3 MPEG-2 Audio Encoder - Lower Sampling Frequency Extension
  4.  *
  5.  * $Id: encode.c,v 1.1 1996/02/14 04:04:23 rowlands Exp $
  6.  *
  7.  * $Log: encode.c,v $
  8.  * Revision 1.1  1996/02/14 04:04:23  rowlands
  9.  * Initial revision
  10.  *
  11.  * Received from Mike Coleman
  12.  **********************************************************************/
  13. /**********************************************************************
  14.  *   date   programmers         comment                               *
  15.  * 3/01/91  Douglas Wong,       start of version 1.1 records          *
  16.  *          Davis Pan                                                 *
  17.  * 3/06/91  Douglas Wong        rename: setup.h to endef.h            *
  18.  *                                      efilter to enfilter           *
  19.  *                                      ewindow to enwindow           *
  20.  *                              integrated "quantizer", "scalefactor",*
  21.  *                              and "transmission" files              *
  22.  *                              update routine "window_subband"       *
  23.  * 3/31/91  Bill Aspromonte     replaced read_filter by               *
  24.  *                              create_an_filter                      *
  25.  * 5/10/91  W. Joseph Carter    Ported to Macintosh and Unix.         *
  26.  *                              Incorporated Jean-Georges Fritsch's   *
  27.  *                              "bitstream.c" package.                *
  28.  *                              Incorporated Bill Aspromonte's        *
  29.  *                              filterbank coefficient matrix         *
  30.  *                              calculation routines and added        *
  31.  *                              roundoff to coincide with specs.      *
  32.  *                              Modified to strictly adhere to        *
  33.  *                              encoded bitstream specs, including    *
  34.  *                              "Berlin changes".                     *
  35.  *                              Modified PCM sound file handling to   *
  36.  *                              process all incoming samples and fill *
  37.  *                              out last encoded frame with zeros     *
  38.  *                              (silence) if needed.                  *
  39.  *                              Located and fixed numerous software   *
  40.  *                              bugs and table data errors.           *
  41.  * 19jun91  dpwe (Aware)        moved "alloc_*" reader to common.c    *
  42.  *                              Globals sblimit, alloc replaced by new*
  43.  *                              struct 'frame_params' passed as arg.  *
  44.  *                              Added JOINT STEREO coding, layers I,II*
  45.  *                              Affects: *_bit_allocation,            *
  46.  *                              subband_quantization, encode_bit_alloc*
  47.  *                              sample_encoding                       *
  48.  * 6/10/91  Earle Jennings      modified II_subband_quantization to   *
  49.  *                              resolve type cast problem for MS_DOS  *
  50.  * 6/11/91  Earle Jennings      modified to avoid overflow on MS_DOS  *
  51.  *                              in routine filter_subband             *
  52.  * 7/10/91  Earle Jennings      port to MsDos from MacIntosh version  *
  53.  * 8/ 8/91  Jens Spille         Change for MS-C6.00                   *
  54.  *10/ 1/91  S.I. Sudharsanan,   Ported to IBM AIX platform.           *
  55.  *          Don H. Lee,                                               *
  56.  *          Peter W. Farrett                                          *
  57.  *10/ 3/91  Don H. Lee          implemented CRC-16 error protection   *
  58.  *                              newly introduced function encode_CRC  *
  59.  *11/ 8/91  Kathy Wang          Documentation of code                 *
  60.  *                              All variablenames are referred to     *
  61.  *                              with surrounding pound (#) signs      *
  62.  * 2/11/92  W. Joseph Carter    Ported new code to Macintosh.  Most   *
  63.  *                              important fixes involved changing     *
  64.  *                              16-bit ints to long or unsigned in    *
  65.  *                              bit alloc routines for quant of 65535 *
  66.  *                              and passing proper function args.     *
  67.  *                              Removed "Other Joint Stereo" option   *
  68.  *                              and made bitrate be total channel     *
  69.  *                              bitrate, irrespective of the mode.    *
  70.  *                              Fixed many small bugs & reorganized.  *
  71.  * 6/16/92  Shaun Astarabadi    Changed I_scale_factor_calc() and     *
  72.  *                              II_scale_factor_calc() to use scale   *
  73.  *                              factor 0 thru 62 only and not to      *
  74.  *                              encode index 63 into the bit stream.  *
  75.  * 7/27/92  Mike Li             (re-)Port to MS-DOS                   *
  76.  * 9/22/92  jddevine@aware.com  Fixed _scale_factor_calc() defs       *
  77.  * 3/31/93  Giogio Dimino       changed II_a_bit_allocation() from:   *
  78.  *                              if( ad > ...) to if(ad >= ...)        *
  79.  * 8/05/93  TEST                changed I_a_bit_allocation() from:    *
  80.  *                              if( ad > ...) to if(ad >= ...)        *
  81.  * 8/02/95  mc@fivebats.com     Changed audio file reading code to    *
  82.  *                              read samples big-endian               *
  83.  *10/15/95  mc@fivebats.com     Modified get_audio() for layer3-LSF   *
  84.  **********************************************************************/
  85.  
  86. #include "common.h"
  87. #include "encoder.h"
  88. #ifdef MS_DOS
  89. extern unsigned _stklen = 16384;
  90. #endif
  91. /*=======================================================================
  92. |                                                                       |
  93. | This segment contains all the core routines of the encoder,           |
  94. | except for the psychoacoustic models.                                 |
  95. |                                                                       |
  96. | The user can select either one of the two psychoacoustic              |
  97. | models. Model I is a simple tonal and noise masking threshold         |
  98. | generator, and Model II is a more sophisticated cochlear masking      |
  99. | threshold generator. Model I is recommended for lower complexity      |
  100. | applications whereas Model II gives better subjective quality at low  |
  101. | bit rates.                                                            |
  102. |                                                                       |
  103. | Layers I and II of mono, stereo, and joint stereo modes are supported.|
  104. | Routines associated with a given layer are prefixed by "I_" for layer |
  105. | 1 and "II_" for layer 2.                                              |
  106. =======================================================================*/
  107.  
  108. /************************************************************************
  109. *
  110. * read_samples()
  111. *
  112. * PURPOSE:  reads the PCM samples from a file to the buffer
  113. *
  114. *  SEMANTICS:
  115. * Reads #samples_read# number of shorts from #musicin# filepointer
  116. * into #sample_buffer[]#.  Returns the number of samples read.
  117. *
  118. ************************************************************************/
  119. unsigned long read_samples(musicin, sample_buffer, num_samples, frame_size)
  120. FILE *musicin;
  121. short sample_buffer[2304];
  122. unsigned long num_samples, frame_size;
  123. {
  124.     unsigned long samples_read;
  125.     static unsigned long samples_to_read;
  126.     static char init = TRUE;
  127.     if (init) {
  128.         samples_to_read = num_samples;
  129.         init = FALSE;
  130.     }
  131.     if (samples_to_read >= frame_size)
  132.         samples_read = frame_size;
  133.     else
  134.         samples_read = samples_to_read;
  135.     if ((samples_read =
  136.          fread(sample_buffer, sizeof(short), (int)samples_read, musicin)) == 0)
  137.         printf("Hit end of audio datan");
  138.     /*
  139.        Samples are big-endian. If this is a little-endian machine
  140.        we must swap
  141.      */
  142.     if ( NativeByteOrder == order_unknown )
  143.       {
  144. NativeByteOrder = DetermineByteOrder();
  145. if ( NativeByteOrder == order_unknown )
  146.   {
  147.     fprintf( stderr, "byte order not determinedn" );
  148.     exit( 1 );
  149.   }
  150.       }
  151.     if ( NativeByteOrder == order_littleEndian )
  152.       SwapBytesInWords( sample_buffer, samples_read );
  153.     samples_to_read -= samples_read;
  154.     if (samples_read < frame_size && samples_read > 0) {
  155.         printf("Insufficient PCM input for one frame - fillout with zerosn");
  156.         for (; samples_read < frame_size; sample_buffer[samples_read++] = 0);
  157.         samples_to_read = 0;
  158.     }
  159.     return(samples_read);
  160. }
  161. /************************************************************************
  162. *
  163. * get_audio()
  164. *
  165. * PURPOSE:  reads a frame of audio data from a file to the buffer,
  166. *   aligns the data for future processing, and separates the
  167. *   left and right channels
  168. *
  169. *
  170. ************************************************************************/
  171.  
  172. unsigned long get_audio( musicin, buffer, num_samples, stereo, info )
  173. FILE *musicin;
  174. short FAR buffer[2][1152];
  175. unsigned long num_samples;
  176. int stereo;
  177. layer *info;
  178. {
  179.     int j;
  180.     short insamp[2304];
  181.     unsigned long samples_read;
  182.     int lay;
  183.     lay = info->lay;
  184.     if ( (lay == 3) && (info->version == 0) )
  185.     {
  186. if ( stereo == 2 )
  187. {
  188.     samples_read = read_samples( musicin, insamp, num_samples,
  189.  (unsigned long) 1152 );
  190.     for ( j = 0; j < 576; j++ )
  191.     {
  192. buffer[0][j] = insamp[2 * j];
  193. buffer[1][j] = insamp[2 * j + 1];
  194.     }
  195. }
  196. else
  197. {
  198.     samples_read = read_samples( musicin, insamp, num_samples,
  199.  (unsigned long) 576 );
  200.     for ( j = 0; j < 576; j++ )
  201.     {
  202. buffer[0][j] = insamp[j];
  203. buffer[1][j] = 0;
  204.     }
  205. }
  206.     }
  207.     else
  208.     {
  209. if (lay == 1){
  210.     if(stereo == 2){ /* layer 1, stereo */
  211. samples_read = read_samples(musicin, insamp, num_samples,
  212.     (unsigned long) 768);
  213. for(j=0;j<448;j++) {
  214.     if(j<64) {
  215. buffer[0][j] = buffer[0][j+384];
  216. buffer[1][j] = buffer[1][j+384];
  217.     }
  218.     else {
  219. buffer[0][j] = insamp[2*j-128];
  220. buffer[1][j] = insamp[2*j-127];
  221.     }
  222. }
  223.     }
  224.     else { /* layer 1, mono */
  225. samples_read = read_samples(musicin, insamp, num_samples,
  226.     (unsigned long) 384);
  227. for(j=0;j<448;j++){
  228.     if(j<64) {
  229. buffer[0][j] = buffer[0][j+384];
  230. buffer[1][j] = 0;
  231.     }
  232.     else {
  233. buffer[0][j] = insamp[j-64];
  234. buffer[1][j] = 0;
  235.     }
  236. }
  237.     }
  238. }
  239. else {
  240.     if(stereo == 2){ /* layer 2 (or 3), stereo */
  241. samples_read = read_samples(musicin, insamp, num_samples,
  242.     (unsigned long) 2304);
  243. for(j=0;j<1152;j++) {
  244.     buffer[0][j] = insamp[2*j];
  245.     buffer[1][j] = insamp[2*j+1];
  246. }
  247.     }
  248.     else { /* layer 2 (or 3), mono */
  249. samples_read = read_samples(musicin, insamp, num_samples,
  250.     (unsigned long) 1152);
  251. for(j=0;j<1152;j++){
  252.     buffer[0][j] = insamp[j];
  253.     buffer[1][j] = 0;
  254. }
  255.     }
  256. }
  257.     }
  258.     return(samples_read);
  259. }
  260.  
  261. /************************************************************************
  262. *
  263. * read_ana_window()
  264. *
  265. * PURPOSE:  Reads encoder window file "enwindow" into array #ana_win#
  266. *
  267. ************************************************************************/
  268.  
  269. void read_ana_window(ana_win)
  270. double FAR ana_win[HAN_SIZE];
  271. {
  272.     int i,j[4];
  273.     FILE *fp;
  274.     double f[4];
  275.     char t[150];
  276.  
  277.     if (!(fp = OpenTableFile("enwindow") ) ) {
  278.        printf("Please check analysis window table 'enwindow'n");
  279.        exit(1);
  280.     }
  281.     for (i=0;i<512;i+=4) {
  282.        fgets(t, 150, fp);
  283.        sscanf(t,"C[%d] = %lf C[%d] = %lf C[%d] = %lf C[%d] = %lfn",
  284.               j, f,j+1,f+1,j+2,f+2,j+3,f+3);
  285.        if (i==j[0]) {
  286.           ana_win[i] = f[0];
  287.           ana_win[i+1] = f[1];
  288.           ana_win[i+2] = f[2];
  289.           ana_win[i+3] = f[3];
  290.        }
  291.        else {
  292.           printf("Check index in analysis window tablen");
  293.           exit(1);
  294.        }
  295.        fgets(t,150,fp);
  296.     }
  297.     fclose(fp);
  298. }
  299. /************************************************************************
  300. *
  301. * window_subband()
  302. *
  303. * PURPOSE:  Overlapping window on PCM samples
  304. *
  305. * SEMANTICS:
  306. * 32 16-bit pcm samples are scaled to fractional 2's complement and
  307. * concatenated to the end of the window buffer #x#. The updated window
  308. * buffer #x# is then windowed by the analysis window #c# to produce the
  309. * windowed sample #z#
  310. *
  311. ************************************************************************/
  312.  
  313. void window_subband(buffer, z, k)
  314. short FAR **buffer;
  315. double FAR z[HAN_SIZE];
  316. int k;
  317. {
  318.     typedef double FAR XX[2][HAN_SIZE];
  319.     static XX FAR *x;
  320.     int i, j;
  321.     static off[2] = {0,0};
  322.     static char init = 0;
  323.     static double FAR *c;
  324.     if (!init) {
  325.         c = (double FAR *) mem_alloc(sizeof(double) * HAN_SIZE, "window");
  326.         read_ana_window(c);
  327.         x = (XX FAR *) mem_alloc(sizeof(XX),"x");
  328.         for (i=0;i<2;i++)
  329.             for (j=0;j<HAN_SIZE;j++)
  330.                 (*x)[i][j] = 0;
  331.         init = 1;
  332.     }
  333.     /* replace 32 oldest samples with 32 new samples */
  334.     for (i=0;i<32;i++) (*x)[k][31-i+off[k]] = (double) *(*buffer)++/SCALE;
  335.     /* shift samples into proper window positions */
  336.     for (i=0;i<HAN_SIZE;i++) z[i] = (*x)[k][(i+off[k])&HAN_SIZE-1] * c[i];
  337.     off[k] += 480;              /*offset is modulo (HAN_SIZE-1)*/
  338.     off[k] &= HAN_SIZE-1;
  339. }
  340.  
  341. /************************************************************************
  342. *
  343. * create_ana_filter()
  344. *
  345. * PURPOSE:  Calculates the analysis filter bank coefficients
  346. *
  347. * SEMANTICS:
  348. * Calculates the analysis filterbank coefficients and rounds to the
  349. * 9th decimal place accuracy of the filterbank tables in the ISO
  350. * document.  The coefficients are stored in #filter#
  351. ************************************************************************/
  352.  
  353. void create_ana_filter(filter)
  354. double FAR filter[SBLIMIT][64];
  355. {
  356.    register int i,k;
  357.  
  358.    for (i=0; i<32; i++)
  359.       for (k=0; k<64; k++) {
  360.           if ((filter[i][k] = 1e9*cos((double)((2*i+1)*(16-k)*PI64))) >= 0)
  361.              modf(filter[i][k]+0.5, &filter[i][k]);
  362.           else
  363.              modf(filter[i][k]-0.5, &filter[i][k]);
  364.           filter[i][k] *= 1e-9;
  365.    }
  366. }
  367. /************************************************************************
  368. *
  369. * filter_subband()
  370. *
  371. * PURPOSE:  Calculates the analysis filter bank coefficients
  372. *
  373. * SEMANTICS:
  374. *      The windowed samples #z# is filtered by the digital filter matrix #m#
  375. * to produce the subband samples #s#. This done by first selectively
  376. * picking out values from the windowed samples, and then multiplying
  377. * them by the filter matrix, producing 32 subband samples.
  378. *
  379. ************************************************************************/
  380.  
  381. void filter_subband(z,s)
  382. double FAR z[HAN_SIZE], s[SBLIMIT];
  383. {
  384.    double y[64];
  385.    int i,j;
  386. static char init = 0;
  387.    typedef double MM[SBLIMIT][64];
  388. static MM FAR *m;
  389. #ifdef MS_DOS
  390.    long    SIZE_OF_MM;
  391.    SIZE_OF_MM      = SBLIMIT*64;
  392.    SIZE_OF_MM      *= 8;
  393.    if (!init) {
  394.        m = (MM FAR *) mem_alloc(SIZE_OF_MM, "filter");
  395.        create_ana_filter(*m);
  396.        init = 1;
  397.    }
  398. #else
  399.    if (!init) {
  400.        m = (MM FAR *) mem_alloc(sizeof(MM), "filter");
  401.        create_ana_filter(*m);
  402.        init = 1;
  403.    }
  404. #endif
  405.    for (i=0;i<64;i++) for (j=0, y[i] = 0;j<8;j++) y[i] += z[i+64*j];
  406.    for (i=0;i<SBLIMIT;i++)
  407.        for (j=0, s[i]= 0;j<64;j++) s[i] += (*m)[i][j] * y[j];
  408. }
  409. /************************************************************************
  410. * encode_info()
  411. *
  412. * PURPOSE:  Puts the syncword and header information on the output
  413. * bitstream.
  414. *
  415. ************************************************************************/
  416.  
  417. void encode_info(fr_ps,bs)
  418. frame_params *fr_ps;
  419. Bit_stream_struc *bs;
  420. {
  421.         layer *info = fr_ps->header;
  422.  
  423.         putbits(bs,0xfff,12);                    /* syncword 12 bits */
  424.         put1bit(bs,info->version);               /* ID        1 bit  */
  425.         putbits(bs,4-info->lay,2);               /* layer     2 bits */
  426.         put1bit(bs,!info->error_protection);     /* bit set => no err prot */
  427.         putbits(bs,info->bitrate_index,4);
  428.         putbits(bs,info->sampling_frequency,2);
  429.         put1bit(bs,info->padding);
  430.         put1bit(bs,info->extension);             /* private_bit */
  431.         putbits(bs,info->mode,2);
  432.         putbits(bs,info->mode_ext,2);
  433.         put1bit(bs,info->copyright);
  434.         put1bit(bs,info->original);
  435.         putbits(bs,info->emphasis,2);
  436. }
  437.  
  438. /************************************************************************
  439. *
  440. * mod()
  441. *
  442. * PURPOSE:  Returns the absolute value of its argument
  443. *
  444. ************************************************************************/
  445.  
  446. double mod(a)
  447. double a;
  448. {
  449.     return (a > 0) ? a : -a;
  450. }
  451.  
  452. /************************************************************************
  453. *
  454. * I_combine_LR    (Layer I)
  455. * II_combine_LR   (Layer II)
  456. *
  457. * PURPOSE:Combines left and right channels into a mono channel
  458. *
  459. * SEMANTICS:  The average of left and right subband samples is put into
  460. * #joint_sample#
  461. *
  462. * Layer I and II differ in frame length and # subbands used
  463. *
  464. ************************************************************************/
  465.  
  466. void I_combine_LR(sb_sample, joint_sample)
  467. double FAR sb_sample[2][3][SCALE_BLOCK][SBLIMIT];
  468. double FAR joint_sample[3][SCALE_BLOCK][SBLIMIT];
  469. {   /* make a filtered mono for joint stereo */
  470.     int sb, smp;
  471.  
  472.    for(sb = 0; sb<SBLIMIT; ++sb)
  473.       for(smp = 0; smp<SCALE_BLOCK; ++smp)
  474.         joint_sample[0][smp][sb] = .5 *
  475.                     (sb_sample[0][0][smp][sb] + sb_sample[1][0][smp][sb]);
  476. }
  477.  
  478. void II_combine_LR(sb_sample, joint_sample, sblimit)
  479. double FAR sb_sample[2][3][SCALE_BLOCK][SBLIMIT];
  480. double FAR joint_sample[3][SCALE_BLOCK][SBLIMIT];
  481. int sblimit;
  482. {  /* make a filtered mono for joint stereo */
  483.    int sb, smp, sufr;
  484.  
  485.    for(sb = 0; sb<sblimit; ++sb)
  486.       for(smp = 0; smp<SCALE_BLOCK; ++smp)
  487.          for(sufr = 0; sufr<3; ++sufr)
  488.             joint_sample[sufr][smp][sb] = .5 * (sb_sample[0][sufr][smp][sb]
  489.                                            + sb_sample[1][sufr][smp][sb]);
  490. }
  491.  
  492. /************************************************************************
  493. *
  494. * I_scale_factor_calc     (Layer I)
  495. * II_scale_factor_calc    (Layer II)
  496. *
  497. * PURPOSE:For each subband, calculate the scale factor for each set
  498. * of the 12 subband samples
  499. *
  500. * SEMANTICS:  Pick the scalefactor #multiple[]# just larger than the
  501. * absolute value of the peak subband sample of 12 samples,
  502. * and store the corresponding scalefactor index in #scalar#.
  503. *
  504. * Layer II has three sets of 12-subband samples for a given
  505. * subband.
  506. *
  507. ************************************************************************/
  508.  
  509. void I_scale_factor_calc(sb_sample,scalar,stereo)
  510. double FAR sb_sample[][3][SCALE_BLOCK][SBLIMIT];
  511. unsigned int scalar[][3][SBLIMIT];
  512. int stereo;
  513. {
  514.    int i,j, k;
  515.    double s[SBLIMIT];
  516.  
  517.    for (k=0;k<stereo;k++) {
  518.      for (i=0;i<SBLIMIT;i++)
  519.        for (j=1, s[i] = mod(sb_sample[k][0][0][i]);j<SCALE_BLOCK;j++)
  520.          if (mod(sb_sample[k][0][j][i]) > s[i])
  521.             s[i] = mod(sb_sample[k][0][j][i]);
  522.  
  523.      for (i=0;i<SBLIMIT;i++)
  524.        for (j=SCALE_RANGE-2,scalar[k][0][i]=0;j>=0;j--) /* $A 6/16/92 */
  525.          if (s[i] <= multiple[j]) {
  526.             scalar[k][0][i] = j;
  527.             break;
  528.          }
  529.    }
  530. }
  531. /******************************** Layer II ******************************/
  532.  
  533. void II_scale_factor_calc(sb_sample,scalar,stereo,sblimit)
  534. double FAR sb_sample[][3][SCALE_BLOCK][SBLIMIT];
  535. unsigned int scalar[][3][SBLIMIT];
  536. int stereo,sblimit;
  537. {
  538.   int i,j, k,t;
  539.   double s[SBLIMIT];
  540.  
  541.   for (k=0;k<stereo;k++) for (t=0;t<3;t++) {
  542.     for (i=0;i<sblimit;i++)
  543.       for (j=1, s[i] = mod(sb_sample[k][t][0][i]);j<SCALE_BLOCK;j++)
  544.         if (mod(sb_sample[k][t][j][i]) > s[i])
  545.              s[i] = mod(sb_sample[k][t][j][i]);
  546.  
  547.   for (i=0;i<sblimit;i++)
  548.     for (j=SCALE_RANGE-2,scalar[k][t][i]=0;j>=0;j--)    /* $A 6/16/92 */
  549.       if (s[i] <= multiple[j]) {
  550.          scalar[k][t][i] = j;
  551.          break;
  552.       }
  553.       for (i=sblimit;i<SBLIMIT;i++) scalar[k][t][i] = SCALE_RANGE-1;
  554.     }
  555. }
  556. /************************************************************************
  557. *
  558. * pick_scale  (Layer II)
  559. *
  560. * PURPOSE:For each subband, puts the smallest scalefactor of the 3
  561. * associated with a frame into #max_sc#.  This is used
  562. * used by Psychoacoustic Model I.
  563. * (I would recommend changin max_sc to min_sc)
  564. *
  565. ************************************************************************/
  566.  
  567. void pick_scale(scalar, fr_ps, max_sc)
  568. unsigned int scalar[2][3][SBLIMIT];
  569. frame_params *fr_ps;
  570. double FAR max_sc[2][SBLIMIT];
  571. {
  572.   int i,j,k,max;
  573.   int stereo  = fr_ps->stereo;
  574.   int sblimit = fr_ps->sblimit;
  575.  
  576.   for (k=0;k<stereo;k++)
  577.     for (i=0;i<sblimit;max_sc[k][i] = multiple[max],i++)
  578.       for (j=1, max = scalar[k][0][i];j<3;j++)
  579.          if (max > scalar[k][j][i]) max = scalar[k][j][i];
  580.   for (i=sblimit;i<SBLIMIT;i++) max_sc[0][i] = max_sc[1][i] = 1E-20;
  581. }
  582. /************************************************************************
  583. *
  584. * put_scale   (Layer I)
  585. *
  586. * PURPOSE:Sets #max_sc# to the scalefactor index in #scalar.
  587. * This is used by Psychoacoustic Model I
  588. *
  589. ************************************************************************/
  590.  
  591. void put_scale(scalar, fr_ps, max_sc)
  592. unsigned int scalar[2][3][SBLIMIT];
  593. frame_params *fr_ps;
  594. double FAR max_sc[2][SBLIMIT];
  595. {
  596.    int i,j,k, max;
  597.    int stereo  = fr_ps->stereo;
  598.    int sblimit = fr_ps->sblimit;
  599.  
  600.    for (k=0;k<stereo;k++) for (i=0;i<SBLIMIT;i++)
  601.         max_sc[k][i] = multiple[scalar[k][0][i]];
  602. }
  603.  
  604. /************************************************************************
  605. *
  606. * II_transmission_pattern (Layer II only)
  607. *
  608. * PURPOSE:For a given subband, determines whether to send 1, 2, or
  609. * all 3 of the scalefactors, and fills in the scalefactor
  610. * select information accordingly
  611. *
  612. * SEMANTICS:  The subbands and channels are classified based on how much
  613. * the scalefactors changes over its three values (corresponding
  614. * to the 3 sets of 12 samples per subband).  The classification
  615. * will send 1 or 2 scalefactors instead of three if the scalefactors
  616. * do not change much.  The scalefactor select information,
  617. * #scfsi#, is filled in accordingly.
  618. *
  619. ************************************************************************/
  620.  
  621. void II_transmission_pattern(scalar, scfsi, fr_ps)
  622. unsigned int scalar[2][3][SBLIMIT];
  623. unsigned int scfsi[2][SBLIMIT];
  624. frame_params *fr_ps;
  625. {
  626.    int stereo  = fr_ps->stereo;
  627.    int sblimit = fr_ps->sblimit;
  628.    int dscf[2];
  629.    int class[2],i,j,k;
  630. static int pattern[5][5] = {0x123, 0x122, 0x122, 0x133, 0x123,
  631.                             0x113, 0x111, 0x111, 0x444, 0x113,
  632.                             0x111, 0x111, 0x111, 0x333, 0x113,
  633.                             0x222, 0x222, 0x222, 0x333, 0x123,
  634.                             0x123, 0x122, 0x122, 0x133, 0x123};
  635.  
  636.    for (k=0;k<stereo;k++)
  637.      for (i=0;i<sblimit;i++) {
  638.        dscf[0] =  (scalar[k][0][i]-scalar[k][1][i]);
  639.        dscf[1] =  (scalar[k][1][i]-scalar[k][2][i]);
  640.        for (j=0;j<2;j++) {
  641.          if (dscf[j]<=-3) class[j] = 0;
  642.          else if (dscf[j] > -3 && dscf[j] <0) class[j] = 1;
  643.               else if (dscf[j] == 0) class[j] = 2;
  644.                    else if (dscf[j] > 0 && dscf[j] < 3) class[j] = 3;
  645.                         else class[j] = 4;
  646.        }
  647.        switch (pattern[class[0]][class[1]]) {
  648.          case 0x123 :    scfsi[k][i] = 0;
  649.                          break;
  650.          case 0x122 :    scfsi[k][i] = 3;
  651.                          scalar[k][2][i] = scalar[k][1][i];
  652.                          break;
  653.          case 0x133 :    scfsi[k][i] = 3;
  654.                          scalar[k][1][i] = scalar[k][2][i];
  655.                          break;
  656.          case 0x113 :    scfsi[k][i] = 1;
  657.                          scalar[k][1][i] = scalar[k][0][i];
  658.                          break;
  659.          case 0x111 :    scfsi[k][i] = 2;
  660.                          scalar[k][1][i] = scalar[k][2][i] = scalar[k][0][i];
  661.                          break;
  662.          case 0x222 :    scfsi[k][i] = 2;
  663.                          scalar[k][0][i] = scalar[k][2][i] = scalar[k][1][i];
  664.                          break;
  665.          case 0x333 :    scfsi[k][i] = 2;
  666.                          scalar[k][0][i] = scalar[k][1][i] = scalar[k][2][i];
  667.                          break;
  668.          case 0x444 :    scfsi[k][i] = 2;
  669.                          if (scalar[k][0][i] > scalar[k][2][i])
  670.                               scalar[k][0][i] = scalar[k][2][i];
  671.                          scalar[k][1][i] = scalar[k][2][i] = scalar[k][0][i];
  672.       }
  673.    }
  674. }
  675.  
  676. /************************************************************************
  677. *
  678. * I_encode_scale  (Layer I)
  679. * II_encode_scale (Layer II)
  680. *
  681. * PURPOSE:The encoded scalar factor information is arranged and
  682. * queued into the output fifo to be transmitted.
  683. *
  684. * For Layer II, the three scale factors associated with
  685. * a given subband and channel are transmitted in accordance
  686. * with the scfsi, which is transmitted first.
  687. *
  688. ************************************************************************/
  689.  
  690. void I_encode_scale(scalar, bit_alloc, fr_ps, bs)
  691. unsigned int scalar[2][3][SBLIMIT];
  692. unsigned int bit_alloc[2][SBLIMIT];
  693. frame_params *fr_ps;
  694. Bit_stream_struc *bs;
  695. {
  696.    int stereo  = fr_ps->stereo;
  697.    int sblimit = fr_ps->sblimit;
  698.    int i,j;
  699.  
  700.    for (i=0;i<SBLIMIT;i++) for (j=0;j<stereo;j++)
  701.       if (bit_alloc[j][i]) putbits(bs,scalar[j][0][i],6);
  702. }
  703.  
  704. /***************************** Layer II  ********************************/
  705.  
  706. void II_encode_scale(bit_alloc, scfsi, scalar, fr_ps, bs)
  707. unsigned int bit_alloc[2][SBLIMIT], scfsi[2][SBLIMIT];
  708. unsigned int scalar[2][3][SBLIMIT];
  709. frame_params *fr_ps;
  710. Bit_stream_struc *bs;
  711. {
  712.    int stereo  = fr_ps->stereo;
  713.    int sblimit = fr_ps->sblimit;
  714.    int jsbound = fr_ps->jsbound;
  715.    int i,j,k;
  716.  
  717.    for (i=0;i<sblimit;i++) for (k=0;k<stereo;k++)
  718.      if (bit_alloc[k][i])  putbits(bs,scfsi[k][i],2);
  719.  
  720.    for (i=0;i<sblimit;i++) for (k=0;k<stereo;k++)
  721.      if (bit_alloc[k][i])  /* above jsbound, bit_alloc[0][i] == ba[1][i] */
  722.         switch (scfsi[k][i]) {
  723.            case 0: for (j=0;j<3;j++)
  724.                      putbits(bs,scalar[k][j][i],6);
  725.                    break;
  726.            case 1:
  727.            case 3: putbits(bs,scalar[k][0][i],6);
  728.                    putbits(bs,scalar[k][2][i],6);
  729.                    break;
  730.            case 2: putbits(bs,scalar[k][0][i],6);
  731.         }
  732. }
  733.  
  734. /*=======================================================================
  735. |                                                                        |
  736. |      The following routines are done after the masking threshold       |
  737. | has been calculated by the fft analysis routines in the Psychoacoustic |
  738. | model. Using the MNR calculated, the actual number of bits allocated   |
  739. | to each subband is found iteratively.                                  |
  740. |                                                                        |
  741. =======================================================================*/
  742.  
  743. /************************************************************************
  744. *
  745. * I_bits_for_nonoise  (Layer I)
  746. * II_bits_for_nonoise (Layer II)
  747. *
  748. * PURPOSE:Returns the number of bits required to produce a
  749. * mask-to-noise ratio better or equal to the noise/no_noise threshold.
  750. *
  751. * SEMANTICS:
  752. * bbal = # bits needed for encoding bit allocation
  753. * bsel = # bits needed for encoding scalefactor select information
  754. * banc = # bits needed for ancillary data (header info included)
  755. *
  756. * For each subband and channel, will add bits until one of the
  757. * following occurs:
  758. * - Hit maximum number of bits we can allocate for that subband
  759. * - MNR is better than or equal to the minimum masking level
  760. *   (NOISY_MIN_MNR)
  761. * Then the bits required for scalefactors, scfsi, bit allocation,
  762. * and the subband samples are tallied (#req_bits#) and returned.
  763. *
  764. * (NOISY_MIN_MNR) is the smallest MNR a subband can have before it is
  765. * counted as 'noisy' by the logic which chooses the number of JS
  766. * subbands.
  767. *
  768. * Joint stereo is supported.
  769. *
  770. ************************************************************************/
  771. static double snr[18] = {0.00, 7.00, 11.00, 16.00, 20.84,
  772.                          25.28, 31.59, 37.75, 43.84,
  773.                          49.89, 55.93, 61.96, 67.98, 74.01,
  774.                          80.03, 86.05, 92.01, 98.01};
  775. int I_bits_for_nonoise(perm_smr, fr_ps)
  776. double FAR perm_smr[2][SBLIMIT];
  777. frame_params *fr_ps;
  778. {
  779.    int i,j,k;
  780.    int stereo  = fr_ps->stereo;
  781.    int sblimit = fr_ps->sblimit;
  782.    int jsbound = fr_ps->jsbound;
  783.    int req_bits = 0;
  784.  
  785.    /* initial b_anc (header) allocation bits */
  786.    req_bits = 32 + 4 * ( (jsbound * stereo) + (SBLIMIT-jsbound) );
  787.  
  788.    for(i=0; i<SBLIMIT; ++i)
  789.      for(j=0; j<((i<jsbound)?stereo:1); ++j) {
  790.        for(k=0;k<14; ++k)
  791.          if( (-perm_smr[j][i] + snr[k]) >= NOISY_MIN_MNR)
  792.            break; /* we found enough bits */
  793.          if(stereo == 2 && i >= jsbound)     /* check other JS channel */
  794.            for(;k<14; ++k)
  795.              if( (-perm_smr[1-j][i] + snr[k]) >= NOISY_MIN_MNR) break;
  796.          if(k>0) req_bits += (k+1)*SCALE_BLOCK + 6*((i>=jsbound)?stereo:1);
  797.    }
  798.    return req_bits;
  799. }
  800.  
  801. /***************************** Layer II  ********************************/
  802.  
  803. int II_bits_for_nonoise(perm_smr, scfsi, fr_ps)
  804. double FAR perm_smr[2][SBLIMIT];
  805. unsigned int scfsi[2][SBLIMIT];
  806. frame_params *fr_ps;
  807. {
  808.    int sb,ch,ba;
  809.    int stereo  = fr_ps->stereo;
  810.    int sblimit = fr_ps->sblimit;
  811.    int jsbound = fr_ps->jsbound;
  812.    al_table *alloc = fr_ps->alloc;
  813.    int req_bits = 0, bbal = 0, berr = 0, banc = 32;
  814.    int maxAlloc, sel_bits, sc_bits, smp_bits;
  815. static int sfsPerScfsi[] = { 3,2,1,2 };    /* lookup # sfs per scfsi */
  816.    /* added 92-08-11 shn */
  817.    if (fr_ps->header->error_protection) berr=16; else berr=0; 
  818.  
  819.    for (sb=0; sb<jsbound; ++sb)
  820.      bbal += stereo * (*alloc)[sb][0].bits;
  821.    for (sb=jsbound; sb<sblimit; ++sb)
  822.      bbal += (*alloc)[sb][0].bits;
  823.    req_bits = banc + bbal + berr;
  824.  
  825.    for(sb=0; sb<sblimit; ++sb)
  826.      for(ch=0; ch<((sb<jsbound)?stereo:1); ++ch) {
  827.        maxAlloc = (1<<(*alloc)[sb][0].bits)-1;
  828.        sel_bits = sc_bits = smp_bits = 0;
  829.        for(ba=0;ba<maxAlloc-1; ++ba)
  830.          if( (-perm_smr[ch][sb] + snr[(*alloc)[sb][ba].quant+((ba>0)?1:0)])
  831.              >= NOISY_MIN_MNR)
  832.             break;      /* we found enough bits */
  833.        if(stereo == 2 && sb >= jsbound) /* check other JS channel */
  834.          for(;ba<maxAlloc-1; ++ba)
  835.            if( (-perm_smr[1-ch][sb]+ snr[(*alloc)[sb][ba].quant+((ba>0)?1:0)])
  836.                >= NOISY_MIN_MNR)
  837.              break;
  838.        if(ba>0) {
  839.          smp_bits = SCALE_BLOCK * ((*alloc)[sb][ba].group * (*alloc)[sb][ba].bits);
  840.          /* scale factor bits required for subband */
  841.          sel_bits = 2;
  842.          sc_bits  = 6 * sfsPerScfsi[scfsi[ch][sb]];
  843.          if(stereo == 2 && sb >= jsbound) {
  844.            /* each new js sb has L+R scfsis */
  845.            sel_bits += 2;
  846.            sc_bits  += 6 * sfsPerScfsi[scfsi[1-ch][sb]];
  847.          }
  848.          req_bits += smp_bits+sel_bits+sc_bits;
  849.        }
  850.    }
  851.    return req_bits;
  852. }
  853.  
  854. /************************************************************************
  855. *
  856. * I_main_bit_allocation   (Layer I)
  857. * II_main_bit_allocation  (Layer II)
  858. *
  859. * PURPOSE:For joint stereo mode, determines which of the 4 joint
  860. * stereo modes is needed.  Then calls *_a_bit_allocation(), which
  861. * allocates bits for each of the subbands until there are no more bits
  862. * left, or the MNR is at the noise/no_noise threshold.
  863. *
  864. * SEMANTICS:
  865. *
  866. * For joint stereo mode, joint stereo is changed to stereo if
  867. * there are enough bits to encode stereo at or better than the
  868. * no-noise threshold (NOISY_MIN_MNR).  Otherwise, the system
  869. * iteratively allocates less bits by using joint stereo until one
  870. * of the following occurs:
  871. * - there are no more noisy subbands (MNR >= NOISY_MIN_MNR)
  872. * - mode_ext has been reduced to 0, which means that all but the
  873. *   lowest 4 subbands have been converted from stereo to joint
  874. *   stereo, and no more subbands may be converted
  875. *
  876. *     This function calls *_bits_for_nonoise() and *_a_bit_allocation().
  877. *
  878. ************************************************************************/
  879.  
  880. void I_main_bit_allocation(perm_smr, bit_alloc, adb, fr_ps)
  881. double FAR perm_smr[2][SBLIMIT];
  882. unsigned int bit_alloc[2][SBLIMIT];
  883. int *adb;
  884. frame_params *fr_ps;
  885. {
  886.    int  noisy_sbs;
  887.    int  mode, mode_ext, lay, i;
  888.    int  rq_db, av_db = *adb;
  889. static  int init = 0;
  890.  
  891.    if(init == 0) {
  892.      /* rearrange snr for layer I */
  893.      snr[2] = snr[3];
  894.      for (i=3;i<16;i++) snr[i] = snr[i+2];
  895.      init = 1;
  896.    }
  897.  
  898.    if((mode = fr_ps->actual_mode) == MPG_MD_JOINT_STEREO) {
  899.      fr_ps->header->mode = MPG_MD_STEREO;
  900.      fr_ps->header->mode_ext = 0;
  901.      fr_ps->jsbound = fr_ps->sblimit;
  902.      if(rq_db = I_bits_for_nonoise(perm_smr, fr_ps) > *adb) {
  903.        fr_ps->header->mode = MPG_MD_JOINT_STEREO;
  904.        mode_ext = 4;           /* 3 is least severe reduction */
  905.        lay = fr_ps->header->lay;
  906.        do {
  907.           --mode_ext;
  908.           fr_ps->jsbound = js_bound(lay, mode_ext);
  909.           rq_db = I_bits_for_nonoise(perm_smr, fr_ps);
  910.        } while( (rq_db > *adb) && (mode_ext > 0));
  911.        fr_ps->header->mode_ext = mode_ext;
  912.      }    /* well we either eliminated noisy sbs or mode_ext == 0 */
  913.    }
  914.    noisy_sbs = I_a_bit_allocation(perm_smr, bit_alloc, adb, fr_ps);
  915. }
  916.  
  917. /***************************** Layer II  ********************************/
  918.  
  919. void II_main_bit_allocation(perm_smr, scfsi, bit_alloc, adb, fr_ps)
  920. double FAR perm_smr[2][SBLIMIT];
  921. unsigned int scfsi[2][SBLIMIT];
  922. unsigned int bit_alloc[2][SBLIMIT];
  923. int *adb;
  924. frame_params *fr_ps;
  925. {
  926.    int  noisy_sbs, nn;
  927.    int  mode, mode_ext, lay;
  928.    int  rq_db, av_db = *adb;
  929.  
  930.    if((mode = fr_ps->actual_mode) == MPG_MD_JOINT_STEREO) {
  931.      fr_ps->header->mode = MPG_MD_STEREO;
  932.      fr_ps->header->mode_ext = 0;
  933.      fr_ps->jsbound = fr_ps->sblimit;
  934.      if((rq_db=II_bits_for_nonoise(perm_smr, scfsi, fr_ps)) > *adb) {
  935.        fr_ps->header->mode = MPG_MD_JOINT_STEREO;
  936.        mode_ext = 4;           /* 3 is least severe reduction */
  937.        lay = fr_ps->header->lay;
  938.        do {
  939.          --mode_ext;
  940.          fr_ps->jsbound = js_bound(lay, mode_ext);
  941.          rq_db = II_bits_for_nonoise(perm_smr, scfsi, fr_ps);
  942.        } while( (rq_db > *adb) && (mode_ext > 0));
  943.        fr_ps->header->mode_ext = mode_ext;
  944.      }    /* well we either eliminated noisy sbs or mode_ext == 0 */
  945.    }
  946.    noisy_sbs = II_a_bit_allocation(perm_smr, scfsi, bit_alloc, adb, fr_ps);
  947. }
  948.  
  949. /************************************************************************
  950. *
  951. * I_a_bit_allocation  (Layer I)
  952. * II_a_bit_allocation (Layer II)
  953. *
  954. * PURPOSE:Adds bits to the subbands with the lowest mask-to-noise
  955. * ratios, until the maximum number of bits for the subband has
  956. * been allocated.
  957. *
  958. * SEMANTICS:
  959. * 1. Find the subband and channel with the smallest MNR (#min_sb#,
  960. *    and #min_ch#)
  961. * 2. Calculate the increase in bits needed if we increase the bit
  962. *    allocation to the next higher level
  963. * 3. If there are enough bits available for increasing the resolution
  964. *    in #min_sb#, #min_ch#, and the subband has not yet reached its
  965. *    maximum allocation, update the bit allocation, MNR, and bits
  966.     available accordingly
  967. * 4. Repeat until there are no more bits left, or no more available
  968. *    subbands. (A subband is still available until the maximum
  969. *    number of bits for the subband has been allocated, or there
  970. *    aren't enough bits to go to the next higher resolution in the
  971.     subband.)
  972. *
  973. ************************************************************************/
  974.  
  975. int I_a_bit_allocation(perm_smr, bit_alloc, adb, fr_ps) /* return noisy sbs */
  976. double FAR perm_smr[2][SBLIMIT];
  977. unsigned int bit_alloc[2][SBLIMIT];
  978. int *adb;
  979. frame_params *fr_ps;
  980. {
  981.    int i, k, smpl_bits, scale_bits, min_sb, min_ch, oth_ch;
  982.    int bspl, bscf, ad, noisy_sbs, done = 0, bbal ;
  983.    double mnr[2][SBLIMIT], small;
  984.    char used[2][SBLIMIT];
  985.    int stereo  = fr_ps->stereo;
  986.    int sblimit = fr_ps->sblimit;
  987.    int jsbound = fr_ps->jsbound;
  988.    al_table *alloc = fr_ps->alloc;
  989. static char init= 0;
  990. static int banc=32, berr=0;
  991.  
  992.    if (!init) {
  993.       init = 1;
  994.       if (fr_ps->header->error_protection) berr = 16;  /* added 92-08-11 shn */
  995.    }
  996.    bbal = 4 * ( (jsbound * stereo) + (SBLIMIT-jsbound) );
  997.    *adb -= bbal + berr + banc;
  998.    ad= *adb;
  999.  
  1000.    for (i=0;i<SBLIMIT;i++) for (k=0;k<stereo;k++) {
  1001.      mnr[k][i]=snr[0]-perm_smr[k][i];
  1002.      bit_alloc[k][i] = 0;
  1003.      used[k][i] = 0;
  1004.    }
  1005.    bspl = bscf = 0;
  1006.  
  1007.    do  {
  1008.      /* locate the subband with minimum SMR */
  1009.      small = mnr[0][0]+1;    min_sb = -1; min_ch = -1;
  1010.      for (i=0;i<SBLIMIT;i++) for (k=0;k<stereo;k++)
  1011.        /* go on only if there are bits left */
  1012.        if (used[k][i] != 2 && small > mnr[k][i]) {
  1013.          small = mnr[k][i];
  1014.          min_sb = i;  min_ch = k;
  1015.        }
  1016.      if(min_sb > -1) {   /* there was something to find */
  1017.        /* first step of bit allocation is biggest */
  1018.        if (used[min_ch][min_sb])  { smpl_bits = SCALE_BLOCK; scale_bits = 0; }
  1019.        else                       { smpl_bits = 24; scale_bits = 6; }
  1020.        if(min_sb >= jsbound)        scale_bits *= stereo;
  1021.  
  1022.        /* check to see enough bits were available for */
  1023.        /* increasing resolution in the minimum band */
  1024.  
  1025.        if (ad >= bspl + bscf + scale_bits + smpl_bits) {
  1026.          bspl += smpl_bits; /* bit for subband sample */
  1027.          bscf += scale_bits; /* bit for scale factor */
  1028.          bit_alloc[min_ch][min_sb]++;
  1029.          used[min_ch][min_sb] = 1; /* subband has bits */
  1030.          mnr[min_ch][min_sb] = -perm_smr[min_ch][min_sb]
  1031.                                + snr[bit_alloc[min_ch][min_sb]];
  1032.          /* Check if subband has been fully allocated max bits */
  1033.          if (bit_alloc[min_ch][min_sb] ==  14 ) used[min_ch][min_sb] = 2;
  1034.        }
  1035.        else            /* no room to improve this band */
  1036.          used[min_ch][min_sb] = 2; /*   for allocation anymore */
  1037.        if(stereo == 2 && min_sb >= jsbound) {
  1038.          oth_ch = 1-min_ch;  /* joint-st : fix other ch */
  1039.          bit_alloc[oth_ch][min_sb] = bit_alloc[min_ch][min_sb];
  1040.          used[oth_ch][min_sb] = used[min_ch][min_sb];
  1041.          mnr[oth_ch][min_sb] = -perm_smr[oth_ch][min_sb]
  1042.                                + snr[bit_alloc[oth_ch][min_sb]];
  1043.        }
  1044.      }
  1045.    } while(min_sb>-1);     /* i.e. still some sub-bands to find */
  1046.    /* Calculate the number of bits left, add on to pointed var */
  1047.    ad -= bspl+bscf;
  1048.    *adb = ad;
  1049.    /* see how many channels are noisy */
  1050.    noisy_sbs = 0; small = mnr[0][0];
  1051.    for(k=0; k<stereo; ++k) {
  1052.      for(i = 0; i< SBLIMIT; ++i) {
  1053.        if(mnr[k][i] < NOISY_MIN_MNR)   ++noisy_sbs;
  1054.        if(small > mnr[k][i])           small = mnr[k][i];
  1055.      }
  1056.    }
  1057.    return noisy_sbs;
  1058. }
  1059. /***************************** Layer II  ********************************/
  1060.  
  1061. int II_a_bit_allocation(perm_smr, scfsi, bit_alloc, adb, fr_ps)
  1062. double FAR perm_smr[2][SBLIMIT];
  1063. unsigned int scfsi[2][SBLIMIT];
  1064. unsigned int bit_alloc[2][SBLIMIT];
  1065. int *adb;
  1066. frame_params *fr_ps;
  1067. {
  1068.    int i, min_ch, min_sb, oth_ch, k, increment, scale, seli, ba;
  1069.    int bspl, bscf, bsel, ad, noisy_sbs, bbal=0;
  1070.    double mnr[2][SBLIMIT], small;
  1071.    char used[2][SBLIMIT];
  1072.    int stereo  = fr_ps->stereo;
  1073.    int sblimit = fr_ps->sblimit;
  1074.    int jsbound = fr_ps->jsbound;
  1075.    al_table *alloc = fr_ps->alloc;
  1076. static char init= 0;
  1077. static int banc=32, berr=0;
  1078. static int sfsPerScfsi[] = { 3,2,1,2 };    /* lookup # sfs per scfsi */
  1079.  
  1080.    if (!init) { 
  1081.        init = 1;  
  1082.        if (fr_ps->header->error_protection) berr=16; /* added 92-08-11 shn */
  1083.    }
  1084.    for (i=0; i<jsbound; ++i)
  1085.      bbal += stereo * (*alloc)[i][0].bits;
  1086.    for (i=jsbound; i<sblimit; ++i)
  1087.      bbal += (*alloc)[i][0].bits;
  1088.    *adb -= bbal + berr + banc;
  1089.    ad = *adb;
  1090.  
  1091.    for (i=0;i<sblimit;i++) for (k=0;k<stereo;k++) {
  1092.      mnr[k][i]=snr[0]-perm_smr[k][i];
  1093.      bit_alloc[k][i] = 0;
  1094.      used[k][i] = 0;
  1095.    }
  1096.    bspl = bscf = bsel = 0;
  1097.  
  1098.    do  {
  1099.      /* locate the subband with minimum SMR */
  1100.      small = 999999.0; min_sb = -1; min_ch = -1;
  1101.      for (i=0;i<sblimit;i++) for(k=0;k<stereo;++k)
  1102.        if (used[k][i]  != 2 && small > mnr[k][i]) {
  1103.          small = mnr[k][i];
  1104.          min_sb = i;  min_ch = k;
  1105.      }
  1106.      if(min_sb > -1) {   /* there was something to find */
  1107.        /* find increase in bit allocation in subband [min] */
  1108.        increment = SCALE_BLOCK * ((*alloc)[min_sb][bit_alloc[min_ch][min_sb]+1].group *
  1109.                         (*alloc)[min_sb][bit_alloc[min_ch][min_sb]+1].bits);
  1110.        if (used[min_ch][min_sb])
  1111.          increment -= SCALE_BLOCK * ((*alloc)[min_sb][bit_alloc[min_ch][min_sb]].group*
  1112.                            (*alloc)[min_sb][bit_alloc[min_ch][min_sb]].bits);
  1113.  
  1114.        /* scale factor bits required for subband [min] */
  1115.        oth_ch = 1 - min_ch;    /* above js bound, need both chans */
  1116.        if (used[min_ch][min_sb]) scale = seli = 0;
  1117.        else {          /* this channel had no bits or scfs before */
  1118.          seli = 2;
  1119.          scale = 6 * sfsPerScfsi[scfsi[min_ch][min_sb]];
  1120.          if(stereo == 2 && min_sb >= jsbound) {
  1121.            /* each new js sb has L+R scfsis */
  1122.            seli += 2;
  1123.            scale += 6 * sfsPerScfsi[scfsi[oth_ch][min_sb]];
  1124.          }
  1125.        }
  1126.        /* check to see enough bits were available for */
  1127.        /* increasing resolution in the minimum band */
  1128.        if (ad >= bspl + bscf + bsel + seli + scale + increment) {
  1129.          ba = ++bit_alloc[min_ch][min_sb]; /* next up alloc */
  1130.          bspl += increment;  /* bits for subband sample */
  1131.          bscf += scale;      /* bits for scale factor */
  1132.          bsel += seli;       /* bits for scfsi code */
  1133.          used[min_ch][min_sb] = 1; /* subband has bits */
  1134.          mnr[min_ch][min_sb] = -perm_smr[min_ch][min_sb] +
  1135.                                snr[(*alloc)[min_sb][ba].quant+1];
  1136.          /* Check if subband has been fully allocated max bits */
  1137.          if (ba >= (1<<(*alloc)[min_sb][0].bits)-1) used[min_ch][min_sb] = 2;
  1138.        }
  1139.        else used[min_ch][min_sb] = 2; /* can't increase this alloc */
  1140.        if(min_sb >= jsbound && stereo == 2) {
  1141.          /* above jsbound, alloc applies L+R */
  1142.          ba = bit_alloc[oth_ch][min_sb] = bit_alloc[min_ch][min_sb];
  1143.          used[oth_ch][min_sb] = used[min_ch][min_sb];
  1144.          mnr[oth_ch][min_sb] = -perm_smr[oth_ch][min_sb] +
  1145.                                snr[(*alloc)[min_sb][ba].quant+1];
  1146.        }
  1147.      }
  1148.    } while(min_sb > -1);   /* until could find no channel */
  1149.    /* Calculate the number of bits left */
  1150.    ad -= bspl+bscf+bsel;   *adb = ad;
  1151.    for (i=sblimit;i<SBLIMIT;i++) for (k=0;k<stereo;k++) bit_alloc[k][i]=0;
  1152.  
  1153.    noisy_sbs = 0;  small = mnr[0][0];      /* calc worst noise in case */
  1154.    for(k=0;k<stereo;++k) {
  1155.      for (i=0;i<sblimit;i++) {
  1156.        if (small > mnr[k][i]) small = mnr[k][i];
  1157.        if(mnr[k][i] < NOISY_MIN_MNR) ++noisy_sbs; /* noise is not masked */
  1158.      }
  1159.    }
  1160.    return noisy_sbs;
  1161. }
  1162.  
  1163. /************************************************************************
  1164. *
  1165. * I_subband_quantization  (Layer I)
  1166. * II_subband_quantization (Layer II)
  1167. *
  1168. * PURPOSE:Quantizes subband samples to appropriate number of bits
  1169. *
  1170. * SEMANTICS:  Subband samples are divided by their scalefactors, which
  1171.  makes the quantization more efficient. The scaled samples are
  1172. * quantized by the function a*x+b, where a and b are functions of
  1173. * the number of quantization levels. The result is then truncated
  1174. * to the appropriate number of bits and the MSB is inverted.
  1175. *
  1176. * Note that for fractional 2's complement, inverting the MSB for a
  1177.  negative number x is equivalent to adding 1 to it.
  1178. *
  1179. ************************************************************************/
  1180.  
  1181. static double a[17] = {
  1182.   0.750000000, 0.625000000, 0.875000000, 0.562500000, 0.937500000,
  1183.   0.968750000, 0.984375000, 0.992187500, 0.996093750, 0.998046875,
  1184.   0.999023438, 0.999511719, 0.999755859, 0.999877930, 0.999938965,
  1185.   0.999969482, 0.999984741 };
  1186.  
  1187. static double b[17] = {
  1188.   -0.250000000, -0.375000000, -0.125000000, -0.437500000, -0.062500000,
  1189.   -0.031250000, -0.015625000, -0.007812500, -0.003906250, -0.001953125,
  1190.   -0.000976563, -0.000488281, -0.000244141, -0.000122070, -0.000061035,
  1191.   -0.000030518, -0.000015259 };
  1192.  
  1193. void I_subband_quantization(scalar, sb_samples, j_scale, j_samps,
  1194.                             bit_alloc, sbband, fr_ps)
  1195. unsigned int scalar[2][3][SBLIMIT];
  1196. double FAR sb_samples[2][3][SCALE_BLOCK][SBLIMIT];
  1197. unsigned int j_scale[3][SBLIMIT];
  1198. double FAR j_samps[3][SCALE_BLOCK][SBLIMIT]; /* L+R for j-stereo if necess */
  1199. unsigned int bit_alloc[2][SBLIMIT];
  1200. unsigned int FAR sbband[2][3][SCALE_BLOCK][SBLIMIT];
  1201. frame_params *fr_ps;
  1202. {
  1203.    int i, j, k, n, sig;
  1204.    int stereo  = fr_ps->stereo;
  1205.    int sblimit = fr_ps->sblimit;
  1206.    int jsbound = fr_ps->jsbound;
  1207.    double d;
  1208. static char init = 0;
  1209.    if (!init) {
  1210.      init = 1;
  1211.      /* rearrange quantization coef to correspond to layer I table */
  1212.      a[1] = a[2]; b[1] = b[2];
  1213.      for (i=2;i<15;i++) { a[i] = a[i+2]; b[i] = b[i+2]; }
  1214.    }
  1215.    for (j=0;j<SCALE_BLOCK;j++) for (i=0;i<SBLIMIT;i++)
  1216.      for (k=0;k<((i<jsbound)?stereo:1);k++)
  1217.        if (bit_alloc[k][i]) {
  1218.          /* for joint stereo mode, have to construct a single subband stream
  1219.             for the js channels.  At present, we calculate a set of mono
  1220.             subband samples and pass them through the scaling system to
  1221.             generate an alternate normalised sample stream.
  1222.  
  1223.             Could normalise both streams (divide by their scfs), then average
  1224.             them.  In bad conditions, this could give rise to spurious
  1225.             cancellations.  Instead, we could just select the sb stream from
  1226.             the larger channel (higher scf), in which case _that_ channel
  1227.             would be 'properly' reconstructed, and the mate would just be a
  1228.             scaled version.  Spec recommends averaging the two (unnormalised)
  1229.             subband channels, then normalising this new signal without
  1230.             actually sending this scale factor... This means looking ahead.
  1231.          */
  1232.          if(stereo == 2 && i>=jsbound)
  1233.            /* use the joint data passed in */
  1234.            d = j_samps[0][j][i] / multiple[j_scale[0][i]];
  1235.          else
  1236.            d = sb_samples[k][0][j][i] / multiple[scalar[k][0][i]];
  1237.          /* scale and quantize floating point sample */
  1238.          n = bit_alloc[k][i];
  1239.          d = d * a[n-1] + b[n-1];
  1240.          /* extract MSB N-1 bits from the floating point sample */
  1241.          if (d >= 0) sig = 1;
  1242.          else { sig = 0; d += 1.0; }
  1243.          sbband[k][0][j][i] = (unsigned int) (d * (double) (1L<<n));
  1244.          /* tag the inverted sign bit to sbband at position N */
  1245.          if (sig) sbband[k][0][j][i] |= 1<<n;
  1246.        }
  1247. }
  1248.  
  1249. /***************************** Layer II  ********************************/
  1250.  
  1251. void II_subband_quantization(scalar, sb_samples, j_scale, j_samps,
  1252.                              bit_alloc, sbband, fr_ps)
  1253. unsigned int scalar[2][3][SBLIMIT];
  1254. double FAR sb_samples[2][3][SCALE_BLOCK][SBLIMIT];
  1255. unsigned int j_scale[3][SBLIMIT];
  1256. double FAR j_samps[3][SCALE_BLOCK][SBLIMIT];
  1257. unsigned int bit_alloc[2][SBLIMIT];
  1258. unsigned int FAR sbband[2][3][SCALE_BLOCK][SBLIMIT];
  1259. frame_params *fr_ps;
  1260. {
  1261.    int i, j, k, s, n, qnt, sig;
  1262.    int stereo  = fr_ps->stereo;
  1263.    int sblimit = fr_ps->sblimit;
  1264.    int jsbound = fr_ps->jsbound;
  1265.    unsigned int stps;
  1266.    double d;
  1267.    al_table *alloc = fr_ps->alloc;
  1268.    for (s=0;s<3;s++)
  1269.      for (j=0;j<SCALE_BLOCK;j++)
  1270.        for (i=0;i<sblimit;i++)
  1271.          for (k=0;k<((i<jsbound)?stereo:1);k++)
  1272.            if (bit_alloc[k][i]) {
  1273.              /* scale and quantize floating point sample */
  1274.              if(stereo == 2 && i>=jsbound)       /* use j-stereo samples */
  1275.                d = j_samps[s][j][i] / multiple[j_scale[s][i]];
  1276.              else
  1277.                d = sb_samples[k][s][j][i] / multiple[scalar[k][s][i]];
  1278.              if (mod(d) > 1.0)
  1279.                printf("Not scaled properly %d %d %d %dn",k,s,j,i);
  1280.              qnt = (*alloc)[i][bit_alloc[k][i]].quant;
  1281.              d = d * a[qnt] + b[qnt];
  1282.              /* extract MSB N-1 bits from the floating point sample */
  1283.              if (d >= 0) sig = 1;
  1284.              else { sig = 0; d += 1.0; }
  1285.              n = 0;
  1286. #ifndef MS_DOS
  1287.              stps = (*alloc)[i][bit_alloc[k][i]].steps;
  1288.              while ((1L<<n) < stps) n++;
  1289. #else
  1290.              while  ( ( (unsigned long)(1L<<(long)n) <
  1291.                        ((unsigned long) ((*alloc)[i][bit_alloc[k][i]].steps)
  1292.                         & 0xffff
  1293.                         )
  1294.                        ) && ( n <16)
  1295.                      ) n++;
  1296. #endif
  1297.              n--;
  1298.              sbband[k][s][j][i] = (unsigned int) (d * (double) (1L<<n));
  1299.              /* tag the inverted sign bit to sbband at position N */
  1300.              /* The bit inversion is a must for grouping with 3,5,9 steps
  1301.                 so it is done for all subbands */
  1302.              if (sig) sbband[k][s][j][i] |= 1<<n;
  1303.            }
  1304.            for (s=0;s<3;s++)
  1305.              for (j=sblimit;j<SBLIMIT;j++)
  1306.                for (i=0;i<SCALE_BLOCK;i++) for (k=0;k<stereo;k++) sbband[k][s][i][j] = 0;
  1307. }
  1308.  
  1309. /*************************************************************************
  1310. * I_encode_bit_alloc  (Layer I)
  1311. * II_encode_bit_alloc (Layer II)
  1312. *
  1313. * PURPOSE:Writes bit allocation information onto bitstream
  1314. *
  1315. * Layer I uses 4 bits/subband for bit allocation information,
  1316. * and Layer II uses 4,3,2, or 0 bits depending on the
  1317. * quantization table used.
  1318. *
  1319. ************************************************************************/
  1320.  
  1321. void I_encode_bit_alloc(bit_alloc, fr_ps, bs)
  1322. unsigned int bit_alloc[2][SBLIMIT];
  1323. frame_params *fr_ps;
  1324. Bit_stream_struc *bs;
  1325. {
  1326.    int i,k;
  1327.    int stereo  = fr_ps->stereo;
  1328.    int sblimit = fr_ps->sblimit;
  1329.    int jsbound = fr_ps->jsbound;
  1330.  
  1331.    for (i=0;i<SBLIMIT;i++)
  1332.      for (k=0;k<((i<jsbound)?stereo:1);k++) putbits(bs,bit_alloc[k][i],4);
  1333. }
  1334.  
  1335. /***************************** Layer II  ********************************/
  1336.  
  1337. void II_encode_bit_alloc(bit_alloc, fr_ps, bs)
  1338. unsigned int bit_alloc[2][SBLIMIT];
  1339. frame_params *fr_ps;
  1340. Bit_stream_struc *bs;
  1341. {
  1342.    int i,k;
  1343.    int stereo  = fr_ps->stereo;
  1344.    int sblimit = fr_ps->sblimit;
  1345.    int jsbound = fr_ps->jsbound;
  1346.    al_table *alloc = fr_ps->alloc;
  1347.  
  1348.    for (i=0;i<sblimit;i++)
  1349.      for (k=0;k<((i<jsbound)?stereo:1);k++)
  1350.        putbits(bs,bit_alloc[k][i],(*alloc)[i][0].bits);
  1351. }
  1352.  
  1353. /************************************************************************
  1354. *
  1355. * I_sample_encoding   (Layer I)
  1356. * II_sample_encoding  (Layer II)
  1357. *
  1358. * PURPOSE:Put one frame of subband samples on to the bitstream
  1359. *
  1360. * SEMANTICS:  The number of bits allocated per sample is read from
  1361. * the bit allocation information #bit_alloc#.  Layer 2
  1362. * supports writing grouped samples for quantization steps
  1363. * that are not a power of 2.
  1364. *
  1365. ************************************************************************/
  1366.  
  1367. void I_sample_encoding(sbband, bit_alloc, fr_ps, bs)
  1368. unsigned int FAR sbband[2][3][SCALE_BLOCK][SBLIMIT];
  1369. unsigned int bit_alloc[2][SBLIMIT];
  1370. frame_params *fr_ps;
  1371. Bit_stream_struc *bs;
  1372. {
  1373.    int i,j,k;
  1374.    int stereo  = fr_ps->stereo;
  1375.    int sblimit = fr_ps->sblimit;
  1376.    int jsbound = fr_ps->jsbound;
  1377.  
  1378.    for(j=0;j<SCALE_BLOCK;j++) {
  1379.      for(i=0;i<SBLIMIT;i++)
  1380.        for(k=0;k<((i<jsbound)?stereo:1);k++)
  1381.          if(bit_alloc[k][i]) putbits(bs,sbband[k][0][j][i],bit_alloc[k][i]+1);
  1382.    }
  1383. }
  1384.  
  1385. /***************************** Layer II  ********************************/
  1386.  
  1387. void II_sample_encoding(sbband, bit_alloc, fr_ps, bs)
  1388. unsigned int FAR sbband[2][3][SCALE_BLOCK][SBLIMIT];
  1389. unsigned int bit_alloc[2][SBLIMIT];
  1390. frame_params *fr_ps;
  1391. Bit_stream_struc *bs;
  1392. {
  1393.    unsigned int temp;
  1394.    unsigned int i,j,k,s,x,y;
  1395.    int stereo  = fr_ps->stereo;
  1396.    int sblimit = fr_ps->sblimit;
  1397.    int jsbound = fr_ps->jsbound;
  1398.    al_table *alloc = fr_ps->alloc;
  1399.  
  1400.    for (s=0;s<3;s++)
  1401.      for (j=0;j<SCALE_BLOCK;j+=3)
  1402.        for (i=0;i<sblimit;i++)
  1403.          for (k=0;k<((i<jsbound)?stereo:1);k++)
  1404.            if (bit_alloc[k][i]) {
  1405.              if ((*alloc)[i][bit_alloc[k][i]].group == 3) {
  1406.                for (x=0;x<3;x++) putbits(bs,sbband[k][s][j+x][i],
  1407.                                          (*alloc)[i][bit_alloc[k][i]].bits);
  1408.              }
  1409.              else {
  1410.                y =(*alloc)[i][bit_alloc[k][i]].steps;
  1411.                temp = sbband[k][s][j][i] +
  1412.                       sbband[k][s][j+1][i] * y +
  1413.                       sbband[k][s][j+2][i] * y * y;
  1414.                putbits(bs,temp,(*alloc)[i][bit_alloc[k][i]].bits);
  1415.              }
  1416.            }
  1417. }
  1418.  
  1419. /************************************************************************
  1420. *
  1421. * encode_CRC
  1422. *
  1423. ************************************************************************/
  1424.  
  1425. void encode_CRC(crc, bs)
  1426. unsigned int crc;
  1427. Bit_stream_struc *bs;
  1428. {
  1429.    putbits(bs, crc, 16);
  1430. }