Wavelet Transform for N dimensions

Question

I came across this amazing response Applying MATLAB's idwt2 several times which I executed to understand it myself. However, I am unable to get how to use the same with work with an RGB image. So, I have 3 Questions.

1. How would the code be applied to an RGB image with only the transformed image displayed in the output that is along with the high and low frequency components along row and column,is it possible to view the fusion of all the components as a single image? I am aware that I have to put cat operator, but I cant understand how to go about it.

2. Secondly, I am also getting a mazed image! I am perplexed since I cannot seem to follow the reason. I have also attached the same code with the statement showing how this image has been generated.

   3.What does the term db1 in the function signature of dwt imply?

CODE:

    load woman;             % Load image data
%startImage=imread('pic_rgb.jpg');  % IF I WANT TO WORK WITH RGB IMAGE
    nLevel = 3;             % Number of decompositions
    nColors = size(map,1);  % Number of colors in colormap
    cA = cell(1,nLevel);    % Approximation coefficients
    cH = cell(1,nLevel);    % Horizontal detail coefficients
    cV = cell(1,nLevel);    % Vertical detail coefficients
    cD = cell(1,nLevel);    % Diagonal detail coefficients
    startImage = X;
    for iLevel = 1:nLevel,
      [cA{iLevel},cH{iLevel},cV{iLevel},cD{iLevel}] = dwt2(startImage,'db1');



     startImage = cA{iLevel};
    end

    figure;colormap(map);
    imagesc(dwt2(startImage,'db1')); %THIS GIVES THE MAZED IMAGE INSTEAD OF THE TRANSFORMED IMAGE
    figure;
    tiledImage = wcodemat(cA{nLevel},nColors);
    for iLevel = nLevel:-1:1,
     tiledImage = [tiledImage                   wcodemat(cH{iLevel},nColors); ...
                    wcodemat(cV{iLevel},nColors) wcodemat(cD{iLevel},nColors)];

    end
    figure;

    imshow(tiledImage,map);

    %reconstruct
    fullRecon = cA{nLevel};
    for iLevel = nLevel:-1:1,
      fullRecon = idwt2(fullRecon,cH{iLevel},cV{iLevel},cD{iLevel},'db1');
    end
    partialRecon = cA{nLevel};
    for iLevel = nLevel:-1:1,
      partialRecon = idwt2(partialRecon,[],[],[],'db1');
    end
    figure;
    imshow([X fullRecon; partialRecon zeros(size(X))],map,...
           'InitialMagnification',50);

Concourse

Didn't you already ask this earlier today ? I don't see the original question any more - did you delete it and repost it or something ? – Matty

Yes sure I did and it went unnoticed and unanswered for the entire day. So, I guessed may be to repost it again to highlight it. I It would have been of immense help if atleast you could juggle out some answers instead of re-editing it all over again,though I had kept your previous edited version intact. Thanx for the added effort! – Concourse

All I did was fix the tags - I didn't edit the question itself. For future reference you should try to improve the original question if you don't get any responses, rather than reposting. – Matty

@RiaGeorge I removed the bit of text about the reputation requirement for uploading an image. Now that you have enough reputation, please upload it. If you do not have enough reputation to upload an image at any point, add the link anyway and someone with enough reputation will edit it in for you. – Froze

Answer 1

The sample image used in my answer to that other question was an indexed image, so there are a few changes that need to be made to get that code working for an RGB image.

I'll first address your question about the 'db1' argument passed to DWT2. This specifies the type of wavelet to use for the decomposition (in this case, a Daubechies wavelet). More information about available wavelets can be found in the documentation for the functions WFILTERS and WAVEINFO.

I'll address your first two questions by showing you how to modify the code from my other answer to work for an RGB image. I'll use the sample 'peppers.png' image. You'll first want to load your image and define the number of values each color component has. Since the sample image is an unsigned 8-bit integer type (the most common situation), nColors will be 256:

X = imread('peppers.png');  %# Load sample image
nColors = 256;              %# Number of values per color component

If your images are larger unsigned integer types (e.g. 'uint16'), a general way to find the number of color values is to use the function INTMAX like so:

nColors = double(intmax(class(X)))+1;

For the ensuing code, an image type of 'uint8' is assumed.

Applying the decompositions is no different than in the indexed image case. The coefficient matrices will simply be M-by-N-by-3 matrices instead of M-by-N matrices:

nLevel = 3;             %# Number of decompositions
cA = cell(1,nLevel);    %# Approximation coefficient storage
cH = cell(1,nLevel);    %# Horizontal detail coefficient storage
cV = cell(1,nLevel);    %# Vertical detail coefficient storage
cD = cell(1,nLevel);    %# Diagonal detail coefficient storage
startImage = X;
for iLevel = 1:nLevel,  %# Apply nLevel decompositions
  [cA{iLevel},cH{iLevel},cV{iLevel},cD{iLevel}] = dwt2(startImage,'db1');
  startImage = cA{iLevel};
end

The code to create the tiled image to show the horizontal, vertical, and diagonal components for each decomposition will change due to the fact that we are now working with 3-D matrices and must use the CAT function instead of the concatenation operator []:

tiledImage = wcodemat(cA{nLevel},nColors);
for iLevel = nLevel:-1:1
  tiledImage = cat(1,cat(2,tiledImage,...
                           wcodemat(cH{iLevel},nColors)),...
                     cat(2,wcodemat(cV{iLevel},nColors),...
                           wcodemat(cD{iLevel},nColors)));
end
figure;
imshow(uint8(tiledImage-1));  %# Convert to unsigned 8-bit integer to display

This will give the following image showing the horizontal (top right), vertical (bottom left), and diagonal (bottom right) components for each decomposition step, along with the reduced image (top left):

The reconstruction steps are unchanged from the other answer. Only the code for displaying the final images needs to be modified:

fullRecon = cA{nLevel};
for iLevel = nLevel:-1:1,
  fullRecon = idwt2(fullRecon,cH{iLevel},cV{iLevel},cD{iLevel},'db1');
end
partialRecon = cA{nLevel};
for iLevel = nLevel:-1:1,
  partialRecon = idwt2(partialRecon,[],[],[],'db1');
end
figure;
tiledImage = cat(1,cat(2,X,uint8(fullRecon)),...
                   cat(2,uint8(partialRecon),zeros(size(X),'uint8')));
imshow(tiledImage,'InitialMagnification',50);

And you will get an image showing the original RGB image (top left), the fully-reconstructed image using all of the stored detail coefficient matrices (top right), and the partially-reconstructed image using none of the stored detail coefficient matrices (bottom left):

Navarrete

Thank you,but is it also possible to view only the final decomposed image excluding all the components?Just one tiled image displaying the wavelet transform of the peppers image.Further, the variable X seems to be inbuilt meaning that no other variable works apart from using X. – Concourse

@Ria: You can view just the transformed image like so: imshow(uint8(wcodemat(cA{N},nColors)-1));. The value for N is either 1, 2, or 3 depending on whether you want to view an image that has been decomposed once, twice, or three times. X is just the variable I chose to store the image data in. You can easily just replace X everywhere in the code with whatever you want to call that variable. – Navarrete

:I included the above line after the first for loop. But it displays the original image instead of the transformed image. I tried with different levels of N but all give the original image. :( – Concourse

@Ria: Are you sure it's displaying the original image? You may think it's displaying the original, but each transformed image is actually a smaller version of the original, with the image size getting reduced by 2 every time it is wavelet transformed. – Navarrete

Yes its exactly the original image. Actually, I intend to view an output as it appears when we do Fast Fourier Transform. So,in this case, should the output be any one of the image displayed in the tiled result? – Concourse

:Is it not possible to view all the components together?Kindly revert will be of immense help. – Concourse