home *** CD-ROM | disk | FTP | other *** search
/ modiromppu / modiromppu.iso / PROGRAMS / ORGPACKS / MPG12304.ZIP / LAYER3.C < prev    next >
C/C++ Source or Header  |  1997-04-08  |  43KB  |  1,619 lines
  1. /* 
  2.  * Mpeg Layer-3 audio decoder 
  3.  * --------------------------
  4.  * copyright (c) 1995,1996,1997 by Michael Hipp.
  5.  * All rights reserved. See also 'README'
  6.  *
  7.  * - I'm currently working on that .. needs a few more optimizations,
  8.  *   though the code is now fast enough to run in realtime on a 133Mhz 486
  9.  * - a few personal notes are in german .. 
  10.  *
  11.  * used source: 
  12.  *   mpeg1_iis package
  13.  */ 
  14.  
  15. #include "mpg123.h"
  16. #include "huffman.h"
  17.  
  18. static real ispow[8207];
  19. static real aa_ca[8],aa_cs[8];
  20. static real COS1[12][6];
  21. static real win[4][36];
  22. static real win1[4][36];
  23. static real gainpow2[256+118];
  24. static real COS9[9];
  25. static real COS6_1,COS6_2;
  26. static real tfcos36[9];
  27. static real tfcos12[3];
  28.  
  29. static int slen[2][16] = {
  30.   {0, 0, 0, 0, 3, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4},
  31.   {0, 1, 2, 3, 0, 1, 2, 3, 1, 2, 3, 1, 2, 3, 2, 3}
  32. };
  33.  
  34. struct bandInfoStruct {
  35.   int longIdx[23];
  36.   int longDiff[22];
  37.   int shortIdx[14];
  38.   int shortDiff[13];
  39. };
  40.  
  41. struct bandInfoStruct bandInfo[3] = { 
  42.  { {0,4,8,12,16,20,24,30,36,44,52,62,74, 90,110,134,162,196,238,288,342,418,576},
  43.    {4,4,4,4,4,4,6,6,8, 8,10,12,16,20,24,28,34,42,50,54, 76,158},
  44.    {0,4*3,8*3,12*3,16*3,22*3,30*3,40*3,52*3,66*3, 84*3,106*3,136*3,192*3},
  45.    {4,4,4,4,6,8,10,12,14,18,22,30,56} } ,
  46.  { {0,4,8,12,16,20,24,30,36,42,50,60,72, 88,106,128,156,190,230,276,330,384,576},
  47.    {4,4,4,4,4,4,6,6,6, 8,10,12,16,18,22,28,34,40,46,54, 54,192},
  48.    {0,4*3,8*3,12*3,16*3,22*3,28*3,38*3,50*3,64*3, 80*3,100*3,126*3,192*3},
  49.    {4,4,4,4,6,6,10,12,14,16,20,26,66} } ,
  50.  { {0,4,8,12,16,20,24,30,36,44,54,66,82,102,126,156,194,240,296,364,448,550,576} ,
  51.    {4,4,4,4,4,4,6,6,8,10,12,16,20,24,30,38,46,56,68,84,102, 26} ,
  52.    {0,4*3,8*3,12*3,16*3,22*3,30*3,42*3,58*3,78*3,104*3,138*3,180*3,192*3} ,
  53.    {4,4,4,4,6,8,12,16,20,26,34,42,12} } 
  54. };
  55.  
  56. static int mapbuf0[3][152];
  57. static int mapbuf1[3][156];
  58. static int mapbuf2[3][44];
  59. static int *map[3][3];
  60. static int *mapend[3][3];
  61.  
  62. static real tan1_1[16],tan2_1[16],tan1_2[16],tan2_2[16];
  63.  
  64. /* 
  65.  * init tables for layer-3 
  66.  */
  67. void init_layer3(void)
  68. {
  69.   int i,j,m;
  70.   static int init=0;
  71.   double sq;
  72.  
  73.   if(init)
  74.     return;
  75.   init = 1;
  76.  
  77.   for(i=-256;i<118;i++)
  78.     gainpow2[i+256] = pow((double)2.0,(double) 0.25 * ( (double) - i - 210.0 ));
  79.  
  80.   for(i=0;i<8207;i++)
  81.     ispow[i] = pow((double)i,(double)4.0/3.0);
  82.  
  83.   for (i=0;i<8;i++)
  84.   {
  85.     static double Ci[8]={-0.6,-0.535,-0.33,-0.185,-0.095,-0.041,-0.0142,-0.0037};
  86.     sq=sqrt(1.0+Ci[i]*Ci[i]);
  87.     aa_cs[i] = 1.0/sq;
  88.     aa_ca[i] = Ci[i]/sq;
  89.   }
  90.  
  91.    for(i=0;i<18;i++)
  92.    {
  93.      win[0][i]    = win[1][i]    = sin( M_PI/36.0 * ((double)i+0.5) );
  94.      win[0][i+18] = win[3][i+18] = sin( M_PI/36.0 * ((double)i+18.5) );
  95.    }
  96.    for(i=0;i<6;i++)
  97.    {
  98.      win[1][i+18] = win[3][i+12] = 1.0;
  99.      win[1][i+24] =                sin( M_PI/12.0 *((double)i+6.5) );
  100.      win[1][i+30] = win[3][i]    = 0.0;
  101.                     win[3][i+6 ] = sin( M_PI/12.0 *((double)i+0.5) );
  102.    }
  103.  
  104.    for(j=0;j<4;j++) {
  105.      if(j == 2)
  106.        continue;
  107.      for(i=0;i<9;i++)
  108.        win[j][i] *= 0.5 / ( cos ( ( M_PI * ( (i+9) * 2.0 + 1.0)) / 72.0 ) );
  109.      for(i=9;i<27;i++)
  110.        win[j][i] *= -0.5 / ( cos ( ( M_PI * ( (26-i) * 2.0 + 1.0)) / 72.0 ) );
  111.      for(i=27;i<36;i++)
  112.        win[j][i] *= -0.5 / ( cos ( ( M_PI * ( (i-27) * 2.0 + 1.0)) / 72.0 ) );
  113.    }
  114.  
  115.   for(i=0;i<9;i++)
  116.     COS9[i] = cos( M_PI / 18.0 * i);
  117.  
  118.   for(i=0;i<9;i++)
  119.     tfcos36[i] = 1.0 / ( 2.0 * cos ( ( M_PI * (i * 2.0 + 1.0)) / 36.0 ) );
  120.   for(i=0;i<3;i++)
  121.     tfcos12[i] = 0.5 / ( cos ( ( M_PI * (i * 2.0 + 1.0)) / 12.0 ) );
  122.  
  123.   COS6_1 = cos( M_PI / 6.0 * 1.0);
  124.   COS6_2 = cos( M_PI / 6.0 * 2.0);
  125.  
  126.   for(i=0;i<12;i++)
  127.   {
  128.     win[2][i] = sin( M_PI/12.0*((double)i+0.5) ) ;
  129.     for(m=0;m<6;m++)
  130.       COS1[i][m] = cos( M_PI/24.0*(2.0*(double)i+7.0)*(2.0*(double)m+1.0) );
  131.   }
  132.  
  133.   for(i=0;i<3;i++) {
  134.     win[2][i+0] *=  0.5 / ( cos ( ( M_PI * ((i+3) * 2.0 + 1.0)) / 24.0 ) );
  135.     win[2][i+3] *= -0.5 / ( cos ( ( M_PI * ((5-i) * 2.0 + 1.0)) / 24.0 ) );
  136.     win[2][i+6] *= -0.5 / ( cos ( ( M_PI * ((2-i) * 2.0 + 1.0)) / 24.0 ) );
  137.     win[2][i+9] *= -0.5 / ( cos ( ( M_PI * ((i+0) * 2.0 + 1.0)) / 24.0 ) );
  138.   }
  139.  
  140.   for(j=0;j<4;j++) {
  141.     if(j == 2) {
  142.       for(i=0;i<12;i+=2)
  143.         win1[2][i] = + win[2][i];
  144.       for(i=1;i<12;i+=2)
  145.         win1[2][i] = - win[2][i];
  146.     }
  147.     else {
  148.       for(i=0;i<36;i+=2)
  149.         win1[j][i] = + win[j][i];
  150.       for(i=1;i<36;i+=2)
  151.         win1[j][i] = - win[j][i];
  152.     }
  153.   }
  154.  
  155.   for(i=0;i<15;i++)
  156.   {
  157.     double t = tan((double) i * (M_PI / 12.0));
  158.     tan1_1[i] = t / (1.0+t);
  159.     tan2_1[i] = 1.0 / (1.0 + t);
  160.     tan1_2[i] = M_SQRT2 * t / (1.0+t);
  161.     tan2_2[i] = M_SQRT2 / (1.0 + t);
  162.   }
  163.  
  164.   for(m=0;m<3;m++)
  165.   {
  166.    struct bandInfoStruct *bi = &bandInfo[m];
  167.    int *mp;
  168.    int cb,lwin;
  169.    int *bdf;
  170.  
  171.    mp = map[m][0] = mapbuf0[m];
  172.    bdf = bi->longDiff;
  173.    for(i=0,cb = 0; cb < 8 ; cb++,i+=*bdf++) {
  174.      *mp++ = (*bdf) >> 1;
  175.      *mp++ = i;
  176.      *mp++ = 3;
  177.      *mp++ = cb;
  178.    }
  179.    bdf = bi->shortDiff+3;
  180.    for(cb=3;cb<13;cb++) {
  181.      int l = (*bdf++) >> 1;
  182.      for(lwin=0;lwin<3;lwin++) {
  183.        *mp++ = l;
  184.        *mp++ = i + lwin;
  185.        *mp++ = lwin;
  186.        *mp++ = cb;
  187.      }
  188.      i += 6*l;
  189.    }
  190.    mapend[m][0] = mp;
  191.  
  192.    mp = map[m][1] = mapbuf1[m];
  193.    bdf = bi->shortDiff+0;
  194.    for(i=0,cb=0;cb<13;cb++) {
  195.      int l = (*bdf++) >> 1;
  196.      for(lwin=0;lwin<3;lwin++) {
  197.        *mp++ = l;
  198.        *mp++ = i + lwin;
  199.        *mp++ = lwin;
  200.        *mp++ = cb;
  201.      }
  202.      i += 6*l;
  203.    }
  204.    mapend[m][1] = mp;
  205.  
  206.    mp = map[m][2] = mapbuf2[m];
  207.    bdf = bi->longDiff;
  208.    for(cb = 0; cb < 22 ; cb++) {
  209.      *mp++ = (*bdf++) >> 1;
  210.      *mp++ = cb;
  211.    }
  212.    mapend[m][2] = mp;
  213.  
  214.   }  
  215. }
  216.  
  217. /*
  218.  * read additional side information
  219.  */
  220. static void III_get_side_info(struct III_sideinfo *si,int stereo,int ms_stereo,long sfreq)
  221. {
  222.    int ch, gr;
  223.  
  224.    si->main_data_begin = getbits(9);
  225.    if (stereo == 1)
  226.      si->private_bits = getbits_fast(5);
  227.    else 
  228.      si->private_bits = getbits_fast(3);
  229.  
  230.    for (ch=0; ch<stereo; ch++) {
  231.        si->ch[ch].gr[0].scfsi = -1;
  232.        si->ch[ch].gr[1].scfsi = getbits_fast(4);
  233.    }
  234.  
  235.    for (gr=0; gr<2; gr++) 
  236.    {
  237.      for (ch=0; ch<stereo; ch++) 
  238.      {
  239.        register struct gr_info_s *gr_info = &(si->ch[ch].gr[gr]);
  240.  
  241.        gr_info->part2_3_length = getbits(12);
  242.        gr_info->big_values = getbits_fast(9);
  243.        gr_info->pow2gain = gainpow2+256 - getbits_fast(8);
  244.        if(ms_stereo)
  245.          gr_info->pow2gain += 2;
  246.        gr_info->scalefac_compress = getbits_fast(4);
  247. /* window-switching flag == 1 for block_Type != 0 .. and block-type == 0 -> win-sw-flag = 0 */
  248.        if(get1bit()) 
  249.        {
  250.          int i;
  251.          gr_info->block_type = getbits_fast(2);
  252.          gr_info->mixed_block_flag = get1bit();
  253.          gr_info->table_select[0] = getbits_fast(5);
  254.          gr_info->table_select[1] = getbits_fast(5);
  255.          /*
  256.           * table_select[2] not needed, because there is no region2,
  257.           * but to satisfy some verifications tools we set it either.
  258.           */
  259.          gr_info->table_select[2] = 0;
  260.          for(i=0;i<3;i++)
  261.            gr_info->full_gain[i] = gr_info->pow2gain + (getbits_fast(3)<<3);
  262.  
  263.          if(gr_info->block_type == 0) {
  264.            fprintf(stderr,"Blocktype == 0 and window-switching == 1 not allowed.\n");
  265.            exit(1);
  266.          }
  267.          /* region_count/start parameters are implicit in this case. */       
  268.          gr_info->region1start = 36>>1;
  269.          gr_info->region2start = 576>>1;
  270.        }
  271.        else 
  272.        {
  273.          int i,r0c,r1c;
  274.          for (i=0; i<3; i++)
  275.            gr_info->table_select[i] = getbits_fast(5);
  276.          r0c = getbits_fast(4);
  277.          r1c = getbits_fast(3);
  278.          gr_info->region1start = bandInfo[sfreq].longIdx[r0c+1] >> 1 ;
  279.          gr_info->region2start = bandInfo[sfreq].longIdx[r0c+1+r1c+1] >> 1;
  280.          gr_info->block_type = 0;
  281.          gr_info->mixed_block_flag = 0;
  282.        }
  283.        gr_info->preflag = get1bit();
  284.        gr_info->scalefac_scale = get1bit();
  285.        gr_info->count1table_select = get1bit();
  286.      }
  287.    }
  288. }
  289.  
  290. /*
  291.  * read scalefactors
  292.  */
  293. static void III_get_scale_factors(int *scf,struct gr_info_s *gr_info)
  294. {
  295.     if (gr_info->block_type == 2) 
  296.     {
  297.       int i=18;
  298.       int num = slen[0][gr_info->scalefac_compress];
  299.  
  300.       if (gr_info->mixed_block_flag) {
  301.          for (i=8;i;i--)
  302.            *scf++ = getbits_fast(num);
  303.          i = 9;
  304.       }
  305.  
  306.       for (;i;i--)
  307.         *scf++ = getbits_fast(num);
  308.       num = slen[1][gr_info->scalefac_compress];
  309.       for (i = 18; i; i--)
  310.         *scf++ = getbits_fast(num);
  311.       *scf++ = 0; *scf++ = 0; *scf++ = 0; /* short[13][0..2] = 0 */
  312.     }
  313.     else 
  314.     {
  315.       int i,num=slen[0][gr_info->scalefac_compress];
  316.       int scfsi = gr_info->scfsi;
  317.  
  318.       if(scfsi < 0) { /* scfsi < 0 => granule == 0 */
  319.          for(i=11;i;i--)
  320.            *scf++ = getbits_fast(num);
  321.          num=slen[1][gr_info->scalefac_compress];
  322.          for(i=10;i;i--)
  323.            *scf++ = getbits_fast(num);
  324.       }
  325.       else {
  326.         if(!(scfsi & 0x8))
  327.           for (i=6;i;i--)
  328.             *scf++ = getbits_fast(num);
  329.         else {
  330.           *scf++ = 0; *scf++ = 0; *scf++ = 0;  /* set to ZERO necessary? */
  331.           *scf++ = 0; *scf++ = 0; *scf++ = 0;
  332.         }
  333.  
  334.         if(!(scfsi & 0x4))
  335.           for (i=5;i;i--)
  336.             *scf++ = getbits_fast(num);
  337.         else {
  338.           *scf++ = 0; *scf++ = 0; *scf++ = 0;  /* set to ZERO necessary? */
  339.           *scf++ = 0; *scf++ = 0;
  340.         }
  341.  
  342.         num=slen[1][gr_info->scalefac_compress];
  343.         if(!(scfsi & 0x2))
  344.           for(i=5;i;i--)
  345.             *scf++ = getbits_fast(num);
  346.         else {
  347.           *scf++ = 0; *scf++ = 0; *scf++ = 0;  /* set to ZERO necessary? */
  348.           *scf++ = 0; *scf++ = 0;
  349.         }
  350.  
  351.         if(!(scfsi & 0x1))
  352.           for (i=5;i;i--)
  353.             *scf++ = getbits_fast(num);
  354.         else {
  355.           *scf++ = 0; *scf++ = 0; *scf++ = 0;  /* set to ZERO necessary? */
  356.           *scf++ = 0; *scf++ = 0;
  357.         }
  358.       }
  359.  
  360.       *scf++ = 0;  /* no l[21] in original sources */
  361.     }
  362. }
  363.  
  364. static int pretab1[22] = {0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,2,2,3,3,3,2,0};
  365. static int pretab2[22] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
  366.  
  367. /*
  368.  * don't forget to apply the same changes to III_dequantize_sample_ms() !!! 
  369.  * (note: maxband stuff would only be necessary for second channel and I-stereo)
  370.  */
  371. static void III_dequantize_sample(real xr[SBLIMIT][SSLIMIT],int *scf,
  372.    struct gr_info_s *gr_info,int sfreq,int part2_start)
  373. {
  374.   int shift = 1 + gr_info->scalefac_scale;
  375.   real *xrpnt = (real *) xr;
  376.   int l[3],l3;
  377.   int part2end = gr_info->part2_3_length + part2_start;
  378.   int *me;
  379.  
  380.   {
  381.     int bv       = gr_info->big_values;
  382.     int region1  = gr_info->region1start;
  383.     int region2  = gr_info->region2start;
  384.  
  385.     l3 = ((576>>1)-bv)>>1;   
  386. /*
  387.  * we may lose the 'odd' bit here !! 
  388.  * check this later gain 
  389.  */
  390.     if(bv <= region1) {
  391.       l[0] = bv; l[1] = 0; l[2] = 0;
  392.     }
  393.     else {
  394.       l[0] = region1;
  395.       if(bv <= region2) {
  396.         l[1] = bv - l[0];  l[2] = 0;
  397.       }
  398.       else {
  399.         l[1] = region2 - l[0]; l[2] = bv - region2;
  400.       }
  401.     }
  402.   }
  403.  
  404.   if(gr_info->block_type == 2) {
  405.     int i,max[4];
  406.     int step=0,lwin=0,cb=0;
  407.     register real v = 0.0;
  408.     register int *m,mc;
  409.  
  410.     if(gr_info->mixed_block_flag) {
  411.       max[3] = -1;
  412.       max[0] = max[1] = max[2] = 2;
  413.       m = map[sfreq][0];
  414.       me = mapend[sfreq][0];
  415.     }
  416.     else {
  417.       max[0] = max[1] = max[2] = max[3] = -1;
  418.       /* max[3] not really needed in this case */
  419.       m = map[sfreq][1];
  420.       me = mapend[sfreq][1];
  421.     }
  422.  
  423.     mc = 0;
  424.     for(i=0;i<2;i++) {
  425.       int lp = l[i];
  426.       struct newhuff *h = ht+gr_info->table_select[i];
  427.       for(;lp;lp--,mc--) {
  428.         register int x,y;
  429.         if( (!mc) ) {
  430.           mc = *m++;
  431.           xrpnt = ((real *) xr) + (*m++);
  432.           lwin = *m++;
  433.           cb = *m++;
  434.           if(lwin == 3) {
  435.             v = gr_info->pow2gain[(*scf++) << shift];
  436.             step = 1;
  437.           }
  438.           else {
  439.             v = gr_info->full_gain[lwin][(*scf++) << shift];
  440.             step = 3;
  441.           }
  442.         }
  443.         {
  444.           register short *val = h->table;
  445.           while((y=*val++)<0)
  446.             if (get1bit())
  447.               val -= y;
  448.           x = y >> 4;
  449.           y &= 0xf;
  450.         }
  451.         if(x == 15) {
  452.           max[lwin] = cb;
  453.           x += getbits(h->linbits);
  454.           if(get1bit())
  455.             *xrpnt = -ispow[x] * v;
  456.           else
  457.             *xrpnt =  ispow[x] * v;
  458.         }
  459.         else if(x) {
  460.           max[lwin] = cb;
  461.           if(get1bit())
  462.             *xrpnt = -ispow[x] * v;
  463.           else
  464.             *xrpnt =  ispow[x] * v;
  465.         }
  466.         else
  467.           *xrpnt = 0.0;
  468.         xrpnt += step;
  469.         if(y == 15) {
  470.           max[lwin] = cb;
  471.           y += getbits(h->linbits);
  472.           if(get1bit())
  473.             *xrpnt = -ispow[y] * v;
  474.           else
  475.             *xrpnt =  ispow[y] * v;
  476.         }
  477.         else if(y) {
  478.           max[lwin] = cb;
  479.           if(get1bit())
  480.             *xrpnt = -ispow[y] * v;
  481.           else
  482.             *xrpnt =  ispow[y] * v;
  483.         }
  484.         else
  485.           *xrpnt = 0.0;
  486.         xrpnt += step;
  487.       }
  488.     }
  489.  
  490.     for(;l3 && hsstell() < part2end;l3--) {
  491.       struct newhuff *h = htc+gr_info->count1table_select;
  492.       register short *val = h->table,a;
  493.  
  494.       while((a=*val++)<0)
  495.         if (get1bit())
  496.           val -= a;
  497.  
  498.       for(i=0;i<4;i++) {
  499.         if(!(i & 1)) {
  500.           if(!mc) {
  501.             mc = *m++;
  502.             xrpnt = ((real *) xr) + (*m++);
  503.             lwin = *m++;
  504.             cb = *m++;
  505.             if(lwin == 3) {
  506.               v = gr_info->pow2gain[(*scf++) << shift];
  507.               step = 1;
  508.             }
  509.             else {
  510.               v = gr_info->full_gain[lwin][(*scf++) << shift];
  511.               step = 3;
  512.             }
  513.           }
  514.           mc--;
  515.         }
  516.         if( (a & (0x8>>i)) ) {
  517.           max[lwin] = cb;
  518.           if(get1bit()) 
  519.             *xrpnt = -v;
  520.           else
  521.             *xrpnt = v;
  522.         }
  523.         else
  524.           *xrpnt = 0.0;
  525.         xrpnt += step;
  526.       }
  527.     }
  528.  
  529.     while( m < me ) {
  530.       if(!mc) {
  531.         mc = *m++;
  532.         xrpnt = ((real *) xr) + *m++;
  533.         if( (*m++) == 3)
  534.           step = 1;
  535.         else
  536.           step = 3;
  537.         m++; /* cb */
  538.       }
  539.       mc--;
  540.       *xrpnt = 0.0;
  541.       xrpnt += step;
  542.       *xrpnt = 0.0;
  543.       xrpnt += step;
  544.     }
  545.  
  546.     gr_info->maxband[0] = max[0]+1;
  547.     gr_info->maxband[1] = max[1]+1;
  548.     gr_info->maxband[2] = max[2]+1;
  549.     gr_info->maxbandl = max[3]+1;
  550.   }
  551.   else {
  552.     int *pretab = gr_info->preflag ? pretab1 : pretab2;
  553.     int i,max = -1;
  554.     int cb = 0;
  555.     register int *m = map[sfreq][2];
  556.     register real v = 0.0;
  557.     register int mc = 0;
  558.     me = mapend[sfreq][2];
  559.  
  560.     for(i=0;i<3;i++) {
  561.       int lp = l[i];
  562.       struct newhuff *h = ht+gr_info->table_select[i];
  563.  
  564.       for(;lp;lp--,mc--) {
  565.         int x,y;
  566.         if(!mc) {
  567.           mc = *m++;
  568.           v = gr_info->pow2gain[((*scf++) + (*pretab++)) << shift];
  569.           cb = *m++;
  570.         }
  571.         {
  572.           register short *val = h->table;
  573.           while((y=*val++)<0)
  574.             if (get1bit())
  575.               val -= y;
  576.           x = y >> 4;
  577.           y &= 0xf;
  578.         }
  579.         if (x == 15) {
  580.           max = cb;
  581.           x += getbits(h->linbits);
  582.           if(get1bit())
  583.             *xrpnt++ = -ispow[x] * v;
  584.           else
  585.             *xrpnt++ =  ispow[x] * v;
  586.         }
  587.         else if(x) {
  588.           max = cb;
  589.           if(get1bit())
  590.             *xrpnt++ = -ispow[x] * v;
  591.           else
  592.             *xrpnt++ =  ispow[x] * v;
  593.         }
  594.         else
  595.           *xrpnt++ = 0.0;
  596.  
  597.         if (y == 15) {
  598.           max = cb;
  599.           y += getbits(h->linbits);
  600.           if(get1bit())
  601.             *xrpnt++ = -ispow[y] * v;
  602.           else
  603.             *xrpnt++ =  ispow[y] * v;
  604.         }
  605.         else if(y) {
  606.           max = cb;
  607.           if(get1bit())
  608.             *xrpnt++ = -ispow[y] * v;
  609.           else
  610.             *xrpnt++ =  ispow[y] * v;
  611.         }
  612.         else
  613.           *xrpnt++ = 0.0;
  614.       }
  615.     }
  616.     for(;l3 && hsstell() < part2end;l3--) {
  617.       struct newhuff *h = htc+gr_info->count1table_select;
  618.       register short *val = h->table,a;
  619.  
  620.       while((a=*val++)<0)
  621.         if (get1bit())
  622.           val -= a;
  623.  
  624.       for(i=0;i<4;i++) {
  625.         if(!(i & 1)) {
  626.           if(!mc) {
  627.             mc = *m++;
  628.             cb = *m++;
  629.             v = gr_info->pow2gain[((*scf++) + (*pretab++)) << shift];
  630.           }
  631.           mc--;
  632.         }
  633.         if ( (a & (0x8>>i)) ) {
  634.           max = cb;
  635.           if(get1bit())
  636.             *xrpnt++ = -v;
  637.           else
  638.             *xrpnt++ = v;
  639.         }
  640.         else
  641.           *xrpnt++ = 0.0;
  642.       }
  643.     }
  644.     for(i=(&xr[SBLIMIT][SSLIMIT]-xrpnt)>>1;i;i--) {
  645.       *xrpnt++ = 0.0;
  646.       *xrpnt++ = 0.0;
  647.     }
  648.  
  649.     gr_info->maxbandl = max+1;
  650.   }
  651.  
  652.   {
  653.     int h;
  654.     if( (h = (part2end-hsstell()) ) ) {
  655.       while ( h > 16 ) {
  656.          getbits(16); /* Dismiss stuffing Bits */
  657.          h -= 16;
  658.       }
  659.       if(h >= 0 )
  660.          getbits(h);
  661.       else {
  662.          fprintf(stderr,"mpg123: Can't rewind stream by %d bits!\n",-h);
  663.          exit(1);
  664.       }
  665.     }
  666.   }
  667. }
  668.  
  669. static void III_dequantize_sample_ms(real xr[2][SBLIMIT][SSLIMIT],int *scf,
  670.    struct gr_info_s *gr_info,int sfreq,int part2_start)
  671. {
  672.   int shift = 1 + gr_info->scalefac_scale;
  673.   real *xrpnt = (real *) xr[1];
  674.   real *xr0pnt = (real *) xr[0];
  675.   int l[3],l3;
  676.   int part2end = gr_info->part2_3_length + part2_start;
  677.   int *me;
  678.  
  679.   {
  680.     int bv       = gr_info->big_values;
  681.     int region1  = gr_info->region1start;
  682.     int region2  = gr_info->region2start;
  683.  
  684.     l3 = ((576>>1)-bv)>>1;   
  685. /*
  686.  * we may lose the 'odd' bit here !! 
  687.  * check this later gain 
  688.  */
  689.     if(bv <= region1) {
  690.       l[0] = bv; l[1] = 0; l[2] = 0;
  691.     }
  692.     else {
  693.       l[0] = region1;
  694.       if(bv <= region2) {
  695.         l[1] = bv - l[0];  l[2] = 0;
  696.       }
  697.       else {
  698.         l[1] = region2 - l[0]; l[2] = bv - region2;
  699.       }
  700.     }
  701.   }
  702.  
  703.   if(gr_info->block_type == 2) {
  704.     int i,max[4];
  705.     int step=0,lwin=0,cb=0;
  706.     register real v = 0.0;
  707.     register int *m,mc = 0;
  708.  
  709.     if(gr_info->mixed_block_flag) {
  710.       max[3] = -1;
  711.       max[0] = max[1] = max[2] = 2;
  712.       m = map[sfreq][0];
  713.       me = mapend[sfreq][0];
  714.     }
  715.     else {
  716.       max[0] = max[1] = max[2] = max[3] = -1;
  717.       /* max[3] not really needed in this case */
  718.       m = map[sfreq][1];
  719.       me = mapend[sfreq][1];
  720.     }
  721.  
  722.     for(i=0;i<2;i++) {
  723.       int lp = l[i];
  724.       struct newhuff *h = ht+gr_info->table_select[i];
  725.       for(;lp;lp--,mc--) {
  726.         int x,y;
  727.         if(!mc) {
  728.           mc = *m++;
  729.           xrpnt = ((real *) xr[1]) + *m;
  730.           xr0pnt = ((real *) xr[0]) + *m++;
  731.           lwin = *m++;
  732.           cb = *m++;
  733.           if(lwin == 3) {
  734.             v = gr_info->pow2gain[(*scf++) << shift];
  735.             step = 1;
  736.           }
  737.           else {
  738.             v = gr_info->full_gain[lwin][(*scf++) << shift];
  739.             step = 3;
  740.           }
  741.         }
  742.         {
  743.           register short *val = h->table;
  744.           while((y=*val++)<0)
  745.             if (get1bit())
  746.               val -= y;
  747.           x = y >> 4;
  748.           y &= 0xf;
  749.         }
  750.         if(x == 15) {
  751.           max[lwin] = cb;
  752.           x += getbits(h->linbits);
  753.           if(get1bit()) {
  754.             real a = ispow[x] * v;
  755.             *xrpnt = *xr0pnt + a;
  756.             *xr0pnt -= a;
  757.           }
  758.           else {
  759.             real a = ispow[x] * v;
  760.             *xrpnt = *xr0pnt - a;
  761.             *xr0pnt += a;
  762.           }
  763.         }
  764.         else if(x) {
  765.           max[lwin] = cb;
  766.           if(get1bit()) {
  767.             real a = ispow[x] * v;
  768.             *xrpnt = *xr0pnt + a;
  769.             *xr0pnt -= a;
  770.           }
  771.           else {
  772.             real a = ispow[x] * v;
  773.             *xrpnt = *xr0pnt - a;
  774.             *xr0pnt += a;
  775.           }
  776.         }
  777.         else
  778.           *xrpnt = *xr0pnt;
  779.         xrpnt += step;
  780.         xr0pnt += step;
  781.  
  782.         if(y == 15) {
  783.           max[lwin] = cb;
  784.           y += getbits(h->linbits);
  785.           if(get1bit()) {
  786.             real a = ispow[y] * v;
  787.             *xrpnt = *xr0pnt + a;
  788.             *xr0pnt -= a;
  789.           }
  790.           else {
  791.             real a = ispow[y] * v;
  792.             *xrpnt = *xr0pnt - a;
  793.             *xr0pnt += a;
  794.           }
  795.         }
  796.         else if(y) {
  797.           max[lwin] = cb;
  798.           if(get1bit()) {
  799.             real a = ispow[y] * v;
  800.             *xrpnt = *xr0pnt + a;
  801.             *xr0pnt -= a;
  802.           }
  803.           else {
  804.             real a = ispow[y] * v;
  805.             *xrpnt = *xr0pnt - a;
  806.             *xr0pnt += a;
  807.           }
  808.         }
  809.         else
  810.           *xrpnt = *xr0pnt;
  811.         xrpnt += step;
  812.         xr0pnt += step;
  813.       }
  814.     }
  815.  
  816.     for(;l3 && hsstell() < part2end;l3--) {
  817.       struct newhuff *h = htc+gr_info->count1table_select;
  818.       register short *val = h->table,a;
  819.  
  820.       while((a=*val++)<0)
  821.         if (get1bit())
  822.           val -= a;
  823.  
  824.       for(i=0;i<4;i++) {
  825.         if(!(i & 1)) {
  826.           if(!mc) {
  827.             mc = *m++;
  828.             xrpnt = ((real *) xr[1]) + *m;
  829.             xr0pnt = ((real *) xr[0]) + *m++;
  830.             lwin = *m++;
  831.             cb = *m++;
  832.             if(lwin == 3) {
  833.               v = gr_info->pow2gain[(*scf++) << shift];
  834.               step = 1;
  835.             }
  836.             else {
  837.               v = gr_info->full_gain[lwin][(*scf++) << shift];
  838.               step = 3;
  839.             }
  840.           }
  841.           mc--;
  842.         }
  843.         if( (a & (0x8>>i)) ) {
  844.           max[lwin] = cb;
  845.           if(get1bit()) {
  846.             *xrpnt = *xr0pnt + v;
  847.             *xr0pnt -= v;
  848.           }
  849.           else {
  850.             *xrpnt = *xr0pnt - v;
  851.             *xr0pnt += v;
  852.           }
  853.         }
  854.         else
  855.           *xrpnt = *xr0pnt;
  856.         xrpnt += step;
  857.         xr0pnt += step;
  858.       }
  859.     }
  860.  
  861.     while( m < me ) {
  862.       if(!mc) {
  863.         mc = *m++;
  864.         xrpnt = ((real *) xr) + *m;
  865.         xr0pnt = ((real *) xr) + *m++;
  866.         if(*m++ == 3)
  867.           step = 1;
  868.         else
  869.           step = 3;
  870.         m++; /* cb */
  871.       }
  872.       mc--;
  873.       *xrpnt = *xr0pnt;
  874.       xrpnt += step;
  875.       xr0pnt += step;
  876.       *xrpnt = *xr0pnt;
  877.       xrpnt += step;
  878.       xr0pnt += step;
  879.     }
  880.  
  881.     gr_info->maxband[0] = max[0]+1;
  882.     gr_info->maxband[1] = max[1]+1;
  883.     gr_info->maxband[2] = max[2]+1;
  884.     gr_info->maxbandl = max[3]+1;
  885.   }
  886.   else {
  887.     int *pretab = gr_info->preflag ? pretab1 : pretab2;
  888.     int i,max = -1;
  889.     int cb = 0;
  890.     register int mc=0,*m = map[sfreq][2];
  891.     register real v = 0.0;
  892.     me = mapend[sfreq][2];
  893.  
  894.     for(i=0;i<3;i++) {
  895.       int lp = l[i];
  896.       struct newhuff *h = ht+gr_info->table_select[i];
  897.  
  898.       for(;lp;lp--,mc--) {
  899.         int x,y;
  900.         if(!mc) {
  901.           mc = *m++;
  902.           cb = *m++;
  903.           v = gr_info->pow2gain[((*scf++) + (*pretab++)) << shift];
  904.         }
  905.         {
  906.           register short *val = h->table;
  907.           while((y=*val++)<0)
  908.             if (get1bit())
  909.               val -= y;
  910.           x = y >> 4;
  911.           y &= 0xf;
  912.         }
  913.         if (x == 15) {
  914.           max = cb;
  915.           x += getbits(h->linbits);
  916.           if(get1bit()) {
  917.             real a = ispow[x] * v;
  918.             *xrpnt++ = *xr0pnt + a;
  919.             *xr0pnt++ -= a;
  920.           }
  921.           else {
  922.             real a = ispow[x] * v;
  923.             *xrpnt++ = *xr0pnt - a;
  924.             *xr0pnt++ += a;
  925.           }
  926.         }
  927.         else if(x) {
  928.           max = cb;
  929.           if(get1bit()) {
  930.             real a = ispow[x] * v;
  931.             *xrpnt++ = *xr0pnt + a;
  932.             *xr0pnt++ -= a;
  933.           }
  934.           else {
  935.             real a = ispow[x] * v;
  936.             *xrpnt++ = *xr0pnt - a;
  937.             *xr0pnt++ += a;
  938.           }
  939.         }
  940.         else
  941.           *xrpnt++ = *xr0pnt++;
  942.  
  943.         if (y == 15) {
  944.           max = cb;
  945.           y += getbits(h->linbits);
  946.           if(get1bit()) {
  947.             real a = ispow[y] * v;
  948.             *xrpnt++ = *xr0pnt + a;
  949.             *xr0pnt++ -= a;
  950.           }
  951.           else {
  952.             real a = ispow[y] * v;
  953.             *xrpnt++ = *xr0pnt - a;
  954.             *xr0pnt++ += a;
  955.           }
  956.         }
  957.         else if(y) {
  958.           max = cb;
  959.           if(get1bit()) {
  960.             real a = ispow[y] * v;
  961.             *xrpnt++ = *xr0pnt + a;
  962.             *xr0pnt++ -= a;
  963.           }
  964.           else {
  965.             real a = ispow[y] * v;
  966.             *xrpnt++ = *xr0pnt - a;
  967.             *xr0pnt++ += a;
  968.           }
  969.         }
  970.         else
  971.           *xrpnt++ = *xr0pnt++;
  972.       }
  973.     }
  974.     for(;l3 && hsstell() < part2end;l3--) {
  975.       struct newhuff *h = htc+gr_info->count1table_select;
  976.       register short *val = h->table,a;
  977.  
  978.       while((a=*val++)<0)
  979.         if (get1bit())
  980.           val -= a;
  981.  
  982.       for(i=0;i<4;i++) {
  983.         if(!(i & 1)) {
  984.           if(!mc) {
  985.             mc = *m++;
  986.             cb = *m++;
  987.             v = gr_info->pow2gain[((*scf++) + (*pretab++)) << shift];
  988.           }
  989.           mc--;
  990.         }
  991.         if ( (a & (0x8>>i)) ) {
  992.           max = cb;
  993.           if(get1bit()) {
  994.             *xrpnt++ = *xr0pnt + v;
  995.             *xr0pnt++ -= v;
  996.           }
  997.           else {
  998.             *xrpnt++ = *xr0pnt - v;
  999.             *xr0pnt++ += v;
  1000.           }
  1001.         }
  1002.         else
  1003.           *xrpnt++ = *xr0pnt++;
  1004.       }
  1005.     }
  1006.     for(i=(&xr[1][SBLIMIT][SSLIMIT]-xrpnt)>>1;i;i--) {
  1007.       *xrpnt++ = *xr0pnt++;
  1008.       *xrpnt++ = *xr0pnt++;
  1009.     }
  1010.  
  1011.     gr_info->maxbandl = max+1;
  1012.   }
  1013.  
  1014.   {
  1015.     int h;
  1016.     if( (h = (part2end-hsstell()) ) ) {
  1017.       while ( h > 16 ) {
  1018.          getbits(16); /* Dismiss stuffing Bits */
  1019.          h -= 16;
  1020.       }
  1021.       if(h >= 0 )
  1022.          getbits(h);
  1023.       else {
  1024.          fprintf(stderr,"mpg123: Can't rewind stream by %d bits!\n",-h);
  1025.          exit(1);
  1026.       }
  1027.     }
  1028.   }
  1029. }
  1030.  
  1031. /* 
  1032.  * III_stereo: calculate real channel values for Joint-I-Stereo-mode
  1033.  */
  1034. static void III_stereo(real xr_buf[2][SBLIMIT][SSLIMIT],int *scalefac,
  1035.    struct gr_info_s *gr_info,int sfreq,int ms_stereo)
  1036. {
  1037.       real (*xr)[SBLIMIT*SSLIMIT] = (real (*)[SBLIMIT*SSLIMIT] ) xr_buf;
  1038.       struct bandInfoStruct *bi = &bandInfo[sfreq];
  1039.       real *tan1,*tan2;
  1040.  
  1041.       if(ms_stereo) {
  1042.         tan1 = tan1_2; tan2 = tan2_2;
  1043.       }
  1044.       else {
  1045.         tan1 = tan1_1; tan2 = tan2_1;
  1046.       }
  1047.  
  1048.       if (gr_info->block_type == 2)
  1049.       {
  1050.          int lwin,do_l = 0;
  1051.          if( gr_info->mixed_block_flag )
  1052.            do_l = 1;
  1053.  
  1054.          for (lwin=0;lwin<3;lwin++) /* process each window */
  1055.          {
  1056.              /* get first band with zero values */
  1057.            int is_p,sb,idx,sfb = gr_info->maxband[lwin];  /* sfb is minimal 3 for mixed mode */
  1058.            if(sfb > 3)
  1059.              do_l = 0;
  1060.  
  1061.            for(;sfb<12;sfb++)
  1062.            {
  1063.              is_p = scalefac[sfb*3+lwin-gr_info->mixed_block_flag]; /* scale: 0-15 */ 
  1064.              if(is_p != 7) {
  1065.                real t1,t2;
  1066.                sb = bi->shortDiff[sfb];
  1067.                idx = bi->shortIdx[sfb] + lwin;
  1068.                t1 = tan1[is_p]; t2 = tan2[is_p];
  1069.                for (; sb > 0; sb--,idx+=3)
  1070.                {
  1071.                  real v = xr[0][idx];
  1072.                  xr[0][idx] = v * t1;
  1073.                  xr[1][idx] = v * t2;
  1074.                }
  1075.              }
  1076.            }
  1077.  
  1078. #if 1
  1079. /* in the original: copy 10 to 11 , here: copy 11 to 12 
  1080. maybe still wrong??? (copy 12 to 13?) */
  1081.            is_p = scalefac[11*3+lwin-gr_info->mixed_block_flag]; /* scale: 0-15 */
  1082.            sb = bi->shortDiff[12];
  1083.            idx = bi->shortIdx[12] + lwin;
  1084. #else
  1085.            is_p = scalefac[10*3+lwin-gr_info->mixed_block_flag]; /* scale: 0-15 */
  1086.            sb = bi->shortDiff[11];
  1087.            idx = bi->shortIdx[11] + lwin;
  1088. #endif
  1089.            if(is_p != 7)
  1090.            {
  1091.              real t1,t2;
  1092.              t1 = tan1[is_p]; t2 = tan2[is_p];
  1093.              for ( ; sb > 0; sb--,idx+=3 )
  1094.              {  
  1095.                real v = xr[0][idx];
  1096.                xr[0][idx] = v * t1;
  1097.                xr[1][idx] = v * t2;
  1098.              }
  1099.            }
  1100.          } /* end for(lwin; .. ; . ) */
  1101.  
  1102.          if (do_l)
  1103.          {
  1104. /* also check l-part, if ALL bands in the three windows are 'empty'
  1105.  * and mode = mixed_mode 
  1106.  */
  1107.            int sfb = gr_info->maxbandl;
  1108.            int idx = bi->longIdx[sfb];
  1109.  
  1110.            for ( ; sfb<8; sfb++ )
  1111.            {
  1112.              int sb = bi->longDiff[sfb];
  1113.              int is_p = scalefac[sfb]; /* scale: 0-15 */
  1114.              if(is_p != 7) {
  1115.                real t1,t2;
  1116.                t1 = tan1[is_p]; t2 = tan2[is_p];
  1117.                for ( ; sb > 0; sb--,idx++)
  1118.                {
  1119.                  real v = xr[0][idx];
  1120.                  xr[0][idx] = v * t1;
  1121.                  xr[1][idx] = v * t2;
  1122.                }
  1123.              }
  1124.              else 
  1125.                idx += sb;
  1126.            }
  1127.          }     
  1128.       } 
  1129.       else /* ((gr_info->block_type != 2)) */
  1130.       {
  1131.         int sfb = gr_info->maxbandl;
  1132.         int is_p,idx = bi->longIdx[sfb];
  1133.         for ( ; sfb<21; sfb++)
  1134.         {
  1135.           int sb = bi->longDiff[sfb];
  1136.           is_p = scalefac[sfb]; /* scale: 0-15 */
  1137.           if(is_p != 7) {
  1138.             real t1,t2;
  1139.             t1 = tan1[is_p]; t2 = tan2[is_p];
  1140.             for ( ; sb > 0; sb--,idx++)
  1141.             {
  1142.                real v = xr[0][idx];
  1143.                xr[0][idx] = v * t1;
  1144.                xr[1][idx] = v * t2;
  1145.             }
  1146.           }
  1147.           else
  1148.             idx += sb;
  1149.         }
  1150.  
  1151.         is_p = scalefac[20]; /* copy l-band 20 to l-band 21 */
  1152.         if(is_p != 7)
  1153.         {
  1154.           int sb;
  1155.           real t1 = tan1[is_p],t2 = tan2[is_p]; 
  1156.  
  1157.           for ( sb = bi->longDiff[21]; sb > 0; sb--,idx++ )
  1158.           {
  1159.             real v = xr[0][idx];
  1160.             xr[0][idx] = v * t1;
  1161.             xr[1][idx] = v * t2;
  1162.           }
  1163.         }
  1164.       } /* ... */
  1165. }
  1166.  
  1167. static void III_antialias(real xr[SBLIMIT][SSLIMIT],struct gr_info_s *gr_info)
  1168. {
  1169.    int sblim;
  1170.  
  1171.    if(gr_info->block_type == 2)
  1172.    {
  1173.       if(!gr_info->mixed_block_flag) 
  1174.         return;
  1175.       else
  1176.        sblim = 1; 
  1177.    }
  1178.    else
  1179.      sblim = SBLIMIT-1;
  1180.  
  1181.    /* 31 alias-reduction operations between each pair of sub-bands */
  1182.    /* with 8 butterflies between each pair                         */
  1183.  
  1184.    {
  1185.      int sb;
  1186.      real *xr1=(real *) xr[1];
  1187.  
  1188.      for(sb=sblim;sb;sb--,xr1+=10)
  1189.      {
  1190.        int ss;
  1191.        real *cs=aa_cs,*ca=aa_ca;
  1192.        real *xr2 = xr1;
  1193.        for(ss=7;ss>=0;ss--)
  1194.        {       /* upper and lower butterfly inputs */
  1195.          register real bu = *--xr2,bd = *xr1;
  1196.          *xr2   = (bu * (*cs)   ) - (bd * (*ca)   );
  1197.          *xr1++ = (bd * (*cs++) ) + (bu * (*ca++) );
  1198.        }
  1199.      }
  1200.   }
  1201. }
  1202.  
  1203. /* This is an optimized DCT from Jeff Tsay's maplay 1.2+ package.
  1204.  I've done the 'manual unroll ;)' and saved one multiplication
  1205.  by doing it together with the window mul. (MH)
  1206.  
  1207. This uses Byeong Gi Lee's Fast Cosine Transform algorithm, but the
  1208. 9 point IDCT needs to be reduced further. Unfortunately, I don't
  1209. know how to do that, because 9 is not an even number. - Jeff.*/
  1210.  
  1211. /*------------------------------------------------------------------*/
  1212. /*                                                                  */
  1213. /*    Function: Calculation of the inverse MDCT                     */
  1214. /*    In the case of short blocks the 3 output vectors are already  */
  1215. /*    overlapped and added in this modul.                           */
  1216. /*                                                                  */
  1217. /*    New layer3                                                    */
  1218. /*                                                                  */
  1219. /*------------------------------------------------------------------*/
  1220.  
  1221. static void dct36(real *inbuf,real *o1,real *o2,real *wintab,real *tsbuf)
  1222. {
  1223.  {
  1224.   register real *in = inbuf;
  1225.  
  1226.   in[17]+=in[16]; in[16]+=in[15]; in[15]+=in[14];
  1227.   in[14]+=in[13]; in[13]+=in[12]; in[12]+=in[11];
  1228.   in[11]+=in[10]; in[10]+=in[9];  in[9] +=in[8];
  1229.   in[8] +=in[7];  in[7] +=in[6];  in[6] +=in[5];
  1230.   in[5] +=in[4];  in[4] +=in[3];  in[3] +=in[2];
  1231.   in[2] +=in[1];  in[1] +=in[0];
  1232.  
  1233.   in[17]+=in[15]; in[15]+=in[13]; in[13]+=in[11]; in[11]+=in[9];
  1234.   in[9] +=in[7];  in[7] +=in[5];  in[5] +=in[3];  in[3] +=in[1];
  1235.  
  1236.   {
  1237.  
  1238. #define MACRO0(v) { \
  1239.     real tmp; \
  1240.     out2[9+(v)] = (tmp = sum0 + sum1) * w[27+(v)]; \
  1241.     out2[8-(v)] = tmp * w[26-(v)];  } \
  1242.     sum0 -= sum1; \
  1243.     ts[SBLIMIT*(8-(v))] = out1[8-(v)] + sum0 * w[8-(v)]; \
  1244.     ts[SBLIMIT*(9+(v))] = out1[9+(v)] + sum0 * w[9+(v)]; 
  1245. #define MACRO1(v) { \
  1246.     real sum0,sum1; \
  1247.     sum0 = tmp1a + tmp2a; \
  1248.     sum1 = (tmp1b + tmp2b) * tfcos36[(v)]; \
  1249.     MACRO0(v); }
  1250. #define MACRO2(v) { \
  1251.     real sum0,sum1; \
  1252.     sum0 = tmp2a - tmp1a; \
  1253.     sum1 = (tmp2b - tmp1b) * tfcos36[(v)]; \
  1254.     MACRO0(v); }
  1255.  
  1256.     register const real *c = COS9;
  1257.     register real *out2 = o2;
  1258.     register real *w = wintab;
  1259.     register real *out1 = o1;
  1260.     register real *ts = tsbuf;
  1261.  
  1262.     real ta33,ta66,tb33,tb66;
  1263.  
  1264.     ta33 = in[2*3+0] * c[3];
  1265.     ta66 = in[2*6+0] * c[6];
  1266.     tb33 = in[2*3+1] * c[3];
  1267.     tb66 = in[2*6+1] * c[6];
  1268.  
  1269.     { 
  1270.       real tmp1a,tmp2a,tmp1b,tmp2b;
  1271.       tmp1a =             in[2*1+0] * c[1] + ta33 + in[2*5+0] * c[5] + in[2*7+0] * c[7];
  1272.       tmp1b =             in[2*1+1] * c[1] + tb33 + in[2*5+1] * c[5] + in[2*7+1] * c[7];
  1273.       tmp2a = in[2*0+0] + in[2*2+0] * c[2] + in[2*4+0] * c[4] + ta66 + in[2*8+0] * c[8];
  1274.       tmp2b = in[2*0+1] + in[2*2+1] * c[2] + in[2*4+1] * c[4] + tb66 + in[2*8+1] * c[8];
  1275.  
  1276.       MACRO1(0);
  1277.       MACRO2(8);
  1278.     }
  1279.  
  1280.     {
  1281.       real tmp1a,tmp2a,tmp1b,tmp2b;
  1282.       tmp1a = ( in[2*1+0] - in[2*5+0] - in[2*7+0] ) * c[3];
  1283.       tmp1b = ( in[2*1+1] - in[2*5+1] - in[2*7+1] ) * c[3];
  1284.       tmp2a = ( in[2*2+0] - in[2*4+0] - in[2*8+0] ) * c[6] - in[2*6+0] + in[2*0+0];
  1285.       tmp2b = ( in[2*2+1] - in[2*4+1] - in[2*8+1] ) * c[6] - in[2*6+1] + in[2*0+1];
  1286.  
  1287.       MACRO1(1);
  1288.       MACRO2(7);
  1289.     }
  1290.  
  1291.     {
  1292.       real tmp1a,tmp2a,tmp1b,tmp2b;
  1293.       tmp1a =             in[2*1+0] * c[5] - ta33 - in[2*5+0] * c[7] + in[2*7+0] * c[1];
  1294.       tmp1b =             in[2*1+1] * c[5] - tb33 - in[2*5+1] * c[7] + in[2*7+1] * c[1];
  1295.       tmp2a = in[2*0+0] - in[2*2+0] * c[8] - in[2*4+0] * c[2] + ta66 + in[2*8+0] * c[4];
  1296.       tmp2b = in[2*0+1] - in[2*2+1] * c[8] - in[2*4+1] * c[2] + tb66 + in[2*8+1] * c[4];
  1297.  
  1298.       MACRO1(2);
  1299.       MACRO2(6);
  1300.     }
  1301.  
  1302.     {
  1303.       real tmp1a,tmp2a,tmp1b,tmp2b;
  1304.       tmp1a =             in[2*1+0] * c[7] - ta33 + in[2*5+0] * c[1] - in[2*7+0] * c[5];
  1305.       tmp1b =             in[2*1+1] * c[7] - tb33 + in[2*5+1] * c[1] - in[2*7+1] * c[5];
  1306.       tmp2a = in[2*0+0] - in[2*2+0] * c[4] + in[2*4+0] * c[8] + ta66 - in[2*8+0] * c[2];
  1307.       tmp2b = in[2*0+1] - in[2*2+1] * c[4] + in[2*4+1] * c[8] + tb66 - in[2*8+1] * c[2];
  1308.  
  1309.       MACRO1(3);
  1310.       MACRO2(5);
  1311.     }
  1312.  
  1313.     {
  1314.         real sum0,sum1;
  1315.         sum0 =  in[2*0+0] - in[2*2+0] + in[2*4+0] - in[2*6+0] + in[2*8+0];
  1316.         sum1 = (in[2*0+1] - in[2*2+1] + in[2*4+1] - in[2*6+1] + in[2*8+1] ) * tfcos36[4];
  1317.         MACRO0(4);
  1318.     }
  1319.   }
  1320.  }
  1321. }
  1322.  
  1323. /*
  1324.  * new DCT12
  1325.  */
  1326. static void dct12(real *in,real *rawout1,real *rawout2,register real *wi,register real *ts)
  1327. {
  1328. #define DCT12_PART1 \
  1329.              in5 = in[5*3];  \
  1330.      in5 += (in4 = in[4*3]); \
  1331.      in4 += (in3 = in[3*3]); \
  1332.      in3 += (in2 = in[2*3]); \
  1333.      in2 += (in1 = in[1*3]); \
  1334.      in1 += (in0 = in[0*3]); \
  1335.                              \
  1336.      in5 += in3; in3 += in1; \
  1337.                              \
  1338.      in2 *= COS6_1; \
  1339.      in3 *= COS6_1; \
  1340.  
  1341. #define DCT12_PART2 \
  1342.      in0 += in4 * COS6_2; \
  1343.                           \
  1344.      in4 = in0 + in2;     \
  1345.      in0 -= in2;          \
  1346.                           \
  1347.      in1 += in5 * COS6_2; \
  1348.                           \
  1349.      in5 = (in1 + in3) * tfcos12[0]; \
  1350.      in1 = (in1 - in3) * tfcos12[2]; \
  1351.                          \
  1352.      in3 = in4 + in5;    \
  1353.      in4 -= in5;         \
  1354.                          \
  1355.      in2 = in0 + in1;    \
  1356.      in0 -= in1;
  1357.  
  1358.  
  1359.    {
  1360.      real in0,in1,in2,in3,in4,in5;
  1361.      register real *out1 = rawout1;
  1362.      ts[SBLIMIT*0] = out1[0]; ts[SBLIMIT*1] = out1[1]; ts[SBLIMIT*2] = out1[2];
  1363.      ts[SBLIMIT*3] = out1[3]; ts[SBLIMIT*4] = out1[4]; ts[SBLIMIT*5] = out1[5];
  1364.  
  1365.      DCT12_PART1
  1366.  
  1367.      {
  1368.        real tmp0,tmp1 = (in0 - in4);
  1369.        {
  1370.          real tmp2 = (in1 - in5) * tfcos12[1];
  1371.          tmp0 = tmp1 + tmp2;
  1372.          tmp1 -= tmp2;
  1373.        }
  1374.        ts[(17-1)*SBLIMIT] = out1[17-1] + tmp0 * wi[11-1];
  1375.        ts[(12+1)*SBLIMIT] = out1[12+1] + tmp0 * wi[6+1];
  1376.        ts[(6 +1)*SBLIMIT] = out1[6 +1] + tmp1 * wi[1];
  1377.        ts[(11-1)*SBLIMIT] = out1[11-1] + tmp1 * wi[5-1];
  1378.      }
  1379.  
  1380.      DCT12_PART2
  1381.  
  1382.      ts[(17-0)*SBLIMIT] = out1[17-0] + in2 * wi[11-0];
  1383.      ts[(12+0)*SBLIMIT] = out1[12+0] + in2 * wi[6+0];
  1384.      ts[(12+2)*SBLIMIT] = out1[12+2] + in3 * wi[6+2];
  1385.      ts[(17-2)*SBLIMIT] = out1[17-2] + in3 * wi[11-2];
  1386.  
  1387.      ts[(6+0)*SBLIMIT]  = out1[6+0] + in0 * wi[0];
  1388.      ts[(11-0)*SBLIMIT] = out1[11-0] + in0 * wi[5-0];
  1389.      ts[(6+2)*SBLIMIT]  = out1[6+2] + in4 * wi[2];
  1390.      ts[(11-2)*SBLIMIT] = out1[11-2] + in4 * wi[5-2];
  1391.   }
  1392.  
  1393.   in++;
  1394.  
  1395.   {
  1396.      real in0,in1,in2,in3,in4,in5;
  1397.      register real *out2 = rawout2;
  1398.  
  1399.      DCT12_PART1
  1400.  
  1401.      {
  1402.        real tmp0,tmp1 = (in0 - in4);
  1403.        {
  1404.          real tmp2 = (in1 - in5) * tfcos12[1];
  1405.          tmp0 = tmp1 + tmp2;
  1406.          tmp1 -= tmp2;
  1407.        }
  1408.        out2[5-1] = tmp0 * wi[11-1];
  1409.        out2[0+1] = tmp0 * wi[6+1];
  1410.        ts[(12+1)*SBLIMIT] += tmp1 * wi[1];
  1411.        ts[(17-1)*SBLIMIT] += tmp1 * wi[5-1];
  1412.      }
  1413.  
  1414.      DCT12_PART2
  1415.  
  1416.      out2[5-0] = in2 * wi[11-0];
  1417.      out2[0+0] = in2 * wi[6+0];
  1418.      out2[0+2] = in3 * wi[6+2];
  1419.      out2[5-2] = in3 * wi[11-2];
  1420.  
  1421.      ts[(12+0)*SBLIMIT] += in0 * wi[0];
  1422.      ts[(17-0)*SBLIMIT] += in0 * wi[5-0];
  1423.      ts[(12+2)*SBLIMIT] += in4 * wi[2];
  1424.      ts[(17-2)*SBLIMIT] += in4 * wi[5-2];
  1425.   }
  1426.  
  1427.   in++; 
  1428.  
  1429.   {
  1430.      real in0,in1,in2,in3,in4,in5;
  1431.      register real *out2 = rawout2;
  1432.      out2[12]=out2[13]=out2[14]=out2[15]=out2[16]=out2[17]=0.0;
  1433.  
  1434.      DCT12_PART1
  1435.  
  1436.      {
  1437.        real tmp0,tmp1 = (in0 - in4);
  1438.        {
  1439.          real tmp2 = (in1 - in5) * tfcos12[1];
  1440.          tmp0 = tmp1 + tmp2;
  1441.          tmp1 -= tmp2;
  1442.        }
  1443.        out2[11-1] = tmp0 * wi[11-1];
  1444.        out2[6 +1] = tmp0 * wi[6+1];
  1445.        out2[0+1] += tmp1 * wi[1];
  1446.        out2[5-1] += tmp1 * wi[5-1];
  1447.      }
  1448.  
  1449.      DCT12_PART2
  1450.  
  1451.      out2[11-0] = in2 * wi[11-0];
  1452.      out2[6 +0] = in2 * wi[6+0];
  1453.      out2[6 +2] = in3 * wi[6+2];
  1454.      out2[11-2] = in3 * wi[11-2];
  1455.  
  1456.      out2[0+0] += in0 * wi[0];
  1457.      out2[5-0] += in0 * wi[5-0];
  1458.      out2[0+2] += in4 * wi[2];
  1459.      out2[5-2] += in4 * wi[5-2];
  1460.   }
  1461. }
  1462.  
  1463.  
  1464. /*
  1465.  * III_hybrid:
  1466.  *   we still need a faster DCT for the dct36 
  1467.  */
  1468. static void III_hybrid(real fsIn[SBLIMIT][SSLIMIT],real tsOut[SSLIMIT][SBLIMIT],
  1469.    int ch,struct gr_info_s *gr_info)
  1470. {
  1471.    real *tspnt = (real *) tsOut;
  1472.    static real block[2][2][SBLIMIT*SSLIMIT] = { { { 0, } } };
  1473.    static int blc[2]={0,0};
  1474.    real *rawout1,*rawout2;
  1475.    int bt1,bt2;
  1476.  
  1477.    {
  1478.      int b = blc[ch];
  1479.      rawout1=block[b][ch];
  1480.      b=-b+1;
  1481.      rawout2=block[b][ch];
  1482.      blc[ch] = b;
  1483.    }
  1484.  
  1485.    bt1 = gr_info->mixed_block_flag ? 0 : gr_info->block_type;
  1486.    bt2 = gr_info->block_type;
  1487.  
  1488.    if(bt2 == 2) {
  1489.      int sb;
  1490.      if(!bt1) {
  1491.          dct36(fsIn[0],rawout1,rawout2,win[0],tspnt);
  1492.          dct36(fsIn[1],rawout1+18,rawout2+18,win1[0],tspnt+1);
  1493.      } 
  1494.      else {
  1495.          dct12(fsIn[0],rawout1,rawout2,win[2],tspnt); /* only called with bt==2 */
  1496.          dct12(fsIn[1],rawout1+18,rawout2+18,win1[2],tspnt+1);
  1497.      }
  1498.      rawout1 += 36; rawout2 += 36; tspnt += 2;
  1499.      for (sb=2; sb<SBLIMIT; sb+=2,tspnt+=2,rawout1+=36,rawout2+=36) {
  1500.          dct12(fsIn[sb],rawout1,rawout2,win[2],tspnt);
  1501.          dct12(fsIn[sb+1],rawout1+18,rawout2+18,win1[2],tspnt+1);
  1502.      }
  1503.    }
  1504.    else {
  1505.      int sb;
  1506.      dct36(fsIn[0],rawout1,rawout2,win[bt1],tspnt);
  1507.      dct36(fsIn[1],rawout1+18,rawout2+18,win1[bt1],tspnt+1);
  1508.      rawout1 += 36; rawout2 += 36; tspnt += 2;
  1509.      for (sb=2; sb<SBLIMIT; sb+=2,tspnt+=2,rawout1+=36,rawout2+=36) {
  1510.        dct36(fsIn[sb],rawout1,rawout2,win[bt2],tspnt);
  1511.        dct36(fsIn[sb+1],rawout1+18,rawout2+18,win1[bt2],tspnt+1);
  1512.      }
  1513.   }
  1514. }
  1515.  
  1516. int do_layer3(struct frame *fr,int outmode,struct audio_info_struct *ai)
  1517. {
  1518.   int gr, ch, ss,clip=0;
  1519.   int scalefacs[2][39];    /* max 39 for short[13][3] mode, mixed: 38, long: 22 */
  1520.   struct III_sideinfo sideinfo;
  1521.   int stereo = fr->stereo;
  1522.   int single = fr->single;
  1523.   int ms_stereo,i_stereo;
  1524.   int sfreq = fr->sampling_frequency;
  1525.   int stereo1;
  1526.  
  1527.   if(stereo == 1) {
  1528.     stereo1 = 1;
  1529.     single = 0; 
  1530.   }
  1531.   else if(single >= 0)
  1532.     stereo1 = 1;
  1533.   else
  1534.     stereo1 = 2;
  1535.  
  1536.   ms_stereo = (fr->mode == MPG_MD_JOINT_STEREO) && (fr->mode_ext & 0x2);
  1537.   i_stereo = (fr->mode == MPG_MD_JOINT_STEREO) && (fr->mode_ext & 0x1);
  1538.  
  1539.   III_get_side_info(&sideinfo,stereo,ms_stereo,sfreq);
  1540.   set_pointer(sideinfo.main_data_begin);
  1541.  
  1542.   for (gr=0;gr<2;gr++) 
  1543.   {
  1544.     static real hybridIn[2][SBLIMIT][SSLIMIT];
  1545.     static real hybridOut[2][SSLIMIT][SBLIMIT];
  1546.  
  1547.     {
  1548.       struct gr_info_s *gr_info = &(sideinfo.ch[0].gr[gr]);
  1549.       long part2_start = hsstell();
  1550.       III_get_scale_factors(scalefacs[0],gr_info);
  1551.       III_dequantize_sample(hybridIn[0], scalefacs[0],gr_info,sfreq,part2_start);
  1552.     }
  1553.     if(stereo == 2) {
  1554.       struct gr_info_s *gr_info = &(sideinfo.ch[1].gr[gr]);
  1555.       long part2_start = hsstell();
  1556.       III_get_scale_factors(scalefacs[1],gr_info);
  1557.       if(ms_stereo)
  1558.         III_dequantize_sample_ms(hybridIn,scalefacs[1],gr_info,sfreq,part2_start); 
  1559.       else
  1560.         III_dequantize_sample(hybridIn[1],scalefacs[1],gr_info,sfreq,part2_start); 
  1561.  
  1562.       if(i_stereo)
  1563.         III_stereo(hybridIn,scalefacs[1],gr_info,sfreq,ms_stereo);
  1564.  
  1565.       if(single >= 0) {
  1566.         if(single == 3) {
  1567.           register int i;
  1568.           register real *in0 = (real *) hybridIn[0],*in1 = (real *) hybridIn[1];
  1569.           for(i=0;i<SSLIMIT*SBLIMIT;i++,in0++)
  1570.               *in0 = (*in0 + *in1++) * 0.5;
  1571.         }
  1572.         if(single == 1) {
  1573.           register int i;
  1574.           register real *in0 = (real *) hybridIn[0],*in1 = (real *) hybridIn[1];
  1575.           for(i=0;i<SSLIMIT*SBLIMIT;i++)
  1576.             *in0++ = *in1++;
  1577.         }
  1578.       }
  1579.     }
  1580.  
  1581.     for(ch=0;ch<stereo1;ch++) {
  1582.       struct gr_info_s *gr_info = &(sideinfo.ch[ch].gr[gr]);
  1583.       III_antialias(hybridIn[ch],gr_info);
  1584.       III_hybrid(hybridIn[ch], hybridOut[ch], ch,gr_info);
  1585.     }
  1586.  
  1587.     for(ss=0;ss<SSLIMIT;ss++) {
  1588.       if(single >= 0) {
  1589.         int i;
  1590.         short *pcm = pcm_sample+pcm_point;
  1591.         clip += (fr->synth)(hybridOut[0][ss],0,pcm);
  1592.         for(i=0;i<32;i++,pcm+=2)
  1593.           pcm[1] = pcm[0];
  1594.       }
  1595.       else {
  1596.         clip += (fr->synth)(hybridOut[0][ss],0,pcm_sample+pcm_point);
  1597.         clip += (fr->synth)(hybridOut[1][ss],1,pcm_sample+pcm_point);
  1598.       }
  1599.       pcm_point += 64;
  1600.  
  1601. #ifdef VARMODESUPPORT
  1602.       if (playlimit < 128) {
  1603.         pcm_point -= playlimit >> 1;
  1604.         playlimit = 0;
  1605.       }
  1606.       else
  1607.         playlimit -= 128;
  1608. #endif
  1609.  
  1610.       if(pcm_point >= audiobufsize)
  1611.         audio_flush(outmode,ai);
  1612.     }
  1613.   }
  1614.   
  1615.   return clip;
  1616. }
  1617.  
  1618.  
  1619.