# （十三）用JAVA编写MP3解码器——IMDCT快速算法

1.初始化加窗系数floatWinIMDCT[][] 初始化floatWinIMDCT[][]时直接代入相关公式，我这事先另外编程用公式把加窗系数计算出来，供程序中直接查表使用，而不是运行时进行初始化。初始化floatWinIMDCT[][]的代码正好还留着，就顺便贴出来：

`public class InitIMDCT {public static void main(String args[]) {int i, j;final double PI12 = 0.26179938779914943653855361527329;final double PI36 = 0.087266462599716478846184538424431;final double PI72 = 0.043633231299858239423092269212215;double[][] doubleWinIMDCT = new double[4][36];// type 0for (i = 0; i < 36; i++)doubleWinIMDCT[0][i] = Math.sin(PI36 * (i + 0.5));// type 1for (i = 0; i < 18; i++)doubleWinIMDCT[1][i] = Math.sin(PI36 * (i + 0.5));for (i = 18; i < 24; i++)doubleWinIMDCT[1][i] = 1.0;for (i = 24; i < 30; i++)doubleWinIMDCT[1][i] = Math.sin(PI12 * (i + 0.5 - 18));for (i = 30; i < 36; i++)doubleWinIMDCT[1][i] = 0.0;// type 2 (not needed...)// type 3for (i = 0; i < 6; i++)doubleWinIMDCT[3][i] = 0.0;for (i = 6; i < 12; i++)doubleWinIMDCT[3][i] = Math.sin(PI12 * (i + 0.5 - 6.0));for (i = 12; i < 18; i++)doubleWinIMDCT[3][i] = 1.0;for (i = 18; i < 36; i++)doubleWinIMDCT[3][i] = Math.sin(PI36 * (i + 0.5));for (i = 0; i < 4; i++) {doubleWinIMDCT[i][0] /= 2 * Math.cos(19 * PI72);doubleWinIMDCT[i][1] /= 2 * Math.cos(21 * PI72);doubleWinIMDCT[i][2] /= 2 * Math.cos(23 * PI72);doubleWinIMDCT[i][3] /= 2 * Math.cos(25 * PI72);doubleWinIMDCT[i][4] /= 2 * Math.cos(27 * PI72);doubleWinIMDCT[i][5] /= 2 * Math.cos(29 * PI72);doubleWinIMDCT[i][6] /= 2 * Math.cos(31 * PI72);doubleWinIMDCT[i][7] /= 2 * Math.cos(33 * PI72);doubleWinIMDCT[i][8] /= 2 * Math.cos(35 * PI72);doubleWinIMDCT[i][9] /= 2 * Math.cos(37 * PI72);doubleWinIMDCT[i][10] /= 2 * Math.cos(39 * PI72);doubleWinIMDCT[i][11] /= 2 * Math.cos(41 * PI72);doubleWinIMDCT[i][12] /= 2 * Math.cos(43 * PI72);doubleWinIMDCT[i][13] /= 2 * Math.cos(45 * PI72);doubleWinIMDCT[i][14] /= 2 * Math.cos(47 * PI72);doubleWinIMDCT[i][15] /= 2 * Math.cos(49 * PI72);doubleWinIMDCT[i][16] /= 2 * Math.cos(51 * PI72);doubleWinIMDCT[i][17] /= 2 * Math.cos(53 * PI72);doubleWinIMDCT[i][18] /= 2 * Math.cos(55 * PI72);doubleWinIMDCT[i][19] /= 2 * Math.cos(57 * PI72);doubleWinIMDCT[i][20] /= 2 * Math.cos(59 * PI72);doubleWinIMDCT[i][21] /= 2 * Math.cos(61 * PI72);doubleWinIMDCT[i][22] /= 2 * Math.cos(63 * PI72);doubleWinIMDCT[i][23] /= 2 * Math.cos(65 * PI72);doubleWinIMDCT[i][24] /= 2 * Math.cos(67 * PI72);doubleWinIMDCT[i][25] /= 2 * Math.cos(69 * PI72);doubleWinIMDCT[i][26] /= 2 * Math.cos(71 * PI72);doubleWinIMDCT[i][27] /= 2 * Math.cos(71 * PI72);doubleWinIMDCT[i][28] /= 2 * Math.cos(69 * PI72);doubleWinIMDCT[i][29] /= 2 * Math.cos(67 * PI72);doubleWinIMDCT[i][30] /= 2 * Math.cos(65 * PI72);doubleWinIMDCT[i][31] /= 2 * Math.cos(63 * PI72);doubleWinIMDCT[i][32] /= 2 * Math.cos(61 * PI72);doubleWinIMDCT[i][33] /= 2 * Math.cos(59 * PI72);doubleWinIMDCT[i][34] /= 2 * Math.cos(57 * PI72);doubleWinIMDCT[i][35] /= 2 * Math.cos(55 * PI72);}System.out.println("private static float[][] doubleWinIMDCT = {");for(i = 0; i < 4; i++) {System.out.printf("{");for(j = 0; j < 35; j++) {System.out.printf("%gf,", doubleWinIMDCT[i][j]);if(j % 6 == 5)System.out.printf("/n");}System.out.printf("%gf},/n", doubleWinIMDCT[i][j]);}System.out.println("};");}}`

2.IMDCT 这里用到两种反向改进型离散余弦变换（IMDCT）：短块用12点的IMDCT；长块用36点IMDCT。在解码端应用IMDCT完成频域到时域的转换，IMDCT计算的公式是：

3.混合滤波带（hybrid filterbank） 先定义一个先进先出的队列floatPrevBlck[][][]，经IMDCT后产生的36个输出中的前18个值和这个队列中的18个值相加输出到xr[]，再将IMDCT产生的36个输出值中的后18个写入队列；每解码完一个粒度组内的一个声道后要将当前的最大子带之后对应的队列数据清零。

【提示】以下代码是Layer3.java的一部分，应遵守《（一）用JAVA编写MP3解码器——前言》中的许可协议。

class Layer3内的hybrid方法的源码如下：

`//7.//>>>>HYBRID(synthesize via iMDCT)=========================================private static final float[][] floatWinIMDCT = {{0.0322824f,0.1072064f,0.2014143f,0.3256164f,0.5f,0.7677747f,1.2412229f,2.3319514f,7.7441506f,-8.4512568f,-3.0390580f,-1.9483297f,-1.4748814f,-1.2071068f,-1.0327232f,-0.9085211f,-0.8143131f,-0.7393892f,-0.6775254f,-0.6248445f,-0.5787917f,-0.5376016f,-0.5f,-0.4650284f,-0.4319343f,-0.4000996f,-0.3689899f,-0.3381170f,-0.3070072f,-0.2751725f,-0.2420785f,-0.2071068f,-0.1695052f,-0.1283151f,-0.0822624f,-0.0295815f},{0.0322824f,0.1072064f,0.2014143f,0.3256164f,0.5f,0.7677747f,1.2412229f,2.3319514f,7.7441506f,-8.4512568f,-3.0390580f,-1.9483297f,-1.4748814f,-1.2071068f,-1.0327232f,-0.9085211f,-0.8143131f,-0.7393892f,-0.6781709f,-0.6302362f,-0.5928445f,-0.5636910f,-0.5411961f,-0.5242646f,-0.5077583f,-0.4659258f,-0.3970546f,-0.3046707f,-0.1929928f,-0.0668476f,-0.0f,-0.0f,-0.0f,-0.0f,-0.0f,-0.0f},{/* block_type = 2 */},{0.0f,0.0f,0.0f,0.0f,0.0f,0.0f,0.3015303f,1.4659259f,6.9781060f,-9.0940447f,-3.5390582f,-2.2903500f,-1.6627548f,-1.3065630f,-1.0828403f,-0.9305795f,-0.8213398f,-0.7400936f,-0.6775254f,-0.6248445f,-0.5787917f,-0.5376016f,-0.5f,-0.4650284f,-0.4319343f,-0.4000996f,-0.3689899f,-0.3381170f,-0.3070072f,-0.2751725f,-0.2420785f,-0.2071068f,-0.1695052f,-0.1283151f,-0.0822624f,-0.0295815f} };/** inv_mdct -- 12点和36点IMDCT快速算法(展开式)*/private void inv_mdct(final float[] in, final float[] out, final int block_type) {int i, idx6 = 0;float in0,in1,in2,in3,in4,in5,in6,in7,in8,in9,in10,in11,in12,in13,in14,in15,in16,in17;float out0,out1,out2,out3,out4,out5,out6,out7,out8,out9,out10,out11,out12,out13,out14,out15,out16,out17, tmp;if(block_type == 2) {out[0]=out[1]=out[2]=out[3]=out[4]=out[5]=out[6]=out[7]=out[8]=out[9]=out[10]=out[11]=out[12]=out[13]=out[14]=out[15]=out[16]=out[17]=out[18]=out[19]=out[20]=out[21]=out[22]=out[23]=out[24]=out[25]=out[26]=out[27]=out[28]=out[29]=out[30]=out[31]=out[32]=out[33]=out[34]=out[35]=0.0f;for (i = 0; i < 3; i++) {//>>>>>>>>>>>> 12-point IMDCT//>>>>>> 6-point IDCTin[15 + i] += (in[12 + i] += in[9 + i]) + in[6 + i];in[9 + i] += (in[6 + i] += in[3 + i]) + in[i];in[3 + i] += in[i];//>>> 3-point IDCT on evenout1 = (in1 = in[i]) - (in2 = in[12+i]);in3 = in1 + in2 * 0.5f;in4 = in[6+i] * 0.8660254f;out0 = in3 + in4;out2 = in3 - in4;//<<< End 3-point IDCT on even//>>> 3-point IDCT on odd (for 6-point IDCT)out4 = ((in1 = in[3+i]) - (in2 = in[15+i])) * 0.7071068f;in3 = in1 + in2 * 0.5f;in4 = in[9+i] * 0.8660254f;out5 = (in3 + in4) * 0.5176381f;out3 = (in3 - in4) * 1.9318516f;//<<< End 3-point IDCT on odd// Output: butterflies on 2,3-point IDCT's (for 6-point IDCT)tmp = out0; out0 += out5; out5 = tmp - out5;tmp = out1; out1 += out4; out4 = tmp - out4;tmp = out2; out2 += out3; out3 = tmp - out3;//<<<<<< End 6-point IDCT//<<<<<<<<<<<< End 12-point IDCT// Shiftin1 = out3 * 0.1072064f;out[6 + idx6] += in1;out[7 + idx6] += out4 * 0.5f;out[8 + idx6] += out5 * 2.3319512f;out[9 + idx6] -= out5 * 3.0390580f;out[10 + idx6] -= out4 * 1.2071068f;out[11 + idx6] -= in1 * 7.5957541f;out[12 + idx6] -= out2 * 0.6248445f;out[13 + idx6] -= out1 * 0.5f;out[14 + idx6] -= out0 * 0.4000996f;out[15 + idx6] -= out0 * 0.3070072f;out[16 + idx6] -= out1 * 0.2071068f;out[17 + idx6] -= out2 * 0.0822623f;idx6 += 6;}} else {//>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 36-point IDCT//>>>>>>>>>>>>>>>>>> 18-point IDCT for oddin[17] += (in[16] += in[15]) + in[14];in[15] += (in[14] += in[13]) + in[12];in[13] += (in[12] += in[11]) + in[10];in[11] += (in[10] += in[9]) + in[8];in[9] += (in[8] += in[7]) + in[6];in[7] += (in[6] += in[5]) + in[4];in[5] += (in[4] += in[3]) + in[2];in[3] += (in[2] += in[1]) + in[0];in[1] += in[0];//>>>>>>>>> 9-point IDCT on even/* * for(i = 0; i < 9; i++) { * sum = 0; * for(j = 0; j < 18; j += 2) * sum += in[j] * cos(PI36 * (2 * i + 1) * j); * out18[i] = sum; * } */in0 = in[0] + in[12] * 0.5f;in1 = in[0] - in[12];in2 = in[8] + in[16] - in[4];out4 = in1 + in2;in3 = in1 - in2 * 0.5f;in4 = (in[10] + in[14] - in[2]) * 0.8660254f;// cos(PI/6)out1 = in3 - in4;out7 = in3 + in4;in5 = (in[4] + in[8]) * 0.9396926f;//cos( PI/9)in6 = (in[16] - in[8]) * 0.1736482f;//cos(4PI/9)in7 = -(in[4] + in[16]) * 0.7660444f;//cos(2PI/9)in17 = in0 - in5 - in7;in8 = in5 + in0 + in6;in9 = in0 + in7 - in6;in12 = in[6] * 0.8660254f;//cos(PI/6)in10 = (in[2] + in[10]) * 0.9848078f;//cos(PI/18)in11 = (in[14] - in[10]) * 0.3420201f;//cos(7PI/18)in13 = in10 + in11 + in12;out0 = in8 + in13;out8 = in8 - in13;in14 = -(in[2] + in[14]) * 0.6427876f;//cos(5PI/18)in15 = in10 + in14 - in12;in16 = in11 - in14 - in12;out3 = in9 + in15;out5 = in9 - in15;out2 = in17 + in16;out6 = in17 - in16;//<<<<<<<<< End 9-point IDCT on even//>>>>>>>>> 9-point IDCT on odd/* * for(i = 0; i < 9; i++) { * sum = 0; * for(j = 0;j < 18; j += 2) * sum += in[j + 1] * cos(PI36 * (2 * i + 1) * j); * out18[17-i] = sum; * } */in0 = in[1] + in[13] * 0.5f;//cos(PI/3)in1 = in[1] - in[13];in2 = in[9] + in[17] - in[5];out13 = (in1 + in2) * 0.7071068f;//cos(PI/4)in3 = in1 - in2 * 0.5f;in4 = (in[11] + in[15] - in[3]) * 0.8660254f;//cos(PI/6)out16 = (in3 - in4) * 0.5176381f;// 0.5/cos( PI/12)out10 = (in3 + in4) * 1.9318517f;// 0.5/cos(5PI/12)in5 = (in[5] + in[9]) * 0.9396926f;// cos( PI/9)in6 = (in[17] - in[9])* 0.1736482f;// cos(4PI/9)in7 = -(in[5] + in[17])*0.7660444f;// cos(2PI/9)in17 = in0 - in5 - in7;in8 = in5 + in0 + in6;in9 = in0 + in7 - in6;in12 = in[7] * 0.8660254f;// cos(PI/6)in10 = (in[3] + in[11]) * 0.9848078f;// cos(PI/18)in11 = (in[15] - in[11])* 0.3420201f;// cos(7PI/18)in13 = in10 + in11 + in12;out17 = (in8 + in13) * 0.5019099f;// 0.5/cos(PI/36)out9 = (in8 - in13) * 5.7368566f;// 0.5/cos(17PI/36)in14 = -(in[3] + in[15]) * 0.6427876f;// cos(5PI/18)in15 = in10 + in14 - in12;in16 = in11 - in14 - in12;out14 = (in9 + in15) * 0.6103873f;// 0.5/cos(7PI/36)out12 = (in9 - in15) * 0.8717234f;// 0.5/cos(11PI/36)out15 = (in17 + in16) * 0.5516890f;// 0.5/cos(5PI/36)out11 = (in17 - in16) * 1.1831008f;// 0.5/cos(13PI/36)//<<<<<<<<< End 9-point IDCT on odd// Butterflies on 9-point IDCT'stmp = out0; out0 += out17; out17 = tmp - out17;tmp = out1; out1 += out16; out16 = tmp - out16;tmp = out2; out2 += out15; out15 = tmp - out15;tmp = out3; out3 += out14; out14 = tmp - out14;tmp = out4; out4 += out13; out13 = tmp - out13;tmp = out5; out5 += out12; out12 = tmp - out12;tmp = out6; out6 += out11; out11 = tmp - out11;tmp = out7; out7 += out10; out10 = tmp - out10;tmp = out8; out8 += out9; out9 = tmp - out9;//<<<<<<<<<<<<<<<<<< End 18-point IDCT//<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< End 36-point IDCT// Shiftfinal float[] winp = floatWinIMDCT[block_type];out[0] = out9 * winp[0];out[1] = out10 * winp[1];out[2] = out11 * winp[2];out[3] = out12 * winp[3];out[4] = out13 * winp[4];out[5] = out14 * winp[5];out[6] = out15 * winp[6];out[7] = out16 * winp[7];out[8] = out17 * winp[8];out[9] = out17 * winp[9];out[10] = out16 * winp[10];out[11] = out15 * winp[11];out[12] = out14 * winp[12];out[13] = out13 * winp[13];out[14] = out12 * winp[14];out[15] = out11 * winp[15];out[16] = out10 * winp[16];out[17] = out9 * winp[17];out[18] = out8 * winp[18];out[19] = out7 * winp[19];out[20] = out6 * winp[20];out[21] = out5 * winp[21];out[22] = out4 * winp[22];out[23] = out3 * winp[23];out[24] = out2 * winp[24];out[25] = out1 * winp[25];out[26] = out0 * winp[26];out[27] = out0 * winp[27];out[28] = out1 * winp[28];out[29] = out2 * winp[29];out[30] = out3 * winp[30];out[31] = out4 * winp[31];out[32] = out5 * winp[32];out[33] = out6 * winp[33];out[34] = out7 * winp[34];out[35] = out8 * winp[35];}}private static float[] floatRawOut;// [36];private static float[][][] floatPrevBlck;// [2][32][18];private void hybrid(final int ch, final int gr) {GRInfo gr_info = (objSI.ch[ch].gr[gr]);int sb, ss, bt;int max_sb = (17 + rzero_index[ch]) / 18;for (sb = 0; sb < max_sb; sb++) {bt = ((gr_info.window_switching_flag != 0)&& (gr_info.mixed_block_flag != 0) && (sb < 2)) ? 0: gr_info.block_type;inv_mdct(xr[ch][sb], floatRawOut, bt);for(ss = 0; ss < 18; ss++) {xr[ch][sb][ss] = floatRawOut[ss] + floatPrevBlck[ch][sb][ss];floatPrevBlck[ch][sb][ss] = floatRawOut[18 + ss];}}for (; sb < 32; sb++)for(ss = 0; ss < 18; ss++) {xr[ch][sb][ss] = floatPrevBlck[ch][sb][ss];floatPrevBlck[ch][sb][ss] = 0;}}//<<<<HYBRID(synthesize via iMDCT)=========================================`

【下载地址】http://jmp123.sourceforge.net/