您的位置：首页 > 其它

图像旋转的原理，实现与优化

2016-09-11 15:19 429 查看

图像旋转的原理

图像旋转的原理其实很简单，为了简化公式的推导，这里我们假设绕原点(0,0)旋转。

本文规定逆时针旋转为正，当然你也可以规定顺时针为正（此时，正角度就表示顺时针旋转）

很容易得到三个方程：

tan(θ+α)=y′x′(1)

tanα=yx(2)

x2+y2=x′2+y′2(3)

由和角公式可知

tan(θ+α)=tanθ+tanα1−tanθtanα

代入（1）和（2）可得

y′=x′(y+xtanθ)x−ytanθ(4)

(4)代入(3)得到：

x′2=(xcosθ−ysinθ)2(5)

（5）代入（3）

y′2=(xsinθ+ycosθ)2

所以

{x′=xcosθ+y(−sinθ)y′=xsinθ+ycosθ(6)

我们知道仿射变换矩阵可以表示为：

A=[cosθsinθ−sinθcosθ](7)

公式7与公式6在形式上是一致的(这里没有加入平移和缩放)

注意公式6中x,y,x′,y′实际表示x−0,y−0,x′−0,y′−0即与旋转中心横坐标以及纵坐标之间的距离，现在我们将旋转中心变为更具有一般性的(x0,y0),得旋转公式为：

{x′=(x−x0)cosθ+(y−y0)(−sinθ)+x0y′=(x−x0)sinθ+(y−y0)cosθ+y0(8)

对于图像来说，规定顺时针旋转角度为正，则旋转公式与上述公式一致，由于任意角度的旋转后会出现像素坐标为负的情况，如果不将旋转后的图像坐标平移，会缺失部分图像，所以，对于旋转后的图像，我们通常会加入一个平移

{x′=(x−x0)cosθ+(y−y0)(−sinθ)+x0+Δxy′=(x−x0)sinθ+(y−y0)cosθ+y0+Δy(9)

图像中使用后向映射，使用x′,y′表示x,y,则上述公式变为

{x=((x′−x0−Δx)cosθ+(y′−y0−Δy)sinθ)/scale+x0y=((x′−x0−Δx)(−sinθ)+(y′−y0−Δy)cosθ)/scale+y0(10)

这里加入了缩放，其中：scale>1表示原图放大 <1表示原图缩小

有了上面的公式，下面就可以写出相应的代码了

注：为什么公式推导的时候选用右手系？因为右手系与图像坐标系推导出来的公式在形式上是统一的，而右手系又是我们熟悉的坐标系。

图像旋转的实现

双线性插值

下面我们对他进行优化，首先我们将最近邻插值改用更为通用的双线性插值，代码如下

void Rotate_Bilinear(const Mat &srcImage, Mat &dstImage, Point center, Size dstSize, double theta, double scale)
{
CV_Assert(srcImage.depth() == CV_8U);
dstImage.create(dstSize, srcImage.type());
dstImage.setTo(Scalar(0, 0, 0));

int x0 = center.x;
int y0 = center.y;
theta = DEGREE2RADIAN(theta);

//
Mat extendedImage;
copyMakeBorder(srcImage, extendedImage, 1, 1, 1, 1, BORDER_CONSTANT,Scalar(0,0,0)); // 使用0填充边界

// dx,dy就是dst与src图像中心的距离
int dx = dstImage.cols / 2 - srcImage.cols / 2;
int dy = dstImage.rows / 2 - srcImage.rows / 2;
int numberOfChannels = srcImage.channels();

int widthOfDst = dstImage.cols;
int heightOfDst = dstImage.rows;

for (int y = 0; y <= heightOfDst - 1; ++y)
{
for (int x = 0; x <= widthOfDst - 1; ++x)
{
// 按照原来的方式计算原图坐标
float srcX = ((x - x0 - dx)*cos(theta) + (y - y0 - dy)*sin(theta)) / scale + x0;
float srcY = ((x0 + dx - x)*sin(theta) + (y - y0 - dy)*cos(theta)) / scale + y0;

// 加1,得到在extendedImage中的坐标
srcX++;
srcY++;

// get the nearest coordinate of src
int x1 = (int)(srcX);
int y1 = (int)(srcY);

// 浮点转化为整数
int dx1 = (srcX - x1)*(1<< SHIFT);
int dy1 = (srcY - y1)*(1<< SHIFT);

if (numberOfChannels == 1)
{
// !！！注意这里的范围，在extendedImage中，原图的范围就是1~cols - 2了
if ((x1 >= 1 && x1 <= extendedImage.cols - 2) && (y1 >= 1 && y1 <= extendedImage.rows - 2))
{
//双线性插值
//周围4个点
//a就是最近邻像素
//a   b
//  p
//c   d
uchar a = extendedImage.at<uchar>(y1, x1);
uchar b = extendedImage.at<uchar>(y1, x1 + 1);
uchar c = extendedImage.at<uchar>(y1 + 1, x1);
uchar d = extendedImage.at<uchar>(y1 + 1, x1 + 1);

//int p = (a*((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b*dx1*((1 << SHIFT) - dy1) + c*((1 << SHIFT) - dx1)*dy1 + d*dx1*dy1)/(1<<(2* SHIFT));

int p = a*((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b*dx1*((1 << SHIFT) - dy1) + c*((1 << SHIFT) - dx1)*dy1 + d*dx1*dy1;
p = DESCALE(p, 2*SHIFT);

dstImage.at<uchar>(y, x) = p;
}
else
{
// 越界赋值0
dstImage.at<uchar>(y, x) = 0;
}
}
else
{
if ((x1 >= 1 && x1 <= extendedImage.cols - 2) && (y1 >= 1 && y1 <= extendedImage.rows - 2))
{
//双线性插值
//周围4个点
//a就是最近邻像素
//a   b
//  p
//c   d
Vec3b a = extendedImage.at<Vec3b>(y1, x1);
Vec3b b = extendedImage.at<Vec3b>(y1, x1 + 1);
Vec3b c = extendedImage.at<Vec3b>(y1 + 1, x1);
Vec3b d = extendedImage.at<Vec3b>(y1 + 1, x1 + 1);

/*int p1 = (a[0] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[0] * dx1*((1 << SHIFT) - dy1) + c[0] * ((1 << SHIFT) - dx1)*dy1 + d[0] * dx1*dy1)/(1<<(2*SHIFT));
int p2 = (a[1] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[1] * dx1*((1 << SHIFT) - dy1) + c[1] * ((1 << SHIFT) - dx1)*dy1 + d[1] * dx1*dy1)/ (1 << (2 * SHIFT));
int p3 = (a[2] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[2] * dx1*((1 << SHIFT) - dy1) + c[2] * ((1 << SHIFT) - dx1)*dy1 + d[2] * dx1*dy1)/ (1 << (2 * SHIFT));*/

int p1 = a[0]*((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[0]*dx1*((1 << SHIFT) - dy1) + c[0]*((1 << SHIFT) - dx1)*dy1 + d[0]*dx1*dy1;
p1 = DESCALE(p1, 2 * SHIFT);
int p2 = a[1]*((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[1]*dx1*((1 << SHIFT) - dy1) + c[1]*((1 << SHIFT) - dx1)*dy1 + d[1] *dx1*dy1;
p2 = DESCALE(p2, 2 * SHIFT);
int p3 = a[2]*((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[2]*dx1*((1 << SHIFT) - dy1) + c[2]*((1 << SHIFT) - dx1)*dy1 + d[2] *dx1*dy1;
p3 = DESCALE(p3, 2 * SHIFT);

dstImage.at<cv::Vec3b>(y, x) = Vec3b(p1,p2,p3);
}
else
{
dstImage.at<cv::Vec3b>(y, x) = cv::Vec3b(0, 0, 0);
}

}
}
}
}

测试方法与上述相同，测试代码如下：

Rotate_Bilinear(srcImage, dstImage, Point(srcImage.cols / 2, srcImage.rows / 2), Size(2500, 2500),30.0, 2);

基本旋转效果都是一样的，下面我们看局部放大图

这是最近邻的眼部区域

双线性在保持图像细节方面要好于最近邻，而且不容易产生锯齿

双线性的优化

上面的双线性插值速度比较慢，我们对其进行优化，主要优化有：

1.将循环内部不变量提取出来

由于sin和cos的计算是比较慢的函数，所以可以将他们提前计算好，而不是每次循环都要计算

如：

double sinTheta = sin(theta);

double cosTheta = cos(theta);

2.改变循环体内循环变量自增的方式

采用加法自增的方式，而不是每次都要计算乘法，详见代码

注：浮点数转化为整数的优化，上面已经涉及，这里不再说明

优化后的代码如下：

void Rotate_Bilinear2(const Mat &srcImage, Mat &dstImage, Point center, Size dstSize, double theta, double scale)
{
CV_Assert(srcImage.depth() == CV_8U);
dstImage.create(dstSize, srcImage.type());
dstImage.setTo(Scalar(0, 0, 0));

int x0 = center.x;
int y0 = center.y;
theta = DEGREE2RADIAN(theta);

//
Mat extendedImage;
copyMakeBorder(srcImage, extendedImage, 1, 1, 1, 1, BORDER_CONSTANT, Scalar(0, 0, 0)); // 使用0填充边界

// dx,dy就是dst与src图像中心的距离
int dx = dstImage.cols / 2 - srcImage.cols / 2;
int dy = dstImage.rows / 2 - srcImage.rows / 2;
int numberOfChannels = srcImage.channels();

int widthOfDst = dstImage.cols;
int heightOfDst = dstImage.rows;

////////////////////////////////// 优化部分/////////////////////////////
// 将循环内的不变量提取出来
double sinTheta = sin(theta);
double cosTheta = cos(theta);
scale = 1.0 / scale;

// 改变了循环内部增量的方式
double temp1= (0 - y0 - dy)*sinTheta;
double temp2 = (0 - y0 - dy)*cosTheta;
double dtemp1 = sinTheta;
double dtemp2 = cosTheta;

for (int y = 0; y <= heightOfDst - 1; ++y,temp1+=dtemp1,temp2+=dtemp2)
{
// 改变了循环内部增量的方式
double temp3= ((0 - x0 - dx)*cosTheta + temp1)*scale + x0;
double temp4= (-(0 - x0 - dx)*sinTheta + temp2)*scale + y0;
double dtemp3 = (cosTheta)*scale;
double dtemp4= (-sinTheta)*scale;
for (int x = 0; x <= widthOfDst - 1; ++x,temp3+=dtemp3,temp4+=dtemp4)
{
// 计算原图坐标
double srcX = temp3;
double srcY = temp4;

// 加1,得到在extendedImage中的坐标
srcX++;
srcY++;

// get the nearest coordinate of src
int x1 = (int)(srcX);
int y1 = (int)(srcY);

// 浮点转化为整数
int dx1 = (srcX - x1)*(1 << SHIFT);
int dy1 = (srcY - y1)*(1 << SHIFT);

if (numberOfChannels == 1)
{
// !！！注意这里的范围，在extendedImage中，原图的范围就是1~cols - 2了
if ((x1 >= 1 && x1 <= extendedImage.cols - 2) && (y1 >= 1 && y1 <= extendedImage.rows - 2))
{
//双线性插值
//周围4个点
//a就是最近邻像素
//a   b
//  p
//c   d
uchar a = extendedImage.at<uchar>(y1, x1);
uchar b = extendedImage.at<uchar>(y1, x1 + 1);
uchar c = extendedImage.at<uchar>(y1 + 1, x1);
uchar d = extendedImage.at<uchar>(y1 + 1, x1 + 1);

//int p = (a*((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b*dx1*((1 << SHIFT) - dy1) + c*((1 << SHIFT) - dx1)*dy1 + d*dx1*dy1)/(1<<(2* SHIFT));

int p = a*((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b*dx1*((1 << SHIFT) - dy1) + c*((1 << SHIFT) - dx1)*dy1 + d*dx1*dy1;
p = DESCALE(p, 2 * SHIFT);

dstImage.at<uchar>(y, x) = p;
}
else
{
// 越界赋值0
dstImage.at<uchar>(y, x) = 0;
}
}
else
{
if ((x1 >= 1 && x1 <= extendedImage.cols - 2) && (y1 >= 1 && y1 <= extendedImage.rows - 2))
{
//双线性插值
//周围4个点
//a就是最近邻像素
//a   b
//  p
//c   d
Vec3b a = extendedImage.at<Vec3b>(y1, x1);
Vec3b b = extendedImage.at<Vec3b>(y1, x1 + 1);
Vec3b c = extendedImage.at<Vec3b>(y1 + 1, x1);
Vec3b d = extendedImage.at<Vec3b>(y1 + 1, x1 + 1);

/*int p1 = (a[0] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[0] * dx1*((1 << SHIFT) - dy1) + c[0] * ((1 << SHIFT) - dx1)*dy1 + d[0] * dx1*dy1)/(1<<(2*SHIFT));
int p2 = (a[1] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[1] * dx1*((1 << SHIFT) - dy1) + c[1] * ((1 << SHIFT) - dx1)*dy1 + d[1] * dx1*dy1)/ (1 << (2 * SHIFT));
int p3 = (a[2] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[2] * dx1*((1 << SHIFT) - dy1) + c[2] * ((1 << SHIFT) - dx1)*dy1 + d[2] * dx1*dy1)/ (1 << (2 * SHIFT));*/

int p1 = a[0] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[0] * dx1*((1 << SHIFT) - dy1) + c[0] * ((1 << SHIFT) - dx1)*dy1 + d[0] * dx1*dy1;
p1 = DESCALE(p1, 2 * SHIFT);
int p2 = a[1] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[1] * dx1*((1 << SHIFT) - dy1) + c[1] * ((1 << SHIFT) - dx1)*dy1 + d[1] * dx1*dy1;
p2 = DESCALE(p2, 2 * SHIFT);
int p3 = a[2] * ((1 << SHIFT) - dx1)*((1 << SHIFT) - dy1) + b[2] * dx1*((1 << SHIFT) - dy1) + c[2] * ((1 << SHIFT) - dx1)*dy1 + d[2] * dx1*dy1;
p3 = DESCALE(p3, 2 * SHIFT);

dstImage.at<cv::Vec3b>(y, x) = Vec3b(p1, p2, p3);
}
else
{
dstImage.at<cv::Vec3b>(y, x) = cv::Vec3b(0, 0, 0);
}

}
}
}
}

下面我们测试一下他们的速度

测试环境：Intel core i5-6200U,12G,测试图像都是1000*580到2500*2500（三通道彩色图），放大2倍


测试方法	Rotate_Nearest	Rotate_Bilinear	Rotate_Bilinear2
速度(单位ms)	58.2	144.9	78.7

可以看出，优化后的代码速度提升还是很明显的。

双线性代码还可以进一步的优化，有兴趣的朋友可以自己实现。

2016-9-11 15:12:24

Last Updated:2017-3-11 23:22:16

非常感谢您的阅读，如果您觉得这篇文章对您有帮助，欢迎扫码进行赞赏。

内容来自用户分享和网络整理，不保证内容的准确性，如有侵权内容，可联系管理员处理

标签：

相关文章推荐

新的分享

章节导航