fastidct.c

快速离散余玄变换idct

#define W1 2841 
#define W2 2676 
#define W3 2408 
#define W5 1609 
#define W6 1108 
#define W7 565  

#define W1D7    2276    /*W1-W7*/
#define W1A7    3405    /*W1+W7*/
#define W2D6    1568    /*W2-W6*/
#define W2A6    3516    /*W2+W6*/
#define W3D5    799     /*W3-W5*/
#define W3A5    4017    /*W3+W5*/

extern   long    *iclp;

void Fast_IDCT(int * block)     /* fast IDCT algorithm*/
{
    int i;
    int x0, x1, x2, x3, x4, x5, x6, x7, x8;
    
    for(i=0; i<8; i++){
        if (!((x1 = block[8*i+4]<<11) | (x2 = block[8*i+6]) | (x3 = block[8*i+2]) |
            (x4 = block[8*i+1]) | (x5 = block[8*i+7]) | (x6 = block[8*i+5]) | (x7 = block[8*i+3])))
        {
            block[8*i]=block[8*i+1]=block[8*i+2]=block[8*i+3]=block[8*i+4]=block[8*i+5]=block[8*i+6]=block[8*i+7]=block[8*i]<<3;
            continue;
        }
        x0 = (block[8*i]<<11) + 128;  
        /*first stage*/
        x8 = W7*(x4+x5);
        x4 = x8 + W1D7*x4;
        x5 = x8 - W1A7*x5;
        x8 = W3*(x6+x7);
        x6 = x8 - W3D5*x6;
        x7 = x8 - W3A5*x7;
        /*second stage*/
        x8 = x0 + x1;
        x0 -= x1;
        x1 = W6*(x3+x2);
        x2 = x1 - W2A6*x2;
        x3 = x1 + W2D6*x3;
        x1 = x4 + x6;
        x4 -= x6;
        x6 = x5 + x7;
        x5 -= x7;
        /*third stage*/
        x7 = x8 + x3;
        x8 -= x3;
        x3 = x0 + x2;
        x0 -= x2;
        x2 = (181*(x4+x5)+128)>>8;
        x4 = (181*(x4-x5)+128)>>8;
        /*fourth stage*/
        block[8*i] = (x7+x1)>>8;
        block[8*i+1] = (x3+x2)>>8;
        block[8*i+2] = (x0+x4)>>8;
        block[8*i+3] = (x8+x6)>>8;
        block[8*i+4] = (x8-x6)>>8;
        block[8*i+5] = (x0-x4)>>8;
        block[8*i+6] = (x3-x2)>>8;
        block[8*i+7] = (x7-x1)>>8;
        
    }
    for(i=0; i<8; i++){
        if (!((x1 = (block[i+32]<<8)) | (x2 = block[i+48]) | (x3 = block[i+16]) |
            (x4 = block[i+8]) | (x5 = block[i+56]) | (x6 = block[i+40]) | (x7 = block[i+24])))
        {
            block[i]=block[i+8]=block[i+16]=block[i+24]=block[i+32]=block[i+40]
                =block[i+48]=block[i+56]=iclp[(block[i]+32)>>6];
            continue;
        }
        x0 = (block[i]<<8) + 8192;
        /*first step*/
        x8 = W7*(x4+x5) + 4;
        x4 = (x8+W1D7*x4)>>3;
        x5 = (x8-W1A7*x5)>>3;
        x8 = W3*(x6+x7) + 4;
        x6 = (x8-W3D5*x6)>>3;
        x7 = (x8-W3A5*x7)>>3;
        /*second step*/
        x8 = x0 + x1;
        x0 -= x1;
        x1 = W6*(x3+x2) + 4;
        x2 = (x1-W2A6*x2)>>3;
        x3 = (x1+W2D6*x3)>>3;
        x1 = x4 + x6;
        x4 -= x6;
        x6 = x5 + x7;
        x5 -= x7;
        /*third step*/
        x7 = x8 + x3;
        x8 -= x3;
        x3 = x0 + x2;
        x0 -= x2;
        x2 = (181*(x4+x5)+128)>>8;
        x4 = (181*(x4-x5)+128)>>8;
        /*fourth step*/
        block[i] = iclp[(x7+x1)>>14];
        block[i+8] = iclp[(x3+x2)>>14];
        block[i+16] = iclp[(x0+x4)>>14];
        block[i+24] = iclp[(x8+x6)>>14];
        block[i+32] = iclp[(x8-x6)>>14];
        block[i+40] = iclp[(x0-x4)>>14];
        block[i+48] = iclp[(x3-x2)>>14];
        block[i+56] = iclp[(x7-x1)>>14];
    }
}