clipConvexPolygon3dHP

CLIPCONVEXPOLYGON3DHP Clip a convex 3D polygon with Half-space.

   POLY2 = clipConvexPolygon3dHP(POLY, PLANE)
   POLY is a N-by-3 array of points, and PLANE is given as:
   [x0 y0 z0 dx1 dy1 dz1 dx2 dy2 dz2].
   The result POLY2 is also an array of 3d points, sometimes smaller than
   poly, and that can be 0-by-3 (empty polygon).

   POLY2 = clipConvexPolygon3dHP(POLY, PT0, NORMAL)
   uses plane with normal NORMAL and containing point PT0.


   See also:
   polygons3d, polyhedra

 ------
 Author: David Legland
 e-mail: david.legland@grignon.inra.fr
 Created: 2007-01-05
 Copyright 2007 INRA - BIA PV Nantes - MIAJ Jouy-en-Josas.

0001 function poly2 = clipConvexPolygon3dHP(poly, plane)
0002 %CLIPCONVEXPOLYGON3DHP Clip a convex 3D polygon with Half-space.
0003 %
0004 %   POLY2 = clipConvexPolygon3dHP(POLY, PLANE)
0005 %   POLY is a N-by-3 array of points, and PLANE is given as:
0006 %   [x0 y0 z0 dx1 dy1 dz1 dx2 dy2 dz2].
0007 %   The result POLY2 is also an array of 3d points, sometimes smaller than
0008 %   poly, and that can be 0-by-3 (empty polygon).
0009 %
0010 %   POLY2 = clipConvexPolygon3dHP(POLY, PT0, NORMAL)
0011 %   uses plane with normal NORMAL and containing point PT0.
0012 %
0013 %
0014 %   See also:
0015 %   polygons3d, polyhedra
0016 %
0017 % ------
0018 % Author: David Legland
0019 % e-mail: david.legland@grignon.inra.fr
0020 % Created: 2007-01-05
0021 % Copyright 2007 INRA - BIA PV Nantes - MIAJ Jouy-en-Josas.
0022 
0023 %   HISTORY
0024 %   2007/14/09 fix postprocessing of last point
0025 
0026 % ensure last point is the same as the first one
0027 if sum(poly(end, :) == poly(1,:)) ~= 3
0028     poly = [poly; poly(1,:)];
0029 end
0030 
0031 % initialize empty polygon
0032 poly2 = zeros(0, 2);
0033 
0034 % compute visible points
0035 below = isBelowPlane(poly, plane);
0036 
0037 % case of empty polygon
0038 if sum(below) == 0
0039     return;
0040 end
0041 
0042 % case of totally clipped polygon
0043 if sum(below) == length(below)
0044     poly2 = poly;
0045     return;
0046 end
0047 
0048 % indices of edges intersecting the plane
0049 ind = find(below ~= below([2:end 1]));
0050 
0051 % compute intersection points: they are 2 for a convex polygon
0052 lines = createLine3d(poly(ind, :), poly(ind+1, :));
0053 pInt = intersectLinePlane(lines, plane);
0054 
0055 % insert intersection points and remove invisible points
0056 if below(1)
0057     poly2 = [poly(1:ind(1), :); pInt; poly(ind(2)+1:end, :)];
0058 else
0059     poly2 = [pInt(1, :); poly(ind(1)+1:ind(2), :); pInt(2, :)];
0060 end
0061 
0062 % remove last point if it is the same as the first one
0063 if sum(poly2(end, :) == poly2(1,:)) == 3
0064     poly2(end, :) = [];
0065 end
0066