本文共 2582 字，大约阅读时间需要 8 分钟。

function [points, pos, faceInds] = intersectLineMesh3d(line, vertices, faces, varargin)%INTERSECTLINEMESH3D Intersection points of a 3D line with a mesh%%   INTERS = intersectLineMesh3d(LINE, VERTICES, FACES)%   Compute the intersection points between a 3D line and a 3D mesh defined%   by vertices and faces.%%   [INTERS POS INDS] = intersectLineMesh3d(LINE, VERTICES, FACES)%   Also returns the position of each intersection point on the input line,%   and the index of the intersected faces.%   If POS > 0, the point is also on the ray corresponding to the line. %   %   Example%   intersectLineMesh3d%%   See also%   meshes3d, triangulateFaces, intersectLineTriangle3d%% ------% Author: David Legland% e-mail: david.legland@grignon.inra.fr% Created: 2011-12-20,    using Matlab 7.9.0.529 (R2009b)% Copyright 2011 INRA - Cepia Software Platform.tol = 1e-12;if ~isempty(varargin)    tol = varargin{1};end% ensure the mesh has triangular facestri2Face = [];if iscell(faces) || size(faces, 2) ~= 3    [faces, tri2Face] = triangulateFaces(faces);end% find triangle edge vectorst0  = vertices(faces(:,1), :);u   = vertices(faces(:,2), :) - t0;v   = vertices(faces(:,3), :) - t0;% triangle normaln   = normalizeVector3d(vectorCross3d(u, v));% direction vector of linedir = line(4:6);% vector between triangle origin and line originw0 = bsxfun(@minus, line(1:3), t0);a = -dot(n, w0, 2);b = dot(n, repmat(dir, size(n, 1), 1), 2);valid = abs(b) > tol & vectorNorm3d(n) > tol;% compute intersection point of line with supporting plane% If pos < 0: point before ray% IF pos > |dir|: point after edgepos = a ./ b;% coordinates of intersection pointpoints = bsxfun(@plus, line(1:3), bsxfun(@times, pos, dir));%% test if intersection point is inside triangle% normalize direction vectors of triangle edgesuu  = dot(u, u, 2);uv  = dot(u, v, 2);vv  = dot(v, v, 2);% coordinates of vector v in triangle basisw   = points - t0;wu  = dot(w, u, 2);wv  = dot(w, v, 2);% normalization constantD = uv.^2 - uu .* vv;% test first coordinates = (uv .* wv - vv .* wu) ./ D;% ind1 = s < 0.0 | s > 1.0;ind1 = s < -tol | s > (1.0 + tol);points(ind1, :) = NaN;pos(ind1) = NaN;% test second coordinate, and third triangle edget = (uv .* wu - uu .* wv) ./ D;% ind2 = t < 0.0 | (s + t) > 1.0;ind2 = t < -tol | (s + t) > (1.0 + tol);points(ind2, :) = NaN;pos(ind2) = NaN;% keep only interesting pointsinds = ~ind1 & ~ind2 & valid;points = points(inds, :);pos = pos(inds);faceInds = find(inds);% convert to face indices of original meshif ~isempty(tri2Face)    faceInds = tri2Face(faceInds);end

转载地址：http://dvxdi.baihongyu.com/