转自:http://zeze0556.tk/2010/01/15/i-mode图像描绘的仿射变换的C++实现
曾经做项目的时候,用到了i-mode仿射变换的函数,对其算法寻找花费了不少时间,都没有找到最相像的。后来花费了一个周末,根据api提供的公式,重新复习下矩阵的计算,终于写了一个基本上完全一样效果的。看来,还是老老实实一些比较好写。
<!--more-->
<pre cpp="brush:cpp">
#define EPSILON_E5 (int)(((int64)(1E-5))* (1<<FIX_POINT))
int fabs( int a )
{
if( a < 0)
return -a;
return a;
}
///
void Mat_Add_3X3(MATRIX3X3_PTR ma,
MATRIX3X3_PTR mb,
MATRIX3X3_PTR msum)
{
// this function adds two 3x3 matrices together and
// and stores the result
for (int row=0; row<3; row++)
{
for (int col=0; col<3; col++)
{
// insert resulting row,col element
msum->M[row][col] = ma->M[row][col] + mb->M[row][col];
} // end for col
} // end for row
} // end Mat_Add_3X3
void Mat_Mul_VECTOR3D_3X3(VECTOR3D_PTR va,
MATRIX3X3_PTR mb,
VECTOR3D_PTR vprod)
{
// this function multiplies a VECTOR3D against a
// 3x3 matrix - ma*mb and stores the result in mprod
for (int col=0; col < 3; col++)
{
// compute dot product from row of ma
// and column of mb
int sum = 0; // used to hold result
for (int row=0; row<3; row++)
{
// add in next product pair
sum+=(va->M[row]*mb->M[col][row]/*>>FIX_POINT*/);
} // end for index
// insert resulting col element
vprod->M[col] = sum;
} // end for col
} // end Mat_Mul_VECTOR3D_3X3
void Mat_Init_3X3(MATRIX3X3_PTR ma,
int m00, int m01, int m02,
int m10, int m11, int m12,
int m20, int m21, int m22)
{
// this function fills a 3x3 matrix with the sent data in row major form
ma->M[0][0] = m00; ma->M[0][1] = m01; ma->M[0][2] = m02;
ma->M[1][0] = m10; ma->M[1][1] = m11; ma->M[1][2] = m12;
ma->M[2][0] = m20; ma->M[2][1] = m21; ma->M[2][2] = m22;
} // end Mat_Init_3X3
/
int Mat_Inverse_3X3(MATRIX3X3_PTR m, MATRIX3X3_PTR mi)
{
// this function computes the inverse of a 3x3
// first compute the determinate to see if there is
// an inverse
int det = ((int64)m->M[0][0]*(((int64)m->M[1][1]*m->M[2][2]>>FIX_POINT) - ((int64)m->M[2][1]*m->M[1][2]>>FIX_POINT))>>FIX_POINT) -
((int64)m->M[0][1]*(((int64)m->M[1][0]*m->M[2][2]>>FIX_POINT) - ((int64)m->M[2][0]*m->M[1][2]>>FIX_POINT))>>FIX_POINT) +
((int64)m->M[0][2]*(((int64)m->M[1][0]*m->M[2][1]>>FIX_POINT) - ((int64)m->M[2][0]*m->M[1][1]>>FIX_POINT))>>FIX_POINT);
if (fabs(det) <= EPSILON_E5)
return(0);
// compute inverse to save divides
int det_inv = ((1<<FIX_POINT)*(1<<FIX_POINT)/det);//<<FIX_POINT;
// compute inverse using m-1 = adjoint(m)/det(m)
mi->M[0][0] = ((int64)det_inv*(((int64)m->M[1][1]*m->M[2][2]>>FIX_POINT) - ((int64)m->M[2][1]*m->M[1][2]>>FIX_POINT)))>>FIX_POINT;
mi->M[1][0] = -((int64)det_inv*(((int64)m->M[1][0]*m->M[2][2]>>FIX_POINT) - ((int64)m->M[2][0]*m->M[1][2]>>FIX_POINT)))>>FIX_POINT;
mi->M[2][0] = (int64)det_inv*(((int64)m->M[1][0]*m->M[2][1]>>FIX_POINT) - ((int64)m->M[2][0]*m->M[1][1]>>FIX_POINT))>>FIX_POINT;
mi->M[0][1] = -((int64)det_inv*(((int64)m->M[0][1]*m->M[2][2]>>FIX_POINT) - ((int64)m->M[2][1]*m->M[0][2]>>FIX_POINT)))>>FIX_POINT;
mi->M[1][1] = (int64)det_inv*(((int64)m->M[0][0]*m->M[2][2]>>FIX_POINT) - ((int64)m->M[2][0]*m->M[0][2]>>FIX_POINT))>>FIX_POINT;
mi->M[2][1] = -((int64)det_inv*(((int64)m->M[0][0]*m->M[2][1]>>FIX_POINT) - ((int64)m->M[2][0]*m->M[0][1]>>FIX_POINT)))>>FIX_POINT;
mi->M[0][2] = (int64)det_inv*(((int64)m->M[0][1]*m->M[1][2]>>FIX_POINT) - ((int64)m->M[1][1]*m->M[0][2]>>FIX_POINT))>>FIX_POINT;
mi->M[1][2] = -((int64)det_inv*(((int64)m->M[0][0]*m->M[1][2]>>FIX_POINT) - ((int64)m->M[1][0]*m->M[0][2]>>FIX_POINT)))>>FIX_POINT;
mi->M[2][2] = (int64)det_inv*(((int64)m->M[0][0]*m->M[1][1]>>FIX_POINT) - ((int64)m->M[1][0]*m->M[0][1]>>FIX_POINT))>>FIX_POINT;
// return success
return(1);
} // end Mat_Inverse_3X3
/
/
int Mat_Det_3X3(MATRIX3X3_PTR m)
{
// computes the determinate of a 3x3 matrix using co-factor
// expansion
return((m->M[0][0]*((m->M[1][1]*m->M[2][2]>>FIX_POINT) - (m->M[2][1]*m->M[1][2]>>FIX_POINT))>>FIX_POINT) -
(m->M[0][1]*((m->M[1][0]*m->M[2][2]>>FIX_POINT) - (m->M[2][0]*m->M[1][2]>>FIX_POINT))>>FIX_POINT) +
(m->M[0][2]*((m->M[1][0]*m->M[2][1]>>FIX_POINT) - (m->M[2][0]*m->M[1][1]>>FIX_POINT))>>FIX_POINT) );
} // end Mat_Det_3X3
</pre>
上面其实是要用到的算法,我记得是从《游戏编程大师》中摘抄的吧,中间为了计算正确,有一个最小值的判断,不然,可能会发生除0错误(<strong>1/det</strong>)。
<pre class="brush:cpp">
void drawImage(Image *img, int *buf, int ax, int ay, int aw, int ah)
{
if( img == 0 )
return;
MATRIX3X3 Morg;
MATRIX3X3 Mdes;
Mat_Init_3X3( &Morg, buf[0], buf[1], buf[2] + (transX<<FIX_POINT), buf[3], buf[4], buf[5] + (transY<<FIX_POINT), 0,0, (1<<FIX_POINT) );
Mat_Inverse_3X3(&Morg, &Mdes);
VECTOR3D org_left_top;
org_left_top.M[0] = 0;
org_left_top.M[1] = 0;
org_left_top.M[2] = 1;
VECTOR3D org_left_bottom;
org_left_bottom.M[0] = 0;
org_left_bottom.M[1] = (ah);
org_left_bottom.M[2] = 1;
VECTOR3D org_right_top;
org_right_top.M[0] = aw;
org_right_top.M[1] = 0;
org_right_top.M[2] = 1;
VECTOR3D org_right_bottom;
org_right_bottom.M[0] = aw;
org_right_bottom.M[1] = ah;
org_right_bottom.M[2] = 1;
VECTOR3D des_left_top, des_left_bottom, des_right_top, des_right_bottom;
Mat_Mul_VECTOR3D_3X3( &org_left_top, &Morg, &des_left_top);
Mat_Mul_VECTOR3D_3X3( &org_left_bottom, &Morg, &des_left_bottom);
Mat_Mul_VECTOR3D_3X3( &org_right_top, &Morg, &des_right_top);
Mat_Mul_VECTOR3D_3X3( &org_right_bottom, &Morg, &des_right_bottom);
int minx = ((int)(min(min(des_left_top.M[0], des_right_top.M[0]), min(des_left_bottom.M[0], des_right_bottom.M[0])))>>FIX_POINT);
int miny = ((int)(min(min(des_left_top.M[1], des_right_top.M[1]), min(des_left_bottom.M[1], des_right_bottom.M[1])))>>FIX_POINT);
int maxx = ((int)(max(max(des_left_top.M[0], des_right_top.M[0]), max(des_left_bottom.M[0], des_right_bottom.M[0])))>>FIX_POINT);
int maxy = ((int)(max(max(des_left_top.M[1], des_right_top.M[1]), max(des_left_bottom.M[1], des_right_bottom.M[1])))>>FIX_POINT);
if(minx >= clip[0] + clip[2] || miny >= clip[1] + clip[3] || maxx <= clip[0] || maxy <= clip[1] )
return;
if( minx < clip[0] )
minx = clip[0];
if( maxx > clip[0] + clip[2] )
maxx = clip[0] + clip[2];
if( miny < clip[1] )
miny = clip[1];
if(maxy > clip[1] + clip[3] )
maxy = clip[1] + clip[3];
IBitmap *pback = IDISPLAY_GetDestination(m_pDisplay);
IDIB* pDib = NULL;
IBITMAP_QueryInterface( pback, AEECLSID_DIB, (void**)&pDib);
IBITMAP_Release(pback);
int* tablex_x = (int*)static_tablex_x;
int* tablex_y = (int*)static_tablex_y;
int* tabley_x = (int*)static_tabley_x;
int* tabley_y = (int*)static_tabley_y;
for( int j = minx; j < maxx; j++)
{
tablex_x[j-minx] = (j*Mdes.M[0][0]+Mdes.M[0][2])>>FIX_POINT;
}
for( int j = minx; j < maxx; j++)
{
tablex_y[j-minx] = (j*Mdes.M[1][0]+Mdes.M[1][2])>>FIX_POINT;
}
for( int j = miny ; j < maxy; j++)
{
tabley_x[j-miny] = (j*Mdes.M[0][1])>>FIX_POINT;
}
for( int j = miny; j < maxy; j++)
{
tabley_y[j-miny] = (j*Mdes.M[1][1])>>FIX_POINT;
}
switch (img->m_pIDIB->nDepth) {
case 16:
{
unsigned short *pdes = (unsigned short*)((unsigned char*)(pDib->pBmp) + pDib->nPitch * miny + (minx << 1));
unsigned short *psrc = (unsigned short*)((unsigned char*)(img->m_pIDIB->pBmp) + img->m_pIDIB->nPitch *(ay)+ ((ax) << 1));
int desstep = (pDib->nPitch >> 1) - (maxx - minx);
int stepsrc = (img->m_pIDIB->nPitch>>1);
int *ptablex_x = tablex_x;
int *ptabley_x = tabley_x;
int *ptablex_y = tablex_y;
int *ptabley_y = tabley_y;
for( int i = miny; i < maxy; i++, pdes += desstep, ptabley_x++, ptabley_y++, ptablex_x = tablex_x, ptablex_y = tablex_y )
for( int j = minx; j < maxx; j++, pdes++, ptablex_x++, ptablex_y++ )
{
int x = (*ptablex_x+*ptabley_x);
int y = (*ptablex_y+*ptabley_y);;
if( x < 0 || y < 0 || x >= aw || y >= ah )
continue;
unsigned short *ptsrc = psrc + y * stepsrc + x;
if( *ptsrc != img->m_pIDIB->ncTransparent )
*pdes = *ptsrc;
}
}
break;
}
IDIB_Release( pDib );
return;
}
</pre>
基本上已经完全优化了,值得注意的是,我是根据公式反过来计算的,即扫描要填充的区域,计算在原图中的位置,如果在原图上,则填充,如果不在原图上,则跳过。这样做的好处就是防止由于误差引起的漏点。旋转也可以采用类似的方法,改天有时间的话,把代码找到放出来。上述函数的使用办法,请参考相关的i-mode的代码。