几何空间变换和图像配准
几何空间变换又称为橡皮膜变换，因为他可以看做是在一幅橡皮膜上印制图像，然后根据一定规则拉伸橡皮膜。由两个基本操作组成：
1）坐标的空间变换
2）灰度内插
最常用的是仿射变换一般形式如下：
[x,y,1] = [v,w,1]*T
[t11 t12 0]
= [v,w,1]*[t21 t22 0]
[t31 t32 0]
这一变换可以根据矩阵T中元素的值，对一组坐标点做尺度变换，旋转，平移，偏移。
我们有两种方法使用该式子，第一种是前向映射，它由扫描输入图像的像素，并在每一个位置（v,w）直接计算输出图像中像素的空间位置（x,y）组成。
第二种方法是反向映射，扫描输出图像的位置，并在每一个位置（x,y）使用(v,w) = T-1(x,y)反向计算输入图像相应的位置，然后进行灰度内插（一般使用最邻近，双线性，双三次内插）。
前向映射一个问题是，输入图像的的两个或多个像素会映射到输出图像的同一个位置，另一种可能是某些输出位置完全没有像素，因此反向映射比前向映射更加有效，被很多的商业公司选用，包括mathworks的 matlab使用的也是反向映射。
一、图像平移
x = v + tx
y = w + ty
原图：
平移后的图像：
二、图像旋转
x = v*cosA - w * sinA
y = v *sinA + w * cosA
进行图像旋转时，需要注意的是当图像的大小发生改变时此时若以左上角点为原点，那么原点会发生改变，一般会选定图像的中心作为原点，这样在图像进行坐标变换时中心不变。但在在变换时要记得加减一个常数，使原点相对应。
旋转后的图像：
观察旋转后的图像，会有很明显的锯齿现象。
若采用双线性插值或者三次立方插值那么锯齿依然难以消除。
三、图像的尺度变换
x = cx * v
y = cy * w
这个就不附上图了，比较简单
四、图像的偏移
垂直偏移变换：
x =v*sv + w
y = w
水平偏移变换：
x = v
y = sh*v+w
水平偏移后的图像：
仿射变换的源代码：
#include<cv.h>
#include<highgui.h>
#include <cmath>
#include <stdlib.h>
const double PI = 3.141592653;
int main(){
IplImage * image,*image2;
image = cvLoadImage("E:\\image\\base.jpg",0);
cvNamedWindow("image",CV_WINDOW_AUTOSIZE);
//cvSaveImage("E:\\image\\pollen.jpg",image,0);
cvShowImage("image",image);
//cvWaitKey(0);
//平移 dx dy
int dx = 10;
int dy = 20;
image2 = cvCreateImage(cvSize(image->width + dx, image->height + dy),image->depth,1);
unsigned char * ptr,*dst;
int i,j ,m,n;//i j为新图中坐标，m n 为原图中坐标
for(i = 0 ; i < image2->height;i++){
for(j = 0; j< image2->width;j++){
dst = (unsigned char *)image2->imageData + i*image2->widthStep + j;
n = j - dx;
m = i - dy;
if(m < 0 || n < 0 || m >= image->height || n >= image->width){
*dst = 0;
}
else {
ptr = (unsigned char * )(image->imageData + m * image->width + n);
*dst = *ptr;
}
}
}
cvNamedWindow("image2",1);
cvShowImage("image2",image2);
//cvWaitKey(0);
cvSaveImage("E:\\image\\basepingyi.jpg",image,0);
/**************************************************************************************/
//旋转图像 这里以逆时针为正 插值选用双线性插值 围绕左上角点进行旋转 双线性插值
IplImage* image3;
double angle = 15;
double Hangle = angle*PI/180*(-1);
int Width,Height;
//double r,s;
double r,s;
int temp[4];
double z1,z2;
//原点会发生移动
double fx,fy;
double Hcos = cos(Hangle);
double Hsin = sin(Hangle);
//x y 存放原图中四个点的位置，以中心为原点 左上角开始 顺时针数
int x[4];
int y[4];
int x1[4];
int y1[4];
x[0] = -(image->width-1)/2;
x[1] = -x[0];
x[2] = -x[0];
x[3] = x[0];
y[0] = -(image->height -1)/2;
y[1] = y[0];
y[2] = -y[0];
y[3] = -y[0];
//x1 y1 分别存放新图中图像的四个角点的位置 以中心为原点 左上角点开始 顺时针数
for (i = 0; i < 4 ; i++)
{
x1[i] = (int)(x[i]*Hcos - y[i]*Hsin + 0.5);
y1[i] = (int)(x[i]*Hsin + y[i]*Hcos + 0.5);
printf("x: %d \n",x[i] );
printf("y: %d\n",y[i]);
printf("x1: %d \n",x1[i] );
printf("y1: %d\n",y1[i]);
}
//确定新图像的长宽
if (abs(y1[2] - y1[0]) > abs(y1[3] - y1[1]))
{
Height = abs(y1[2] - y1[0]);
Width = abs(x1[3] - x1[1]);
}
else{
Height = abs(y1[3] - y1[1]);
Width = abs( x1[2] - x1[0]);
}
image3 = cvCreateImage(cvSize(Width,Height),image->depth,1);
//两个偏移常量
fx = -1*(Width -1)* Hcos*0.5 -(Height -1)*Hsin*0.5 + (image->width-1)/2;
fy = (Width -1)*Hsin*0.5 -(Height -1)*Hcos*0.5 + (image->height -1)/2;
for (i = 0 ; i< Height ; i++)
{
for (j =0 ; j< Width; j++)
{
dst = (unsigned char *)(image3->imageData + i*image3->widthStep + j);
r = (int)(i*Hcos- j * Hsin + fy + 0.5);
s = (int)(i*Hsin + j*Hcos + fx + 0.5 );
if (r <0 || s < 0 || r >= image->height || s >= image->width)
{
*dst = 0;
}
else{
ptr = (unsigned char *)(image->imageData + image->widthStep*(int)r + (int)s );
temp[0] = *ptr;
temp[1] = *(ptr+1);
ptr = (unsigned char*)(image->imageData + image->widthStep *(int)(r+1) + int (s));
temp[2] = *ptr;
temp[3] = *(ptr+1);
z1 = (temp[1] - temp[0])*(s - int(s)) + temp[0];
z2 = (temp[3] -temp[2])*(s- int(s)) + temp[2];
*dst = (int)(z1 + (z2 -z1)*(r- (int)r));
}
}
}
cvSaveImage("bldpingyi.bmp",image2);
/**************************************************************************************/
/**************************************************************************************/
/* //这里是采用最简单的取离映射点最近的那个点 另外 的方法是采用双线性插值的方式
//旋转图像 这里以逆时针为正 插值选用双线性插值 围绕左上角点进行旋转
IplImage* image3;
double angle = 15;
double Hangle = angle*PI/180*(-1);
int Width,Height;
//double r,s;
int r,s;
//原点会发生移动
double fx,fy;
double Hcos = cos(Hangle);
double Hsin = sin(Hangle);
//x y 存放原图中四个点的位置，以中心为原点 左上角开始 顺时针数
int x[4];
int y[4];
int x1[4];
int y1[4];
x[0] = -(image->width-1)/2;
x[1] = -x[0];
x[2] = -x[0];
x[3] = x[0];
y[0] = -(image->height -1)/2;
y[1] = y[0];
y[2] = -y[0];
y[3] = -y[0];
//x1 y1 分别存放新图中图像的四个角点的位置 以中心为原点 左上角点开始 顺时针数
for (i = 0; i < 4 ; i++)
{
x1[i] = (int)(x[i]*Hcos - y[i]*Hsin + 0.5);
y1[i] = (int)(x[i]*Hsin + y[i]*Hcos + 0.5);
printf("x: %d \n",x[i] );
printf("y: %d\n",y[i]);
printf("x1: %d \n",x1[i] );
printf("y1: %d\n",y1[i]);
}
//确定新图像的长宽
if (abs(y1[2] - y1[0]) > abs(y1[3] - y1[1]))
{
Height = abs(y1[2] - y1[0]);
Width = abs(x1[3] - x1[1]);
}
else{
Height = abs(y1[3] - y1[1]);
Width = abs( x1[2] - x1[0]);
}
image3 = cvCreateImage(cvSize(Width,Height),image->depth,1);
//两个偏移常量
fx = -1*(Width -1)* Hcos*0.5 -(Height -1)*Hsin*0.5 + (image->width-1)/2;
fy = (Width -1)*Hsin*0.5 -(Height -1)*Hcos*0.5 + (image->height -1)/2;
for (i = 0 ; i< Height ; i++)
{
for (j =0 ; j< Width; j++)
{
dst = (unsigned char *)(image3->imageData + i*image3->widthStep + j);
r = (int)(i*Hcos- j * Hsin + fy + 0.5);
s = (int)(i*Hsin + j*Hcos + fx + 0.5 );
if (r <0 || s < 0 || r >= image->height || s >= image->width)
{
*dst = 0;
}
else{
ptr = (unsigned char *)(image->imageData + image->widthStep*r + s );
*dst = *ptr;
}
}
}
*/
/**************************************************************************************/
cvNamedWindow("image3",1);
cvShowImage("image3",image3);
cvSaveImage("bldxuanzhuan2.bmp",image3);
//cvWaitKey(0);
//垂直偏移 x = v*Sv + w y = w
//水平偏移 x = v y = Sh*v + w
//这里实现水平偏移
double sh = 0.8;
//确定新图的长和宽
Width = (int)(image->width + sh* image->height);
Height = image->height;
int t;
IplImage * image4 = cvCreateImage(cvSize(Width , Height),image->depth,1);
for (i = 0 ; i< Height ; i++)
{
for (j =0 ; j< Width; j++)
{
dst = (unsigned char *)(image4->imageData + i*image4->widthStep + j);
//r = (int)(i*Hcos- j * Hsin + fy + 0.5);
//s = (int)(i*Hsin + j*Hcos + fx + 0.5 );
t = i;
s = j - sh*i;
if (t <0 || s < 0 || t >= image->height || s >= image->width)
{
*dst = 0;
}
else{
ptr = (unsigned char *)(image->imageData + image->widthStep*t + (int)s );
z1 = *ptr;
z2 = *(ptr+1);
*dst = (int)(z1 + (z2 -z1)*(s- (int)s));
}
}
}
cvNamedWindow("image4",1);
cvShowImage("image4",image4);
cvSaveImage("bldpianyi.bmp",image4);
cvWaitKey(0);
cvReleaseImage(&image);
cvReleaseImage(&image2);
cvReleaseImage(&image3);
cvDestroyAllWindows();
return 0;
}
//图像的仿射变换和图像配准 //平移、缩放、旋转、偏移 //若想变换之后的图像时可以恢复的，那么在进行变换的时候图像的大小应该可以变化，是为了保留更多的信息 //使用反向映射的方式
五、图片配准
图像配准，一般用于对齐两幅或多幅相同场景的图像。前面已经讲了所有的仿射变换的函数形式，更为复杂的是多种仿射变化的连续作用，显然可以得出，多种仿射变换的连续作用可以用如下式子概括：
x =c1*v + c2*w + c3*v*w + c4
y =c5*v + c6*w + c7*v*w + c8
方程组含有8个未知量，假设我们已知一对图像中相对应的4个点的坐标，那个既可以解出8个未知量，从而实现两张图片的配准。
使用photoshop打开一种图片，放上至少4个小方块，能够肉眼识别，但又不要太大，太大会使误差增大，对图像进行一定的仿射变换，并记下四对匹配点的坐标。
这里
image1中取四个点
47 44
433 39
47 430
462 444
image2中取四个点
140 40
461 106
55 305
395 390
采用matlab解矩阵，求出8个参数，然后根据该参数进行图像配准：
matlab程序：************************************
% 图像配准 首先是解方程组获取参数
d = ones(4,4);
d(1:2,:) = [47 433 47 462;
44, 39 ,430 ,444];
d(3,:) = d(1,:).*d(2,:)
x= [140 461 55 395]
C1 = x/d
%得 C1 = 0.8290 -0.2200 -0.0000 110.7284
%事实上这个图像当时因为只是进行了缩放和旋转，因此不存在偏移，C（3） = 0是合理的
y = [40 ;106 ;305 ;390];
C2 = y/d
%C2 = 0.1797 0.6863 0.0000 1.3466
%理由同上，事实上这个图像当时因为只是进行了缩放和旋转，因此不存在偏移，C（3） = 0是合理的
图像配准程序，也是基于反向映射
/*******************************************************************
已知两幅图像上的四个对应点，很据变换式，采用matlab求出二者之间的变换方程
根据变换式，假如需要将图片2变换到图片1，那么需要求出图片1到2 的变换式，反向映射，
假如映射区域在2区域之内，那么赋予该处的值，否者需要赋予0
************************************************************************/
//中值滤波和均值滤波
#include<cv.h>
#include<highgui.h>
int main(){
IplImage * image,*image2,*image3;
image = cvLoadImage("E:\\image\\base.jpg",0);//以灰度图像的形式读入图片
image2 = cvLoadImage("E:\\image\\base2.jpg",0);
cvNamedWindow("image",CV_WINDOW_AUTOSIZE);
cvNamedWindow("image2",CV_WINDOW_AUTOSIZE);
cvNamedWindow("image3",CV_WINDOW_AUTOSIZE);
cvShowImage("image",image);
image3 = cvCreateImage(cvGetSize(image),image->depth,1);
//cvWaitKey(0);
unsigned char * ptr,*dst;
int i,j;
double m, n;
int temp[4];
double C[8] = {0.8290, -0.2200, -0.0000, 110.7284,0.1797,0.6863,0.0000,1.3466};//图片1到2 的变换矩阵
double Z1,Z2;
for( i = 0 ; i < image3->height;i++){
for( j = 0; j< image3->width;j++){
dst = (unsigned char *)image3->imageData+ i*image3->widthStep+ j;
n = C[0]*j + C[1]*i + C[2]*i*j + C[3];
m = C[4]*j + C[5]*i + C[6]*i*j + C[7];
if (m < 0.0 || n < 0.0 || m >= image2->height || n >= image2->width)
{
*dst = 0;
}
else{
ptr = (unsigned char *)(image2->imageData + image2->widthStep * (int)m + (int)n);
temp[0] = *ptr;
temp[1] = *(ptr+1);
ptr = (unsigned char *)(image2->imageData + image2->widthStep * (int)(m +1) + (int)n);
temp[2] = *ptr;
temp[3] = *(ptr+1);
Z1 = (temp[1] - temp[0])*(m - int(m)) + temp[0];
Z2 = (temp[3] -temp[2])*(m- int(m)) + temp[2];
*dst = (int)(Z1 + (Z2 -Z1)*(n- (int)n));
}
}
}
cvShowImage("image2",image2);
cvShowImage("image3",image3);
cvSaveImage("bldpeizhun.bmp",image3);
//比较配准效果
IplImage* image4;
image4 = cvCreateImage(cvSize(image->width , image->height),image->depth , 1);
cvSub(image , image3 , image4,0);
cvNamedWindow("sub");
cvShowImage("sub", image4);
cvSaveImage("bldsub.bmp",image4);
cvWaitKey(0);
return 0;
}配准的两个源图像：
图1：
图2
配准之后的图像：
配准图像和原图像的差值：
