function dual = computeDualPolytope(poly, rep)
% calculate the dual polytope
%
% Syntax: dual = computeDualPolytope(polytope)
%   or    dual = computeDualPolytope(polytope, rep)
%
% The dual of a point is a hyperplane and vice versa. In that way the 
%   extreme points of a V-representation can be mapped to the halfspaces of
%   a dual polytope. For the polytope and the dual polytope there exists a
%   one-to-one relationship of extreme points and hyperplanes. Especially
%   the following terms are equivalent by duality:
%     * a hyperplane is uniquely determined by n linearly independent
%       points
%     * n linearly independent hyperplanes intersect in a single point
%
% Usually, if a H-representation is provided, a V-representation is
%   returned and vice versa. Using rep, you can specify the output
%   representation directly. Use [] or '' to maintain the default behavior.
%
% The dual mapping is connected to the unit sphere. Hyperplanes that 
%   intersect the sphere (do not intersect; touch) are mapped to points 
%   outside (inside; on the circumfence) of the sphere and vice versa. For 
%   this mapping the following must be considered:
%     * no point of a V-representation is equal to the origin. This would
%       map to an hyperplane at infinity.
%     * no facet distance equals zero. This would map to a point at
%       infinity.
%
% See also: scaleVrep, scaleHrep, normalizeHrep

% 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
infoPolytope = identifyPolytope(poly);
if ~exist('rep', 'var') || isempty(rep)
    switch infoPolytope.rep
        case 'Hrep'
            rep = 'Vrep';
        case 'Vrep'
            rep = 'Hrep';
    end
end

switch infoPolytope.rep
    case 'Hrep'
        for iFacet = 1:infoPolytope.m
            poly.A(iFacet,:) = poly.A(iFacet,:) / ...
                                poly.h(iFacet);
            poly.h(iFacet) = 1;
        end
        dual.V = poly.A;
    case 'Vrep'
        dual.A = poly.V;
        dual.h = ones(infoPolytope.m, 1);
        dual = normalizeHrep(dual);
end

dual = convertPolytope(dual, rep);

end
