1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
|
function [isinside, pt, coord] = linextriangle(p0, p1, plane)
% [isinside,pt,coord]=linextriangle(p0,p1,plane)
%
% calculate the intersection of a 3d line (passing two points)
% with a plane (determined by 3 points)
%
% author: Qianqian Fang <q.fang at neu.edu>
% date: 12/12/2008
%
% parameters:
% p0: a 3d point in form of (x,y,z)
% p1: another 3d point in form of (x,y,z), p0 and p1 determins the line
% plane: a 3x3 matrix, each row is a 3d point in form of (x,y,z)
% this is used to define a plane
% outputs:
% isinside: a boolean variable, 1 for the intersection is within the
% 3d triangle determined by the 3 points in plane; 0 is outside
% pt: the coordinates of the intersection pint
% coord: 1x3 vector, if isinside=1, coord will record the barycentric
% coordinate for the intersection point within the triangle;
% otherwise it will be all zeros.
%
% for degenerated lines or triangles, this will stop
%
% Please find more information at http://iso2mesh.sf.net/cgi-bin/index.cgi?metch
%
% this function is part of "metch" toobox, see COPYING for license
[a, b, c, d] = getplanefrom3pt(plane);
if (a * a + b * b + c * c == 0.0)
error('degenerated plane');
end
dl_n = sum([a b c] .* (p1 - p0));
if (dl_n == 0.0)
error('degenerated line');
end
% solve for the intersection point
t = -(a * p0(1) + b * p0(2) + c * p0(3) + d) / dl_n;
pt = p0 + (p1 - p0) * t;
dist = sum(abs(diff(plane)));
[md, imax] = sort(dist);
if (md(2) == 0.0)
error('degenerated triangle');
end
goodidx = imax(2:end);
ptproj = pt(goodidx);
mat0 = [plane(:, goodidx)', ptproj'; 1 1 1 1];
isinside = 0;
coord = [0 0 0];
det1 = det(mat0(:, [4 2 3], :));
det2 = det(mat0(:, [1 4 3], :));
if (det1 * det2 < 0)
return
end
det3 = det(mat0(:, [1 2 4], :));
if (det2 * det3 < 0)
return
end
if (det1 * det3 < 0)
return
end
isinside = 1;
det0 = det(mat0(:, 1:3));
coord = [det1 det2 det3] / det0;
|
