目录

  • 图像平移、旋转、缩放、镜像的MATLAB实现(仿照MATLAB内置函数实现)
  • 原理
  • 代码
  • 平移
  • 镜像
  • 缩放与旋转
  • 辅助工具
  • 图像预处理
  • 最近邻元法
  • 双线性内插法
  • 三次内插法
  • 旋转
  • 缩放
  • 提升运行效率的一些尝试
  • 实验结果



文中附上的代码可能有缺失,建议通过该GitHub仓库访问代码。

图像平移、旋转、缩放、镜像的MATLAB实现(仿照MATLAB内置函数实现)

(《数字图像处理》课程实验3)
参考MATLAB内置im系列函数实现。本次实验结果上没有做到与内置函数100%一致,但相差无几,应该是取整原则带来的差异。旋转与缩放运行效率相对较低,待改进。

原理

(待填坑)

代码

平移

function RESULT=myshift(IMAGE,X,Y)
%MYSHIFT Shift image.
%   RESULT = myshift(IMAGE, X, Y) shifts image IMAGE by X pixels
%   horizontally to the right and Y pixels vertically down. To shift image
%   IMAGE horizontally to the left or vertically up, specify a negative
%   value for X or Y.
%
%   Example
%   -------
%   Shift 100 pixels horizontally to the right and 100 pixels vertically
%   down.
%
%       filename='color.jpg';
%       img=imread(filename);
%       x=100;y=100;
%       res=myshift(img,x,y);
%       imshow(res);
%       storagepath='results/shift';
%       mkdir(storagepath);
%       imwrite(res,[storagepath,'/',filename]);
[rows,cols,channels]=size(IMAGE);
RESULT=uint8(zeros(rows,cols,channels));
 
if X>0
    xstart=1;
    xend=cols-X;
else
    xstart=-X+1;
    xend=cols;
end
 
if Y>0
    ystart=1;
    yend=rows-Y;
else
    ystart=-Y+1;
    yend=rows;
end
 
RESULT(ystart+Y:yend+Y,xstart+X:xend+X,:)=IMAGE(ystart:yend,xstart:xend,:);

镜像

function RESULT=myflip(IMAGE,DIM)
%MYFLIP Flip image.
%   RESULT = myflip(IMAGE, DIM) flips image IMAGE horizontally by
%   specifying 2 fod DIM or vertically by specifying 1 for DIM.
%
%   Example
%   -------
%   Flip vertically.
%
%       filename='color.jpg';
%       img=imread(filename);
%       dim=1;
%       res=myflip(img,dim);
%       imshow(res);
%       mode={'vertical','horizontal'};
%       storagepath=['results/flip/',mode{dim}];
%       mkdir(storagepath);
%       imwrite(res,[storagepath,'/',filename]);
 
[rows,cols,channels]=size(IMAGE);
RESULT=uint8(zeros(rows,cols,channels));
switch DIM
    case 1
        for i=1:rows
            RESULT(i,:,:)=IMAGE(rows-i+1,:,:);
        end
    case 2
        for j=1:cols
            RESULT(:,j,:)=IMAGE(:,cols-j+1,:);
        end
    otherwise
        warning("Only horizontal flip or vertical flip is supported.");
end

缩放与旋转

辅助工具
图像预处理
function RESULT=preprocessimage(IMAGE,N)
%PREPROCESSIMAGE Enlarge the image by N-circle pixels.
%   RESULT=preprocessimage(IMAGE,N) enlarges the IMAGE by N-circle pixels.
%   For example, if the size of IMAGE is (98,98) and N is 1, the size of
%   RESULT will be (100, 100). The new added pixel will use the pixel value
%   of the pixel of the original image closest to it as its own value.
%
%   Note
%   ----
%   PREPROCESSIMAGE is only a precautionary measure to prevent
%   cross-borders and should not be used for output images.
[rows,cols,channels]=size(IMAGE);
RESULT=uint8(zeros(rows+2*N,cols+2*N,channels));
for k=1:channels
    % process original image
    RESULT(N+1:N+rows,N+1:N+cols,k)=IMAGE(:,:,k);
    % process four lines
    RESULT(1:N,N+1:N+cols,k)=repmat(IMAGE(1,:,k),N,1);
    RESULT(N+rows+1:2*N+rows,N+1:N+cols,k)=repmat(IMAGE(rows,:,k),N,1);
    RESULT(N+1:N+rows,1:N,k)=repmat(IMAGE(:,1,k),1,N);
    RESULT(N+1:N+rows,N+cols+1:2*N+cols,k)=repmat(IMAGE(:,cols,k),1,N);
    % process four corners
    RESULT(1:N,1:N,k)=IMAGE(1,1,k);
    RESULT(1:N,N+cols+1:2*N+cols,k)=IMAGE(1,cols,k);
    RESULT(N+rows+1:2*N+rows,1:N,k)=IMAGE(rows,1,k);
    RESULT(N+rows+1:2*N+rows,N+cols+1:2*N+cols,k)=IMAGE(rows,cols,k);
End
最近邻元法
function PIXELVALUE=nearestinterpolation(IMAGE,X,Y,K)
%NEARESTINTERPOLATION Implement the nearest-neighbor interpolation.
%   PIXELVALUE=nearestinterpolation(IMAGE,X,Y,K) calculates the
%   interpolation value of point(X, Y) from IMAGE on the K channel, using 
%   the nearest-neighbor interpolation.
%
%   Note
%   ----
%   Refer to function MYRESIZE or function MYROTATE to learn how to use
%   this function.
PIXELVALUE=IMAGE(round(Y),round(X),K);
双线性内插法
function PIXELVALUE=bilinearinterpolation(IMAGE,X,Y,K)
%BILINEARINTERPOLATION Implement the bilinear interpolation.
%   PIXELVALUE = bilinearinterpolation(IMAGE, X, Y, K) calculates the
%   interpolation value of point(X, Y) from IMAGE on the K channel, using 
%   the bilinear interpolation.
%
%   Note
%   ----
%   Refer to function MYRESIZE or function MYROTATE to learn how to use
%   this function.
ii=floor(X);
jj=floor(Y);
f1=((ii+1)-X)*IMAGE(jj,ii,K)+(X-ii)*IMAGE(jj,ii+1,K);
f2=((ii+1)-X)*IMAGE(jj+1,ii,K)+(X-ii)*IMAGE(jj+1,ii+1,K);
f12=((jj+1)-Y)*f1+(Y-jj)*f2;
PIXELVALUE=f12;
三次内插法
function PIXELVALUE=bicubicinterpolation(IMAGE,X,Y,K)
%BICUBICINTERPOLATION Implement the bicubic interpolation.
%   PIXELVALUE = bicubicinterpolation(IMAGE, X, Y, K) calculates the
%   interpolation value of point(X, Y) from IMAGE on the K channel, using 
%   the bicubicinterpolation interpolation.
%
%   Note
%   ----
%   Refer to function MYRESIZE or function MYROTATE to learn how to use
%   this function.
i=floor(X);
j=floor(Y);
u=X-i;
v=Y-j;
A=[S(1+v),S(v),S(1-v),S(2-v)];
B=double([
    IMAGE(j-1,i-1,K),IMAGE(j-1,i,K),IMAGE(j-1,i+1,K),IMAGE(j-1,i+2,K);
    IMAGE(j,i-1,K),IMAGE(j,i,K),IMAGE(j,i+1,K),IMAGE(j,i+2,K);
    IMAGE(j+1,i-1,K),IMAGE(j+1,i,K),IMAGE(j+1,i+1,K),IMAGE(j+1,i+2,K);
    IMAGE(j+2,i-1,K),IMAGE(j+2,i,K),IMAGE(j+2,i+1,K),IMAGE(j+2,i+2,K);
    ]);
C=[S(1+u),S(u),S(1-u),S(2-u)]';
PIXELVALUE=uint8(A*B*C);
 
% S(x) for bicubic interpolation
function res=S(x)
x=abs(x);
if 0<=x && x<1
    res=1-2*x^2+x^3;
elseif 1<=x && x<2
    res=4-8*x+5*x^2-x^3;
elseif x>=2
    res=0;
end
旋转
function RESULT=myrotate(IMAGE,ANGLE,METHOD)
%MYROTATE Rotate image.
%   RESULT = myrotate(IMAGE, ANGLE, METHOD) rotates image IMAGE by ANGLE 
%   degrees in a clockwise direction around its center point.  To rotate 
%   the image counterclockwise, specify a negative value for ANGLE.  
%   MYROTATE makes the output image RESULT large enough to contain the 
%   entire rotated image. MYROTATE sets the values of pixels in RESULT that 
%   are outside the rotated image to 0 (zero).
%
%   MYROTATE rotates image IMAGE, using the interpolation method specified
%   by METHOD. METHOD is a string that can have one of the following
%   values. 
%
%        'nearest'    Nearest neighbor interpolation
%
%        'bilinear'   Bilinear interpolation
%
%        'bicubic'    Bicubic interpolation.
%
%   Note that you must specify METHOD argument, because there is
%   no default value for METHOD.
%
%   Example
%   -------
%   Rotates the color.jpg by 45 degrees in a clockwise direction.
%
%       filename='color.jpg';
%       img=imread(filename);
%       res=myrotate(img,45,'bilinear');
%       imshow(res);
%       storagepath='results/rotate';
%       mkdir(storagepath);
%       imwrite(res,[storagepath,'/',filename]);
[M,N,channels]=size(IMAGE);
theta=mod(ANGLE,360)/180*pi;
% preprocess the input image according to the used method
% and get the required method
switch METHOD
    case 'nearest'
        n=1;
        interpolation=@nearestinterpolation;
    case 'bilinear'
        n=1;
        interpolation=@bilinearinterpolation;
    otherwise 
        n=2;
        interpolation=@bicubicinterpolation;
end
IMAGE=preprocessimage(IMAGE,n);
 
rotationmatrix=[cos(theta),-sin(theta);sin(theta),cos(theta)];
 
northwest=[1,1]*rotationmatrix;
northeast=[1,N]*rotationmatrix;
southwest=[M,1]*rotationmatrix;
southeast=[M,N]*rotationmatrix;
 
M1=ceil(max([abs(northwest(1)-southeast(1)) abs(northeast(1)-southwest(1))]));    
N1=ceil(max([abs(northwest(2)-southeast(2)) abs(northeast(2)-southwest(2))]));   
 
RESULT=uint8(zeros(M1,N1,channels));
 
DY=abs(min([northwest(1) northeast(1) southwest(1) southeast(1)])); 
DX=abs(min([northwest(2) northeast(2) southwest(2) southeast(2)])); 
 
 
for i=1-DY:M1-DY
    for j=1-DX:N1-DX
        originalpoint=[i,j]/rotationmatrix;
        y=originalpoint(1);
        x=originalpoint(2);
        if 1<=y && y<=M && 1<=x && x<=N
            for k=1:channels
                RESULT(i+DY,j+DX,k)=interpolation(IMAGE,x+n,y+n,k);
            end
        end
    end
end
缩放
function RESULT=myresize(IMAGE,XSCALE,YSCALE,METHOD)
%MYRESIZE Resize image.
%   RESULT = myresize(IMAGE, XSCALE, YSCALE, METHOD) returns an image
%   that is XSCALE times the column size and YSCALE times the row size of
%   IMAGE, which is a grayscale, RGB, or binary image. If A has more than
%   two dimensions, only the first two dimensions are resized.
%
%   To control the interpolation method used by MYRESIZE, use the METHOD
%   argument. METHOD can be a string naming a general interpolation method:
%  
%       'nearest'    - nearest-neighbor interpolation
% 
%       'bilinear'   - bilinear interpolation
% 
%       'bicubic'    - cubic interpolation; the default method
%
%
%   Note that you must specify METHOD argument, because there is
%   no default value for METHOD.
%
%   Examples
%   --------
%   Shrink by factor of two using bilinear interpolation
%
%       filename='color.jpg';
%       img=imread(filename);
%       xscale=0.5; yscale=0.5;
%       res=myresize(img,xscale,yscale,'bilinear');
%       imshow(res);
%       storagepath='results/resize';
%       mkdir(storagepath);
%       imwrite(res,[storagepath,'/',filename]);
 
% get the size of input image and create the result image
[rows,cols,channels]=size(IMAGE);
RESULT=uint8(zeros(round(rows*YSCALE),round(cols*XSCALE),channels));
% preprocess the input image according to the used method
switch METHOD
    case 'nearest'
        n=1;
        interpolation=@nearestinterpolation;
    case 'bilinear'
        n=1;
        interpolation=@bilinearinterpolation;
    otherwise 
        n=2;
        interpolation=@bicubicinterpolation;
end
IMAGE=preprocessimage(IMAGE,n);
% resize the image
for k=1:channels
    for i=1:size(RESULT,1)
        for j=1:size(RESULT,2)
            x=j/XSCALE;
            y=i/YSCALE;
            RESULT(i,j,k)=interpolation(IMAGE,x+n,y+n,k);
        end
    end
end

提升运行效率的一些尝试

(可能有误)
在执行图像旋转函数myrotate.m时发现问题,运行一次(放大1.5倍)往往需要近20s的时间。相对地,内置函数imrotate.m仅需要不到1s;参考另一位同学代码发现也只需要6s左右。

经过检查与分析发现原始代码中主要有三个影响运行效率的地方:

  1. 程序使用了switch语句判断执行哪一种插值算法,但初始时为了写起来简便将switch语句写在了三重for循环内,可能导致了较大分支开销。将这部分移到循环外,使用@实现后,节省了约3s时间;但使用@相对直接调用函数仍然多消耗了3s左右。(即如果完全去掉分支直接使用双线性内插法,会比原来节省6s左右。这部分原因未知,可能是由于@带来了一定的参数传递消耗)。
  2. 由于习惯性将三通道的处理放在最外层循环,导致由放大后图像某个像素求其在原始图像中的位置(浮点坐标)以及其他一些运算重复运算了多次。将三通道的处理由最外层循环改为最内层循环后,程序运行时间近乎变为原来的1/3,达到6s左右。
  3. 试图进一步提高程序运行效率时,搜索了大量资料,得知MATLAB中应尽量使用矩阵运算而不是for循环,使用内置的矩阵运算往往能提高4-5倍运行效率(甚至更多)。经过一些尝试,暂时没有实现将现有的三重循环以矩阵运算形式实现。尽管最内层的处理三通道的循环可以很轻易地使用矩阵运算实现,但其竟比使用for循环还要缓慢,原因未知。

由于处理图像缩放的主体思路与处理图像旋转基本一致,对应修改了myresize的影响执行效率部分,一定程度上提高了执行效率。

实验结果

(图片过多,待补充)