home *** CD-ROM | disk | FTP | other *** search
/ AmigActive 4 / AACD04.ISO / AACD / Sound / LAME / Source / mdct.c < prev    next >
C/C++ Source or Header  |  1999-06-01  |  19KB  |  509 lines
  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: mdct.c,v 1.1 1996/02/14 04:04:23 rowlands Exp $
  6.  *
  7.  * $Log: mdct.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. #include "common.h"
  15. #include "l3side.h"
  16. #include "mdct.h"
  17.  
  18. double ca[8], cs[8];
  19.  
  20. /*
  21.   This is table B.9: coefficients for aliasing reduction
  22. */
  23. static double c[8] = { -0.6,-0.535,-0.33,-0.185,-0.095,-0.041,-0.0142, -0.0037 };
  24.  
  25. void mdct_sub( L3SBS (*sb_sample), double (*mdct_freq)[2][576], int stereo, III_side_info_t *l3_side, int mode_gr )
  26. {
  27.     gr_info *cod_info;
  28.     double mdct_in[36];
  29.     int ch,gr,band,k,j;
  30.     double bu,bd;
  31.     static int init = 0;
  32.     int    block_type;
  33.     double (*mdct_enc)[2][32][18] = (double (*)[2][32][18]) mdct_freq;
  34.     
  35.     if ( init == 0 )
  36.     {
  37.     /* prepare the aliasing reduction butterflies */
  38.     for ( k = 0; k < 8; k++ )
  39.     {
  40.         double sq;
  41.         sq = sqrt( 1.0 + c[k] * c[k] );
  42.         ca[k] = c[k] / sq;
  43.         cs[k] = 1.0 / sq;
  44.     }
  45.     init++;
  46.     }
  47.     
  48.     for ( gr = 0; gr < mode_gr; gr++ )
  49.     for ( ch = 0; ch < stereo; ch++ )
  50.     {
  51.         cod_info = (gr_info *) &(l3_side->gr[gr].ch[ch]) ;
  52.         block_type = cod_info->block_type;
  53.         
  54.         /*
  55.           Compensate for inversion in the analysis filter
  56.         */
  57.         for ( band = 1; band < 32; band+=2 )
  58.         for ( k = 1; k < 18; k+=2 )
  59.             (*sb_sample)[ch][gr+1][k][band] *= -1.0;
  60.         
  61.         /*
  62.           Perform imdct of 18 previous subband samples
  63.           + 18 current subband samples
  64.         */
  65.         for ( band = 0; band < 32; band++ )
  66.         {
  67.         for ( k = 0; k < 18; k++ )
  68.         {
  69.             mdct_in[k]    = (*sb_sample)[ch][ gr ][k][band];
  70.             mdct_in[k+18] = (*sb_sample)[ch][gr+1][k][band];
  71.         }
  72. #ifdef ALLOW_MIXED
  73.         if ( cod_info->mixed_block_flag && (band < 2) )
  74.             block_type = 0;
  75. #endif
  76.         mdct( mdct_in, &mdct_enc[gr][ch][band][0], block_type );
  77.         }
  78.         
  79.         /*
  80.           Perform aliasing reduction butterfly
  81.           on long blocks
  82.         */
  83.         if ( block_type != 2 )
  84.         for ( band = 0; band < 31; band++ )
  85.             for ( k = 0; k < 8; k++ )
  86.             {
  87.             bu = mdct_enc[gr][ch][band][17-k] * cs[k] + mdct_enc[gr][ch][band+1][k] * ca[k];
  88.             bd = mdct_enc[gr][ch][band+1][k] * cs[k] - mdct_enc[gr][ch][band][17-k] * ca[k];
  89.             mdct_enc[gr][ch][band][17-k] = bu;
  90.             mdct_enc[gr][ch][band+1][k]  = bd;
  91.             }
  92.         
  93.     }
  94.     
  95.     /*
  96.       Save latest granule's subband samples to be used in
  97.       the next mdct call
  98.     */
  99.     for ( ch = 0; ch < stereo; ch++ )
  100.     for ( j = 0; j < 18; j++ )
  101.         for ( band = 0; band < 32; band++ )
  102.         (*sb_sample)[ch][0][j][band] = (*sb_sample)[ch][mode_gr][j][band];
  103. }
  104.  
  105. void mdct( double *in, double *out, int block_type )
  106. {
  107. /*-------------------------------------------------------------------*/
  108. /*                                                                   */
  109. /*   Function: Calculation of the MDCT                               */
  110. /*   In the case of long blocks ( block_type 0,1,3 ) there are       */
  111. /*   36 coefficents in the time domain and 18 in the frequency       */
  112. /*   domain.                                                         */
  113. /*   In the case of short blocks (block_type 2 ) there are 3         */
  114. /*   transformations with short length. This leads to 12 coefficents */
  115. /*   in the time and 6 in the frequency domain. In this case the     */
  116. /*   results are stored side by side in the vector out[].            */
  117. /*                                                                   */
  118. /*   New layer3                                                      */
  119. /*                                                                   */
  120. /*-------------------------------------------------------------------*/
  121.  
  122.   int l,k,i,m,N;
  123.   double sum;
  124.   static double win[4][36];
  125.   static int init = 0;
  126.   static double cos_s[6][12], cos_l[18][18];
  127.   static double fin[18];
  128.   static double *cos_l0;
  129.  
  130.   if ( init == 0 )
  131.   {
  132.     /* type 0 */
  133.     for ( i = 0; i < 36; i++ )
  134.       win[0][i] = sin( PI/36 * (i + 0.5) );
  135.     /* type 1*/
  136.     for ( i = 0; i < 18; i++ ) 
  137.       win[1][i] = sin( PI/36 * (i + 0.5) );
  138.     for ( i = 18; i < 24; i++ )
  139.       win[1][i] = 1.0;
  140.     for ( i = 24; i < 30; i++ )
  141.       win[1][i] = sin( PI/12 * ( i + 0.5 - 18) );
  142.     for ( i = 30; i < 36; i++ )
  143.       win[1][i] = 0.0;
  144.     /* type 3*/
  145.     for ( i = 0; i < 6; i++ )
  146.       win[3][i] = 0.0;
  147.     for ( i = 6; i < 12; i++ ) 
  148.       win[3][i] = sin( PI/12 * (i + 0.5 - 6) );
  149.     for ( i = 12; i < 18; i++ )
  150.       win[3][i] = 1.0;
  151.     for ( i = 18; i < 36; i++ )
  152.       win[3][i] = sin( PI/36 * (i + 0.5) );
  153.  
  154.     N = 12;
  155.     for ( m = 0; m < N / 2; m++ )
  156.       for ( k = 0; k < N; k++ )
  157.         cos_s[m][k] = cos( (PI /(2 * N)) * (2 * k + 1 + N / 2) *
  158.                      (2 * m + 1) ) / (N / 4) * sin( PI/12 * (k + 0.5) );
  159.  
  160.     N = 36;
  161.     for ( m = 0; m < N / 2; m++ )
  162.     {
  163.       for ( k = 0; k < 9; k++ )
  164.       {
  165.         cos_l[m][k] = cos( (PI / (2 * N)) * (2 * k + 1 + N / 2) *
  166.                      (2 * m + 1) ) / (N / 4);
  167.         cos_l[m][9+k] = cos( (PI / (2 * N)) * (2 * (18 + k) + 1 + N / 2) *
  168.                      (2 * m + 1) ) / (N / 4);
  169.       }
  170.     }
  171.     init++;
  172.   }
  173.  
  174.   if ( block_type == 2)
  175.   {
  176.     N = 12;
  177.     for ( l = 0; l < 3; l++ )
  178.     {
  179.       for ( m = 0; m < N / 2; m++ )
  180.       {
  181.         for ( sum = 0.0, k = 0; k < N; k++ )
  182.           sum += in[k + 6 * l + 6] * cos_s[m][k];
  183.         out[ 3 * m + l] = sum;
  184.       }
  185.     }
  186.   }
  187.   else 
  188.   {
  189.       for (k=0;k<9;k++)
  190.       {
  191.      fin[k] = win[block_type][k] * in[k] - win[block_type][17-k] * in[17-k];
  192.      fin[9+k] = win[block_type][18+k] * in[18+k] + win[block_type][35-k] * in[35-k];
  193.       }
  194.  
  195.       /* 0 */  
  196.       cos_l0= &cos_l[0][0];
  197.       sum = ( fin[0] ) * *cos_l0; /* 17 */
  198.       sum += ( fin[1] ) * *++cos_l0; /* 15 */
  199.       sum += ( fin[2] ) * *++cos_l0; /* 13 */
  200.       sum += ( fin[3] ) * *++cos_l0; /* 11 */
  201.       sum += ( fin[4] ) * *++cos_l0; /* 9 */
  202.       sum += ( fin[5] ) * *++cos_l0; /* 7  */
  203.       sum += ( fin[6] ) * *++cos_l0; /* 5  */
  204.       sum += ( fin[7] ) * *++cos_l0; /* 3 */
  205.       sum += ( fin[8] ) * *++cos_l0; /* 1  */
  206.       sum += ( fin[9] ) * *++cos_l0; /* 19*/
  207.       sum += ( fin[10] ) * *++cos_l0; /* 21 */
  208.       sum += ( fin[11] ) * *++cos_l0; /* 23 */
  209.       sum += ( fin[12] ) * *++cos_l0; /*25 */
  210.       sum += ( fin[13] ) * *++cos_l0; /* 27 */
  211.       sum += ( fin[14] ) * *++cos_l0; /* 29 */
  212.       sum += ( fin[15] ) * *++cos_l0; /* 31*/
  213.       sum += ( fin[16] ) * *++cos_l0; /* 33 */
  214.       sum += ( fin[17] ) * *++cos_l0; /* 35 */    
  215.       out[0]=sum; 
  216.  
  217.  
  218.       /* 1 */
  219.       sum = ( fin[2]+fin[3]+fin[17] ) * cos_l[1][2];  /* mfc=3 2*/
  220.       sum += ( fin[1]+fin[4]+fin[16] ) * cos_l[1][1];  /* mfc=9 1*/
  221.       sum += ( fin[0]+fin[5]+fin[15] ) * *++cos_l0;  /* mfc=15 0*/
  222.       sum += ( fin[6]-fin[9]+fin[14] ) * *(cos_l0+=6); /* mfc=21 6*/
  223.       sum += ( fin[7]-fin[10]+fin[13] ) * *++cos_l0; /* mfc=27 7*/
  224.       sum += ( fin[8]-fin[11]+fin[12] ) * *++cos_l0; /* mfc = 28 8*/
  225.       out[1]=sum;
  226.  
  227.       /* 2 */
  228.       cos_l0= &cos_l[2][0];
  229.       sum = ( fin[0] ) * *cos_l0;  /* mfc 23 */
  230.       sum += ( fin[1] ) * *++cos_l0;  /* mfc 33 */
  231.       sum += ( fin[2] ) * *++cos_l0;   /* mfc 29 */
  232.       sum += ( fin[3] ) * *++cos_l0;   /* mfc 19 */
  233.       sum += ( fin[4] ) * *++cos_l0;   /* mfc 9  */
  234.       sum += ( fin[5] ) * *++cos_l0;   /* mfc 1  */
  235.       sum += ( fin[6] ) * *++cos_l0;   /* mfc 11 */
  236.       sum += ( fin[7] ) * *++cos_l0;   /* mfc 21 */
  237.       sum += ( fin[8] ) * *++cos_l0;    /* mfc 31 */
  238.       sum += ( fin[9] ) * *++cos_l0; /* mfc 13 */
  239.       sum += ( fin[10] ) * *++cos_l0; /* mfc 3  */
  240.       sum += ( fin[11] ) * *++cos_l0; /* mfc 7  */
  241.       sum += ( fin[12] ) * *++cos_l0; /* mfc 17 */
  242.       sum += ( fin[13] ) * *++cos_l0; /* mfc 27 */
  243.       sum += ( fin[14] ) * *++cos_l0; /* mfc 35 */
  244.       sum += ( fin[15] ) * *++cos_l0; /* mfc 25 */
  245.       sum += ( fin[16] ) * *++cos_l0; /* mfc 15 */
  246.       sum += ( fin[17] ) * *++cos_l0; /* mfc 5 */
  247.       out[2]=sum;
  248.  
  249.       /* 3 */;
  250.       sum = ( fin[0] ) * *++cos_l0; /* mfc 11 */
  251.       sum += ( fin[1] ) * *++cos_l0; /* mfc 3 */
  252.       sum += ( fin[2] ) * *++cos_l0; /* 17 */
  253.       sum += ( fin[3] ) * *++cos_l0; /* 31 */
  254.       sum += ( fin[4] ) * *++cos_l0; /* 27 */
  255.       sum += ( fin[5] ) * *++cos_l0; /* 13 */
  256.       sum += ( fin[6] ) * *++cos_l0; /* 1 */
  257.       sum += ( fin[7] ) * *++cos_l0; /* 15 */
  258.       sum += ( fin[8] ) * *++cos_l0; /* 29 */
  259.       sum += ( fin[9] ) * *++cos_l0; /* 25 */
  260.       sum += ( fin[10] ) * *++cos_l0; /* 33 */
  261.       sum += ( fin[11] ) * *++cos_l0; /* 19 */
  262.       sum += ( fin[12] ) * *++cos_l0; /* 5 */
  263.       sum += ( fin[13] ) * *++cos_l0; /* 9 */
  264.       sum += ( fin[14] ) * *++cos_l0; /* 23 */
  265.       sum += ( fin[15] ) * *++cos_l0; /* 35 */
  266.       sum += ( fin[16] ) * *++cos_l0; /* 21 */
  267.       sum += ( fin[17] ) * *++cos_l0; /* 7 */
  268.       out[3]=sum;
  269.  
  270.       /* 4 */
  271.       /* mfc 27 */
  272.       sum = ( fin[0]-fin[1]-fin[4]+fin[5]+fin[8] \
  273.            -fin[11]+fin[12]+fin[15]-fin[16] ) * *++cos_l0;
  274.       /* mfc 9 */
  275.       sum += ( -fin[2]-fin[3]+fin[6]+fin[7]-fin[9]-fin[10] \
  276.            +fin[13]+fin[14]-fin[17] ) * cos_l[4][6];
  277.       out[4]=sum;
  278.  
  279.       /* 5 */
  280.       cos_l0= &cos_l[5][0];
  281.       sum = ( fin[0] ) * *cos_l0; /* 7 */
  282.       sum += ( fin[1] ) * *++cos_l0; /* 15 */
  283.       sum += ( fin[2] ) * *++cos_l0; /* 35 */
  284.       sum += ( fin[3] ) * *++cos_l0; /* 13 */
  285.       sum += ( fin[4] ) * *++cos_l0; /* 9 */
  286.       sum += ( fin[5] ) * *++cos_l0; /* 31 */
  287.       sum += ( fin[6] ) * *++cos_l0; /* 19 */
  288.       sum += ( fin[7] ) * *++cos_l0; /* 3 */
  289.       sum += ( fin[8] ) * *++cos_l0; /* 8 */
  290.       sum += ( fin[9] ) * *++cos_l0; /* 29 */
  291.       sum += ( fin[10] ) * *++cos_l0; /* 21 */
  292.       sum += ( fin[11] ) * *++cos_l0; /* 1 */
  293.       sum += ( fin[12] ) * *++cos_l0; /* 23 */
  294.       sum += ( fin[13] ) * *++cos_l0; /* 27 */
  295.       sum += ( fin[14] ) * *++cos_l0; /* 5 */
  296.       sum += ( fin[15] ) * *++cos_l0; /* 17 */
  297.       sum += ( fin[16] ) * *++cos_l0; /* 33 */
  298.       sum += ( fin[17] ) * *++cos_l0; /* 11 */
  299.       out[5]=sum;
  300.  
  301.       /* 6 */
  302.       sum = ( fin[0] ) * *++cos_l0; /* 31 */
  303.       sum += ( fin[1] ) * *++cos_l0; /* 15 */
  304.       sum += ( fin[2] ) * *++cos_l0; /* 11 */
  305.       sum += ( fin[3] ) * *++cos_l0; /* 35 */
  306.       sum += ( fin[4] ) * *++cos_l0; /* 9 */
  307.       sum += ( fin[5] ) * *++cos_l0; /* 17 */
  308.       sum += ( fin[6] ) * *++cos_l0; /* 29 */
  309.       sum += ( fin[7] ) * *++cos_l0; /* 3 */
  310.       sum += ( fin[8] ) * *++cos_l0; /* 23 */
  311.       sum += ( fin[9] ) * *++cos_l0; /* 5 */
  312.       sum += ( fin[10] ) * *++cos_l0; /* 21 */
  313.       sum += ( fin[11] ) * *++cos_l0; /* 25 */
  314.       sum += ( fin[12] ) * *++cos_l0; /* 1 */
  315.       sum += ( fin[13] ) * *++cos_l0; /* 27 */
  316.       sum += ( fin[14] ) * *++cos_l0; /* 19 */
  317.       sum += ( fin[15] ) * *++cos_l0; /* 7 */
  318.       sum += ( fin[16] ) * *++cos_l0; /* 33 */
  319.       sum += ( fin[17] ) * *++cos_l0; /* 13 */    
  320.       out[6]=sum;
  321.  
  322.       /* 7 */
  323.       sum = ( fin[0]+fin[5]+fin[15] ) * *++cos_l0; /* 3 */
  324.       sum += ( fin[1]+fin[4]+fin[16] ) * *++cos_l0; /* 27*/
  325.       sum += ( fin[2]+fin[3]+fin[17] ) * *++cos_l0; /* 15 */
  326.       sum += ( fin[6]-fin[9]+fin[14] ) * *(cos_l0+=4); /* 33 */
  327.       sum += ( fin[7]-fin[10]+fin[13] ) * *++cos_l0; /* 9 */
  328.       sum += ( fin[8]-fin[11]+fin[12] ) * *++cos_l0; /* 21 */
  329.       out[7]=sum;
  330.  
  331.       /* 8 */
  332.       cos_l0= &cos_l[8][0];
  333.       sum = (  fin[0] ) * *cos_l0; /* 35 */
  334.       sum += ( fin[1] ) * *++cos_l0; /* 3  */
  335.       sum += ( fin[2] ) * *++cos_l0; /* 31 */
  336.       sum += ( fin[3] ) * *++cos_l0; /*  7 */
  337.       sum += ( fin[4] ) * *++cos_l0; /* 27*/
  338.       sum += ( fin[5] ) * *++cos_l0; /* 11 */
  339.       sum += ( fin[6] ) * *++cos_l0; /* 23 */
  340.       sum += ( fin[7] ) * *++cos_l0; /* 15*/
  341.       sum += ( fin[8] ) * *++cos_l0; /* 19 */
  342.       sum += ( fin[9] ) * *++cos_l0; /* 1 */
  343.       sum += ( fin[10] ) * *++cos_l0; /* 33 */
  344.       sum += ( fin[11] ) * *++cos_l0; /*  5 */
  345.       sum += ( fin[12] ) * *++cos_l0; /* 29*/
  346.       sum += ( fin[13] ) * *++cos_l0; /*  9 */
  347.       sum += ( fin[14] ) * *++cos_l0; /* 25 */
  348.       sum += ( fin[15] ) * *++cos_l0; /* 13*/
  349.       sum += ( fin[16] ) * *++cos_l0; /* 21 */
  350.       sum += ( fin[17] ) * *++cos_l0; /* 17 */    
  351.       out[8]=sum;
  352.  
  353.       /* 9 */
  354.       sum = ( fin[0] ) * *++cos_l0; /* 1  */
  355.       sum += ( fin[1] ) * *++cos_l0; /* 33 */
  356.       sum += ( fin[2] ) * *++cos_l0; /*  5 */
  357.       sum += ( fin[3] ) * *++cos_l0; /*  29*/
  358.       sum += ( fin[4] ) * *++cos_l0; /* 9 */
  359.       sum += ( fin[5] ) * *++cos_l0; /* 25 */
  360.       sum += ( fin[6] ) * *++cos_l0; /* 13 */
  361.       sum += ( fin[7] ) * *++cos_l0; /* 21*/
  362.       sum += ( fin[8] ) * *++cos_l0; /* 17 */
  363.       sum += ( fin[9] ) * *++cos_l0; /* 35*/
  364.       sum += ( fin[10] ) * *++cos_l0; /* 3  */
  365.       sum += ( fin[11] ) * *++cos_l0; /* 31 */
  366.       sum += ( fin[12] ) * *++cos_l0; /* 7 */
  367.       sum += ( fin[13] ) * *++cos_l0; /* 27 */
  368.       sum += ( fin[14] ) * *++cos_l0; /* 11 */
  369.       sum += ( fin[15] ) * *++cos_l0; /* 23*/
  370.       sum += ( fin[16] ) * *++cos_l0; /* 15 */
  371.       sum += ( fin[17] ) * *++cos_l0; /* 19 */    
  372.       out[9]=sum; 
  373.  
  374.       /* 10 */
  375.       sum = ( fin[0]+fin[5]+fin[15] ) * *++cos_l0; /* 33 */
  376.       sum += ( fin[1]+fin[4]+fin[16] ) * *++cos_l0; /* 9 */
  377.       sum += ( fin[2]+fin[3]+fin[17] ) * *++cos_l0; /* 21 */
  378.       sum += ( fin[6]-fin[9]+fin[14] ) * *(cos_l0+=4); /* 3 */
  379.       sum += ( fin[7]-fin[10]+fin[13] ) * *++cos_l0; /* 27 */
  380.       sum += ( fin[8]-fin[11]+fin[12] ) * *++cos_l0; /* 15 */
  381.       out[10]=sum;
  382.  
  383.       /* 11 */
  384.       cos_l0= &cos_l[11][0];
  385.       sum = ( fin[0] ) * *cos_l0; /* 5  */
  386.       sum += ( fin[1] ) * *++cos_l0; /* 21 */
  387.       sum += ( fin[2] ) * *++cos_l0; /* 25 */
  388.       sum += ( fin[3] ) * *++cos_l0; /*  1 */
  389.       sum += ( fin[4] ) * *++cos_l0; /*27 */
  390.       sum += ( fin[5] ) * *++cos_l0; /* 19 */
  391.       sum += ( fin[6] ) * *++cos_l0; /* 7  */
  392.       sum += ( fin[7] ) * *++cos_l0; /* 33*/
  393.       sum += ( fin[8] ) * *++cos_l0; /* 13 */
  394.       sum += ( fin[9] ) * *++cos_l0; /* 31*/
  395.       sum += ( fin[10] ) * *++cos_l0; /* 15 */
  396.       sum += ( fin[11] ) * *++cos_l0; /* 11 */
  397.       sum += ( fin[12] ) * *++cos_l0; /*35 */
  398.       sum += ( fin[13] ) * *++cos_l0; /*  9 */
  399.       sum += ( fin[14] ) * *++cos_l0; /* 17 */
  400.       sum += ( fin[15] ) * *++cos_l0; /* 29*/
  401.       sum += ( fin[16] ) * *++cos_l0; /* 3  */
  402.       sum += ( fin[17] ) * *++cos_l0; /* 23 */    
  403.       out[11]=sum; 
  404.  
  405.       /* 12 */
  406.       sum = ( fin[0] ) * *++cos_l0; /* 29 */
  407.       sum += ( fin[1] ) * *++cos_l0; /* 21 */
  408.       sum += ( fin[2] ) * *++cos_l0; /* 1  */
  409.       sum += ( fin[3] ) * *++cos_l0; /* 23 */
  410.       sum += ( fin[4] ) * *++cos_l0; /*27 */
  411.       sum += ( fin[5] ) * *++cos_l0; /* 5  */
  412.       sum += ( fin[6] ) * *++cos_l0; /* 17 */
  413.       sum += ( fin[7] ) * *++cos_l0; /* 33*/
  414.       sum += ( fin[8] ) * *++cos_l0; /* 11 */
  415.       sum += ( fin[9] ) * *++cos_l0; /* 7 */
  416.       sum += ( fin[10] ) * *++cos_l0; /* 15 */
  417.       sum += ( fin[11] ) * *++cos_l0; /* 35 */
  418.       sum += ( fin[12] ) * *++cos_l0; /*13 */
  419.       sum += ( fin[13] ) * *++cos_l0; /*  9 */
  420.       sum += ( fin[14] ) * *++cos_l0; /* 31 */
  421.       sum += ( fin[15] ) * *++cos_l0; /* 19*/
  422.       sum += ( fin[16] ) * *++cos_l0; /* 3  */
  423.       sum += ( fin[17] ) * *++cos_l0; /* 25 */    
  424.       out[12]=sum; 
  425.  
  426.       /* 13 */
  427.       sum = ( -fin[0]+fin[1]+fin[4]-fin[5]-fin[8]+fin[11]-fin[12] \
  428.           -fin[15]+fin[16] ) * -*++cos_l0; /* 9 */
  429.       sum += ( -fin[2]-fin[3]+fin[6]+fin[7]-fin[9]-fin[10]+fin[13] \
  430.            +fin[14]-fin[17] ) * -cos_l[13][2]; /* 27 */
  431.       out[13]=sum;
  432.  
  433.       /* 14 */
  434.       cos_l0= &cos_l[14][0];
  435.       sum = ( fin[0] ) * *cos_l0; /* 25 */
  436.       sum += ( fin[1] ) * *++cos_l0; /* 33 */
  437.       sum += ( fin[2] ) * *++cos_l0; /* 19 */
  438.       sum += ( fin[3] ) * *++cos_l0; /* 5  */
  439.       sum += ( fin[4] ) * *++cos_l0; /* 9 */
  440.       sum += ( fin[5] ) * *++cos_l0; /* 23 */
  441.       sum += ( fin[6] ) * *++cos_l0; /* 35 */
  442.       sum += ( fin[7] ) * *++cos_l0; /* 21*/
  443.       sum += ( fin[8] ) * *++cos_l0; /* 7  */
  444.       sum += ( fin[9] ) * *++cos_l0; /* 11*/
  445.       sum += ( fin[10] ) * *++cos_l0; /* 3  */
  446.       sum += ( fin[11] ) * *++cos_l0; /* 17 */
  447.       sum += ( fin[12] ) * *++cos_l0; /*31 */
  448.       sum += ( fin[13] ) * *++cos_l0; /* 27 */
  449.       sum += ( fin[14] ) * *++cos_l0; /* 13 */
  450.       sum += ( fin[15] ) * *++cos_l0; /* 1 */
  451.       sum += ( fin[16] ) * *++cos_l0; /* 15 */
  452.       sum += ( fin[17] ) * *++cos_l0; /* 29 */    
  453.       out[14]=sum; 
  454.  
  455.       /* 15 */
  456.       sum = ( fin[0] ) * *++cos_l0; /* 13 */
  457.       sum += ( fin[1] ) * *++cos_l0; /* 3  */
  458.       sum += ( fin[2] ) * *++cos_l0; /* 7  */
  459.       sum += ( fin[3] ) * *++cos_l0; /* 17 */
  460.       sum += ( fin[4] ) * *++cos_l0; /* 27*/
  461.       sum += ( fin[5] ) * *++cos_l0; /* 35 */
  462.       sum += ( fin[6] ) * *++cos_l0; /* 25 */
  463.       sum += ( fin[7] ) * *++cos_l0; /* 15*/
  464.       sum += ( fin[8] ) * *++cos_l0; /* 5  */
  465.       sum += ( fin[9] ) * *++cos_l0; /* 23*/
  466.       sum += ( fin[10] ) * *++cos_l0; /* 33 */
  467.       sum += ( fin[11] ) * *++cos_l0; /* 29 */
  468.       sum += ( fin[12] ) * *++cos_l0; /*19 */
  469.       sum += ( fin[13] ) * *++cos_l0; /* 9  */
  470.       sum += ( fin[14] ) * *++cos_l0; /* 1  */
  471.       sum += ( fin[15] ) * *++cos_l0; /* 11*/
  472.       sum += ( fin[16] ) * *++cos_l0; /* 21 */
  473.       sum += ( fin[17] ) * *++cos_l0; /* 31 */    
  474.       out[15]=sum; 
  475.  
  476.       /* 16 */
  477.       sum = ( fin[0]+fin[5]+fin[15] ) * *++cos_l0; /* 21 */
  478.       sum += ( fin[1]+fin[4]+fin[16] ) * *++cos_l0;  /* 27 */
  479.       sum += ( fin[2]+fin[3]+fin[17] ) * *++cos_l0; /* 33 */
  480.       sum += ( fin[6]-fin[9]+fin[14] ) * *(cos_l0+=4); /* 15 */
  481.       sum += ( fin[7]-fin[10]+fin[13] ) * *++cos_l0; /* 9 */
  482.       sum += ( fin[8]-fin[11]+fin[12] ) * *++cos_l0; /* 3 */
  483.       out[16]=sum;
  484.  
  485.       /* 17 */
  486.       cos_l0= &cos_l[17][0];
  487.       sum = ( fin[0] ) * *cos_l0; /* 17 */
  488.       sum += ( fin[1] ) * *++cos_l0; /* 15 */
  489.       sum += ( fin[2] ) * *++cos_l0; /* 13 */
  490.       sum += ( fin[3] ) * *++cos_l0; /* 11 */
  491.       sum += ( fin[4] ) * *++cos_l0; /* 9 */
  492.       sum += ( fin[5] ) * *++cos_l0; /* 7  */
  493.       sum += ( fin[6] ) * *++cos_l0; /* 5  */
  494.       sum += ( fin[7] ) * *++cos_l0; /* 3 */
  495.       sum += ( fin[8] ) * *++cos_l0; /* 1  */
  496.       sum += ( fin[9] ) * *++cos_l0; /* 19*/
  497.       sum += ( fin[10] ) * *++cos_l0; /* 21 */
  498.       sum += ( fin[11] ) * *++cos_l0; /* 23 */
  499.       sum += ( fin[12] ) * *++cos_l0; /*25 */
  500.       sum += ( fin[13] ) * *++cos_l0; /* 27 */
  501.       sum += ( fin[14] ) * *++cos_l0; /* 29 */
  502.       sum += ( fin[15] ) * *++cos_l0; /* 31*/
  503.       sum += ( fin[16] ) * *++cos_l0; /* 33 */
  504.       sum += ( fin[17] ) * *++cos_l0; /* 35 */    
  505.       out[17]=sum; 
  506.    }
  507. }
  508.  
  509.