function [vrep, scale, x0] = scaleVrep(vrep, zerotol, flagIsometric)
% scale polytope to reasonable size and position
%
% Syntax: [vrep, scale, x0] = scaleVrep(vrep)
%   or    [vrep, scale, x0] = scaleVrep(vrep, zerotol)
%
% The bounding box of the polytope is calculated. The polytope is then 
%   first moved and secondly stretched so that the bounding box of this
%   polytop equals the cube [-1, 1]^n. In a second step the polytope is
%   moved so that the mean fits the origin.
%
% Stretching values become large, when a polytope is flat and well aligned 
%   with one of the coordinate axes. Zerotol relates to the flatness a
%   polytope can have before stretching in that direction is disabled.
%   Default: see elkZerotol.
%
% With scaleVrep(.., flagIsometric), the polytope can be scaled
%   isometrically (equally in all dimensions). This is usually required for
%   measure calculations. Default: 0
%
% The polytope is scaled back to original by:
%   originalVrep = stretchVrep(scaledVrep, scale);
%   originalVrep = moveVrep(originalVrep, x0);
%
% See also: convertVrepToHrep, scaleHrep, reduceVrepDimension, reduceVrep

% The elk-library: convex geometry applied to crystallization modeling.
%   Copyright (C) 2012 Alexander Reinhold
%
% This program is free software: you can redistribute it and/or modify it 
%   under the terms of the GNU General Public License as published by the 
%   Free Software Foundation, either version 3 of the License, or (at your 
%   option) any later version.
%
% This program is distributed in the hope that it will be useful, but 
%   WITHOUT ANY WARRANTY; without even the implied warranty of 
%   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU 
%   General Public License for more details.
%
% You should have received a copy of the GNU General Public License along 
%   with this program.  If not, see <http://www.gnu.org/licenses/>.

%% Input
if ~exist('zerotol', 'var') || isempty(zerotol)
    zerotol = elkZerotol;
end
if ~exist('flagIsometric', 'var')
    flagIsometric = 0;
end
if isempty(vrep.V)
    scale = 1;
    x0 = ones(size(vrep.V, 2), 1);
    return
end

%% Gather bounding box
dim = size(vrep.V, 2);
minPos = ones(1, dim);
maxPos = ones(1, dim);
for iDim = 1:dim
    minPos(iDim) = min(vrep.V(:, iDim));
    maxPos(iDim) = max(vrep.V(:, iDim));
end

%% Determine scaling coefficients
if flagIsometric
    % the center of the bounding box
    x0 = (maxPos + minPos)/2;
    % the size of the bounding box (but as it is scaled to the box [-1, 1]^dim
    %   the scaling factor is only half of the original size
    scale = 0.5*(maxPos - minPos);
    scale = max(scale);
else
    % the center of the bounding box
    x0 = mean(vrep.V, 1);
    % the size of the bounding box (but as it is scaled to the box [-1, 1]^dim
    %   the scaling factor is only half of the original size
    scale = 0.5*(maxPos - minPos);
    scale(abs(scale) < 0.5*zerotol) = 1;
end
%% Scale polytope
vrep = moveVrep(vrep, -1*x0);
vrep = stretchVrep(vrep, 1./scale, zeros(1, dim));