如果按照你的程序的话,根据提示内容,在计算

i1=round(i*cos(a) - j*sin(a)+ n * sin(a))+1;

j1=round(i*sin(a) + j*cos(a));

这两个时会出现零值,那么,在MATLAB中索引F矩阵就是错误的了,你可以设置断点,然后单步运行一下看看在哪一步出现的零值,你根据旋转矩阵的计算应该是没有问题的,关键是灰度映射时可能出现零位置,最好在其中加入判断的语句,j1计算会出现0值.

你这样写,是沿坐标轴原点旋转,你说的我原来写的那个是沿固定点旋转,实现都是使用旋转矩阵,沿固定点只是先平移到原点,然后再平移回去,你可以看看,具体的我已经记不太清楚了,呵呵,已经不做图像了.

我写的那个的程序：

function im_final = imrotate_my(im_path,theta,options)

% IM_ROTATE 两维图像旋转以及双线性灰度插值算法的实现

% im_path 图像存储路径

% theta 旋转角度,正数表示顺时针旋转

% options 可以为circular（超出范围部分,按照周期形式扩展）

% crop（超出部分置零,即全黑）

% Ref. 章毓晋. 图像工程（上册）——图像处理. 清华大学出版社

% Author: lskyp Date: 2009.08.12

% Version: V1.2 Original Version: V1.0 im_bilinear.m;im_rotate.m

% with the parameter options added

error(nargchk(2,3,nargin,'string'))

if nargin == 2

options = 'circular';

else

if ~ (strcmp(options,'circular') || strcmp(options,'crop'))

error('错误的输出方法')

end

end

im_init = imread(im_path);

im_init = double(im_init);

im_height = size(im_init,1);

im_width = size(im_init,2);

% 分别处理灰度图像和RGB图像

if ndims(im_init) == 3

im_final = zeros(im_height,im_width,3);

R = im_init(:,:,1);

G = im_init(:,:,2);

B = im_init(:,:,3);

R_final = im_final(:,:,1);

G_final = im_final(:,:,2);

B_final = im_final(:,:,3);

else

im_final = zeros(im_height,im_width);

end

rot_matrix = [cos(theta) -sin(theta);sin(theta) cos(theta)];

orig_h = (im_height + 1)/2;

orig_w = (im_width + 1)/2;

for h = 1:im_height

for w = 1:im_width

% 平移至原点,旋转,然后再平移回去

new_position = rot_matrix*[h - orig_h;w - orig_w] + [orig_h;orig_w];

% 超出范围按周期扩展控制,options参数控制

if strcmp(options,'circular')

new_position(1) = mod(new_position(1),im_height);

new_position(2) = mod(new_position(2),im_width);

if new_position(1) == 0

new_position(1) = im_height;

end

if new_position(2) == 0

new_position(2) = im_width;

end

end

% 如果新位置为整数,那么直接赋予灰度值或者RGB值,否则,按照双线性插值计算.

% 使用后向映射

if new_position == round(new_position)

if new_position(1) == 0

new_position(1) = 1;

end

if new_position(2) == 0

new_position(2) = 1;

end

% 超出范围控制,options为crop选项,超出范围置零

if strcmp(options,'crop') && (new_position(1) >= im_height || ...

new_position(2) >= im_width || new_position(1) < 0 || ...

new_position(2) < 0)

if ndims(im_init) == 3

R_final(h,w) = 0;

G_final(h,w) = 0;

B_final(h,w) = 0;

else

im_final(h,w) = 0;

end

else

if ndims(im_init) == 3

R_final(h,w) = R(new_position(1),new_position(2));

G_final(h,w) = G(new_position(1),new_position(2));

B_final(h,w) = B(new_position(1),new_position(2));

else

im_final(h,w) = im_init(new_position(1),new_position(2));

end

end

else

h_new = floor(new_position(1));

w_new = floor(new_position(2));

if h_new == 0

h_new = 1;

end

if w_new == 0

w_new = 1;

end

% 超出范围控制,options为crop选项,超出范围置零

if strcmp(options,'crop') && (h_new >= im_height || ...

w_new >= im_width || h_new < 0 || ...

w_new < 0)

if ndims(im_init) == 3

R_final(h,w) = 0;

G_final(h,w) = 0;

B_final(h,w) = 0;

else

im_final(h,w) = 0;

end

else

% 边界控制

h1 = h_new + 1;

w1 = w_new + 1;

if h1 >= im_height + 1

h1 = mod(h1,im_height);

end

if w1 >= im_width + 1

w1 = mod(w1,im_width);

end

if ndims(im_init) == 3

% 双线性插值的实现过程

% Ref. 章毓晋. 图像工程（上册）——图像处理. 清华大学出版社

R_temp1 = R(h1,w_new)*(new_position(1) - h_new) + ...

R(h_new,w_new)*(h_new + 1 - new_position(1));

R_temp2 = R(h1,w1)*(new_position(1) - h_new) + ...

R(h_new,w1)*(h_new + 1 - new_position(1));

R_final(h,w) = R_temp1*(w_new + 1 - new_position(2)) + ...

R_temp2*(new_position(2) - w_new);

G_temp1 = G(h1,w_new)*(new_position(1) - h_new) + ...

G(h_new,w_new)*(h_new + 1 - new_position(1));

G_temp2 = G(h1,w1)*(new_position(1) - h_new) + ...

G(h_new,w1)*(h_new + 1 - new_position(1));

G_final(h,w) = G_temp1*(w_new + 1 - new_position(2)) + ...

G_temp2*(new_position(2) - w_new);

B_temp1 = B(h1,w_new)*(new_position(1) - h_new) + ...

B(h_new,w_new)*(h_new + 1 - new_position(1));

B_temp2 = B(h1,w1)*(new_position(1) - h_new) + ...

B(h_new,w1)*(h_new + 1 - new_position(1));

B_final(h,w) = B_temp1*(w_new + 1 - new_position(2)) + ...

B_temp2*(new_position(2) - w_new);

else

gray_temp1 = im_init(h1,w_new)*(new_position(1) - h_new) + ...

im_init(h_new,w_new)*(h_new + 1 - new_position(1));

gray_temp2 = im_init(h1,w1)*(new_position(1) - h_new) + ...

im_init(h_new,w1)*(h_new + 1 - new_position(1));

im_final(h,w) = gray_temp1*(w_new + 1 - new_position(2)) + ...

gray_temp2*(new_position(2) - w_new);

end

end

end

end

end

if ndims(im_init) == 3

im_final(:,:,1) = R_final;

im_final(:,:,2) = G_final;

im_final(:,:,3) = B_final;

end

im_final = im2uint8(mat2gray(im_final));

实现部分主要就在控制边界那个地方,我用的周期扩展,如果为零则认为是边界地方,实现部分你可以看看.

具体的就只有这么多了,呵呵,多谢信任.

声明一下,程序内容如有朋友转,请注明.

这个函数只能帮助理解算法,如果真正用起来,还是使用IPT的imrotate函数.