Jpeg2000 2-level wavelet transform-lichtenstein

The picture shows a 2 level wavelet transform. For the low-pass and high-pass filters, I have used the same of the JPEG2000, lossy transform, so that this is the wavelet that is used in the JPEG2000 format. In general, it is quite hard to show a wavelet transform properly, because of the different ranges of the different parts. I have made the following choices.

First of all, I wanted to use a picture with both flat surfaces, sharp edges and details. The standard test image Lenna is commonly used for this purpose, but I didn't want to use it because it is not in the public domain. I took a look at the big featured picture gallery con Commons and I have chosen Image:Lichtenstein img processing test.png because it has lots of details in the bottom left, while the top right has almost no details because of the sky. Moreover, the edges of the castle are really sharp, so they can be clearly seen in the wavelet transform.

I have saved the original image in a file called lichtenstein.png, to make it simpler.

Then I have used the following Matlab script:

clear all

%input must be 512x512
img=imread('lichtenstein.png');

img=im2double(img);
img=rgb2gray(img);

% filters in JPEG2000
Lo_D=[0.0267 -0.0168 -0.0782 0.2668 0.6029 0.2668 -0.0782 -0.0168 0.0267];
Hi_D=[0.0912 -0.0575 -0.5912 1.1150 -0.5912 -0.0575 0.0912];

%calculate the 2-level Wavelet transform
[scaled, vertical, horizontal, diagonal]=dwt2(img,Lo_D,Hi_D);
[scaled2, vertical2, horizontal2, diagonal2]=dwt2(scaled,Lo_D,Hi_D);

%crop pictures (for displaying purpose)
scaled2=imcrop(scaled2,[4 4 127 127]);
horizontal2=imcrop(horizontal2,[4 4 127 127]);
vertical2=imcrop(vertical2,[4 4 127 127]);
diagonal2=imcrop(diagonal2,[4 4 127 127]);
horizontal=imcrop(horizontal,[4 4 255 255]);
vertical=imcrop(vertical,[4 4 255 255]);
diagonal=imcrop(diagonal,[4 4 255 255]);

%increase contrast (for displaying purpose)
vertical=imadjust(abs(vertical));
horizontal=imadjust(abs(horizontal));
diagonal=imadjust(abs(diagonal));
vertical2=imadjust(abs(vertical2));
horizontal2=imadjust(abs(horizontal2));
diagonal2=imadjust(abs(diagonal2));

%put them all together
transform=zeros(512,512);

transform(1:128,1:128)=scaled2;
transform(1:128,129:256)=horizontal2;
transform(129:256,1:128)=vertical2;
transform(129:256,129:256)=diagonal2;
transform(1:256,257:512)=horizontal;
transform(257:512,1:256)=vertical;
transform(257:512,257:512)=diagonal;

imwrite(transform,'Jpeg2000_2-level_wavelet_transform-lichtenstein.png');

If I had placed all the transforms one next to the other, I would have got a picture wider than 512x512, with some odd behaviour close to the edges. In order to avoid those problems, I have cropped the transforms so that I can easily put one next to the other to get a 512x512 picture. In general, it is not a good idea to discard information, but here I just wanted to create a simple and clear example.