examp_compression.m

linear time－frequency toolbox

%EXAMP_COMPRESSION  Image compression example
%
%   This examples shows how to do image compression using a
%   two-dimensional Wilson transform.
%
%   The image is transformed using an orthonormal DWILT2 transform
%   followed by N-term approximation and quantization. Then the image is
%   reconstructed and compared with the original.
%
%   This example demonstrates both lossless and lossy compression. An image
%   compression is lossless if the truncation error is below the quantization
%   error, as the truncation error is not visible in that case.
%
%   FIGURE 1 The cameraman.
%
%     This figure to the left shows lossless compression of the
%     cameraman. The figure to the right shows a lossy compression.
%

disp('Type "help examp_compression" to see a description of how this example works.');

% -------- defining constants -------------


% Image has dimension 256x256.
L=256;

% The original signal uses 8 bits:
origbits=8;

% Number of channels.
M=32;

% Bitrate
nbits=8;

% Number of coefficients to keep for lossy compression.
Nkeep=2000;

% Constant to use for mulaw encoding / decoding.
mu=2;

% -------- initial setup ------------------

% Load the cameraman picture
f=cameraman();
f = double(f);

% Create window for Wilson orhonormal basis.
g=wilorth(M,L);

% Get coefficients.
c=dwilt2(f,g,M);


disp('----- Lossless compression -----');

% Find max-value of image
max_c=max(abs(f(:)));

% Determine a threshold
thr_val=.5*max_c/(2^origbits);

% Do threshholding
[cN,N_thresh]=thresh('hard',c,thr_val);

% Do quantization.
cQ=uquant('unsigned',cN,nbits);

% Reconstruct.
f_lossless=idwilt2(cQ,g);

% Display it
figure(1);

subplot(1,2,1);
colormap(gray);
imagesc(f_lossless);

disp('');
disp('Compression ratio:');
nbits*N_thresh/(origbits*L*L)

disp('')
disp('relative l^2 error:'); 
norm(f_lossless(:)-double(f(:)))/norm(double(f(:)))

disp('----- Lossy compression ------');

% Do N-term approximation.
cN=largestn(c,Nkeep);

% Do mulaw encoding
[c_mulaw, sigweight] = mulawencode(cN,mu);

% Do quantization.
c_mulaw_Q=uquant('unsigned',c_mulaw,nbits);

% Do mulaw deconding
c_Q=mulawdecode(c_mulaw_Q,mu,sigweight);

% Reconstruct.
f_lossy=idwilt2(c_Q,g);

% Display it
subplot(1,2,2);
colormap(gray);
imagesc(f_lossy);

disp('');
disp('Compression ratio:');
nbits*Nkeep/(origbits*L*L)

disp('')
disp('relative l^2 error:'); 
norm(f_lossy(:)-double(f(:)))/norm(double(f(:)))