clear all;
close all;
im=imread('m640.tif','tiff');

[xs,ys]=size(im(:,:,1));
xss=64;
fal=zeros(xs,ys,3);
for color=1:3;
fal(xs/2+1:xs,ys/2+1:ys,color)=foveate3(im(xs/2+1:xs,ys/2+1:ys,color),xss)/255;
color
fal(1:xs/2,ys/2+1:ys,color)=rot90(foveate3(rot90(im(1:xs/2,ys/2+1:ys,color),3),xss))/255;
fal(1:xs/2,1:ys/2,color)=rot90(foveate3(rot90(im(1:xs/2,1:ys/2,color),2),xss),2)/255;
fal(xs/2+1:xs,1:ys/2,color)=rot90(foveate3(rot90(im(xs/2+1:xs,1:ys/2,color)),xss),3)/255;
end;

%imShow(fal);
imwrite(fal,'out.tif','tiff');

%===================================================

function out=foveate(in,max);
[xs,ys]=size(in);
out=zeros(xs,ys);
xsize=1;
ysize=1;
nx=max;
ny=max;
x=1;
y=1;
while y<=ys;

while x<=xs;
%[x,y,nx,ny, xsize, ysize]
for xx=1:nx
 for yy=1:ny
%  [x,x+xsize-1,y,y+ysize-1]
  out(x:x+xsize-1,y:y+ysize-1)=mean(mean(in(x:x+xsize-1,y:y+ysize-1)));
  y=y+ysize;
 end;
 x=x+xsize;
 y=y-ny*ysize;
end;
xsize=xsize*2;
nx=nx/2;
end;

%'x>32'
x=1;
xsize=1;
nx=max;
ysize=ysize*2;
ny=ny/2;
if nx==0.5
 x=1;
 nx=max;
 xsize=1;
end;
y=y+ny*ysize;
%ysize=1;
%ny=max;
end;
