pkg load image;

%% render a Droste image by 'ray tracing'

img = im2double(imread('test.jpg'));

xs = linspace(0, 1, size(img, 2));
ys = linspace(0, 1, size(img, 1));

[X, Y] = meshgrid(xs, ys);

%% Droste regions (y, x, shrink)
p = [.4  .325  .6; ...
     .2  .8    .1; ...
     .75  .825  .3];


tic
[X, Y, hops] = bkup(X, Y, p);
toc


% now that we know where every pixel comes from, render image
V = nan(length(ys), length(xs), 3);
tic
for l = 1:3
  V(:,:,l) = interp2(img(:,:,l), 1+(size(img,2)-1)*X, 1+(size(img,1)-1)*Y, 'cubic');
end
toc

figure(1);
imshow(V);
axis equal

figure(2);
imagesc(hops);
colorbar;
axis equal

imwrite(V, 'exact.png');
imwrite(V, 'exact.jpg');
