%% show level 1 decomposition of an image
% using the wavelet tool box  
% S. Allen Broughton  7 Nov 2008

%% set dwt parameters
dwtmode('per')
wn = 'db2'
% load image
load wbarb
finegray = (0:255)'*[1,1,1]/256;

%% compute DWT
[CA,CH,CV,CD] = dwt2(X,wn);

% compute max values
ma = max(abs(CA(:))); 
mh = max(abs(CH(:)));
mv = max(abs(CV(:)));
md = max(abs(CD(:)));

figure(1)
colormap(finegray);
imagesc(X);
axis equal
axis tight
title('original')
figure(2) 
colormap(finegray);
imagesc([CA,CV; CH,CD]);
title('level 1 transform')
axis equal
axis tight
figure(3)
colormap(finegray);
imagesc([CA/ma,CV/mv;CH/mh,CD/md]);
axis equal
axis tight
title('normalized level 1 transform')

%% compute and show approximation and details
Z = zeros(size(CA));
XA = idwt2(CA,Z,Z,Z,wn);
XH = idwt2(Z,CH,Z,Z,wn);
XV = idwt2(Z,Z,CV,Z,wn);
XD = idwt2(Z,Z,Z,CD,wn);

figure(5)
colormap(finegray);
subplot(2,2,1)
imagesc(XA);
title('approximation')
axis equal
axis tight
subplot(2,2,2);
imagesc(XV);
title('vertical details')
axis equal
axis tight
subplot(2,2,3);
imagesc(XH);
title('horizontal details')
axis equal
axis tight
subplot(2,2,4);
imagesc(XD);
title('diagonal details')
axis equal
axis tight
figure(6)
colormap(finegray);
imagesc([XA,XH+XV+XD]);
title('approximation and details')
axis equal
axis tight


%% Energy calculations
% recall energy = X.X/(#elements) 
% and that XA has 4 times as many elments as CA 

wae = energyCA); 
whe = energy(CH);
wve = energy(CV);
wde = energy(CD);
xae = energy(XA); 
xhe = energy(XH); 
xve = energy(XV); 
xde = energy(XD); 

before = [wae,whe,wve,wde]
after = [xae,xhe,xve,xde]
enery_diff = before/4-after