steinerPolytope

STEINERPOLYTOPE Create a steiner polytope from a set of vectors.

   [VERTICES FACES] = steinerPolygon(VECTORS)
   Creates the Steiner polytope defined by the set of vectors VECTORS.

   Example
     % Creates and display a planar Steiner polytope (ie, a polygon)
     [v f] = steinerPolytope([1 0;0 1;1 1]);
     fillPolygon(v);

     % Creates and display a 3D Steiner polytope 
     [v f] = steinerPolytope([1 0 0;0 1 0;0 0 1;1 1 1]);
     drawMesh(v, f);
     view(3); axis vis3d

   See also
   meshes3d, drawMesh, steinerPolygon, mergeCoplanarFaces

0001 function [vertices, faces] = steinerPolytope(vectors)
0002 %STEINERPOLYTOPE Create a steiner polytope from a set of vectors.
0003 %
0004 %   [VERTICES FACES] = steinerPolygon(VECTORS)
0005 %   Creates the Steiner polytope defined by the set of vectors VECTORS.
0006 %
0007 %   Example
0008 %     % Creates and display a planar Steiner polytope (ie, a polygon)
0009 %     [v f] = steinerPolytope([1 0;0 1;1 1]);
0010 %     fillPolygon(v);
0011 %
0012 %     % Creates and display a 3D Steiner polytope
0013 %     [v f] = steinerPolytope([1 0 0;0 1 0;0 0 1;1 1 1]);
0014 %     drawMesh(v, f);
0015 %     view(3); axis vis3d
0016 %
0017 %   See also
0018 %   meshes3d, drawMesh, steinerPolygon, mergeCoplanarFaces
0019 %
0020 
0021 % ------
0022 % Author: David Legland
0023 % e-mail: david.legland@grignon.inra.fr
0024 % Created: 2006-04-28
0025 % Copyright 2006 INRA - CEPIA Nantes - MIAJ (Jouy-en-Josas).
0026 
0027 % History
0028 % 2013-02-22 merge coplanar faces, add management of 2D case, update doc
0029 
0030 
0031 % compute vectors dimension
0032 nd = size(vectors, 2);
0033 
0034 % create candidate vertices
0035 vertices = zeros(1, size(vectors, 2));
0036 for i = 1:length(vectors)
0037     nv = size(vertices, 1);
0038     vertices = [vertices; vertices+repmat(vectors(i,:), [nv 1])]; %#ok<AGROW>
0039 end
0040 
0041 if nd == 2
0042     % for planar case, use specific function convhull
0043     K = convhull(vertices(:,1), vertices(:,2));
0044     vertices = vertices(K, :);
0045     faces = 1:length(K);
0046     
0047 else 
0048     % Process the general case (tested only for nd==3)
0049     
0050     % compute convex hull
0051     K = convhulln(vertices);
0052     
0053     % keep only relevant points, and update faces indices
0054     ind = unique(K);
0055     for i = 1:length(ind)
0056         K(K==ind(i)) = i;
0057     end
0058     
0059     % return results
0060     vertices = vertices(ind, :);
0061     faces = K;
0062     
0063     % in case of 3D meshes, merge coplanar faces
0064     if nd == 3
0065         faces = mergeCoplanarFaces(vertices, faces);
0066     end
0067 end