% file:      	wavelet_box.m, (c) Matthew Roughan, Mon Aug 14 2006
%
N = 256;
Image = MakeImage('Square',N);
figure(1)
imagesc(Image);
colormap(gray);
opts = struct('height',8, 'Color', 'gray');
axis image
axis('off')
set(gca,'DataAspectRatio', [1 1 1]);
exportfig(gcf,sprintf('Plots/wavelet_box.eps'), opts, 'format', 'eps');

[n,J] = quadlength(Image);
qmf = MakeONFilter('Daubechies',8);
L = 5;
wc = FWT2_PO(Image,L,qmf);

figure(3)
Display2dProjV(wc,L,qmf);
axis image
axis('off')
set(gca,'DataAspectRatio', [1 1 1]);
exportfig(gcf,sprintf('Plots/wavelet_box_approx.eps'), opts, 'format', 'eps');

% taken from wt07fig26 (Mallat Fig 7.26)
wc2 = wc;
avg = wc(1:2^L,1:2^L);
wc2(1:2^L,1:2^L) = 1-(avg ./max(max(abs(avg))));
wc2(1:2^L,2^L) = zeros(2^L,1);
wc2(2^L,1:2^L) = zeros(1,2^L);
for j = L:J-1,
	low = 1:2^j;
	high = 2^j+1:2^(j+1);
	max1 = max(max(max(abs(wc2(low,high)))));
	wc2(low,high) = (wc2(low,high)+max1) ./ (2*max1);
	max2 = max(max(max(abs(wc2(high,low)))));
	wc2(high,low) = (wc2(high,low)+max2) ./(2*max2);
	max3 = max(max(max(abs(wc2(high,high)))));
	wc2(high,high) = (wc2(high,high)+max3) ./(2*max3);
	wc2(2^(j),2^j+1:2^(j+1)) = zeros(1,2^(j));
	wc2(2^(j+1),1:2^j) = zeros(1,2^(j));
	wc2(2^j+1:2^(j+1),2^(j)) = zeros(2^(j),1);
	wc2(1:2^j,2^(j+1)) = zeros(2^(j),1);
	wc2(2^(j+1),2^j+1:2^(j+1)) = zeros(1,2^(j));
	wc2(2^j+1:2^(j+1),2^(j+1)) = zeros(2^(j),1);
end
figure(2);clf
set(gcf, 'NumberTitle','off', 'Name','Window 2')
imagesc(wc2)
axis image
axis('off')
colormap(1-gray(256))
exportfig(gcf,sprintf('Plots/wavelet_box_wt.eps'), opts, 'format', 'eps');
