ff=imread('lena256.bmp');           %讀取圖像到ff
[mm,nn]=size(ff);              %將圖像隔行隔列抽取元素，得到縮小一半的圖像f
m=mm/2;n=nn/2;
for i=1:m
  for j=1:n
    f(i,j)=ff(2*i,2*j);
  end
end
k=3;                       %設置放大倍數
%將待插值圖像矩陣前后各擴展兩行兩列,共擴展四行四列
a=f(1,:);		%取矩陣f的第一行的全部元素
c=f(m,:);             %取矩陣f的第m行的全部元素
b=[f(1,1),f(1,1),f(:,1)',f(m,1),f(m,1)];d=[f(1,n),f(1,n),f(:,n)',f(m,n),f(m,n)];
a1=[a;a;f;c;c];		%a1為f上擴展兩行，下擴展兩列
a1';			%a1轉置
b1=[b;b;a1';d;d];	%b1為a1轉置后上下各擴展兩行
f=b1';f1=double(f);	%b1再轉置，成為f的周圍擴展兩圈
for i=1:k*m;                 %利用雙三次插值公式對新圖象所有像素賦值
  u=rem(i,k)/k; i1=floor(i/k)+2;		%i能否被k整除，整除時u=0，floor指朝負無窮方向舍入
  A=[sw(1+u) sw(u) sw(1-u) sw(2-u)]; 	%矩陣A  
  for j=1:k*n;
    v=rem(j,k)/k;j1=floor(j/k)+2;
    C=[sw(1+v);sw(v);sw(1-v);sw(2-v)];	%矩陣C
    B=[f1(i1-1,j1-1) f1(i1-1,j1) f1(i1-1,j1+1) f1(i1-1,j1+2)
      f1(i1,j1-1)   f1(i1,j1)   f1(i1,j1+1)   f1(i1,j1+2)
      f1(i1,j1-1)   f1(i1+1,j1) f1(i1+1,j1+1) f1(i1+1,j1+2)
      f1(i1+2,j1-1) f1(i1+2,j1) f1(i1+2,j1+1) f1(i1+2,j1+2)];	%矩陣B
    g1(i,j)=(A*B*C);
  end
end
g=uint8(g1);
imshow(f); title('縮小一半的圖像');             %顯示縮小的圖像
figure,imshow(ff);title('原圖像');               %顯示原圖像
figure,imshow(g);title('雙三次插值放大后的圖像');     %顯示插值后的圖像