l3psy.c
上传用户:bjsgzm
上传日期:2007-01-08
资源大小:256k
文件大小:25k
源码类别:

mpeg/mp3

开发平台:

Visual C++

  1. /*
  2. (c) Copyright 1998, 1999 - Tord Jansson
  3. =======================================
  4. This file is part of the BladeEnc MP3 Encoder, based on
  5. ISO's reference code for MPEG Layer 3 compression, and might
  6. contain smaller or larger sections that are directly taken
  7. from ISO's reference code.
  8. All changes to the ISO reference code herein are either
  9. copyrighted by Tord Jansson (tord.jansson@swipnet.se)
  10. or sublicensed to Tord Jansson by a third party.
  11. BladeEnc is free software; you can redistribute this file
  12. and/or modify it under the terms of the GNU Lesser General Public
  13. License as published by the Free Software Foundation; either
  14. version 2.1 of the License, or (at your option) any later version.
  15. */
  16. #include <stdlib.h>
  17. #include "common.h"
  18. #include "encoder.h"
  19. #include "l3psy.h"
  20. #include "l3side.h"
  21. #define maximum(x,y) ( (x>y) ? x : y )
  22. #define minimum(x,y) ( (x<y) ? x : y )
  23. /*____ Global Static Variables ______________________________________________*/
  24. static double ratio[2][21];
  25. static double ratio_s[2][12][3];
  26. /* The static variables "r", "phi_sav", "new_", "old" and "oldest" have    */
  27. /* to be remembered for the unpredictability measure.  For "r" and        */
  28. /* "phi_sav", the first index from the left is the channel select and     */
  29. /* the second index is the "age" of the data.                             */
  30.    static FLOAT window_s[BLKSIZE_s] ;
  31. static int     new_ = 0, old = 1, oldest = 0;
  32. static int     flush, sync_flush, syncsize, sfreq_idx;
  33. static double  cw[HBLKSIZE], eb[CBANDS];
  34. static double  ctb[CBANDS];
  35. static double SNR_l[CBANDS], SNR_s[CBANDS_s];
  36. static double minval[CBANDS],qthr_l[CBANDS],norm_l[CBANDS];
  37. static double qthr_s[CBANDS_s],norm_s[CBANDS_s];
  38. static double nb_1[2][CBANDS], nb_2[2][CBANDS];
  39. static double s3_l[CBANDS][CBANDS]; /* s3_s[CBANDS_s][CBANDS_s]; */
  40. /* Scale Factor Bands */
  41. static int cbw_l[SBMAX_l],bu_l[SBMAX_l],bo_l[SBMAX_l] ;
  42. static int cbw_s[SBMAX_s],bu_s[SBMAX_s],bo_s[SBMAX_s] ;
  43. static double w1_l[SBMAX_l], w2_l[SBMAX_l];
  44. static double w1_s[SBMAX_s], w2_s[SBMAX_s];
  45. static double en[SBMAX_l],   thm[SBMAX_l] ;
  46. static int blocktype_old[2] ;
  47. static int partition_l[HBLKSIZE],partition_s[HBLKSIZE_s];
  48. /* The following static variables are constants.                           */
  49. static float  crit_band[27] = {0,  100,  200, 300, 400, 510, 630,  770,
  50.                                920, 1080, 1270,1480,1720,2000,2320, 2700,
  51.                               3150, 3700, 4400,5300,6400,7700,9500,12000,
  52.                              15500,25000,30000};
  53. /* The following pointer variables point to large areas of memory         */
  54. /* dynamically allocated by the mem_alloc() function.  Dynamic memory     */
  55. /* allocation is used in order to avoid stack frame or data area          */
  56. /* overflow errors that otherwise would have occurred at compile time     */
  57. /* on the Macintosh computer.                                             */
  58. static FLOAT energy_s[3][256];
  59. static FLOAT phi_s[3][256] ; /* 256 samples not 129 */
  60. static int *numlines ;
  61. static int     *partition;
  62. static FLOAT   *cbval, *rnorm;
  63. static FLOAT   *window;
  64. static FLOAT   *absthr;
  65. static double  *tmn;
  66. static FCB     *s;
  67. static FHBLK   *lthr;
  68. static F2HBLK  *r, *phi_sav;
  69. /*____ read_absthr() ________________________________________________________*/
  70. void read_absthr(float * absthr, int table)
  71. {
  72. float * pSource;
  73. int j;
  74. switch(table)
  75. {
  76.      case 0 :
  77. pSource = absthr_0;
  78.        break;
  79.      case 1 :
  80. pSource = absthr_1;
  81.      break;
  82.      case 2 :
  83. pSource = absthr_2;
  84.        break;
  85.      default :
  86. printf("absthr table: Not valid table numbern");
  87. return;
  88. }
  89.   for(j=0; j<HBLKSIZE; j++)
  90.      absthr[j] =  pSource[j];
  91. }
  92. void L3para_read( int sfreq, int numlines[CBANDS], int partition_l[HBLKSIZE],
  93.   double minval[CBANDS], double qthr_l[CBANDS], double norm_l[CBANDS],
  94.   double s3_l[CBANDS][CBANDS], int partition_s[HBLKSIZE_s], double qthr_s[CBANDS_s],
  95.   double norm_s[CBANDS_s], double SNR_s[CBANDS_s],
  96.   int cbw_l[SBMAX_l], int bu_l[SBMAX_l], int bo_l[SBMAX_l],
  97.   double w1_l[SBMAX_l], double w2_l[SBMAX_l],
  98.   int cbw_s[SBMAX_s], int bu_s[SBMAX_s], int bo_s[SBMAX_s],
  99.   double w1_s[SBMAX_s], double w2_s[SBMAX_s] );
  100. /*____ psycho_anal_init() ___________________________________________________*/
  101. void psycho_anal_init( double sfreq )
  102. {
  103. unsigned int   i, j;
  104. float           freq_mult, bval_lo;
  105. double          temp1,temp2,temp3;
  106.    float  fthr[sizeof(FHBLK)/4];
  107. new_ = 0;
  108. old = 1;
  109. oldest = 0;
  110. for( i = 0 ; i < 21 ; i++ )
  111. {
  112. ratio[0][i] = 0.0;
  113. ratio[1][i] = 0.0;
  114. }
  115. for( i = 0 ; i < 2 ; i++ )
  116. for( j = 0 ; j < 12 ; j++ )
  117. {
  118. ratio_s[i][j][0] = 0.0;
  119. ratio_s[i][j][0] = 0.0;
  120. ratio_s[i][j][0] = 0.0;
  121. }
  122. /* These dynamic memory allocations simulate "static" variables placed    */
  123. /* in the data space.  Each mem_alloc() call here occurs only once at     */
  124. /* initialization time.  The mem_free() function must not be called.      */
  125.     numlines = (int *) mem_alloc(sizeof(ICB), "numlines");
  126.     partition  = (int *) mem_alloc(sizeof(IHBLK), "partition");
  127.     cbval  = (FLOAT *) mem_alloc(sizeof(FCB), "cbval");
  128.     rnorm  = (FLOAT *) mem_alloc(sizeof(FCB), "rnorm");
  129.     window  = (FLOAT *) mem_alloc(sizeof(FBLK), "window");
  130.     absthr  = (FLOAT *) mem_alloc(sizeof(FHBLK), "absthr");
  131.     tmn  = (double *) mem_alloc(sizeof(DCB), "tmn");
  132.     s  = (FCB *) mem_alloc(sizeof(FCBCB), "s");
  133.     lthr  = (FHBLK *) mem_alloc(sizeof(F2HBLK), "lthr");
  134.     r  = (F2HBLK *) mem_alloc(sizeof(F22HBLK), "r");
  135.     phi_sav  = (F2HBLK *) mem_alloc(sizeof(F22HBLK), "phi_sav");
  136. /*#if 0 */
  137.     i = sfreq + 0.5;
  138.     switch(i)
  139. {
  140.      case 32000:
  141. sfreq_idx = 0;
  142. break;
  143.         case 44100:
  144. sfreq_idx = 1;
  145. break;
  146.         case 48000:
  147. sfreq_idx = 2;
  148. break;
  149. /*      default:    printf("error, invalid sampling frequency: %d Hzn",i);
  150.          exit(-1); */
  151.     }
  152.     read_absthr(absthr, sfreq_idx);
  153. sync_flush=768;
  154. flush=576;
  155. syncsize=1344;
  156. /* calculate HANN window coefficients */
  157.     for(i=0;i<BLKSIZE;i++)
  158. window[i]  =0.5*(1-cos(2.0*PI*(i-0.5)/BLKSIZE));
  159.     for(i=0;i<BLKSIZE_s;i++)
  160. window_s[i]=0.5*(1-cos(2.0*PI*(i-0.5)/BLKSIZE_s));
  161. /* reset states used in unpredictability measure */
  162. for(i=0;i<HBLKSIZE;i++)
  163. {
  164.      r[0][0][i]=r[1][0][i]=r[0][1][i]=r[1][1][i]=0;
  165.         phi_sav[0][0][i]=phi_sav[1][0][i]=0;
  166.         phi_sav[0][1][i]=phi_sav[1][1][i]=0;
  167.         lthr[0][i] = 60802371420160.0;
  168.         lthr[1][i] = 60802371420160.0;
  169. }
  170. /*****************************************************************************
  171.  * Initialization: Compute the following constants for use later             *
  172.  *    partition[HBLKSIZE] = the partition number associated with each        *
  173.  *                          frequency line                                   *
  174.  *    cbval[CBANDS]       = the center (average) bark value of each          *
  175.  *                          partition                                        *
  176.  *    numlines[CBANDS]    = the number of frequency lines in each partition  *
  177.  *    tmn[CBANDS]         = tone masking noise                               *
  178.  *****************************************************************************/
  179. /* compute fft frequency multiplicand */
  180. freq_mult = sfreq/BLKSIZE;
  181. /* calculate fft frequency, then bval of each line (use fthr[] as tmp storage)*/
  182.     for(i=0;i<HBLKSIZE;i++)
  183. {
  184.         temp1 = i*freq_mult;
  185.         j = 1;
  186.         while(temp1>crit_band[j])
  187. j++;
  188.      fthr[i]=j-1+(temp1-crit_band[j-1])/(crit_band[j]-crit_band[j-1]);
  189.     }
  190.     partition[0] = 0;
  191. /* temp2 is the counter of the number of frequency lines in each partition */
  192. temp2 = 1;
  193.     cbval[0]=fthr[0];
  194.     bval_lo=fthr[0];
  195.     for(i=1;i<HBLKSIZE;i++)
  196. {
  197.         if((fthr[i]-bval_lo)>0.33)
  198. {
  199.          partition[i]=partition[i-1]+1;
  200.          cbval[partition[i-1]] = cbval[partition[i-1]]/temp2;
  201.          cbval[partition[i]] = fthr[i];
  202.          bval_lo = fthr[i];
  203.          numlines[partition[i-1]] = temp2;
  204.          temp2 = 1;
  205.         }
  206.         else
  207. {
  208.          partition[i]=partition[i-1];
  209.          cbval[partition[i]] += fthr[i];
  210.          temp2++;
  211. }
  212.     }
  213.     numlines[partition[i-1]] = temp2;
  214.     cbval[partition[i-1]] = cbval[partition[i-1]]/temp2;
  215. /************************************************************************
  216.  * Now compute the spreading function, s[j][i], the value of the spread-*
  217.  * ing function, centered at band j, for band i, store for later use    *
  218.  ************************************************************************/
  219. for(j=0;j<CBANDS;j++)
  220. {
  221.      for(i=0;i<CBANDS;i++)
  222. {
  223.          temp1 = (cbval[i] - cbval[j])*1.05;
  224. if(temp1>=0.5 && temp1<=2.5)
  225. {
  226.              temp2 = temp1 - 0.5;
  227.              temp2 = 8.0 * (temp2*temp2 - 2.0 * temp2);
  228. }
  229. else
  230. temp2 = 0;
  231. temp1 += 0.474;
  232. temp3 = 15.811389+7.5*temp1-17.5*sqrt((double) (1.0+temp1*temp1));
  233. if(temp3 <= -100)
  234. s[i][j] = 0;
  235. else
  236. {
  237.              temp3 = (temp2 + temp3)*LN_TO_LOG10;
  238.              s[i][j] = exp(temp3);
  239. }
  240. }
  241.     }
  242. /* Calculate Tone Masking Noise values */
  243.     for(j=0;j<CBANDS;j++)
  244. {
  245.         temp1 = 15.5 + cbval[j];
  246.         tmn[j] = (temp1>24.5) ? temp1 : 24.5;
  247.    /* Calculate normalization factors for the net spreading functions */
  248.         rnorm[j] = 0;
  249.         for(i=0;i<CBANDS;i++)
  250. {
  251.            rnorm[j] += s[j][i];
  252.         }
  253. }
  254.     L3para_read( (int)sfreq,numlines,partition_l,minval,qthr_l,norm_l,s3_l,
  255.    partition_s,qthr_s,norm_s,SNR_s,
  256.    cbw_l,bu_l,bo_l,w1_l,w2_l, cbw_s,bu_s,bo_s,w1_s,w2_s );
  257. }
  258. /*____ psycho_anal_exit() ___________________________________________________*/
  259. void psycho_anal_exit( void )
  260. {
  261.     free( numlines );
  262.     free( partition );
  263.     free( cbval );
  264.     free( rnorm );
  265.     free( window );
  266.     free( absthr );
  267.     free( tmn );
  268.     free( s );
  269.     free( lthr );
  270.     free( r );
  271.     free( phi_sav );
  272. }
  273. /*____ psycho_anal() ________________________________________________________*/
  274. void psycho_anal( short int *buffer, short int savebuf[1344], int chn, int lay, float snr32[32],
  275.        double ratio_d[21], double ratio_ds[12][3],
  276.        double *pe, gr_info *cod_info )
  277. {
  278. int  blocktype;
  279. unsigned int   b, i, j, k;
  280. double          r_prime, phi_prime; /* not FLOAT */
  281. double          temp1,temp2,temp3;
  282. double    thr[CBANDS];
  283. int sb,sblock;
  284. float nb[sizeof(FCB)/4];
  285. float  cb[sizeof(FCB)/4];
  286.    float  ecb[sizeof(FCB)/4];
  287.    float  bc[sizeof(FCB)/4];
  288.    float wsamp_r[sizeof(FBLK)/4];
  289.    float  wsamp_i[sizeof(FBLK)/4];
  290.    float  phi[sizeof(FBLK)/4];
  291.    float  energy[sizeof(FBLK)/4];
  292.    float  c[sizeof(FHBLK)/4];
  293.    float  fthr[sizeof(FHBLK)/4];
  294.    F32    snrtmp[2];
  295.    memset(nb, 0, sizeof(FCB));
  296.    memset(cb, 0, sizeof(FCB));
  297.    memset(ecb, 0, sizeof(FCB));
  298.    memset(bc, 0, sizeof(FCB));
  299.    memset(wsamp_r, 0, sizeof(FBLK));
  300.    memset(wsamp_i, 0, sizeof(FBLK));
  301.    memset(phi, 0, sizeof(FBLK));
  302.    memset(energy, 0, sizeof(FBLK));
  303.    memset(c, 0, sizeof(FHBLK));
  304.    memset(fthr, 0, sizeof(FHBLK));
  305.    memset(snrtmp, 0, sizeof(F32)*2);
  306. for ( j = 0; j < 21; j++ )
  307.     ratio_d[j] = ratio[chn][j];
  308. for ( j = 0; j < 12; j++ )
  309. for ( i = 0; i < 3; i++ )
  310. ratio_ds[j][i] = ratio_s[chn][j][i];
  311. if ( chn == 0 )
  312. {
  313. if ( new_ == 0 )
  314.     {
  315. new_ = 1;
  316. old = 0;
  317. oldest = 1;
  318.     }
  319.     else
  320.     {
  321. new_ = 0;
  322. old = 1;
  323. oldest = 0;
  324.     }
  325. }
  326. /**********************************************************************
  327. *  Delay signal by sync_flush=768 samples                             *
  328. **********************************************************************/
  329. for ( j = 0; j < sync_flush; j++ ) /* for long window samples */
  330.     savebuf[j] = savebuf[j+flush];
  331. for ( j = sync_flush; j < syncsize; j++ )
  332.     savebuf[j] = *buffer++;
  333. for ( j = 0; j < BLKSIZE; j++ )
  334. { /**window data with HANN window**/
  335.     wsamp_r[j] = window[j] * savebuf[j];  
  336.     wsamp_i[j] = 0.0;
  337. }
  338. /**********************************************************************
  339. *    compute unpredicatability of first six spectral lines            * 
  340. **********************************************************************/
  341. fft( wsamp_r, wsamp_i, energy, phi, 1024 ); /**long FFT**/
  342. for ( j = 0; j < 6; j++ )
  343. {  /* calculate unpredictability measure cw */
  344.     r_prime = 2.0 * r[chn][old][j] - r[chn][oldest][j];
  345.     phi_prime = 2.0 * phi_sav[chn][old][j]-phi_sav[chn][oldest][j];
  346.     r[chn][new_][j] = sqrt((double) energy[j]);
  347.     phi_sav[chn][new_][j] = phi[j];
  348.     temp1 = r[chn][new_][j] * cos((double) phi[j]) - r_prime * cos(phi_prime);
  349.     temp2 = r[chn][new_][j] * sin((double) phi[j]) - r_prime * sin(phi_prime);
  350.     temp3=r[chn][new_][j] + fabs(r_prime);
  351.     
  352.     if ( temp3 != 0.0 )
  353. cw[j] = sqrt( temp1*temp1+temp2*temp2 ) / temp3;
  354.     else
  355. cw[j] = 0;
  356. }
  357. /**********************************************************************
  358. *     compute unpredicatibility of next 200 spectral lines            *
  359. **********************************************************************/ 
  360. for ( sblock = 0; sblock < 3; sblock++ )
  361. { /**window data with HANN window**/
  362.     for ( j = 0, k = 128 * (2 + sblock); j < 256; j++, k++ )
  363.     {
  364. wsamp_r[j] = window_s[j]* savebuf[k];
  365. wsamp_i[j] = 0.0;
  366.     } /* short FFT*/
  367.     
  368.     fft( wsamp_r, wsamp_i, &energy_s[sblock][0], &phi_s[sblock][0], 256 );
  369. }
  370.  
  371.     sblock = 1;
  372. for ( j = 6; j < 206; j += 4 )
  373. {/* calculate unpredictability measure cw */
  374.     double r2, phi2, temp1, temp2, temp3;
  375.     
  376.     k = (j+2) / 4; 
  377.     r_prime = 2.0 * sqrt((double) energy_s[0][k]) - sqrt((double) energy_s[2][k]);
  378.     phi_prime = 2.0 * phi_s[0][k] - phi_s[2][k];
  379.     r2 = sqrt((double) energy_s[1][k]);
  380.     phi2 = phi_s[1][k];
  381.     temp1 = r2 * cos( phi2 ) - r_prime * cos( phi_prime );
  382.     temp2 = r2 * sin( phi2 ) - r_prime * sin( phi_prime );
  383.     temp3 = r2 + fabs( r_prime );
  384.     if ( temp3 != 0.0 )
  385. cw[j] = sqrt( temp1 * temp1 + temp2 * temp2 ) / temp3;
  386.     else
  387. cw[j] = 0.0;
  388.     cw[j+1] = cw[j+2] = cw[j+3] = cw[j];
  389. }
  390. /**********************************************************************
  391. *    Set unpredicatiblility of remaining spectral lines to 0.4        *
  392. **********************************************************************/
  393. for ( j = 206; j < HBLKSIZE; j++ )
  394.     cw[j] = 0.4;
  395. /**********************************************************************
  396. *    Calculate the energy and the unpredictability in the threshold   *
  397. *    calculation partitions                                           *
  398. **********************************************************************/
  399. for ( b = 0; b < CBANDS; b++ )
  400. {
  401.     eb[b] = 0.0;
  402.     cb[b] = 0.0;
  403. }
  404. for ( j = 0; j < HBLKSIZE; j++ )
  405. {
  406.     int tp = partition_l[j];
  407.     if ( tp >= 0 )
  408.     {
  409. eb[tp] += energy[j];
  410. cb[tp] += cw[j] * energy[j];
  411.     }
  412. }
  413. /**********************************************************************
  414. *      convolve the partitioned energy and unpredictability           *
  415. *      with the spreading function, s3_l[b][k]                        *
  416. ******************************************************************** */
  417. for ( b = 0; b < CBANDS; b++ )
  418. {
  419.     ecb[b] = 0.0;
  420.     ctb[b] = 0.0;
  421. }
  422. for ( b = 0;b < CBANDS; b++ )
  423. {
  424.     for ( k = 0; k < CBANDS; k++ )
  425.     {
  426. ecb[b] += s3_l[b][k] * eb[k]; /* sprdngf for Layer III */
  427. ctb[b] += s3_l[b][k] * cb[k];
  428.     }
  429. }
  430. /* calculate the tonality of each threshold calculation partition */
  431. /* calculate the SNR in each threshhold calculation partition */
  432. for ( b = 0; b < CBANDS; b++ )
  433. {
  434.     double cbb,tbb;
  435.     if (ecb[b] != 0.0 )
  436.         {
  437. cbb = ctb[b]/ecb[b];
  438.          if (cbb <0.01) cbb = 0.01;
  439. cbb = log( cbb);
  440.         }
  441.     else
  442. cbb = 0.0 ;
  443.     tbb = -0.299 - 0.43*cbb;  /* conv1=-0.299, conv2=-0.43 */
  444.     tbb = minimum( 1.0, maximum( 0.0, tbb) ) ;  /* 0<tbb<1 */
  445.     SNR_l[b] = maximum( minval[b], 29.0*tbb+6.0*(1.0-tbb) );
  446. } /* TMN=29.0,NMT=6.0 for all calculation partitions */
  447. for ( b = 0; b < CBANDS; b++ ) /* calculate the threshold for each partition */
  448.     nb[b] = ecb[b] * norm_l[b] * exp( -SNR_l[b] * LN_TO_LOG10 );
  449. for ( b = 0; b < CBANDS; b++ )
  450. { /* pre-echo control */
  451.     double temp_1; /* BUG of IS */
  452.     temp_1 = minimum( nb[b], minimum(2.0*nb_1[chn][b],16.0*nb_2[chn][b]) );
  453.     thr[b] = maximum( qthr_l[b], temp_1 );/* rpelev=2.0, rpelev2=16.0 */
  454.     nb_2[chn][b] = nb_1[chn][b];
  455.     nb_1[chn][b] = nb[b];
  456. }
  457. *pe = 0.0; /*  calculate percetual entropy */
  458. for ( b = 0; b < CBANDS; b++ )
  459. {
  460.     double tp ;
  461.     tp = minimum( 0.0, log((thr[b]+1.0) / (eb[b]+1.0) ) ) ; /*not log*/
  462.     *pe -= numlines[b] * tp ;
  463. } /* thr[b] -> thr[b]+1.0 : for non sound portition */
  464. #define switch_pe  1800
  465.         blocktype = NORM_TYPE;
  466. if ( *pe < switch_pe )
  467. { /* no attack : use long blocks */
  468. switch( blocktype_old[chn] )
  469.     {
  470.      case NORM_TYPE:
  471.      case STOP_TYPE:
  472. blocktype = NORM_TYPE;
  473. break;
  474.     
  475.        case SHORT_TYPE:
  476. blocktype = STOP_TYPE;
  477. break;
  478.     
  479. /*        case START_TYPE:
  480. fprintf( stderr, "Error in block selectingn" );
  481. abort();
  482. break; problem */
  483.     }
  484.     /* threshold calculation (part 2) */
  485.     for ( sb = 0; sb < SBMAX_l; sb++ )
  486.     {
  487. en[sb] = w1_l[sb] * eb[bu_l[sb]] + w2_l[sb] * eb[bo_l[sb]];
  488. thm[sb] = w1_l[sb] *thr[bu_l[sb]] + w2_l[sb] * thr[bo_l[sb]];
  489. for ( b = bu_l[sb]+1; b < bo_l[sb]; b++ )
  490. {
  491.     en[sb]  += eb[b];
  492.     thm[sb] += thr[b];
  493. }
  494. if ( en[sb] != 0.0 )
  495. ratio[chn][sb] = thm[sb]/en[sb];
  496. else
  497.     ratio[chn][sb] = 0.0;
  498.     }
  499. }
  500. else 
  501. {
  502.     /* attack : use short blocks */
  503.     blocktype = SHORT_TYPE;
  504.     
  505.     if ( blocktype_old[chn] == NORM_TYPE ) 
  506. blocktype_old[chn] = START_TYPE;
  507.     if ( blocktype_old[chn] == STOP_TYPE )
  508. blocktype_old[chn] = SHORT_TYPE ;
  509.     
  510.     /* threshold calculation for short blocks */
  511.     
  512.     for ( sblock = 0; sblock < 3; sblock++ )
  513.     {
  514. for ( b = 0; b < CBANDS_s; b++ )
  515. {
  516.     eb[b] = 0.0;
  517.     ecb[b] = 0.0;
  518. }
  519. for ( j = 0; j < HBLKSIZE_s; j++ )
  520.     eb[partition_s[j]] += energy_s[sblock][j];
  521. for ( b = 0; b < CBANDS_s; b++ )
  522.     for ( k = 0; k < CBANDS_s; k++ )
  523. ecb[b] += s3_l[b][k] * eb[k];
  524. for ( b = 0; b < CBANDS_s; b++ )
  525. {
  526.     nb[b] = ecb[b] * norm_l[b] * exp( (double) SNR_s[b] * LN_TO_LOG10 );
  527.     thr[b] = maximum (qthr_s[b],nb[b]);
  528. }
  529. for ( sb = 0; sb < SBMAX_s; sb++ )
  530. {
  531.     en[sb] = w1_s[sb] * eb[bu_s[sb]] + w2_s[sb] * eb[bo_s[sb]];
  532.     thm[sb] = w1_s[sb] *thr[bu_s[sb]] + w2_s[sb] * thr[bo_s[sb]];
  533.     for ( b = bu_s[sb]+1; b < bo_s[sb]; b++ )
  534.     {
  535. en[sb] += eb[b];
  536. thm[sb] += thr[b];
  537.     }
  538.     if ( en[sb] != 0.0 )
  539. ratio_s[chn][sb][sblock] = thm[sb]/en[sb];
  540.     else
  541. ratio_s[chn][sb][sblock] = 0.0;
  542. }
  543.     }
  544. cod_info->block_type = blocktype_old[chn];
  545. blocktype_old[chn] = blocktype;
  546. if ( cod_info->block_type == NORM_TYPE )
  547.     cod_info->window_switching_flag = 0;
  548. else
  549.     cod_info->window_switching_flag = 1;
  550. cod_info->mixed_block_flag = 0;
  551. }
  552. /*____ L3para_read() __________________________________________________________*/
  553. void L3para_read( int sfreq, int *numlines, int *partition_l, double *minval,
  554.   double *qthr_l, double *norm_l, double (*s3_l)[63], int *partition_s,
  555.   double *qthr_s, double *norm_s, double *SNR, int *cbw_l, int *bu_l,
  556.   int *bo_l, double *w1_l, double *w2_l, int *cbw_s, int *bu_s,
  557.   int *bo_s, double *w1_s, double *w2_s )
  558. {
  559.    static double bval_l[CBANDS], bval_s[CBANDS];
  560.    int   cbmax, cbmax_tp;
  561.    static double s3_s[CBANDS][CBANDS];
  562.    int  sbmax ;
  563.    int  i,j,k,k2, part_max ;
  564.  psyDataElem * rpa1;
  565.  psyDataElem2 * rpa2;
  566.  psyDataElem3 * rpa3;
  567.  /* Read long block data */
  568. switch( sfreq )
  569. {
  570. case 32000:
  571. rpa1 = psy_longBlock__32000_58;
  572. cbmax_tp = 59;
  573. break;
  574. case 44100:
  575. rpa1 = psy_longBlock_44100_62;
  576. cbmax_tp = 63;
  577. break;
  578. case 48000:
  579. rpa1 = psy_longBlock_48000_61;
  580. cbmax_tp = 62;
  581. break;
  582. default:
  583. return; /* Just to avoid compiler warnings */
  584. }
  585. cbmax = cbmax_tp;
  586. for(i=0,k2=0;i<cbmax_tp;i++)
  587. {
  588. numlines[i] = rpa1->lines;
  589. minval[i] = rpa1->minVal;
  590. qthr_l[i] = rpa1->qthr;
  591. norm_l[i] = rpa1->norm;
  592. bval_l[i] = rpa1->bVal;
  593. rpa1++;
  594. for(k=0;k<numlines[i];k++)
  595. partition_l[k2++] = i ;
  596. }
  597. /************************************************************************
  598.  * Now compute the spreading function, s[j][i], the value of the spread-*
  599.  * ing function, centered at band j, for band i, store for later use    *
  600.  ************************************************************************/
  601.   part_max = cbmax ;
  602.       for(i=0;i<part_max;i++)
  603.   {
  604. double tempx,x,tempy,temp;
  605.         for(j=0;j<part_max;j++)
  606.     {
  607. /*            tempx = (bval_l[i] - bval_l[j])*1.05; */
  608. if (j>=i)
  609. tempx = (bval_l[i] - bval_l[j])*3.0;
  610.             else
  611. tempx = (bval_l[i] - bval_l[j])*1.5;
  612. /*             if (j>=i) tempx = (bval_l[j] - bval_l[i])*3.0;
  613.                else    tempx = (bval_l[j] - bval_l[i])*1.5; */
  614.              if(tempx>=0.5 && tempx<=2.5)
  615.      {
  616.                temp = tempx - 0.5;
  617.                x = 8.0 * (temp*temp - 2.0 * temp);
  618.              }
  619.              else x = 0.0;
  620.              tempx += 0.474;
  621.              tempy = 15.811389 + 7.5*tempx - 17.5*sqrt(1.0+tempx*tempx);
  622.              if (tempy <= -60.0) s3_l[i][j] = 0.0;
  623.              else                s3_l[i][j] = exp( (x + tempy)*LN_TO_LOG10 );
  624.             }
  625.           }
  626. /* Read short block data */
  627. switch( sfreq )
  628. {
  629. case 32000:
  630. rpa2 = psy_shortBlock_32000_41;
  631. cbmax_tp = 42;
  632. break;
  633. case 44100:
  634. rpa2 = psy_shortBlock_44100_38;
  635. cbmax_tp = 39;
  636. break;
  637. case 48000:
  638. rpa2 = psy_shortBlock_48000_37;
  639. cbmax_tp = 38;
  640. break;
  641. default:
  642. return; /* Just to avoid compiler warnings */
  643. }
  644. cbmax = cbmax_tp;
  645. for(i=0,k2=0;i<cbmax_tp;i++)
  646. {
  647. numlines[i] = rpa2->lines;
  648. qthr_s[i] = rpa2->qthr;
  649. norm_s[i] = rpa2->norm;
  650. SNR[i] = rpa2->snr;
  651. bval_s[i] = rpa2->bVal;
  652. rpa2++;
  653. for(k=0;k<numlines[i];k++)
  654. partition_s[k2++] = i ;
  655. }
  656. /************************************************************************
  657.  * Now compute the spreading function, s[j][i], the value of the spread-*
  658.  * ing function, centered at band j, for band i, store for later use    *
  659.  ************************************************************************/
  660.   part_max = cbmax ;
  661.           for(i=0;i<part_max;i++)
  662.   {
  663.   double tempx,x,tempy,temp;
  664.             for(j=0;j<part_max;j++)
  665.     {
  666.              tempx = (bval_s[i] - bval_s[j])*1.05;
  667.              if (j>=i) tempx = (bval_s[i] - bval_s[j])*3.0;
  668.                else    tempx = (bval_s[i] - bval_s[j])*1.5;
  669.              if(tempx>=0.5 && tempx<=2.5)
  670.      {
  671.                temp = tempx - 0.5;
  672.                x = 8.0 * (temp*temp - 2.0 * temp);
  673.              }
  674.              else x = 0.0;
  675.              tempx += 0.474;
  676.              tempy = 15.811389 + 7.5*tempx - 17.5*sqrt(1.0+tempx*tempx);
  677.              if (tempy <= -60.0) s3_s[i][j] = 0.0;
  678.              else                s3_s[i][j] = exp( (x + tempy)*LN_TO_LOG10 );
  679.             }
  680.           }
  681. /* Read long block data for converting threshold calculation
  682.    partitions to scale factor bands */
  683. switch( sfreq )
  684. {
  685. case 32000:
  686. rpa3 = psy_data3_32000_20;
  687. sbmax = 21;
  688. break;
  689. case 44100:
  690. rpa3 = psy_data3_44100_20;
  691. sbmax = 21;
  692. break;
  693. case 48000:
  694. rpa3 = psy_data3_48000_20;
  695. sbmax = 21;
  696. break;
  697. default:
  698. return; /* Just to avoid compiler warnings */
  699. }
  700.   for(i=0;i<sbmax;i++)
  701.   {
  702. cbw_l[i] = rpa3->cbw;
  703. bu_l[i] = rpa3->bu;
  704. bo_l[i] = rpa3->bo;
  705. w1_l[i] = rpa3->w1;
  706. w2_l[i] = rpa3->w2;
  707. rpa3++;
  708. }
  709. /* Read short block data for converting threshold calculation
  710.    partitions to scale factor bands */
  711. switch( sfreq )
  712. {
  713. case 32000:
  714. rpa3 = psy_data4_32000_11;
  715. sbmax = 12;
  716. break;
  717. case 44100:
  718. rpa3 = psy_data4_44100_11;
  719. sbmax = 12;
  720. break;
  721. case 48000:
  722. rpa3 = psy_data4_48000_11;
  723. sbmax = 12;
  724. break;
  725. }
  726.   for(i=0;i<sbmax;i++)
  727.   {
  728. cbw_s[i] = rpa3->cbw;
  729. bu_s[i] = rpa3->bu;
  730. bo_s[i] = rpa3->bo;
  731. w1_s[i] = rpa3->w1;
  732. w2_s[i] = rpa3->w2;
  733. rpa3++;
  734. }
  735. }