(转)图形图像处理-之-高质量的快速的图像缩放 中篇 二次线性插值和三次卷积插值

来源:互联网 发布:网络精灵国语版 编辑:程序博客网 时间:2024/05/01 05:39

        图形图像处理-之-高质量的快速的图像缩放 中篇 二次线性插值和三次卷积插值

                     HouSisong@GMail.com   2006255.12.13

(2009.03.07  可以到这里下载缩放算法的完整的可以编译的项目源代码: http://blog.csdn.net/housisong/archive/2009/03/07/3967270.aspx  )

(2007.11.12 替换了二次线性插值的实现(以前偷懒使用了一个近似公式),改进后在图片边缘的插值效果更好(包括三次卷积插值的边界也更精确); 
(2007.09.14 修正三次卷积的MMX版本中表的精度太低(7bit),造成卷积结果误差较大的问题,该版本提高了插值质量,并且速度加快12-25%)
(2007.09.07 PicZoom_ThreeOrder2和PicZoom_ThreeOrder_MMX在缩放的图片宽或高
小于3个像素的时候有一个Bug(边界计算错误);将unsigned long xrIntFloat_16, 
yrIntFloat_16的定义改成long xrIntFloat_16,yrIntFloat_16就可以了)
(2007.07.02 ThreeOrder2_Fast一点小的改进,加快14%)
(2007.06.18 优化PicZoom_BilInear_MMX的实现(由138.5fps提高到147.9fps),
            并添加更快的两路展开的实现版本BilInear_MMX_expand2函数;
            补充新的SSE2的实现PicZoom_BilInear_SSE2函数)
(2007.06.06 更新测试数据,编译器由vc6改为vc2005,CPU由赛扬2G改为AMD64x2 4200+(2.1G) )
(2007.03.06 更新)


tag:图像缩放,速度优化,定点数优化,近邻取样插值,二次线性插值,三次线性插值,
   MipMap链,三次卷积插值,MMX,SSE,SSE2,CPU缓存优化

摘要:首先给出一个基本的图像缩放算法,然后一步一步的优化其速度和缩放质量;

高质量的快速的图像缩放 全文 分为:
     上篇 近邻取样插值和其速度优化
     中篇 二次线性插值和三次卷积插值
     下篇 三次线性插值和MipMap链

     补充 使用SSE2优化

正文:
  为了便于讨论,这里只处理32bit的ARGB颜色;
  代码使用C++;涉及到汇编优化的时候假定为x86平台;使用的编译器为vc2005;
  为了代码的可读性,没有加入异常处理代码;
   测试使用的CPU为AMD64x2 4200+(2.37G)  和 Intel Core2 4400(2.00G);


速度测试说明:
  只测试内存数据到内存数据的缩放
  测试图片都是800*600缩放到1024*768; fps表示每秒钟的帧数,值越大表示函数越快


A:近邻取样插值、二次线性插值、三次卷积插值 缩放效果对比

                                   
       原图         近邻取样缩放到0.6倍     近邻取样缩放到1.6倍

                                         
                二次线性插值缩放到0.6倍   二次线性插值缩放到1.6倍

                                        
               三次卷积插值缩放到0.6倍   三次卷积插值缩放到1.6倍

      
 原图 近邻取样缩放到8倍 二次线性插值缩放到8倍 三次卷积插值缩放到8倍 二次线性插值(近似公式)


     近邻取样插值缩放简单、速度快,但很多时候缩放出的图片质量比较差(特别是对于人物、景色等),
图片的缩放有比较明显的锯齿;使用二次或更高次插值有利于改善缩放效果;


B: 首先定义图像数据结构:

#define asm __asm

typedef unsigned 
char TUInt8; // [0..255]
struct TARGB32      //32 bit color
{
    TUInt8  b,g,r,a;          
//a is alpha
};

struct TPicRegion  //一块颜色数据区的描述,便于参数传递
{
    TARGB32
*    pdata;         //颜色数据首地址
    long        byte_width;    //一行数据的物理宽度(字节宽度);
                
//abs(byte_width)有可能大于等于width*sizeof(TARGB32);
    long        width;         //像素宽度
    long        height;        //像素高度
};

//那么访问一个点的函数可以写为:
inline TARGB32& Pixels(const TPicRegion& pic,const long x,const long y)
{
    
return ( (TARGB32*)((TUInt8*)pic.pdata+pic.byte_width*y) )[x];
}



二次线性插值缩放:

C: 二次线性插值缩放原理和公式图示:

        

             缩放后图片                 原图片
            (宽DW,高DH)              (宽SW,高SH)

  缩放映射原理:
  (Sx-0)/(SW-0)=(Dx-0)/(DW-0)   (Sy-0)/(SH-0)=(Dy-0)/(DH-0)
 =>   Sx=Dx*SW/DW                    Sy=Dy*SH/DH

  聚焦看看(Sx,Sy)坐标点(Sx,Sy为浮点数)附近的情况;


         


  对于近邻取样插值的缩放算法,直接取Color0颜色作为缩放后点的颜色;
二次线性插值需要考虑(Sx,Sy)坐标点周围的4个颜色值Color0/Color1/Color2/Color3,
把(Sx,Sy)到A/B/C/D坐标点的距离作为系数来把4个颜色混合出缩放后点的颜色;
( u=Sx-floor(Sx); v=Sy-floor(Sy); 说明:floor函数的返回值为小于等于参数的最大整数 )  
  二次线性插值公式为:
 tmpColor0=Color0*(1-u) + Color2*u;
 tmpColor1=Color1*(1-u) + Color3*u;
        DstColor =tmpColor0*(1-v) + tmpColor2*v;

  展开公式为:
        pm0=(1-u)*(1-v);
        pm1=v*(1-u);
        pm2=u*(1-v);
        pm3=u*v;
  则颜色混合公式为:
        DstColor = Color0*pm0 + Color1*pm1 + Color2*pm2 + Color3*pm3;

参数函数图示:

      

                                              二次线性插值函数图示

对于上面的公式,它将图片向右下各移动了半个像素,需要对此做一个修正;
  =>   Sx=(Dx+0.5)*SW/DW-0.5; Sy=(Dy+0.5)*SH/DH-0.5;
而实际的程序,还需要考虑到边界(访问源图片可能超界)对于算法的影响,边界的处理可能有各种
方案(不处理边界或边界回绕或边界饱和或边界映射或用背景颜色混合等;文章中默认使用边界饱和来处理超界);
比如:边界饱和函数: 

//访问一个点的函数,(x,y)坐标可能超出图片边界; //边界处理模式:边界饱和
inline TARGB32 Pixels_Bound(const TPicRegion& pic,long x,long y)
{
    
//assert((pic.width>0)&&(pic.height>0));
    bool IsInPic=true;
    
if (x<0) {x=0; IsInPic=false; } else if (x>=pic.width ) {x=pic.width -1; IsInPic=false; }
    
if (y<0) {y=0; IsInPic=false; } else if (y>=pic.height) {y=pic.height-1; IsInPic=false; }
    TARGB32 result
=Pixels(pic,x,y);
    
if (!IsInPic) result.a=0;
    
return result;
}


D: 二次线性插值缩放算法的一个参考实现:PicZoom_BilInear0
  该函数并没有做什么优化,只是一个简单的浮点实现版本;


    inline void Bilinear0(const TPicRegion& pic,float fx,float fy,TARGB32* result)
    {
        
long x=(long)fx; if (x>fx) --x; //x=floor(fx);    
        long y=(long)fy; if (y>fy) --y; //y=floor(fy);
        
        TARGB32 Color0
=Pixels_Bound(pic,x,y);
        TARGB32 Color2
=Pixels_Bound(pic,x+1,y);
        TARGB32 Color1
=Pixels_Bound(pic,x,y+1);
        TARGB32 Color3
=Pixels_Bound(pic,x+1,y+1);

        
float u=fx-x;
        
float v=fy-y;
        
float pm3=u*v;
        
float pm2=u*(1-v);
        
float pm1=v*(1-u);
        
float pm0=(1-u)*(1-v);

        result
->a=(pm0*Color0.a+pm1*Color1.a+pm2*Color2.a+pm3*Color3.a);
        result
->r=(pm0*Color0.r+pm1*Color1.r+pm2*Color2.r+pm3*Color3.r);
        result
->g=(pm0*Color0.g+pm1*Color1.g+pm2*Color2.g+pm3*Color3.g);
        result
->b=(pm0*Color0.b+pm1*Color1.b+pm2*Color2.b+pm3*Color3.b);
    }

void PicZoom_Bilinear0(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    unsigned 
long dst_width=Dst.width;
    TARGB32
* pDstLine=Dst.pdata;
    
for (unsigned long y=0;y<Dst.height;++y)
    {
        
float srcy=(y+0.4999999)*Src.height/Dst.height-0.5;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            
float srcx=(x+0.4999999)*Src.width/Dst.width-0.5;
            Bilinear0(Src,srcx,srcy,
&pDstLine[x]);
        }
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_BilInear0      8.3 fps
////////////////////////////////////////////////////////////////////////////////

E: 把PicZoom_BilInear0的浮点计算改写为定点数实现:PicZoom_BilInear1

    inline void Bilinear1(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        
long x=x_16>>16;
        
long y=y_16>>16;
        TARGB32 Color0
=Pixels_Bound(pic,x,y);
        TARGB32 Color2
=Pixels_Bound(pic,x+1,y);
        TARGB32 Color1
=Pixels_Bound(pic,x,y+1);
        TARGB32 Color3
=Pixels_Bound(pic,x+1,y+1);

        unsigned 
long u_8=(x_16 & 0xFFFF)>>8;
        unsigned 
long v_8=(y_16 & 0xFFFF)>>8;
        unsigned 
long pm3_16=(u_8*v_8);
        unsigned 
long pm2_16=(u_8*(unsigned long)(256-v_8));
        unsigned 
long pm1_16=(v_8*(unsigned long)(256-u_8));
        unsigned 
long pm0_16=((256-u_8)*(256-v_8));

        result
->a=((pm0_16*Color0.a+pm1_16*Color1.a+pm2_16*Color2.a+pm3_16*Color3.a)>>16);
        result
->r=((pm0_16*Color0.r+pm1_16*Color1.r+pm2_16*Color2.r+pm3_16*Color3.r)>>16);
        result
->g=((pm0_16*Color0.g+pm1_16*Color1.g+pm2_16*Color2.g+pm3_16*Color3.g)>>16);
        result
->b=((pm0_16*Color0.b+pm1_16*Color1.b+pm2_16*Color2.b+pm3_16*Color3.b)>>16);
    }

void PicZoom_Bilinear1(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    
long xrIntFloat_16=((Src.width)<<16)/Dst.width+1
    
long yrIntFloat_16=((Src.height)<<16)/Dst.height+1;
    
const long csDErrorX=-(1<<15)+(xrIntFloat_16>>1);
    
const long csDErrorY=-(1<<15)+(yrIntFloat_16>>1);

    unsigned 
long dst_width=Dst.width;

    TARGB32
* pDstLine=Dst.pdata;
    
long srcy_16=csDErrorY;
    
long y;
    
for (y=0;y<Dst.height;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear1(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }


////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_BilInear1     17.7 fps
////////////////////////////////////////////////////////////////////////////////

F: 二次线性插值需要考略边界访问超界的问题,我们可以将边界区域和内部区域分开处理,这样就可以优化内部的插值实现函数了:比如不需要判断访问超界、减少颜色数据复制、减少一些不必要的重复坐标计算等等



    inline void Bilinear2_Fast(TARGB32* PColor0,TARGB32* PColor1,unsigned long u_8,unsigned long v_8,TARGB32* result)
    {
        unsigned 
long pm3_16=u_8*v_8;
        unsigned 
long pm2_16=(u_8<<8)-pm3_16;
        unsigned 
long pm1_16=(v_8<<8)-pm3_16;
        unsigned 
long pm0_16=(1<<16)-pm1_16-pm2_16-pm3_16;
   
        result
->a=((pm0_16*PColor0[0].a+pm2_16*PColor0[1].a+pm1_16*PColor1[0].a+pm3_16*PColor1[1].a)>>16);
        result
->r=((pm0_16*PColor0[0].r+pm2_16*PColor0[1].r+pm1_16*PColor1[0].r+pm3_16*PColor1[1].r)>>16);
        result
->g=((pm0_16*PColor0[0].g+pm2_16*PColor0[1].g+pm1_16*PColor1[0].g+pm3_16*PColor1[1].g)>>16);
        result
->b=((pm0_16*PColor0[0].b+pm2_16*PColor0[1].b+pm1_16*PColor1[0].b+pm3_16*PColor1[1].b)>>16);
    }

    inline 
void Bilinear2_Border(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        
long x=(x_16>>16);
        
long y=(y_16>>16);
        unsigned 
long u_16=((unsigned short)(x_16));
        unsigned 
long v_16=((unsigned short)(y_16));

        TARGB32 pixel[
4];
        pixel[
0]=Pixels_Bound(pic,x,y);
        pixel[
1]=Pixels_Bound(pic,x+1,y);
        pixel[
2]=Pixels_Bound(pic,x,y+1);
        pixel[
3]=Pixels_Bound(pic,x+1,y+1);
        
        Bilinear2_Fast(
&pixel[0],&pixel[2],u_16>>8,v_16>>8,result);
    }

void PicZoom_Bilinear2(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    
long xrIntFloat_16=((Src.width)<<16)/Dst.width+1
    
long yrIntFloat_16=((Src.height)<<16)/Dst.height+1;
    
const long csDErrorX=-(1<<15)+(xrIntFloat_16>>1);
    
const long csDErrorY=-(1<<15)+(yrIntFloat_16>>1);

    unsigned 
long dst_width=Dst.width;

    
//计算出需要特殊处理的边界
    long border_y0=-csDErrorY/yrIntFloat_16+1;              //y0+y*yr>=0; y0=csDErrorY => y>=-csDErrorY/yr
    if (border_y0>=Dst.height) border_y0=Dst.height;
    
long border_x0=-csDErrorX/xrIntFloat_16+1;     
    
if (border_x0>=Dst.width ) border_x0=Dst.width; 
    
long border_y1=(((Src.height-2)<<16)-csDErrorY)/yrIntFloat_16+1//y0+y*yr<=(height-2) => y<=(height-2-csDErrorY)/yr
    if (border_y1<border_y0) border_y1=border_y0;
    
long border_x1=(((Src.width-2)<<16)-csDErrorX)/xrIntFloat_16+1
    
if (border_x1<border_x0) border_x1=border_x0;

    TARGB32
* pDstLine=Dst.pdata;
    
long Src_byte_width=Src.byte_width;
    
long srcy_16=csDErrorY;
    
long y;
    
for (y=0;y<border_y0;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear2_Border(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y0;y<border_y1;++y)
    {
        
long srcx_16=csDErrorX;
        
long x;
        
for (x=0;x<border_x0;++x)
        {
            Bilinear2_Border(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }

        {
            unsigned 
long v_8=(srcy_16 & 0xFFFF)>>8;
            TARGB32
* PSrcLineColor= (TARGB32*)((TUInt8*)(Src.pdata)+Src_byte_width*(srcy_16>>16)) ;
            
for (unsigned long x=border_x0;x<border_x1;++x)
            {
                TARGB32
* PColor0=&PSrcLineColor[srcx_16>>16];
                TARGB32
* PColor1=(TARGB32*)((TUInt8*)(PColor0)+Src_byte_width);
                Bilinear2_Fast(PColor0,PColor1,(srcx_16 
& 0xFFFF)>>8,v_8,&pDstLine[x]);
                srcx_16
+=xrIntFloat_16;
            }
        }

        
for (x=border_x1;x<dst_width;++x)
        {
            Bilinear2_Border(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y1;y<Dst.height;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear2_Border(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_BilInear2     43.4 fps
////////////////////////////////////////////////////////////////////////////////

(F'补充:
  如果不想处理边界访问超界问题,可以考虑扩大源图片的尺寸,加一个边框 (“哨兵”优化);
这样插值算法就不用考虑边界问题了,程序写起来也简单很多! 
  如果对缩放结果的边界像素级精度要求不是太高,我还有一个方案,一个稍微改变的缩放公式:
  Sx=Dx*(SW-1)/DW; Sy=Dy*(SH-1)/DH;  (源图片宽和高:SW>=2;SH>=2)
  证明这个公式不会造成内存访问超界: 
   要求Dx=DW-1时: sx+1=int( (dw-1)/dw*(dw-1) ) +1 <= (sw-1)
      有:  int( (sw-1)*(dw-1)/dw ) <=sw-2
             (sw-1)*(dw-1)/dw <(sw-1)
             (dw-1) /dw<1
             (dw-1) <dw
  比如,按这个公式的一个简单实现: (缩放效果见前面的"二次线性插值(近似公式)"图示)

void PicZoom_ftBilinear_Common(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(2>Src.width)||(2>Src.height)) return;

    
long xrIntFloat_16=((Src.width-1)<<16)/Dst.width; 
    
long yrIntFloat_16=((Src.height-1)<<16)/Dst.height;

    unsigned 
long dst_width=Dst.width;
    
long Src_byte_width=Src.byte_width;
    TARGB32
* pDstLine=Dst.pdata;
    
long srcy_16=0;
    
for (unsigned long y=0;y<Dst.height;++y)
    {
        unsigned 
long v_8=(srcy_16 & 0xFFFF)>>8;
        TARGB32
* PSrcLineColor= (TARGB32*)((TUInt8*)(Src.pdata)+Src_byte_width*(srcy_16>>16)) ;
        
long srcx_16=0;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            TARGB32
* PColor0=&PSrcLineColor[srcx_16>>16];
            Bilinear_Fast_Common(PColor0,(TARGB32
*)((TUInt8*)(PColor0)+Src_byte_width),(srcx_16 & 0xFFFF)>>8,v_8,&pDstLine[x]);
            srcx_16
+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }


)

G:利用单指令多数据处理的MMX指令一般都可以加快颜色的运算;在使用MMX改写之前,利用
32bit寄存器(或变量)来模拟单指令多数据处理;
数据储存原理:一个颜色数据分量只有一个字节,用2个字节来储存单个颜色分量的计算结果,
对于很多颜色计算来说精度就够了;那么一个32bit寄存器(或变量)就可以储存2个计算出的
临时颜色分量;从而达到了单个指令两路数据处理的目的;
单个指令两路数据处理的计算: 
  乘法: ((0x00AA*a)<<16) | (0x00BB*a) = 0x00AA00BB * a 
    可见只要保证0x00AA*a和0x00BB*a都小于(1<<16)那么乘法可以直接使用无符号数乘法了
  加法: ((0x00AA+0x00CC)<<16) | (0x00BB+0x00DD) = 0x00AA00BB + 0x00CC00DD 
    可见只要0x00AA+0x00CC和0x00BB+0x00DD小于(1<<16)那么加法可以直接使用无符号数加法了
  (移位、减法等稍微复杂一点,因为这里没有用到就不推倒运算公式了)


    inline void Bilinear_Fast_Common(TARGB32* PColor0,TARGB32* PColor1,unsigned long u_8,unsigned long v_8,TARGB32* result)
    {
        unsigned 
long pm3_8=(u_8*v_8)>>8;
        unsigned 
long pm2_8=u_8-pm3_8;
        unsigned 
long pm1_8=v_8-pm3_8;
        unsigned 
long pm0_8=256-pm1_8-pm2_8-pm3_8;

        unsigned 
long Color=*(unsigned long*)(PColor0);
        unsigned 
long BR=(Color & 0x00FF00FF)*pm0_8;
        unsigned 
long GA=((Color & 0xFF00FF00)>>8)*pm0_8;
                      Color
=((unsigned long*)(PColor0))[1];
                      GA
+=((Color & 0xFF00FF00)>>8)*pm2_8;
                      BR
+=(Color & 0x00FF00FF)*pm2_8;
                      Color
=*(unsigned long*)(PColor1);
                      GA
+=((Color & 0xFF00FF00)>>8)*pm1_8;
                      BR
+=(Color & 0x00FF00FF)*pm1_8;
                      Color
=((unsigned long*)(PColor1))[1];
                      GA
+=((Color & 0xFF00FF00)>>8)*pm3_8;
                      BR
+=(Color & 0x00FF00FF)*pm3_8;

        
*(unsigned long*)(result)=(GA & 0xFF00FF00)|((BR & 0xFF00FF00)>>8);
    }

    inline 
void Bilinear_Border_Common(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        
long x=(x_16>>16);
        
long y=(y_16>>16);
        unsigned 
long u_16=((unsigned short)(x_16));
        unsigned 
long v_16=((unsigned short)(y_16));

        TARGB32 pixel[
4];
        pixel[
0]=Pixels_Bound(pic,x,y);
        pixel[
1]=Pixels_Bound(pic,x+1,y);
        pixel[
2]=Pixels_Bound(pic,x,y+1);
        pixel[
3]=Pixels_Bound(pic,x+1,y+1);
        
        Bilinear_Fast_Common(
&pixel[0],&pixel[2],u_16>>8,v_16>>8,result);
    }

void PicZoom_Bilinear_Common(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    
long xrIntFloat_16=((Src.width)<<16)/Dst.width+1
    
long yrIntFloat_16=((Src.height)<<16)/Dst.height+1;
    
const long csDErrorX=-(1<<15)+(xrIntFloat_16>>1);
    
const long csDErrorY=-(1<<15)+(yrIntFloat_16>>1);

    unsigned 
long dst_width=Dst.width;

    
//计算出需要特殊处理的边界
    long border_y0=-csDErrorY/yrIntFloat_16+1;              //y0+y*yr>=0; y0=csDErrorY => y>=-csDErrorY/yr
    if (border_y0>=Dst.height) border_y0=Dst.height;
    
long border_x0=-csDErrorX/xrIntFloat_16+1;     
    
if (border_x0>=Dst.width ) border_x0=Dst.width; 
    
long border_y1=(((Src.height-2)<<16)-csDErrorY)/yrIntFloat_16+1//y0+y*yr<=(height-2) => y<=(height-2-csDErrorY)/yr
    if (border_y1<border_y0) border_y1=border_y0;
    
long border_x1=(((Src.width-2)<<16)-csDErrorX)/xrIntFloat_16+1
    
if (border_x1<border_x0) border_x1=border_x0;

    TARGB32
* pDstLine=Dst.pdata;
    
long Src_byte_width=Src.byte_width;
    
long srcy_16=csDErrorY;
    
long y;
    
for (y=0;y<border_y0;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y0;y<border_y1;++y)
    {
        
long srcx_16=csDErrorX;
        
long x;
        
for (x=0;x<border_x0;++x)
        {
            Bilinear_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }

        {
            unsigned 
long v_8=(srcy_16 & 0xFFFF)>>8;
            TARGB32
* PSrcLineColor= (TARGB32*)((TUInt8*)(Src.pdata)+Src_byte_width*(srcy_16>>16)) ;
            
for (unsigned long x=border_x0;x<border_x1;++x)
            {
                TARGB32
* PColor0=&PSrcLineColor[srcx_16>>16];
                TARGB32
* PColor1=(TARGB32*)((TUInt8*)(PColor0)+Src_byte_width);
                Bilinear_Fast_Common(PColor0,PColor1,(srcx_16 
& 0xFFFF)>>8,v_8,&pDstLine[x]);
                srcx_16
+=xrIntFloat_16;
            }
        }

        
for (x=border_x1;x<dst_width;++x)
        {
            Bilinear_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y1;y<Dst.height;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_BilInear_Common   65.3 fps
////////////////////////////////////////////////////////////////////////////////

H:使用MMX指令改写:PicZoom_Bilinear_MMX

    inline void  Bilinear_Fast_MMX(TARGB32* PColor0,TARGB32* PColor1,unsigned long u_8,unsigned long v_8,TARGB32* result)
    {
        asm
        {    
              MOVD      MM6,v_8
              MOVD      MM5,u_8
              mov       edx,PColor0
              mov       eax,PColor1
              PXOR      mm7,mm7

              MOVD         MM2,dword ptr [eax]  
              MOVD         MM0,dword ptr [eax
+4]
              PUNPCKLWD    MM5,MM5
              PUNPCKLWD    MM6,MM6
              MOVD         MM3,dword ptr [edx]  
              MOVD         MM1,dword ptr [edx
+4]
              PUNPCKLDQ    MM5,MM5 
              PUNPCKLBW    MM0,MM7
              PUNPCKLBW    MM1,MM7
              PUNPCKLBW    MM2,MM7
              PUNPCKLBW    MM3,MM7
              PSUBw        MM0,MM2
              PSUBw        MM1,MM3
              PSLLw        MM2,
8
              PSLLw        MM3,
8
              PMULlw       MM0,MM5
              PMULlw       MM1,MM5
              PUNPCKLDQ    MM6,MM6 
              PADDw        MM0,MM2
              PADDw        MM1,MM3

              PSRLw        MM0,
8
              PSRLw        MM1,
8
              PSUBw        MM0,MM1
              PSLLw        MM1,
8
              PMULlw       MM0,MM6
              mov       eax,result
              PADDw        MM0,MM1

              PSRLw        MM0,
8
              PACKUSwb     MM0,MM7
              movd      [eax],MM0 
              
//emms
        }
    }

    
void Bilinear_Border_MMX(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        
long x=(x_16>>16);
        
long y=(y_16>>16);
        unsigned 
long u_16=((unsigned short)(x_16));
        unsigned 
long v_16=((unsigned short)(y_16));

        TARGB32 pixel[
4];
        pixel[
0]=Pixels_Bound(pic,x,y);
        pixel[
1]=Pixels_Bound(pic,x+1,y);
        pixel[
2]=Pixels_Bound(pic,x,y+1);
        pixel[
3]=Pixels_Bound(pic,x+1,y+1);
        
        Bilinear_Fast_MMX(
&pixel[0],&pixel[2],u_16>>8,v_16>>8,result);
    }

void PicZoom_Bilinear_MMX(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    
long xrIntFloat_16=((Src.width)<<16)/Dst.width+1
    
long yrIntFloat_16=((Src.height)<<16)/Dst.height+1;
    
const long csDErrorX=-(1<<15)+(xrIntFloat_16>>1);
    
const long csDErrorY=-(1<<15)+(yrIntFloat_16>>1);

    unsigned 
long dst_width=Dst.width;

    
//计算出需要特殊处理的边界
    long border_y0=-csDErrorY/yrIntFloat_16+1;              //y0+y*yr>=0; y0=csDErrorY => y>=-csDErrorY/yr
    if (border_y0>=Dst.height) border_y0=Dst.height;
    
long border_x0=-csDErrorX/xrIntFloat_16+1;     
    
if (border_x0>=Dst.width ) border_x0=Dst.width; 
    
long border_y1=(((Src.height-2)<<16)-csDErrorY)/yrIntFloat_16+1//y0+y*yr<=(height-2) => y<=(height-2-csDErrorY)/yr
    if (border_y1<border_y0) border_y1=border_y0;
    
long border_x1=(((Src.width-2)<<16)-csDErrorX)/xrIntFloat_16+1
    
if (border_x1<border_x0) border_x1=border_x0;

    TARGB32
* pDstLine=Dst.pdata;
    
long Src_byte_width=Src.byte_width;
    
long srcy_16=csDErrorY;
    
long y;
    
for (y=0;y<border_y0;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y0;y<border_y1;++y)
    {
        
long srcx_16=csDErrorX;
        
long x;
        
for (x=0;x<border_x0;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }

        {
            unsigned 
long v_8=(srcy_16 & 0xFFFF)>>8;
            TARGB32
* PSrcLineColor= (TARGB32*)((TUInt8*)(Src.pdata)+Src_byte_width*(srcy_16>>16)) ;
            
for (unsigned long x=border_x0;x<border_x1;++x)
            {
                TARGB32
* PColor0=&PSrcLineColor[srcx_16>>16];
                TARGB32
* PColor1=(TARGB32*)((TUInt8*)(PColor0)+Src_byte_width);
                Bilinear_Fast_MMX(PColor0,PColor1,(srcx_16 
& 0xFFFF)>>8,v_8,&pDstLine[x]);
                srcx_16
+=xrIntFloat_16;
            }
        }

        
for (x=border_x1;x<dst_width;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y1;y<Dst.height;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    asm emms
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_BilInear_MMX 132.9 fps
////////////////////////////////////////////////////////////////////////////////

 

H' 对BilInear_MMX简单改进:PicZoom_Bilinear_MMX_Ex


void PicZoom_Bilinear_MMX_Ex(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    
long xrIntFloat_16=((Src.width)<<16)/Dst.width+1
    
long yrIntFloat_16=((Src.height)<<16)/Dst.height+1;
    
const long csDErrorX=-(1<<15)+(xrIntFloat_16>>1);
    
const long csDErrorY=-(1<<15)+(yrIntFloat_16>>1);

    unsigned 
long dst_width=Dst.width;

    
//计算出需要特殊处理的边界
    long border_y0=-csDErrorY/yrIntFloat_16+1;              //y0+y*yr>=0; y0=csDErrorY => y>=-csDErrorY/yr
    if (border_y0>=Dst.height) border_y0=Dst.height;
    
long border_x0=-csDErrorX/xrIntFloat_16+1;     
    
if (border_x0>=Dst.width ) border_x0=Dst.width; 
    
long border_y1=(((Src.height-2)<<16)-csDErrorY)/yrIntFloat_16+1//y0+y*yr<=(height-2) => y<=(height-2-csDErrorY)/yr
    if (border_y1<border_y0) border_y1=border_y0;
    
long border_x1=(((Src.width-2)<<16)-csDErrorX)/xrIntFloat_16+1
    
if (border_x1<border_x0) border_x1=border_x0;

    TARGB32
* pDstLine=Dst.pdata;
    
long Src_byte_width=Src.byte_width;
    
long srcy_16=csDErrorY;
    
long y;
    
for (y=0;y<border_y0;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }

    
for (y=border_y0;y<border_y1;++y)
    {
        
long srcx_16=csDErrorX;
        
long x;
        
for (x=0;x<border_x0;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }

        {
            
long dst_width_fast=border_x1-border_x0;
            
if (dst_width_fast>0)
            {
                unsigned 
long v_8=(srcy_16 & 0xFFFF)>>8;
                TARGB32
* PSrcLineColor= (TARGB32*)((TUInt8*)(Src.pdata)+Src_byte_width*(srcy_16>>16)) ;
                TARGB32
* PSrcLineColorNext= (TARGB32*)((TUInt8*)(PSrcLineColor)+Src_byte_width) ;
                TARGB32
* pDstLine_Fast=&pDstLine[border_x0];
                asm
                {
                      movd         mm6,v_8
                      pxor         mm7,mm7 
//mm7=0
                      PUNPCKLWD    MM6,MM6
                      PUNPCKLDQ    MM6,MM6
//mm6=v_8
                    
                      mov       esi,PSrcLineColor
                      mov       ecx,PSrcLineColorNext
                      mov       edx,srcx_16
                      mov       ebx,dst_width_fast
                      mov       edi,pDstLine_Fast
                      lea       edi,[edi
+ebx*4]
                      push      ebp
                      mov       ebp,xrIntFloat_16
                      neg       ebx

                loop_start:

                          mov       eax,edx
                          shl       eax,
16
                          shr       eax,
24
                          
//== movzx       eax,dh  //eax=u_8
                          MOVD      MM5,eax
                          mov       eax,edx
                          shr       eax,
16     //srcx_16>>16

                          MOVD         MM2,dword ptr [ecx
+eax*4]  
                          MOVD         MM0,dword ptr [ecx
+eax*4+4]
                          PUNPCKLWD    MM5,MM5
                          MOVD         MM3,dword ptr [esi
+eax*4]  
                          MOVD         MM1,dword ptr [esi
+eax*4+4]
                          PUNPCKLDQ    MM5,MM5 
//mm5=u_8
                          PUNPCKLBW    MM0,MM7
                          PUNPCKLBW    MM1,MM7
                          PUNPCKLBW    MM2,MM7
                          PUNPCKLBW    MM3,MM7
                          PSUBw        MM0,MM2
                          PSUBw        MM1,MM3
                          PSLLw        MM2,
8
                          PSLLw        MM3,
8
                          PMULlw       MM0,MM5
                          PMULlw       MM1,MM5
                          PADDw        MM0,MM2
                          PADDw        MM1,MM3

                          PSRLw        MM0,
8
                          PSRLw        MM1,
8
                          PSUBw        MM0,MM1
                          PSLLw        MM1,
8
                          PMULlw       MM0,MM6
                          PADDw        MM0,MM1

                          PSRLw     MM0,
8
                          PACKUSwb  MM0,MM7
                          MOVd   dword ptr    [edi
+ebx*4],MM0 //write DstColor
                                      
                          add       edx,ebp 
//srcx_16+=xrIntFloat_16
                          inc       ebx
                          jnz       loop_start

                      pop       ebp
                      mov       srcx_16,edx
                }
            }
        }

        
for (x=border_x1;x<dst_width;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y1;y<Dst.height;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            Bilinear_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    asm emms
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_Bilinear_MMX_Ex 157.0 fps
////////////////////////////////////////////////////////////////////////////////

I: 把测试成绩放在一起:

////////////////////////////////////////////////////////////////////////////////
//CPU: AMD64x2 4200+(2.37G)  zoom 800*600 to 1024*768
//==============================================================================
// StretchBlt                   232.7 fps   
// PicZoom3_SSE                 711.7 fps 
// 
// PicZoom_BilInear0              8.3 fps
// PicZoom_BilInear1             17.7 fps
// PicZoom_BilInear2             43.4 fps
// PicZoom_BilInear_Common       65.3 fps
// PicZoom_BilInear_MMX         132.9 fps
// PicZoom_BilInear_MMX_Ex      157.0 fps
////////////////////////////////////////////////////////////////////////////////

补充Intel Core2 4400上的测试成绩:

////////////////////////////////////////////////////////////////////////////////
//CPU: Intel Core2 4400(2.00G)  zoom 800*600 to 1024*768
//==============================================================================
// PicZoom3_SSE                1099.7 fps  
// 
// PicZoom_BilInear0             10.7 fps
// PicZoom_BilInear1             24.2 fps
// PicZoom_BilInear2             54.3 fps
// PicZoom_BilInear_Common       59.8 fps
// PicZoom_BilInear_MMX         118.4 fps
// PicZoom_BilInear_MMX_Ex      142.9 fps
////////////////////////////////////////////////////////////////////////////////

 

 

三次卷积插值:

 


J: 二次线性插值缩放出的图片很多时候让人感觉变得模糊(术语叫低通滤波),特别是在放大
的时候;使用三次卷积插值来改善插值结果;三次卷积插值考虑映射点周围16个点(4x4)的颜色来
计算最终的混合颜色,如图;

          
         P(0,0)所在像素为映射的点,加上它周围的15个点,按一定系数混合得到最终输出结果;

         混合公式参见PicZoom_ThreeOrder0的实现;

    插值曲线公式sin(x*PI)/(x*PI),如图:


                             三次卷积插值曲线sin(x*PI)/(x*PI) (其中PI=3.1415926...)

 

K:三次卷积插值缩放算法的一个参考实现:PicZoom_ThreeOrder0
  该函数并没有做过多的优化,只是一个简单的浮点实现版本;



        inline 
double SinXDivX(double x) 
        {
            
//该函数计算插值曲线sin(x*PI)/(x*PI)的值 //PI=3.1415926535897932385; 
            
//下面是它的近似拟合表达式
            const float a = -1//a还可以取 a=-2,-1,-0.75,-0.5等等,起到调节锐化或模糊程度的作用

            
if (x<0) x=-x; //x=abs(x);
            double x2=x*x;
            
double x3=x2*x;
            
if (x<=1)
              
return (a+2)*x3 - (a+3)*x2 + 1;
            
else if (x<=2
              
return a*x3 - (5*a)*x2 + (8*a)*- (4*a);
            
else
              
return 0;
        } 

        inline TUInt8 border_color(
long Color)
        {
            
if (Color<=0)
                
return 0;
            
else if (Color>=255)
                
return 255;
            
else
                
return Color;
        }
        
    
void ThreeOrder0(const TPicRegion& pic,const float fx,const float fy,TARGB32* result)
    {
        
long x0=(long)fx; if (x0>fx) --x0; //x0=floor(fx);    
        long y0=(long)fy; if (y0>fy) --y0; //y0=floor(fy);
        float fu=fx-x0;
        
float fv=fy-y0;

        TARGB32 pixel[
16];
        
long i,j;

        
for (i=0;i<4;++i)
        {
            
for (j=0;j<4;++j)
            {
                
long x=x0-1+j;
                
long y=y0-1+i;
                pixel[i
*4+j]=Pixels_Bound(pic,x,y);
            }
        }

        
float afu[4],afv[4];
        
//
        afu[0]=SinXDivX(1+fu);
        afu[
1]=SinXDivX(fu);
        afu[
2]=SinXDivX(1-fu);
        afu[
3]=SinXDivX(2-fu);
        afv[
0]=SinXDivX(1+fv);
        afv[
1]=SinXDivX(fv);
        afv[
2]=SinXDivX(1-fv);
        afv[
3]=SinXDivX(2-fv);

        
float sR=0,sG=0,sB=0,sA=0;
        
for (i=0;i<4;++i)
        {
            
float aR=0,aG=0,aB=0,aA=0;
            
for (long j=0;j<4;++j)
            {
                aA
+=afu[j]*pixel[i*4+j].a;
                aR
+=afu[j]*pixel[i*4+j].r;
                aG
+=afu[j]*pixel[i*4+j].g;
                aB
+=afu[j]*pixel[i*4+j].b;
            }
            sA
+=aA*afv[i];
            sR
+=aR*afv[i];
            sG
+=aG*afv[i];
            sB
+=aB*afv[i];
        }

        result
->a=border_color((long)(sA+0.5));
        result
->r=border_color((long)(sR+0.5));
        result
->g=border_color((long)(sG+0.5));
        result
->b=border_color((long)(sB+0.5));
    }

void PicZoom_ThreeOrder0(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;


    unsigned 
long dst_width=Dst.width;
    TARGB32
* pDstLine=Dst.pdata;
    
for (unsigned long y=0;y<Dst.height;++y)
    {
        
float srcy=(y+0.4999999)*Src.height/Dst.height-0.5;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            
float srcx=(x+0.4999999)*Src.width/Dst.width-0.5;
            ThreeOrder0(Src,srcx,srcy,
&pDstLine[x]);
        }
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_ThreeOrder0    3.6 fps
////////////////////////////////////////////////////////////////////////////////

L: 使用定点数来优化缩放函数;边界和内部分开处理;对SinXDivX做一个查找表;对border_color做一个查找表;


    static long SinXDivX_Table_8[(2<<8)+1];
    
class _CAutoInti_SinXDivX_Table {
    
private
        
void _Inti_SinXDivX_Table()
        {
            
for (long i=0;i<=(2<<8);++i)
                SinXDivX_Table_8[i]
=long(0.5+256*SinXDivX(i*(1.0/(256))))*1;
        };
    
public:
        _CAutoInti_SinXDivX_Table() { _Inti_SinXDivX_Table(); }
    };
    
static _CAutoInti_SinXDivX_Table __tmp_CAutoInti_SinXDivX_Table;


    
//颜色查表
    static TUInt8 _color_table[256*3];
    
static const TUInt8* color_table=&_color_table[256];
    
class _CAuto_inti_color_table
    {
    
public:
        _CAuto_inti_color_table() {
            
for (int i=0;i<256*3;++i)
                _color_table[i]
=border_color(i-256);
        }
    };
    
static _CAuto_inti_color_table _Auto_inti_color_table;

    
void ThreeOrder_Fast_Common(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        unsigned long u_8=(unsigned char)((x_16)>>8);
        unsigned 
long v_8=(unsigned char)((y_16)>>8);
        
const TARGB32* pixel=&Pixels(pic,(x_16>>16)-1,(y_16>>16)-1);
        
long pic_byte_width=pic.byte_width;

        
long au_8[4],av_8[4];
        
//
        au_8[0]=SinXDivX_Table_8[(1<<8)+u_8];
        au_8[
1]=SinXDivX_Table_8[u_8];
        au_8[
2]=SinXDivX_Table_8[(1<<8)-u_8];
        au_8[
3]=SinXDivX_Table_8[(2<<8)-u_8];
        av_8[
0]=SinXDivX_Table_8[(1<<8)+v_8];
        av_8[
1]=SinXDivX_Table_8[v_8];
        av_8[
2]=SinXDivX_Table_8[(1<<8)-v_8];
        av_8[
3]=SinXDivX_Table_8[(2<<8)-v_8];

        
long sR=0,sG=0,sB=0,sA=0;
        
for (long i=0;i<4;++i)
        {
            
long aA=au_8[0]*pixel[0].a + au_8[1]*pixel[1].a + au_8[2]*pixel[2].a + au_8[3]*pixel[3].a;
            
long aR=au_8[0]*pixel[0].r + au_8[1]*pixel[1].r + au_8[2]*pixel[2].r + au_8[3]*pixel[3].r;
            
long aG=au_8[0]*pixel[0].g + au_8[1]*pixel[1].g + au_8[2]*pixel[2].g + au_8[3]*pixel[3].g;
            
long aB=au_8[0]*pixel[0].b + au_8[1]*pixel[1].b + au_8[2]*pixel[2].b + au_8[3]*pixel[3].b;
            sA
+=aA*av_8[i];
            sR
+=aR*av_8[i];
            sG
+=aG*av_8[i];
            sB
+=aB*av_8[i];
            ((TUInt8
*&)pixel)+=pic_byte_width;
        }

        result
->a=color_table[sA>>16];
        result
->r=color_table[sR>>16];
        result
->g=color_table[sG>>16];
        result
->b=color_table[sB>>16];
    }

    
void ThreeOrder_Border_Common(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        
long x0_sub1=(x_16>>16)-1;
        
long y0_sub1=(y_16>>16)-1;
        unsigned 
long u_16_add1=((unsigned short)(x_16))+(1<<16);
        unsigned 
long v_16_add1=((unsigned short)(y_16))+(1<<16);

        TARGB32 pixel[
16];
        
long i;

        
for (i=0;i<4;++i)
        {
            
long y=y0_sub1+i;
            pixel[i
*4+0]=Pixels_Bound(pic,x0_sub1+0,y);
            pixel[i
*4+1]=Pixels_Bound(pic,x0_sub1+1,y);
            pixel[i
*4+2]=Pixels_Bound(pic,x0_sub1+2,y);
            pixel[i
*4+3]=Pixels_Bound(pic,x0_sub1+3,y);
        }
        
        TPicRegion npic;
        npic.pdata     
=&pixel[0];
        npic.byte_width
=4*sizeof(TARGB32);
        
//npic.width     =4;
        
//npic.height    =4;
        ThreeOrder_Fast_Common(npic,u_16_add1,v_16_add1,result);
    }

void PicZoom_ThreeOrder_Common(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    
long xrIntFloat_16=((Src.width)<<16)/Dst.width+1
    
long yrIntFloat_16=((Src.height)<<16)/Dst.height+1;
    
const long csDErrorX=-(1<<15)+(xrIntFloat_16>>1);
    
const long csDErrorY=-(1<<15)+(yrIntFloat_16>>1);

    unsigned 
long dst_width=Dst.width;

    
//计算出需要特殊处理的边界
    long border_y0=((1<<16)-csDErrorY)/yrIntFloat_16+1;//y0+y*yr>=1; y0=csDErrorY => y>=(1-csDErrorY)/yr
    if (border_y0>=Dst.height) border_y0=Dst.height;
    
long border_x0=((1<<16)-csDErrorX)/xrIntFloat_16+1;
    
if (border_x0>=Dst.width ) border_x0=Dst.width;
    
long border_y1=(((Src.height-3)<<16)-csDErrorY)/yrIntFloat_16+1//y0+y*yr<=(height-3) => y<=(height-3-csDErrorY)/yr
    if (border_y1<border_y0) border_y1=border_y0;
    
long border_x1=(((Src.width-3)<<16)-csDErrorX)/xrIntFloat_16+1;; 
    
if (border_x1<border_x0) border_x1=border_x0;

    TARGB32
* pDstLine=Dst.pdata;
    
long srcy_16=csDErrorY;
    
long y;
    
for (y=0;y<border_y0;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            ThreeOrder_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y0;y<border_y1;++y)
    {
        
long srcx_16=csDErrorX;
        
long x;
        
for (x=0;x<border_x0;++x)
        {
            ThreeOrder_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        
for (x=border_x0;x<border_x1;++x)
        {
            ThreeOrder_Fast_Common(Src,srcx_16,srcy_16,
&pDstLine[x]);//fast  !
            srcx_16+=xrIntFloat_16;
        }
        
for (x=border_x1;x<dst_width;++x)
        {
            ThreeOrder_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y1;y<Dst.height;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            ThreeOrder_Border_Common(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_ThreeOrder_Common    16.9 fps
////////////////////////////////////////////////////////////////////////////////

 M: 用MMX来优化ThreeOrder_Common函数:ThreeOrder_MMX


    typedef   unsigned long TMMXData32;
    
static TMMXData32 SinXDivX_Table_MMX[(2<<8)+1];
    
class _CAutoInti_SinXDivX_Table_MMX {
    
private
        
void _Inti_SinXDivX_Table_MMX()
        {
            
for (long i=0;i<=(2<<8);++i)
            {
                unsigned 
short t=long(0.5+(1<<14)*SinXDivX(i*(1.0/(256))));
                unsigned 
long tl=| (((unsigned long)t)<<16);
                SinXDivX_Table_MMX[i]
=tl;
            }
        };
    
public:
        _CAutoInti_SinXDivX_Table_MMX() { _Inti_SinXDivX_Table_MMX(); }
    };
    
static _CAutoInti_SinXDivX_Table_MMX __tmp_CAutoInti_SinXDivX_Table_MMX;


    
void __declspec(naked) _private_ThreeOrder_Fast_MMX()
    {
        asm
        {
            movd        mm1,dword ptr [edx]
            movd        mm2,dword ptr [edx
+4]
            movd        mm3,dword ptr [edx
+8]
            movd        mm4,dword ptr [edx
+12]
            movd        mm5,dword ptr [(offset SinXDivX_Table_MMX)
+256*4+eax*4]
            movd        mm6,dword ptr [(offset SinXDivX_Table_MMX)
+eax*4]
            punpcklbw   mm1,mm7
            punpcklbw   mm2,mm7
            punpcklwd   mm5,mm5
            punpcklwd   mm6,mm6
            psllw       mm1,
7
            psllw       mm2,
7
            pmulhw      mm1,mm5
            pmulhw      mm2,mm6
            punpcklbw   mm3,mm7
            punpcklbw   mm4,mm7
            movd        mm5,dword ptr [(offset SinXDivX_Table_MMX)
+256*4+ecx*4]
            movd        mm6,dword ptr [(offset SinXDivX_Table_MMX)
+512*4+ecx*4]
            punpcklwd   mm5,mm5
            punpcklwd   mm6,mm6
            psllw       mm3,
7
            psllw       mm4,
7
            pmulhw      mm3,mm5
            pmulhw      mm4,mm6
            paddsw      mm1,mm2
            paddsw      mm3,mm4
            movd        mm6,dword ptr [ebx] 
//v
            paddsw      mm1,mm3
            punpcklwd   mm6,mm6

            pmulhw      mm1,mm6
            add     edx,esi  
//+pic.byte_width
            paddsw      mm0,mm1

            ret
        }
    }

    inline 
void ThreeOrder_Fast_MMX(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        asm
        {
            mov     ecx,pic
            mov     eax,y_16
            mov     ebx,x_16
            movzx   edi,ah 
//v_8
            mov     edx,[ecx+TPicRegion::pdata]
            shr     eax,
16
            mov     esi,[ecx
+TPicRegion::byte_width]
            dec     eax
            movzx   ecx,bh 
//u_8
            shr     ebx,16
            imul    eax,esi
            lea     edx,[edx
+ebx*4-4]
            add     edx,eax 
//pixel

            mov     eax,ecx
            neg     ecx

            pxor    mm7,mm7  
//0
            
//mov     edx,pixel
            pxor    mm0,mm0  //result=0
            
//lea     eax,auv_7

            lea    ebx,[(offset SinXDivX_Table_MMX)
+256*4+edi*4]
            call  _private_ThreeOrder_Fast_MMX
            lea    ebx,[(offset SinXDivX_Table_MMX)
+edi*4]
            call  _private_ThreeOrder_Fast_MMX
            neg    edi
            lea    ebx,[(offset SinXDivX_Table_MMX)
+256*4+edi*4]
            call  _private_ThreeOrder_Fast_MMX
            lea    ebx,[(offset SinXDivX_Table_MMX)
+512*4+edi*4]
            call  _private_ThreeOrder_Fast_MMX

            psraw     mm0,
3
            mov       eax,result
            packuswb  mm0,mm7
            movd      [eax],mm0
            
//emms
        }
    }

    
void ThreeOrder_Border_MMX(const TPicRegion& pic,const long x_16,const long y_16,TARGB32* result)
    {
        unsigned 
long x0_sub1=(x_16>>16)-1;
        unsigned 
long y0_sub1=(y_16>>16)-1;
        
long u_16_add1=((unsigned short)(x_16))+(1<<16);
        
long v_16_add1=((unsigned short)(y_16))+(1<<16);

        TARGB32 pixel[
16];

        
for (long i=0;i<4;++i)
        {
            
long y=y0_sub1+i;
            pixel[i
*4+0]=Pixels_Bound(pic,x0_sub1  ,y);
            pixel[i
*4+1]=Pixels_Bound(pic,x0_sub1+1,y);
            pixel[i
*4+2]=Pixels_Bound(pic,x0_sub1+2,y);
            pixel[i
*4+3]=Pixels_Bound(pic,x0_sub1+3,y);
        }
        
        TPicRegion npic;
        npic.pdata     
=&pixel[0];
        npic.byte_width
=4*sizeof(TARGB32);
        
//npic.width     =4;
        
//npic.height    =4;
        ThreeOrder_Fast_MMX(npic,u_16_add1,v_16_add1,result);
    }

void PicZoom_ThreeOrder_MMX(const TPicRegion& Dst,const TPicRegion& Src)
{
    
if (  (0==Dst.width)||(0==Dst.height)
        
||(0==Src.width)||(0==Src.height)) return;

    
long xrIntFloat_16=((Src.width)<<16)/Dst.width+1
    
long yrIntFloat_16=((Src.height)<<16)/Dst.height+1;
    
const long csDErrorX=-(1<<15)+(xrIntFloat_16>>1);
    
const long csDErrorY=-(1<<15)+(yrIntFloat_16>>1);

    unsigned 
long dst_width=Dst.width;

    
//计算出需要特殊处理的边界
    long border_y0=((1<<16)-csDErrorY)/yrIntFloat_16+1;//y0+y*yr>=1; y0=csDErrorY => y>=(1-csDErrorY)/yr
    if (border_y0>=Dst.height) border_y0=Dst.height;
    
long border_x0=((1<<16)-csDErrorX)/xrIntFloat_16+1;
    
if (border_x0>=Dst.width ) border_x0=Dst.width;
    
long border_y1=(((Src.height-3)<<16)-csDErrorY)/yrIntFloat_16+1//y0+y*yr<=(height-3) => y<=(height-3-csDErrorY)/yr
    if (border_y1<border_y0) border_y1=border_y0;
    
long border_x1=(((Src.width-3)<<16)-csDErrorX)/xrIntFloat_16+1;; 
    
if (border_x1<border_x0) border_x1=border_x0;

    TARGB32
* pDstLine=Dst.pdata;
    
long srcy_16=csDErrorY;
    
long y;
    
for (y=0;y<border_y0;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            ThreeOrder_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y0;y<border_y1;++y)
    {
        
long srcx_16=csDErrorX;
        
long x;
        
for (x=0;x<border_x0;++x)
        {
            ThreeOrder_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        
for (x=border_x0;x<border_x1;++x)
        {
            ThreeOrder_Fast_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]);//fast MMX !
            srcx_16+=xrIntFloat_16;
        }
        
for (x=border_x1;x<dst_width;++x)
        {
            ThreeOrder_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]);//border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }
    
for (y=border_y1;y<Dst.height;++y)
    {
        
long srcx_16=csDErrorX;
        
for (unsigned long x=0;x<dst_width;++x)
        {
            ThreeOrder_Border_MMX(Src,srcx_16,srcy_16,
&pDstLine[x]); //border
            srcx_16+=xrIntFloat_16;
        }
        srcy_16
+=yrIntFloat_16;
        ((TUInt8
*&)pDstLine)+=Dst.byte_width;
    }

    asm emms
}

////////////////////////////////////////////////////////////////////////////////
//速度测试:
//==============================================================================
// PicZoom_ThreeOrder_MMX   34.3 fps
////////////////////////////////////////////////////////////////////////////////

 

 N:将测试结果放到一起:

////////////////////////////////////////////////////////////////////////////////
//CPU: AMD64x2 4200+(2.37G)  zoom 800*600 to 1024*768
//==============================================================================
// StretchBlt                   232.7 fps   
// PicZoom3_SSE                 711.7 fps  
// PicZoom_BilInear_MMX_Ex      157.0 fps
// 
// PicZoom_ThreeOrder0            3.6 fps
// PicZoom_ThreeOrder_Common     16.9 fps
// PicZoom_ThreeOrder_MMX        34.3 fps
////////////////////////////////////////////////////////////////////////////////

补充Intel Core2 4400上的测试成绩:

////////////////////////////////////////////////////////////////////////////////
//CPU: Intel Core2 4400(2.00G)  zoom 800*600 to 1024*768
//==============================================================================
// PicZoom3_SSE                1099.7 fps  
// PicZoom_BilInear_MMX_Ex      142.9 fps
// 
// PicZoom_ThreeOrder0            4.2 fps
// PicZoom_ThreeOrder_Common     17.6 fps
// PicZoom_ThreeOrder_MMX        34.4 fps
////////////////////////////////////////////////////////////////////////////////


转自:http://blog.csdn.net/housisong/article/details/1452249

推荐看原文链接,内中的评论很有价值