clear;
M=50;
N=50;
L=2;
load JuleszOddDipole;
% since dipole calculation is quite time consuming, it is done separately 
% and stored in JuleszOddDipole.mat
% the program is very simple
%
% load Julesz ( load the odd, even and coin-flip texture images)
% d=dipole(iodd);
% save JuleszOddDipole
%
% reconstruct the original image by full positive d of dipole histogram
ihat=reconstr(d);
% shrink d to dx<=ceil((M+1)/2) & dy<=ceil((N+1)/2)
ddx=ceil((M+1)/2);
ddy=ceil((N+1)/2);
D=reshape(d(M-ddx:M+ddx,N-ddy:N+ddy,:,:),[ddx*2+1 ddy*2+1 L L]);
% reconstruct the original image by half distance dipole histogram
iihat=reconstrHalf(D,M,N);
if ihat==iodd
  ihatEQiodd=1
else
  ihatEQiodd=0
end;
if iihat==iodd
  iihatEQiodd=1
else
  iihatEQiodd=0
end;
figure(1);
AutoImage(iodd);
figure(2);
AutoImage(ihat);
figure(3);
AutoImage(iihat);
