图形图像处理-之-任意角度的高质量的快速的图像旋转 上篇 纯软件的任意角度的快速旋转

FW; http://blog.csdn.net/housisong/archive/2007/04/27/1586717.aspx

图形图像处理-之-任意角度的高质量的快速的图像旋转 上篇 纯软件的任意角度的快速旋转
                               
HouSisong@GMail.com   2007.04.26

tag:图像旋转、任意角度、图像缩放、速度优化、定点数优化、近邻取样插值、二次线性插值、
   三次卷积插值、MipMap链、三次线性插值、MMX/SSE优化、CPU缓存优化

摘要:首先给出一个基本的图像旋转算法,然后一步一步的优化其速度和旋转质量,打破不能软件旋转的神话!

任意角度的高质量的快速的图像旋转 全文 分为:
     上篇 纯软件的任意角度的快速旋转
     下篇 高质量的旋转

(2007.04.29修正一个TRotaryClipData.find_begin的bug)

正文:
  为了便于讨论,这里只处理32bit的ARGB颜色;
  代码使用C++;涉及到汇编优化的时候假定为x86平台;使用的编译器为vc6;
  为了代码的可读性,没有加入异常处理代码;
   测试使用的CPU为赛扬2G;
  (一些基础代码和插值原理的详细说明参见作者的《图形图像处理-之-高质量的快速的图像缩放》系列文章)


速度测试说明:
  只测试内存数据到内存数据的缩放
  测试图片都是800*600旋转到1004*1004,测试成绩取各个旋转角度的平均速度值; fps表示每秒钟的帧数,值越大表示函数越快


A:旋转原理和旋转公式:
  推导旋转公式:
             
                       旋转示意图
   有:  tg(b)=y/x                             ----(1)
         tg(a+b)=y'/x'                         ----(2)
         x*x + y*y = x'*x' + y'*y'             ----(3)
   有公式:tg(a+b) = ( tg(a)+tg(b) ) / ( 1-tg(a)*tg(b) )  ----(4)
     把(1)代入(4)从而消除参数b;
     tg(a)+y/x = y'/x'*( 1-tg(a)*y/x )                ----(5)
     由(5)可以得x'=y'*(x-y*tg(a))/( x*tg(a)+y )       ----(6)
   把(6)代入(3)从而消除参数x',化简后求得:
     y'=x*sin(a)+y*cos(a);                     ----(7)
   把(7)代入(6),有:
     x'=x*cos(a)-y*sin(a);                     ----(8)

  OK,旋转公式有了,那么来看看在图片旋转中的应用;
  假设对图片上任意点(x,y),绕一个坐标点(rx0,ry0)逆时针旋转RotaryAngle角度后的新的坐标设为(x', y'),有公式:
  (x平移rx0,y平移ry0,角度a对应-RotaryAngle , 带入方程(7)、(8)后有: ) 
  x'= (x - rx0)*cos(RotaryAngle) + (y - ry0)*sin(RotaryAngle) + rx0 ;
  y'=-(x - rx0)*sin(RotaryAngle) + (y - ry0)*cos(RotaryAngle) + ry0 ;

那么,根据新的坐标点求源坐标点的公式为:
  x=(x'- rx0)*cos(RotaryAngle) - (y'- ry0)*sin(RotaryAngle) + rx0 ;
  y=(x'- rx0)*sin(RotaryAngle) + (y'- ry0)*cos(RotaryAngle) + ry0 ;

旋转的时候还可以顺便加入x轴和y轴的缩放和平移,而不影响速度,那么完整的公式为:           
  x=(x'- move_x-rx0)/ZoomX*cos(RotaryAngle) - (y'- move_y-ry0)/ZoomY*sin(RotaryAngle) + rx0 ;
  y=(x'- move_x-rx0)/ZoomX*sin(RotaryAngle) + (y'- move_y-ry0)/ZoomY*cos(RotaryAngle) + ry0 ;
  其中: RotaryAngle为逆时针旋转的角度;
         ZoomX,ZoomY为x轴y轴的缩放系数(支持负的系数,相当于图像翻转);
         move_x,move_y为x轴y轴的平移量;

 一些颜色和图片的数据定义:

#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];
}
// 判断一个点是否在图片中
inline  bool  PixelsIsInPic( const  TPicRegion &  pic, const   long  x, const   long  y)
{
    
return  ( (x >= 0 ) && (x < pic.width)  &&  (y >= 0 ) && (y < pic.height) );
}

 

 B:一个简单的浮点实现版本

//
//函数假设以原图片的中心点坐标为旋转和缩放的中心

void PicRotary0(const TPicRegion& Dst,const TPicRegion& Src,double RotaryAngle,double ZoomX,double ZoomY,double move_x,double  move_y)
{
    
if ( (fabs(ZoomX*Src.width)<1.0e-4|| (fabs(ZoomY*Src.height)<1.0e-4) ) return//太小的缩放比例认为已经不可见

    double rx0=Src.width*0.5;  //(rx0,ry0)为旋转中心 
    double ry0=Src.height*0.5
    
for (long y=0;y<Dst.height;++
y)
    {
        
for (long x=0;x<Dst.width;++
x)
        {
            
long srcx=(long)((x- move_x-rx0)/ZoomX*cos(RotaryAngle) - (y- move_y-ry0)/ZoomY*sin(RotaryAngle) +
 rx0) ;
            
long srcy=(long)((x- move_x-rx0)/ZoomX*sin(RotaryAngle) + (y- move_y-ry0)/ZoomY*cos(RotaryAngle) +
 ry0) ;
            
if
 (PixelsIsInPic(Src,srcx,srcy))
                Pixels(Dst,x,y)
=
Pixels(Src,srcx,srcy);
        }
    }
}

(调用方法比如:
   PicRotary0(ppicDst,ppicSrc,PI/6,0.9,0.9,(dst_wh-ppicSrc.width)*0.5,(dst_wh-ppicSrc.height)*0.5);
   //作用:将图片ppicSrc按0.9的缩放比例旋转PI/6幅度后绘制到图片ppicDst的中心
)

//注:测试图片都是800*600的图片旋转到1004*1004的图片中心 测试成绩取各个旋转角度的平均速度值

//速度测试:
//==============================================================================
// PicRotary0              12.6 fps
 

旋转结果图示(小图):

   
                
 30度                              60度                          90度

   
                 120度                             150度                         180度

  
                 210度                             240度                         270度

   
                 300度                             330度                         360度


C:优化循环内部,化简系数
  1.sin和cos函数是很慢的计算函数,可以在循环前预先计算好sin(RotaryAngle)和cos(RotaryAngle)的值:
    double sinA=sin(RotaryAngle);
    double cosA=cos(RotaryAngle);
  2.可以将除以ZoomX、ZoomY改成乘法,预先计算出倒数:
    double rZoomX=1.0/ZoomX;
    double rZoomY=1.0/ZoomY;
  3.优化内部的旋转公式,将能够预先计算的部分提到循环外(即:拆解公式):
     原:  long srcx=(long)((x- move_x-rx0)/ZoomX*cos(RotaryAngle) - (y- move_y-ry0)/ZoomY*sin(RotaryAngle) + rx0) ;
           long srcy=(long)((x- move_x-rx0)/ZoomX*sin(RotaryAngle) + (y- move_y-ry0)/ZoomY*cos(RotaryAngle) + ry0) ;
     变形为:
           long srcx=(long)( Ax*x + Bx*y +Cx ) ;
           long srcy=(long)( Ay*x + By*y +Cy ) ;
      其中: Ax=(rZoomX*cosA); Bx=(-rZoomY*sinA); Cx=(-(rx0+move_x)*rZoomX*cosA+(ry0+move_y)*rZoomY*sinA+rx0);
            Ay=(rZoomX*sinA); By=(rZoomY*cosA);  Cy=(-(rx0+move_x)*rZoomX*sinA-(ry0+move_y)*rZoomY*cosA+ry0);
      (提示: Ax,Bx,Cx,Ay,By,Cy都可以在旋转之前预先计算出来)
  改进后的函数为:

void PicRotary1(const TPicRegion& Dst,const TPicRegion& Src,double RotaryAngle,double ZoomX,double ZoomY,double move_x,double  move_y)
{
    
if ( (fabs(ZoomX*Src.width)<1.0e-4|| (fabs(ZoomY*Src.height)<1.0e-4) ) return//太小的缩放比例认为已经不可见

    double rZoomX=1.0/ ZoomX;
    
double rZoomY=1.0/
ZoomY;
    
double sinA=
sin(RotaryAngle);
    
double cosA=
cos(RotaryAngle);
    
double Ax=(rZoomX*
cosA); 
    
double Ay=(rZoomX*
sinA); 
    
double Bx=(-rZoomY*
sinA); 
    
double By=(rZoomY*
cosA); 
    
double rx0=Src.width*0.5;  //(rx0,ry0)为旋转中心 

    double ry0=Src.height*0.5
    
double Cx=(-(rx0+move_x)*rZoomX*cosA+(ry0+move_y)*rZoomY*sinA+
rx0);
    
double Cy=(-(rx0+move_x)*rZoomX*sinA-(ry0+move_y)*rZoomY*cosA+
ry0); 

    TARGB32
* pDstLine=
Dst.pdata;
    
double srcx0_f=
(Cx);
    
double srcy0_f=
(Cy);
    
for (long y=0;y<Dst.height;++
y)
    {
        
double srcx_f=
srcx0_f;
        
double srcy_f=
srcy0_f;
        
for (long x=0;x<Dst.width;++
x)
        {
            
long srcx=(long
)(srcx_f);
            
long srcy=(long
)(srcy_f);
            
if
 (PixelsIsInPic(Src,srcx,srcy))
                pDstLine[x]
=
Pixels(Src,srcx,srcy);
            srcx_f
+=
Ax;
            srcy_f
+=
Ay;
        }
        srcx0_f
+=
Bx;
        srcy0_f
+=
By;
        ((TUInt8
*&)pDstLine)+=
Dst.byte_width;
    }
}


//速度测试:
//==============================================================================
// PicRotary1               8.4 fps
 

( PicRotary1优化后速度反而更慢的说明: vc6编译器将PicRotary0优化的很好,sin/cos等重复计算都自动提到了循环外;对于PicRotary1,编译器的优化反而降低了,可能是因为有了太多的变量;经过测试用汇编重写该函数(并优化掉取整)得到的速度只比PicRotary2慢些而已 ) 

D:更深入的优化、定点数优化
  (浮点数到整数的转化也是应该优化的一个地方,这里不再处理,可以参见
    <图形图像处理-之-高质量的快速的图像缩放 上篇 近邻取样插值和其速度优化>中的PicZoom3_float函数)
  1.优化除法:
    原: double rZoomX=1.0/ZoomX;
         double rZoomY=1.0/ZoomY;
    改写为(优化掉了一次除法):
         double tmprZoomXY=1.0/(ZoomX*ZoomY); 
         double rZoomX=tmprZoomXY*ZoomY;
         double rZoomY=tmprZoomXY*ZoomX;
  2.x86的浮点计算单元(FPU)有一条指令"fsincos"可以同时计算出sin和cos值
    //定义SinCos函数: 同时计算sin(Angle)和cos(Angle)的内嵌x86汇编函数
    void __declspec(naked) __stdcall SinCos(const double Angle,double& sina,double& cosa)
    {
        asm
        {
            fld  qword ptr [esp+4]//Angle  
            mov  eax,[esp+12]//&sina
            mov  edx,[esp+16]//&cosa
            fsincos  
            fstp qword ptr [edx]  
            fstp qword ptr [eax] 
            ret 16
        }
    }
  3.用定点数代替浮点计算

void PicRotary2(const TPicRegion& Dst,const TPicRegion& Src,double RotaryAngle,double ZoomX,double ZoomY,double move_x,double  move_y)
{
    
if ( (fabs(ZoomX*Src.width)<1.0e-4|| (fabs(ZoomY*Src.height)<1.0e-4) ) return//太小的缩放比例认为已经不可见

    double tmprZoomXY=1.0/(ZoomX* ZoomY);  
    
double rZoomX=tmprZoomXY*
ZoomY;
    
double rZoomY=tmprZoomXY*
ZoomX;
    
double
 sinA,cosA;
    SinCos(RotaryAngle,sinA,cosA);
    
long Ax_16=(long)(rZoomX*cosA*(1<<16
)); 
    
long Ay_16=(long)(rZoomX*sinA*(1<<16
)); 
    
long Bx_16=(long)(-rZoomY*sinA*(1<<16
)); 
    
long By_16=(long)(rZoomY*cosA*(1<<16
)); 
    
double rx0=Src.width*0.5;  //(rx0,ry0)为旋转中心 

    double ry0=Src.height*0.5
    
long Cx_16=(long)((-(rx0+move_x)*rZoomX*cosA+(ry0+move_y)*rZoomY*sinA+rx0)*(1<<16
));
    
long Cy_16=(long)((-(rx0+move_x)*rZoomX*sinA-(ry0+move_y)*rZoomY*cosA+ry0)*(1<<16
)); 

    TARGB32
* pDstLine=
Dst.pdata;
    
long srcx0_16=
(Cx_16);
    
long srcy0_16=
(Cy_16);
    
for (long y=0;y<Dst.height;++
y)
    {
        
long srcx_16=
srcx0_16;
        
long srcy_16=
srcy0_16;
        
for (long x=0;x<Dst.width;++
x)
        {
            
long srcx=(srcx_16>>16
);
            
long srcy=(srcy_16>>16
);
            
if
 (PixelsIsInPic(Src,srcx,srcy))
                pDstLine[x]
=
Pixels(Src,srcx,srcy);
            srcx_16
+=
Ax_16;
            srcy_16
+=
Ay_16;
        }
        srcx0_16
+=
Bx_16;
        srcy0_16
+=
By_16;
        ((TUInt8
*&)pDstLine)+=
Dst.byte_width;
    }
}


//速度测试:
//==============================================================================
// PicRotary2              55.8 fps
 


E:优化内部循环的判断函数PixelsIsInPic,将判断展开到内部循环之外,跳过不需要处理的目标像素;
   TRotaryClipData类用于寻找旋转需要处理的边界范围;算法思路是首先寻找原图片中心点对应的;
那条
扫描线,然后依次向上和向下寻找边界;  如果想要更快速的边界寻找算法,可以考虑利用像素的直线
  绘制原理来自动生成边界(有机会的时候再来实现它);

                  

                       边界寻找算法图示

 

 

struct  TRotaryClipData{
public :
    
long  src_width;
    
long  src_height;
    
long  dst_width;
    
long  dst_height;
    
long  Ax_16; 
    
long  Ay_16; 
    
long  Bx_16; 
    
long  By_16; 
    
long  Cx_16;
    
long  Cy_16; 
public :
    
long  out_dst_up_y;
    
long  out_dst_down_y;

    
long  out_src_x0_16;
    
long  out_src_y0_16;
private :
    
long  out_dst_up_x0;
    
long  out_dst_up_x1;
    
long  out_dst_down_x0;
    
long  out_dst_down_x1;
public :
    inline 
long  get_up_x0(){  if  (out_dst_up_x0 < 0 return   0 ;   else   return  out_dst_up_x0; }
    inline 
long  get_up_x1(){  if  (out_dst_up_x1 >= dst_width)  return  dst_width;   else   return  out_dst_up_x1; }
    inline 
long  get_down_x0(){  if  (out_dst_down_x0 < 0 return   0 ;   else   return  out_dst_down_x0; }
    inline 
long  get_down_x1(){  if  (out_dst_down_x1 >= dst_width)  return  dst_width;   else   return  out_dst_down_x1; }

    inline 
bool  is_in_src( long  src_x_16, long  src_y_16)
    {
         
return  ( ( (src_x_16 >= 0 ) && ((src_x_16 >> 16 ) < src_width) )
               
&&  ( (src_y_16 >= 0 ) && ((src_y_16 >> 16 ) < src_height) ) );
    }
    
void  find_begin_in( long  dst_y, long &  out_dst_x, long &  src_x_16, long &  src_y_16)
    {
        src_x_16
-= Ax_16;
        src_y_16
-= Ay_16;
        
while  (is_in_src(src_x_16,src_y_16))
        {
            
-- out_dst_x;
            src_x_16
-= Ax_16;
            src_y_16
-= Ay_16;
        }
        src_x_16
+= Ax_16;
        src_y_16
+= Ay_16;
    }
    
bool  find_begin( long  dst_y, long &  out_dst_x0, long  dst_x1)
    {
        
long  test_dst_x0 = out_dst_x0 - 1 ;
        
long  src_x_16 = Ax_16 * test_dst_x0  +  Bx_16 * dst_y  +  Cx_16;
        
long  src_y_16 = Ay_16 * test_dst_x0  +  By_16 * dst_y  +  Cy_16;
        
for  ( long  i = test_dst_x0;i <= dst_x1; ++ i)
        {
            
if  (is_in_src(src_x_16,src_y_16))
            {
                out_dst_x0
= i;
                
if  (i == test_dst_x0)
                    find_begin_in(dst_y,out_dst_x0,src_x_16,src_y_16);
                
if  (out_dst_x0 < 0 )
                {
                    src_x_16
-= (Ax_16 * out_dst_x0);
                    src_y_16
-= (Ay_16 * out_dst_x0);
                }
                out_src_x0_16
= src_x_16;
                out_src_y0_16
= src_y_16;
                
return   true ;
            }
            
else
            {
                src_x_16
+= Ax_16;
                src_y_16
+= Ay_16;
            }
        }
        
return   false ;
    }
    
void  find_end( long  dst_y, long  dst_x0, long &  out_dst_x1)
    {
        
long  test_dst_x1 = out_dst_x1;
        
if  (test_dst_x1 < dst_x0) test_dst_x1 = dst_x0;

        
long  src_x_16 = Ax_16 * test_dst_x1  +  Bx_16 * dst_y  +  Cx_16;
        
long  src_y_16 = Ay_16 * test_dst_x1  +  By_16 * dst_y  +  Cy_16;
        
if  (is_in_src(src_x_16,src_y_16))
        {
            
++ test_dst_x1;
            src_x_16
+= Ax_16;
            src_y_16
+= Ay_16;
            
while  (is_in_src(src_x_16,src_y_16))
            {
                
++ test_dst_x1;
                src_x_16
+= Ax_16;
                src_y_16
+= Ay_16;
            }
            out_dst_x1
= test_dst_x1;
        }
        
else
        {
            src_x_16
-= Ax_16;
            src_y_16
-= Ay_16;
            
while  ( ! is_in_src(src_x_16,src_y_16))
            {
                
-- test_dst_x1;
                src_x_16
-= Ax_16;
                src_y_16
-= Ay_16;
            }
            out_dst_x1
= test_dst_x1;
        }
    }

    
bool  inti_clip( double  move_x, double  move_y)
    {
        
// 计算src中心点映射到dst后的坐标
        out_dst_down_y = ( long )(src_height * 0.5 + move_y);
        out_dst_down_x0
= ( long )(src_width * 0.5 + move_x);
        out_dst_down_x1
= out_dst_down_x0;
        
// 得到初始out_???
         if  (find_begin(out_dst_down_y,out_dst_down_x0,out_dst_down_x1))
            find_end(out_dst_down_y,out_dst_down_x0,out_dst_down_x1);
        out_dst_up_y
= out_dst_down_y;
        out_dst_up_x0
= out_dst_down_x0;
        out_dst_up_x1
= out_dst_down_x1;
        
return  (out_dst_down_x0 < out_dst_down_x1);
    }
    
bool  next_clip_line_down()
    {
        
++ out_dst_down_y;
        
if  ( ! find_begin(out_dst_down_y,out_dst_down_x0,out_dst_down_x1))  return   false ;
        find_end(out_dst_down_y,out_dst_down_x0,out_dst_down_x1);
        
return  (out_dst_down_x0 < out_dst_down_x1);
    }
    
bool  next_clip_line_up()
    {
        
-- out_dst_up_y;
        
if  ( ! find_begin(out_dst_up_y,out_dst_up_x0,out_dst_up_x1))  return   false ;
        find_end(out_dst_up_y,out_dst_up_x0,out_dst_up_x1);
        
return  (out_dst_up_x0 < out_dst_up_x1);
    }
};

void  PicRotary3( const  TPicRegion &  Dst, const  TPicRegion &  Src, double  RotaryAngle, double  ZoomX, double  ZoomY, double  move_x, double  move_y)
{
    
if  ( (fabs(ZoomX * Src.width) < 1.0e-4 ||  (fabs(ZoomY * Src.height) < 1.0e-4 ) )  return // 太小的缩放比例认为已经不可见
     double  tmprZoomXY = 1.0 / (ZoomX * ZoomY);  
    
double  rZoomX = tmprZoomXY * ZoomY;
    
double  rZoomY = tmprZoomXY * ZoomX;
    
double  sinA,cosA;
    SinCos(RotaryAngle,sinA,cosA);
    
long  Ax_16 = ( long )(rZoomX * cosA * ( 1 << 16 )); 
    
long  Ay_16 = ( long )(rZoomX * sinA * ( 1 << 16 )); 
    
long  Bx_16 = ( long )( - rZoomY * sinA * ( 1 << 16 )); 
    
long  By_16 = ( long )(rZoomY * cosA * ( 1 << 16 )); 
    
double  rx0 = Src.width * 0.5 ;   // (rx0,ry0)为旋转中心 
     double  ry0 = Src.height * 0.5
    
long  Cx_16 = ( long )(( - (rx0 + move_x) * rZoomX * cosA + (ry0 + move_y) * rZoomY * sinA + rx0) * ( 1 << 16 ));
    
long  Cy_16 = ( long )(( - (rx0 + move_x) * rZoomX * sinA - (ry0 + move_y) * rZoomY * cosA + ry0) * ( 1 << 16 )); 

    TRotaryClipData rcData;
    rcData.Ax_16
= Ax_16;
    rcData.Bx_16
= Bx_16;
    rcData.Cx_16
= Cx_16;
    rcData.Ay_16
= Ay_16;
    rcData.By_16
= By_16;
    rcData.Cy_16
= Cy_16;
    rcData.dst_width
= Dst.width;
    rcData.dst_height
= Dst.height;
    rcData.src_width
= Src.width;
    rcData.src_height
= Src.height;
    
if  ( ! rcData.inti_clip(move_x,move_y))  return ;

    TARGB32
*  pDstLine = Dst.pdata;
    ((TUInt8
*& )pDstLine) += (Dst.byte_width * rcData.out_dst_down_y);
    
while  ( true // to down
    {
        
long  y = rcData.out_dst_down_y;
        
if  (y >= Dst.height)  break ;
        
if  (y >= 0 )
        {
            
long  srcx_16 = rcData.out_src_x0_16;
            
long  srcy_16 = rcData.out_src_y0_16;
            
long  x1 = rcData.get_down_x1();
            
for  ( long  x = rcData.get_down_x0();x < x1; ++ x)
            {
                pDstLine[x]
= Pixels(Src,(srcx_16 >> 16 ),(srcy_16 >> 16 ));
                srcx_16
+= Ax_16;
                srcy_16
+= Ay_16;
            }
        }
        
if  ( ! rcData.next_clip_line_down())  break ;
        ((TUInt8
*& )pDstLine) += Dst.byte_width;
    }
   
    pDstLine
= Dst.pdata;
    ((TUInt8
*& )pDstLine) += (Dst.byte_width * rcData.out_dst_up_y);
    
while  (rcData.next_clip_line_up())  // to up
    {
        
long  y = rcData.out_dst_up_y;
        
if  (y < 0 break ;
        ((TUInt8
*& )pDstLine) -= Dst.byte_width;
        
if  (y < Dst.height)
        {
            
long  srcx_16 = rcData.out_src_x0_16;
            
long  srcy_16 = rcData.out_src_y0_16;
            
long  x1 = rcData.get_up_x1();
            
for  ( long  x = rcData.get_up_x0();x < x1; ++ x)
            {
                pDstLine[x]
= Pixels(Src,(srcx_16 >> 16 ),(srcy_16 >> 16 ));
                srcx_16
+= Ax_16;
                srcy_16
+= Ay_16;
            }
        }
    }
}

 
//速度测试:
//==============================================================================
// PicRotary3              74.1 fps
 


F:使用SSE的MOVNTQ指令优化CPU写缓冲

 

struct  TRotaryCopyLineData
{
  TARGB32
*  psrc;
  
long      src_byte_width;
  
long      Ax_16;
  
long      Ay_16;
  TARGB32
*  pDstLine;
  
long  xCount;
  
long  srcx0_16;
  
long  srcy0_16;
};

void  PicRotarySSE_copyLine(TRotaryCopyLineData *  rclData)
{
    
/*  //汇编实现的功能等价的对应高级语言代码示例
    for (long x=0;x<rclData->xCount;++x)
    {
        rclData->pDstLine[x]=((TARGB32*)((TUInt8*)rclData->psrc
            +(rclData->srcy0_16>>16)*rclData->src_byte_width))[rclData->srcx0_16>>16];
        rclData->srcx0_16+=rclData->Ax_16;
        rclData->srcy0_16+=rclData->Ay_16;
    }
    //
*/

    
/*  //汇编实现
    asm  
    {
            mov  esi,rclData
            mov  ecx,[esi]TRotaryCopyLineData.xCount
            push ebp
            test   ecx,ecx   
            jle    EndWrite
            mov  edi,[esi]TRotaryCopyLineData.pDstLine
            mov  edx,[esi]TRotaryCopyLineData.srcx0_16
            mov  ebx,[esi]TRotaryCopyLineData.srcy0_16

            lea  edi,[edi+ecx*4]
            neg  ecx

        loop_start:
              mov       eax,ebx
              mov       ebp,edx
              shr       eax,16
              imul      eax,[esi]TRotaryCopyLineData.src_byte_width
              add       ebx,[esi]TRotaryCopyLineData.Ay_16           
              add       edx,[esi]TRotaryCopyLineData.Ax_16           
              shr       ebp,16
              lea       eax,[eax+ebp*4]
              mov       ebp,[esi]TRotaryCopyLineData.psrc
              mov       eax,[eax+ebp]
              mov       [edi+ecx*4],eax

              inc       ecx
              jnz       loop_start


        EndWrite:
            pop ebp
    }
    //
*/

    
// * // 优化写缓冲的汇编实现
    asm
    {
            mov  esi,rclData
            mov  ecx,[esi]TRotaryCopyLineData.xCount
            push ebp
            test   ecx,ecx   
            jle    EndWrite
            mov  edi,[esi]TRotaryCopyLineData.pDstLine
            mov  edx,[esi]TRotaryCopyLineData.srcx0_16
            mov  ebx,[esi]TRotaryCopyLineData.srcy0_16

            push ecx
            and  ecx, (not 
1 // 循环2次展开

            test   ecx,ecx   
// nop
            jle    EndWriteLoop2

            lea  edi,[edi
+ ecx * 4 ]
            neg  ecx

        loop2_start:
              mov       eax,ebx
              mov       ebp,edx
              shr       eax,
16
              imul      eax,[esi]TRotaryCopyLineData.src_byte_width
              add       ebx,[esi]TRotaryCopyLineData.Ay_16           
              add       edx,[esi]TRotaryCopyLineData.Ax_16           
              shr       ebp,
16
              lea       eax,[eax
+ ebp * 4 ]
              mov       ebp,[esi]TRotaryCopyLineData.psrc
              MOVD          mm0,[eax
+ ebp]

              mov       eax,ebx
              mov       ebp,edx
              shr       eax,
16
              imul      eax,[esi]TRotaryCopyLineData.src_byte_width
              add       ebx,[esi]TRotaryCopyLineData.Ay_16           
              add       edx,[esi]TRotaryCopyLineData.Ax_16           
              shr       ebp,
16
              lea       eax,[eax
+ ebp * 4 ]
              mov       ebp,[esi]TRotaryCopyLineData.psrc
              PUNPCKlDQ     mm0,[eax
+ ebp]

              
//  MOVNTQ qword ptr [edi+ecx*4], mm0   // 不使用缓存的写入指令
              asm _emit  0x0F  asm _emit  0xE7  asm _emit  0x04  asm _emit  0x8F   

              add       ecx,
2
              jnz       loop2_start

              
// sfence  // 刷新写入
              
// asm _emit 0x0F asm _emit 0xAE asm _emit 0xF8  
            emms

        EndWriteLoop2:

            pop    ecx
            and    ecx,
1   
            test   ecx,ecx
            jle    EndWrite

              mov       eax,ebx
              shr       eax,
16
              imul      eax,[esi]TRotaryCopyLineData.src_byte_width
              shr       edx,
16
              lea       eax,[eax
+ edx * 4 ]
              mov       edx,[esi]TRotaryCopyLineData.psrc
              mov       eax,[eax
+ edx]
              mov       [edi],eax
        EndWrite:

            pop ebp

    }
    
// */
}


void  PicRotarySSE( const  TPicRegion &  Dst, const  TPicRegion &  Src, double  RotaryAngle, double  ZoomX, double  ZoomY, double  move_x, double  move_y)
{
    
if  ( (fabs(ZoomX * Src.width) < 1.0e-4 ||  (fabs(ZoomY * Src.height) < 1.0e-4 ) )  return // 太小的缩放比例认为已经不可见
     double  tmprZoomXY = 1.0 / (ZoomX * ZoomY);  
    
double  rZoomX = tmprZoomXY * ZoomY;
    
double  rZoomY = tmprZoomXY * ZoomX;
    
double  sinA,cosA;
    SinCos(RotaryAngle,sinA,cosA);
    
long  Ax_16 = ( long )(rZoomX * cosA * ( 1 << 16 )); 
    
long  Ay_16 = ( long )(rZoomX * sinA * ( 1 << 16 )); 
    
long  Bx_16 = ( long )( - rZoomY * sinA * ( 1 << 16 )); 
    
long  By_16 = ( long )(rZoomY * cosA * ( 1 << 16 )); 
    
double  rx0 = Src.width * 0.5 ;   // (rx0,ry0)为旋转中心 
     double  ry0 = Src.height * 0.5
    
long  Cx_16 = ( long )(( - (rx0 + move_x) * rZoomX * cosA + (ry0 + move_y) * rZoomY * sinA + rx0) * ( 1 << 16 ));
    
long  Cy_16 = ( long )(( - (rx0 + move_x) * rZoomX * sinA - (ry0 + move_y) * rZoomY * cosA + ry0) * ( 1 << 16 )); 

    TRotaryClipData rcData;
    rcData.Ax_16
= Ax_16;
    rcData.Bx_16
= Bx_16;
    rcData.Cx_16
= Cx_16;
    rcData.Ay_16
= Ay_16;
    rcData.By_16
= By_16;
    rcData.Cy_16
= Cy_16;
    rcData.dst_width
= Dst.width;
    rcData.dst_height
= Dst.height;
    rcData.src_width
= Src.width;
    rcData.src_height
= Src.height;
    
if  ( ! rcData.inti_clip(move_x,move_y))  return ;

    TRotaryCopyLineData rclData;
    rclData.Ax_16
= rcData.Ax_16;
    rclData.Ay_16
= rcData.Ay_16;
    rclData.psrc
= Src.pdata;
    rclData.src_byte_width
= Src.byte_width;

    TARGB32
*  pDstLine = Dst.pdata;
    ((TUInt8
*& )pDstLine) += (Dst.byte_width * rcData.out_dst_down_y);
    
while  ( true // to down
    {
        
long  y = rcData.out_dst_down_y;
        
if  (y >= Dst.height)  break ;
        
if  (y >= 0 )
        {
            rclData.srcx0_16
= rcData.out_src_x0_16;
            rclData.srcy0_16
= rcData.out_src_y0_16;
            
long  x0 = rcData.get_down_x0();
            rclData.pDstLine
=& pDstLine[x0];
            rclData.xCount
= rcData.get_down_x1() - x0;
            PicRotarySSE_copyLine(
& rclData);
        }
        
if  ( ! rcData.next_clip_line_down())  break ;
        ((TUInt8
*& )pDstLine) += Dst.byte_width;
    }

    pDstLine
= Dst.pdata;
    ((TUInt8
*& )pDstLine) += (Dst.byte_width * rcData.out_dst_up_y);
    
while  (rcData.next_clip_line_up())  // to up
    {
        
long  y = rcData.out_dst_up_y;
        
if  (y < 0 break ;
        ((TUInt8
*& )pDstLine) -= Dst.byte_width;
        
if  (y < Dst.height)
        {
            rclData.srcx0_16
= rcData.out_src_x0_16;
            rclData.srcy0_16
= rcData.out_src_y0_16;
            
long  x0 = rcData.get_up_x0();
            rclData.pDstLine
=& pDstLine[x0];
            rclData.xCount
= rcData.get_up_x1() - x0;
            PicRotarySSE_copyLine(
& rclData);
        }
    }
}

 
//速度测试:
//==============================================================================
// PicRotarySEE            80.8 fps

一张效果图:
 //程序使用的调用参数:
    const long testcount=2000;
    long dst_wh=1004;
    for (int i=0;i<testcount;++i)
    {
        double zoom=rand()*(1.0/RAND_MAX)+0.5;
        PicRotarySSE(ppicDst,ppicSrc,rand()*(PI*2/RAND_MAX),zoom,zoom,((dst_wh+ppicSrc.width)*rand()*(1.0/RAND_MAX)-ppicSrc.width),(dst_wh+ppicSrc.height)*rand()*(1.0/RAND_MAX)-ppicSrc.height);
    }

 //ps:如果很多时候源图片绘制时可能落在目标区域的外面,那么需要写一个剪切算法快速排除不必要的绘制

 

一张测试函数速度的时候生成的图像:


G:旋转测试的结果放到一起:

//注:测试图片都是800*600的图片旋转到1004*1004的图片中心,测试成绩取各个旋转角度的平均速度值 

//速度测试:
//==============================================================================
// PicRotary0              12.6 fps
// PicRotary1               8.4 fps
// PicRotary2              55.8 fps
// PicRotary3              74.1 fps
// PicRotarySEE            80.8 fps

//ps:文章的下篇将进一步优化图片旋转的质量(使用二次线性插值、三次卷积插值和MipMap链),并完美的处理边缘的锯齿,并考虑介绍颜色的Alpha Blend混合


(希望读者能在这一系列的文章中不仅能学到旋转和缩放,还能够学到一些优化的基本技巧和思路;也欢迎指出文章中的错误、我没有做到的优化、改进意见 等)

 
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值