function [xval, yval] = compress2Dh(x, y, width, height, alpha)
% COMPRESS2DH: Rudimental compression algorithm for 2D matlab plots
%
%
% It is a variant of compress2D(x, y, width, height, alpha). compress2Dh is
% MUCH slower and A BIT more accurate than compress2D.
% Use them if you need to plot highly populated vectors when you expect
% high density of points per pixel.
%
%
% [XVAL, YVAL] = COMPRESS2DH(x, y, width, height) reduces the number of
% points (x, y) which is necessary to represent in a figure, eliminating
% the redundancies. It returns two new vectors xval and yval which can be
% plotted instead of inputs x and y.
%
% It is necessary to specify the width and the height of the desired figure
% and a resolution parameter alpha. The plane is initially fragmented in a
% number of cells which is a function of width, height and alpha:
% ceil(alpha * width) * ceil(height / width * ceil(alpha * width))
%
% alpha must be a positive number such that width * alpha < 3000: the
% smaller alpha, the higher the compression will be, the lower quality the
% final image will have. A suitable value for alpha, when the figure is at
% its maximizum size, is 0.27.
%
% ATTENTION: alpha, width and height all contribute to the quality of the
% final, compressed image. Try different combinations in order to obtain
% the desired compression at the desired quality.
%
%
% If uncertain on how to specify width and height, these lines may help:
% h = figure; % open a new figure (manually choose the desired size)
% set(h, 'Units', 'pixel'); % unit must be pixels
% pos = get(h, 'Position'); % get the Position
% close; % close figure
% pos(1:2) = []; % keep only the last two
% width = pos(1); % find the desired value of width
% height = pos(2); % find the desired value of height
% clear h pos; % clear unnecessary variables
%
% I found this function useful because I needed to save my matlab figures
% in eps format. Since my vectors were big, I obtained very heavy images
% which could not be easily imported in my LaTeX documents. I have been
% looking around for a function that could help me to compress the images
% but I could not find it. So I wrote functions compress2D and compress2Dh
% which are very rudimental and probably need important implementation, but
% they did their job well enough for my purposes. I'm happy to share them:
% suggestions on how to implement them more efficiently are all welcome.
%
% % -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % EXAMPLE1: how to use the function
% N = 1000000; % A big number
% x = rand(N, 1); % random x values
% y = x .^ 2 + rand(N, 1) / 10; % y values along the square of x, plus small random noise
% h = plot(x, y, '.'); % heavy, uncompressed image ...
% saveas(h, 'heavyuncompressedimage.eps'); % ... heavy, uncompressed file
% h1 = dir('heavyuncompressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the uncompressed image, is equal to %0.2f KB', h1.bytes / 1024));
%
% h = figure; % open a new figure (manually choose the desired size)
% set(h, 'Units', 'pixel'); % unit must be pixels
% pos = get(h, 'Position'); % get the Position
% pos(1:2) = []; % keep only the last two
% width = pos(1); % find the desired value of width
% height = pos(2); % find the desired value of height
% clear pos; % clear unnecessary variables
% alpha = 0.27; % a resolution coefficient
% [xind, yind] = compress2Dh(x, y, width, height, alpha); % new vectors for compressed image
% h = plot(xind, yind, '.'); % compressed image ...
% saveas(h, 'compressedimage.eps'); % ... small file
% h2 = dir('compressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the compressed image, is equal to %0.2f KB', h2.bytes / 1024));
% disp(sprintf('Congratulations! You have saved %0.2f%% of space', 100 - h2.bytes / h1.bytes * 100));
%
% % -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % EXAMPLE2: an extremely lucky case
% N = 1000000; % A big number
% x = randsample(50, N, 1); % random x values, integers from 1 to 50
% y = x .^ 2 + rand(N, 1); % y values along the square of x, plus small random noise
% h = plot(x, y, '.'); % heavy, uncompressed image ...
% saveas(h, 'heavyuncompressedimage.eps'); % ... heavy, uncompressed file
% h1 = dir('heavyuncompressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the uncompressed image, is equal to %0.2f KB', h1.bytes / 1024));
%
% h = figure; % open a new figure (manually choose the desired size)
% set(h, 'Units', 'pixel'); % unit must be pixels
% pos = get(h, 'Position'); % get the Position
% pos(1:2) = []; % keep only the last two
% width = pos(1); % find the desired value of width
% height = pos(2); % find the desired value of height
% clear pos; % clear unnecessary variables
% alpha = 3000 / width; % a resolution coefficient
% [xind, yind] = compress2Dh(x, y, width, height, alpha); % new vectors for compressed image
% h = plot(xind, yind, '.'); % compressed image ...
% saveas(h, 'compressedimage.eps'); % ... small file
% h2 = dir('compressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the compressed image, is equal to %0.2f KB', h2.bytes / 1024));
% disp(sprintf('Congratulations! You have saved %0.2f%% of space', 100 - h2.bytes / h1.bytes * 100));
%
% % -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % EXAMPLE3: an unlucky case
% N = 1000000; % A big number
% x = rand(N, 1); % random x values
% y = rand(N, 1); % random y values
% h = plot(x, y, '.'); % heavy, uncompressed image ...
% saveas(h, 'heavyuncompressedimage.eps'); % ... heavy, uncompressed file
% h1 = dir('heavyuncompressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the uncompressed image, is equal to %0.2f KB', h1.bytes / 1024));
%
% h = figure; % open a new figure (manually choose the desired size)
% set(h, 'Units', 'pixel'); % unit must be pixels
% pos = get(h, 'Position'); % get the Position
% pos(1:2) = []; % keep only the last two
% width = pos(1); % find the desired value of width
% height = pos(2); % find the desired value of height
% clear pos; % clear unnecessary variables
% alpha = 0.4; % a resolution coefficient
% [xind, yind] = compress2Dh(x, y, width, height, alpha); % new vectors for compressed image
% h = plot(xind, yind, '.'); % compressed image ...
% saveas(h, 'compressedimage.eps'); % ... small file
% h2 = dir('compressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the compressed image, is equal to %0.2f KB', h2.bytes / 1024));
% disp(sprintf('Congratulations! You have saved %0.2f%% of space', 100 - h2.bytes / h1.bytes * 100));
%
% % -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % EXAMPLE4: lucky (obvious case)
% N = 1000000; % A big number
% x1 = rand(N, 1); % random x1 values
% y1 = 0.5 * ones(N, 1); % horizontal line
% x2 = 0.5 * ones(N, 1); % vertical line
% y2 = rand(N, 1); % random y2 values
% h = plot(x1, y1, '.'); % heavy, uncompressed image ...
% hold on
% plot(x2, y2, '.r');
% saveas(h, 'heavyuncompressedimage.eps'); % ... heavy, uncompressed file
% h1 = dir('heavyuncompressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the uncompressed image, is equal to %0.2f KB', h1.bytes / 1024));
%
% h = figure; % open a new figure (manually choose the desired size)
% set(h, 'Units', 'pixel'); % unit must be pixels
% pos = get(h, 'Position'); % get the Position
% pos(1:2) = []; % keep only the last two
% width = pos(1); % find the desired value of width
% height = pos(2); % find the desired value of height
% clear pos; % clear unnecessary variables
% alpha = 1.50; % a resolution coefficient
% [xind1, yind1] = compress2Dh(x1, y1, width, height, alpha); % find indexes for compressed image
% [xind2, yind2] = compress2Dh(x2, y2, width, height, alpha); % find indexes for compressed image
% h = plot(xind1, yind1, '.'); % compressed image ...
% hold on
% plot(xind2, yind2, '.r');
% saveas(h, 'compressedimage.eps'); % ... small file
% h2 = dir('compressedimage.eps'); % retrieve file size and display the result
%
% disp(sprintf('The file size, for the compressed image, is equal to %0.2f KB', h2.bytes / 1024));
% disp(sprintf('Congratulations! You have saved %0.2f%% of space', 100 - h2.bytes / h1.bytes * 100));
%
% % -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% Author: Francesco Pozzi
% E-Mail: francesco.pozzi@anu.edu.au
% Date: 08-07-2008
%
% Check inputs
ctrl = isvector(x) & isnumeric(x) & isreal(x) & isvector(y) & isnumeric(y) & isreal(y);
if ctrl
ctrl = ctrl & (length(x) == length(y));
ctrl = ctrl & all(~isinf(x)) & all(~isinf(y));
end
if ctrl
x = x(:);
y = y(:);
else
error('x and y must be vectors of real numbers and of the same size!')
end
ctrl = isscalar(width) & isscalar(height) & isscalar(alpha);
if ctrl
ctrl = (width > 0) & (height > 0) & (alpha > 0) & (alpha * width <= 3000);
end
if ~ctrl
error('Specify a valid Position and a valid value for parameter alpha!')
end
% Range of x and y:
xmin = min(x);
xmax = max(x);
ymin = min(y);
ymax = max(y);
% resol1 and resol2 are intervals between minimum and maximum for x and y
% In each area of the grid thus obtained, only one point at most will be
% represented. All others will be discarded
if (xmin == xmax) & (ymin == ymax)
xval = xmin;
yval = ymin;
warning('Hey! Your entire data lie on a single 2D point! Are you sure data are correct?');
return;
elseif ymin == ymax % if horizontal line ...
resol1 = ceil(alpha * width); % resolution along x-axis
resol2 = 1; % resolution along y-axis
j = floor((x - xmin) / (xmax - xmin) * (resol1 - 1)) + 1; % discretized x
k = ones(size(j));
elseif xmin == xmax % if vertical line ...
resol1 = ceil(alpha * width); % resolution along x-axis
resol2 = ceil(height / width * resol1); % resolution along y-axis
resol1 = 1; % resolution along x-axis
k = floor((y - ymin) / (ymax - ymin) * (resol2 - 1)) + 1; % discretized y
j = ones(size(k));
else
resol1 = ceil(alpha * width); % resolution along x-axis
resol2 = ceil(height / width * resol1); % resolution along y-axis
j = floor((x - xmin) / (xmax - xmin) * (resol1 - 1)) + 1; % discretized x
k = floor((y - ymin) / (ymax - ymin) * (resol2 - 1)) + 1; % discretized y
end
xval = sparse(resol1, resol2); % new values for x
yval = sparse(resol1, resol2); % new values for y
m = sparse(resol1, resol2); % number of elements in each cell
for i = 1:length(x)
xval(j(i), k(i)) = xval(j(i), k(i)) + x(i); % for each cell, sum x
yval(j(i), k(i)) = yval(j(i), k(i)) + y(i); % for each cell, sum y
m(j(i), k(i)) = m(j(i), k(i)) + 1; % for each cell, remember number of points
end
temp = find(m); % keep only those elements where data exist
m = m(temp);
xval = xval(temp); % keep only those x elements that exist
yval = yval(temp); % keep only those y elements that exist
xval = xval ./ m; % average values for x
yval = yval ./ m; % average values for y
