[数字图像处理]简单的几何学图像变换与图像配准

【修订版】本文已整理到我的个人博客订正了原文的几处错误旋转公式里的减号应为加号原文与代码矛盾坐标记号与倾斜系数的命名前后不一示意图也已重新绘制。CSDN 这一版保留原样之后的勘误与更新只在博客进行 Geometric Transform : 几何学变换与图像配准 | Finns blog1.图像的几何学变换之前的博文里我简单的介绍了图像的放大与缩小。放大与缩小也算是图像的几何学变换本文介绍了其他的几何学变换包括旋转、水平倾斜和垂直倾斜当然还有水平移动与垂直移动。这些变换很简单不需要插值所以这里就不着重介绍了。假设输入图像为G(u,v)其变换后的图像为F(x,y)。其变化的方法如下所示。图像的几何学变换主要有两种向前映射与向后映射。1.1向前映射所谓向前映射就是从输入图像为g(u,v)的(0,0)点开始将g(u,v)遍历一遍依次计算g(u,v)变换后的坐标。当然计算出来的坐标不会是整数很大程度上不会是整数就像下图所示那样。借由上图我们重新理解一些向前映射。假设图像的一部分这5个点经过变换我们得到了右边的5个点。其实这个5个点经过变换之后计算出来的坐标并不是整数。假设现在遍历到这个点图中以蓝色表示通过计算得到的点是。我们将的坐标进行四舍五入将的值赋值给这种处理相当于选择最近的点属于最初级的插值方法。使用上述的向前映射的思想将图像旋转我们能得到如下效果。可以看出来所得到的结果非常“斑驳”。其原因是我们将g(u,v)进行变换计算之后所得到的坐标四舍五入之后有的点没有被赋值而有点被赋值了多次。这就产生了“空穴”所以让变换后的图像非常斑驳。我看了一些文章基本到这里都会说“由于向前映射会产生空穴所以向前映射一般不使用”。诸如此类的话其实这样理解不太对产生空穴的根本原因是由于插值法的不恰当所产生的。使用向前映射也可以不产生“空穴”所使用的方法如下所示。我们可以通过最接近的四个点去确定的值。这样的话就可以不产生空穴。但是这个方法实现起来很困难所以向前映射才一般不被使用。1.2 向后映射向后映射就是将输出图像f(x,y)遍历一遍然后计算输出(x,y)的时候所需要g(x,y)的坐标。当然这个数不一定是整数。为了方便理解还是看下图。假设我们遍历到通过计算我们得到这样一个点。也就意味着如果需要确定的灰度值那么我们需要点处的灰度值。这样几何变换的问题又归结到了插值问题。最方便的办法就是选择最近的点入上图的例子的灰度值等于的灰度值。当然还有更加“高级”点的方法精度也还算可以。的灰度值的确定方法如下所示。如上图我们使用四个点去确定的值。首先与之间使用线性插值确定出的值。同样的方法确定出的值。与之间进行线性插值就可以得到的值了。图画的非常清楚了使用向后映射将图像旋转我们可以得到如下结果。1.3 上述两部分的Matlab代码close all; clear all; clc; %% -----------Geometric_spatial_transform------------------ f imread(./letter_T.tif); f mat2gray(f,[0 255]); [M,N] size(f); %-----by forward mapping-----% g_fm zeros(M,N); seta -pi/8; for v (-M/2):1:(M/2)-1 for w (-N/2):1:(N/2)-1 x round(v*cos(seta) - w*sin(seta)); y round(v*sin(seta) w*cos(seta)); if (((y(-N/2))(y(N/2)))((x(-M/2))(x(M/2)))) g_fm(x(M/2)1,y(N/2)1) f(v(M/2)1,w(N/2)1); end end end figure(); subplot(1,2,1); imshow(f,[0 1]); xlabel(a).Original Image); subplot(1,2,2); imshow(g_fm,[0 1]); xlabel(b).Ruselt of Geometric spatial transform); %-----by inverse mapping---- g_im zeros(M,N); seta -pi/8; for x (-M/2):1:(M/2)-1 for y (-N/2):1:(N/2)-1 v x*cos(-seta) - y*sin(-seta); w x*sin(-seta) y*cos(-seta); if (((w(-N/2)1)(w(N/2)-1))((v(-M/2)1)(v(M/2)-1))) %g_im(x(M/2)1,y(N/2)1) 1;%f(round(v(M/2)1),round(w(N/2)1)); Q_11 f(floor(v)(M/2)1,floor(w)(N/2)1); Q_21 f(floor(v)(M/2)1,ceil(w)(N/2)1); Q_12 f(ceil(v)(M/2)1,floor(w)(N/2)1); Q_22 f(ceil(v)(M/2)1,ceil(w)(N/2)1); R1 (Q_21 - Q_11)*(w-floor(w)) Q_11; R2 (Q_22 - Q_12)*(w-floor(w)) Q_12; g_im(x(M/2)1,y(N/2)1) (R2-R1)*(v-floor(v)) R1; end end end figure(); subplot(1,2,1); imshow(f,[0 1]); xlabel(a).Original Image); subplot(1,2,2); imshow(g_im,[0 1]); xlabel(b).Ruselt of Geometric spatial transform);2.图像的配准图像的配准常常用于超解析领域。当然了作为基础学习博主没学那么深。这次的图像配准试图去还原被倾斜变换的图像。首先先将图像进行两个方向的倾斜。我们这次的目的是将图像b).还原为a).。这次变换相对简单我们使用如下模型去拟合变换关系。我们在原图与变换后的图像上选择出四个标准点然后带入方程并求出系数。使用这个方程将图像b).还原为a)。这个过程比较简单下面是结果。另外所选择的标准点我已经在上图中用红圈标定出来了。可以看到还原的效果非常的好额其实这也是应为畸变比较简单的缘故。为了看出于原图的区别我们做出了差分图像。可以看出来还原不是完美的。至此上个代码作为本文的结尾。额由于特征点的选择是手工的这个代码可能没有多少意义。不要吐槽最后这一句话。一般的在使用特征点做图象配准的时候一般选择在某个实际的物体上放入某个实际的标志然后通过检测这个实际的标志去进行变换。【这里可以参考《Digital Image Processing》 Rafael C. Gonzalez / Richard E. Woods的第二章的相关内容这里有叙述的】close all; clear all; clc; %% --------------image registration--------------------------- f_Original imread(./characters_test_pattern.tif); f_Original mat2gray(f_Original,[0 255]); [M,N] size(f_Original); f zeros(M4,N4); for x 1:M f(x2,:) [0 0 f_Original(x,:) 0 0]; end [M,N] size(f); s_v 0.3; s_w 0.02; P round(Ms_v*N); Q round(Ns_w*M); g_1 zeros(P,N); for x (-P/2):1:(P/2)-1 for y (-N/2):1:(N/2)-1 v x - s_v * y; w y; if (((w(-N/2))(w(N/2)-1))((v(-M/2))(v(M/2)-1))) Q_11 f(floor(v)(M/2)1,floor(w)(N/2)1); Q_21 f(floor(v)(M/2)1,ceil(w)(N/2)1); Q_12 f(ceil(v)(M/2)1,floor(w)(N/2)1); Q_22 f(ceil(v)(M/2)1,ceil(w)(N/2)1); R1 (Q_21 - Q_11)*(w-floor(w)) Q_11; R2 (Q_22 - Q_12)*(w-floor(w)) Q_12; g_1(x(P/2)1,y(N/2)1) (R2-R1)*(v-floor(v)) R1; end end end g zeros(P,Q); for x (-P/2):1:(P/2)-1 for y (-Q/2):1:(Q/2)-1 v x; w y - s_w * x; if (((w(-N/2))(w(N/2)-1))((v(-P/2))(v(P/2)-1))) Q_11 g_1(floor(v)(P/2)1,floor(w)(N/2)1); Q_21 g_1(floor(v)(P/2)1,ceil(w)(N/2)1); Q_12 g_1(ceil(v)(P/2)1,floor(w)(N/2)1); Q_22 g_1(ceil(v)(P/2)1,ceil(w)(N/2)1); R1 (Q_21 - Q_11)*(w-floor(w)) Q_11; R2 (Q_22 - Q_12)*(w-floor(w)) Q_12; g(x(P/2)1,y(Q/2)1) (R2-R1)*(v-floor(v)) R1; end end end figure(); subplot(1,2,1); imshow(f,[0 1]); xlabel(a).Original Image); hold on; plot(621, 78,ro,MarkerSize,7); plot(116,113,ro,MarkerSize,7); plot( 85,649,ro,MarkerSize,7); plot(624,641,ro,MarkerSize,7); subplot(1,2,2); imshow(g,[0 1]); xlabel(b).Ruselt of Geometric spatial transform); hold on; plot(625,264,ro,MarkerSize,7); plot(116,147,ro,MarkerSize,7); plot( 96,674,ro,MarkerSize,7); plot(638,828,ro,MarkerSize,7); %% f_reg zeros(M,N); Original [621-(N/2) 78-(M/2) (621-(N/2))*(78-(M/2)) 1; 116-(N/2) 113-(M/2) (116-(N/2))*(113-(M/2)) 1; 85-(N/2) 649-(M/2) (85-(N/2))*(649-(M/2)) 1; 624-(N/2) 641-(M/2) (624-(N/2))*(641-(M/2)) 1]; Ruselt [625-(Q/2);116-(Q/2);96-(Q/2);638-(Q/2)]; c (Original)\Ruselt; %C1~C4 Ruselt [264-(P/2);147-(P/2);674-(P/2);828-(P/2)]; c [c (Original)\Ruselt]; %C5~C8 wwwwww for y (-M/2):1:(M/2)-1 for x (-N/2):1:(N/2)-1 w c(1,1)*y c(2,1)*x c(3,1)*x*y c(4,1); v c(1,2)*y c(2,2)*x c(3,2)*x*y c(4,2); if (((w(-Q/2))(w(Q/2)-1))((v(-P/2))(v(P/2)-1))) Q_11 g(floor(v)(P/2)1,floor(w)(Q/2)1); Q_21 g(floor(v)(P/2)1,ceil(w)(Q/2)1); Q_12 g(ceil(v)(P/2)1,floor(w)(Q/2)1); Q_22 g(ceil(v)(P/2)1,ceil(w)(Q/2)1); R1 (Q_21 - Q_11)*(w-floor(w)) Q_11; R2 (Q_22 - Q_12)*(w-floor(w)) Q_12; f_reg(x(N/2)1,y(M/2)1) (R2-R1)*(v-floor(v)) R1; %g(round(v(P/2)1),round(w(Q/2)1)) 0.5; end end end g_diff abs(f - f_reg); figure(); subplot(1,2,1); imshow(f_reg,[0 1]); xlabel(c).Ruselt of image registration); subplot(1,2,2); imshow(g_diff,[0 1]); xlabel(d).Difference image);原文发于博客 zhoufan900428-CSDN博客