2
// This unit is part of the GLScene Engine https://github.com/glscene
5
Octree management classes and structures.
7
TODO: move the many public vars/fields to private/protected
10
10/11/14 - PW - Renamed from "Octree.pas" to "GLOctree.pas" and added ListHandle in TOctreeNode
11
10/11/12 - PW - Added CPP compatibility: changed vector arrays to records
12
30/03/07 - DaStr - Added $I GLScene.inc
13
28/03/07 - DaStr - Renamed parameters in some methods
14
(thanks Burkhard Carstens) (Bugtracker ID = 1678658)
15
24/03/07 - DaStr - Added explicit pointer dereferencing
16
(thanks Burkhard Carstens) (Bugtracker ID = 1678644)
17
02/08/04 - LR, YHC - BCB corrections: use record instead array
18
19/06/04 - LucasG - Moved triangleFiler and WalkSphereToLeaf to public
19
31/01/04 - Mathx - Bugfix on DisposeTree (thanks dikoe Kanguru)
20
20/07/03 - DanB - Modified SphereSweepIntersect to deal with embedded spheres better
21
08/05/03 - DanB - name changes + added ClosestPointOnTriangle + fixes
22
08/05/03 - DanB - added AABBIntersect (Matheus Degiovani)
23
22/01/03 - EG - GetTrianglesInCube moved in (Bernd Klaiber)
24
29/11/02 - EG - Added triangleInfo
25
14/07/02 - EG - Dropped GLvectorFileObjects dependency
26
17/03/02 - EG - Added SphereIntersectAABB from Robert Hayes
27
13/03/02 - EG - Made in a standalone unit, based on Robert Hayes code
38
GLVectorTypes, GLVectorGeometry, GLVectorLists, GLGeometryBB, GLContext;
42
TProcInt = procedure(I: Integer) of object;
43
TProcAffineAffineAffine = procedure(V1, V2, V3: TAffineFLTVector) of object;
45
// TOctreeTriangleInfo
47
{ : Stores information about an intersected triangle. }48
TOctreeTriangleInfo = record
50
Vertex: array [0 .. 2] of TAffineVector;
53
POctreeTriangleInfo = ^TOctreeTriangleInfo;
57
POctreeNode = ^TOctreeNode;
60
MinExtent: TAffineFLTVector;
61
MaxExtent: TAffineFLTVector;
62
// TextureHandle:TGLTextureHandle;
63
ListHandle: TGLListHandle;
64
WillDraw: Boolean; // temporary change
66
// Duplicates possible?
67
TriArray: array of Integer; // array of triangle references
69
ChildArray: array [0 .. 7] of POctreeNode; // Octree's 8 children
74
{ : Manages an Octree containing references to triangles. }75
TOctree = class(TObject)
79
Intersections: Integer;
80
// for debugging - number of triangles intersecting an AABB plane
84
// Find the exact centre of an AABB
85
function GetMidPoint(Min, Max: Single): Single;
86
// Check if a point lies within the AABB specified by min and max entents
87
function PointInNode(const Min, Max, APoint: TAffineFLTVector): Boolean;
88
// Check if a triangle (vertices v1, v2, v3) lies within the AABB specified by min and max entents
89
function TriIntersectNode(const MinExtent, MaxExtent, V1, V2,
90
V3: TAffineFLTVector): BOOLEAN;
91
// Check if a sphere (at point C with radius) lies within the AABB specified by min and max entents
92
function SphereInNode(const MinExtent, MaxExtent: TAffineVector;
93
const C: TVector; Radius: Single): Boolean;
95
procedure WalkTriToLeafx(Onode: POctreeNode;
96
const V1, V2, V3: TAffineFLTVector);
97
procedure WalkPointToLeafx(ONode: POctreeNode; const P: TAffineVector);
98
procedure WalkSphereToLeafx(Onode: POctreeNode; const P: TVector;
100
procedure WalkRayToLeafx(Onode: POctreeNode; const P, V: TVector);
102
function GetExtent(const Flags: array of Byte; ParentNode: POctreeNode)
105
{ : Recursive routine to build nodes from parent to max depth level. }106
procedure Refine(ParentNode: POctreeNode; Level: Integer);
108
// Main "walking" routines. Walks the item through the Octree down to a leaf node.
109
procedure WalkPointToLeaf(ONode: POctreeNode; const P: TAffineVector);
110
procedure WalkTriToLeaf(Onode: POctreeNode;
111
const V1, V2, V3: TAffineVector);
112
procedure WalkRayToLeaf(Onode: POctreeNode; const P, V: TVector);
114
// : Example of how to process each node in the tree
115
procedure ConvertR4(ONode: POctreeNode; const Scale: TAffineFLTVector);
117
procedure CreateTree(Depth: Integer);
122
WorldMinExtent, WorldMaxExtent: TAffineFLTVector;
123
RootNode: POctreeNode; // always points to root node
124
MaxOlevel: Integer; // max depth level of TOctreeNode
126
// number of nodes (ex: 8 for level 1, 72 for level 2 etc).
127
TriCountMesh: Integer; // total number of triangles in the mesh
128
TriCountOctree: Integer; // total number of triangles cut into the octree
129
MeshCount: Integer; // number of meshes currently cut into the Octree
131
ResultArray: array of POctreeNode;
132
// holds the result nodes of various calls
134
{ 19/06/2004 - Lucas G. - Needed this change - Used in ECMisc.pas }135
TriangleFiler: TAffineVectorList;
136
procedure WalkSphereToLeaf(Onode: POctreeNode; const P: TVector;
139
{ : Initializes the tree from the triangle list.140
All triangles must be contained in the world extent to be properly
141
taken into account. }
142
procedure InitializeTree(const AWorldMinExtent, AWorldMaxExtent
143
: TAffineVector; const ATriangles: TAffineVectorList;
144
const ATreeDepth: Integer);
145
procedure DisposeTree;
147
destructor Destroy; override;
149
function RayCastIntersect(const RayStart, RayVector: TVector;
150
IntersectPoint: PVector = nil; IntersectNormal: PVector = nil;
151
TriangleInfo: POctreeTriangleInfo = nil): Boolean;
152
function SphereSweepIntersect(const RayStart, RayVector: TVector;
153
const Velocity, Radius: Single; IntersectPoint: PVector = nil;
154
IntersectNormal: PVector = nil): Boolean;
156
function TriangleIntersect(const V1, V2, V3: TAffineVector): Boolean;
157
{ : Returns all triangles in the AABB. }158
function GetTrianglesFromNodesIntersectingAABB(const ObjAABB: TAABB)
160
{ : Returns all triangles in an arbitrarily placed cube }161
function GetTrianglesFromNodesIntersectingCube(const ObjAABB: TAABB;
162
const ObjToSelf, SelfToObj: TMatrix): TAffineVectorList;
163
{ : Checks if an AABB intersects a face on the octree }164
function AABBIntersect(const AABB: TAABB; M1to2, M2to1: TMatrix;
165
Triangles: TAffineVectorList = nil): Boolean;
166
// function SphereIntersect(position:TAffineVector; radius:single);
169
// ------------------------------------------------------------------
170
// ------------------------------------------------------------------
171
// ------------------------------------------------------------------
174
// ------------------------------------------------------------------
175
// ------------------------------------------------------------------
176
// ------------------------------------------------------------------
178
// ----------------------------------------------------------------------
179
// Name : CheckPointInSphere()
180
// Input : point - point we wish to check for inclusion
181
// sO - Origin of sphere
182
// sR - radius of sphere
184
// Return: TRUE if point is in sphere, FALSE if not.
185
// -----------------------------------------------------------------------
187
function CheckPointInSphere(const Point, SO: TVector; const SR: Single)
190
// Allow small margin of error
191
Result := (VectorDistance2(Point, SO) <= Sqr(SR));
194
// ----------------------------------------------------------------------
195
// Name : CheckPointInTriangle()
196
// Input : point - point we wish to check for inclusion
197
// a - first vertex in triangle
198
// b - second vertex in triangle
199
// c - third vertex in triangle
200
// Notes : Triangle should be defined in clockwise order a,b,c
201
// Return: TRUE if point is in triangle, FALSE if not.
202
// -----------------------------------------------------------------------
204
function CheckPointInTriangle(Point, A, B, C: TAffineVector): Boolean;
206
Total_angles: Single;
207
V1, V2, V3: TAffineVector;
211
// make the 3 vectors
212
V1 := VectorSubtract(Point, A);
213
V2 := VectorSubtract(Point, B);
214
V3 := VectorSubtract(Point, C);
220
Total_angles := Total_angles + Arccos(VectorDotProduct(V1, V2));
221
Total_angles := Total_angles + Arccos(VectorDotProduct(V2, V3));
222
Total_angles := Total_angles + Arccos(VectorDotProduct(V3, V1));
224
if (Abs(Total_angles - 2 * PI) <= 0.005) then
230
// ----------------------------------------------------------------------
231
// Name : ClosestPointOnLine()
232
// Input : a - first end of line segment
233
// b - second end of line segment
234
// p - point we wish to find closest point on line from
235
// Notes : Helper function for closestPointOnTriangle()
236
// Return: closest point on line segment
237
// -----------------------------------------------------------------------
239
function ClosestPointOnLine(const A, B, P: TAffineVector): TAffineVector;
242
C, V: TAffineFLTVector;
244
VectorSubtract(P, A, C);
245
VectorSubtract(B, A, V);
247
D := VectorLength(V);
249
T := VectorDotProduct(V, C);
251
// Check to see if t is beyond the extents of the line segment
261
Result := VectorAdd(A, V);
265
// ----------------------------------------------------------------------
266
// Name : ClosestPointOnTriangle()
267
// Input : a - first vertex in triangle
268
// b - second vertex in triangle
269
// c - third vertex in triangle
270
// p - point we wish to find closest point on triangle from
272
// Return: closest point on triangle
273
// -----------------------------------------------------------------------
275
function ClosestPointOnTriangle(const a, b, c, n, p: TAffineVector): TAffineVector;
277
dAB, dBC, dCA : Single;
278
Rab, Rbc, Rca, intPoint : TAffineFLTVector;
281
//this would be faster if RayCastTriangleIntersect detected backwards hits
282
hit:=RayCastTriangleIntersect(VectorMake(p),VectorMake(n),a,b,c,@intPoint) or
283
RayCastTriangleIntersect(VectorMake(p),VectorMake(VectorNegate(n)),a,b,c,@intPoint);
290
Rab:=ClosestPointOnLine(a, b, p);
291
Rbc:=ClosestPointOnLine(b, c, p);
292
Rca:=ClosestPointOnLine(c, a, p);
294
dAB:=VectorDistance2(p, Rab);
295
dBC:=VectorDistance2(p, Rbc);
296
dCA:=VectorDistance2(p, Rca);
309
// ----------------------------------------------------------------------
310
// Name : ClosestPointOnTriangleEdge()
311
// Input : a - first vertex in triangle
312
// b - second vertex in triangle
313
// c - third vertex in triangle
314
// p - point we wish to find closest point on triangle from
316
// Return: closest point on line triangle edge
317
// -----------------------------------------------------------------------
319
function ClosestPointOnTriangleEdge(const A, B, C, P: TAffineVector)
322
DAB, DBC, DCA: Single;
323
Rab, Rbc, Rca: TAffineFLTVector;
325
Rab := ClosestPointOnLine(A, B, P);
326
Rbc := ClosestPointOnLine(B, C, P);
327
Rca := ClosestPointOnLine(C, A, P);
329
DAB := VectorDistance2(P, Rab);
330
DBC := VectorDistance2(P, Rbc);
331
DCA := VectorDistance2(P, Rca);
338
else if DCA < DAB then
346
function HitBoundingBox(const MinB, MaxB: TAffineFLTVector;
347
const Origin, Dir: TVector; var Coord: TVector): BOOLEAN;
354
I, Whichplane: Integer;
356
Quadrant: array [0 .. NUMDIM] of Byte;
357
MaxT: array [0 .. NUMDIM] of Double;
358
CandidatePlane: array [0 .. NUMDIM] of Double;
363
// Find candidate planes; this loop can be avoided if
364
// rays cast all from the eye(assume perpsective view)
365
for I := 0 to NUMDIM do
367
if (Origin.V[I] < MinB.V[I]) then
370
CandidatePlane[I] := MinB.V[I];
373
else if (Origin.V[I] > MaxB.V[I]) then
375
Quadrant[I] := RIGHT;
376
CandidatePlane[I] := MaxB.V[I];
380
Quadrant[I] := MIDDLE;
383
// * Ray origin inside bounding box */
386
SetVector(Coord, Origin);
391
// * Calculate T distances to candidate planes */
392
for I := 0 to NUMDIM do
394
if (Quadrant[I] <> MIDDLE) AND (Dir.V[I] <> 0) then
395
MaxT[I] := (CandidatePlane[I] - Origin.V[I]) / Dir.V[I]
400
// * Get largest of the maxT's for final choice of intersection */
402
for I := 1 to NUMDIM do
403
if (MaxT[WhichPlane] < MaxT[I]) then
406
// * Check final candidate actually inside box */
407
if (MaxT[WhichPlane] < 0) then
413
for I := 0 to NUMDIM do
415
if WhichPlane <> I then
417
Coord.V[I] := Origin.V[I] + MaxT[WhichPlane] * Dir.V[I];
418
if (Coord.V[I] < MinB.V[I]) OR (Coord.V[I] > MaxB.V[I])
426
Coord.V[I] := CandidatePlane[I];
429
Result := TRUE; // * ray hits box */
434
USE_EPSILON_TEST = TRUE;
439
function Coplanar_tri_tri(const N, V0, V1, V2, U0, U1,
440
U2: TAffineFLTVEctor): Integer;
445
function EDGE_AGAINST_TRI_EDGES(const V0, V1, U0, U1,
446
U2: TAffineFLTVector): Integer;
448
Ax, Ay, Bx, By, Cx, Cy, E, D, F: Single;
450
// * this edge to edge test is based on Franlin Antonio's gem:
451
// "Faster Line Segment Intersection", in Graphics Gems III,
453
function EDGE_EDGE_TEST(const V0, U0, U1: TAffineFLTVector): Integer;
456
Bx := U0.V[I0] - U1.V[I0];
457
By := U0.V[I1] - U1.V[I1];
458
Cx := V0.V[I0] - U0.V[I0];
459
Cy := V0.V[I1] - U0.V[I1];
460
F := Ay * Bx - Ax * By;
461
D := By * Cx - Bx * Cy;
462
if ((F > 0) and (D >= 0) and (D <= F)) or
463
((F < 0) and (D <= 0) and (D >= F)) then
465
E := Ax * Cy - Ay * Cx;
468
if (E >= 0) and (E <= F) then
471
else if (E <= 0) and (E >= F) then
477
Ax := V1.V[I0] - V0.V[I0];
478
Ay := V1.V[I1] - V0.V[I1];
479
// * test edge U0,U1 against V0,V1 */
480
Result := EDGE_EDGE_TEST(V0, U0, U1);
483
// * test edge U1,U2 against V0,V1 */
484
Result := EDGE_EDGE_TEST(V0, U1, U2);
487
// * test edge U2,U1 against V0,V1 */
488
Result := EDGE_EDGE_TEST(V0, U2, U0);
491
function POINT_IN_TRI(const V0, U0, U1, U2: TAffineFLTVector): Integer;
493
A, B, C, D0, D1, D2: Single;
496
// * is T1 completly inside T2? */
497
// * check if V0 is inside tri(U0,U1,U2) */
498
A := U1.V[I1] - U0.V[I1];
499
B := -(U1.V[I0] - U0.V[I0]);
500
C := -A * U0.V[I0] - B * U0.V[I1];
501
D0 := A * V0.V[I0] + B * V0.V[I1] + C;
503
A := U2.V[I1] - U1.V[I1];
504
B := -(U2.V[I0] - U1.V[I0]);
505
C := -A * U1.V[I0] - B * U1.V[I1];
506
D1 := A * V0.V[I0] + B * V0.V[I1] + C;
508
A := U0.V[I1] - U2.V[I1];
509
B := -(U0.V[I0] - U2.V[I0]);
510
C := -A * U2.V[I0] - B * U2.V[I1];
511
D2 := A * V0.V[I0] + B * V0.V[I1] + C;
512
if (D0 * D1 > 0.0) then
513
if (D0 * D2 > 0.0) then
517
/// Begin Main logic ///////////////////////////////
519
// * first project onto an axis-aligned plane, that maximizes the area */
520
// * of the triangles, compute indices: i0,i1. */
521
A.V[0] := Abs(N.V[0]);
522
A.V[1] := Abs(N.V[1]);
523
A.V[2] := Abs(N.V[2]);
524
if (A.V[0] > A.V[1]) then
526
if (A.V[0] > A.V[2]) then
528
I0 := 1; // * A[0] is greatest */
533
I0 := 0; // * A[2] is greatest */
538
begin // * A[0]<=A[1] */
539
if (A.V[2] > A.V[1]) then
541
I0 := 0; // * A[2] is greatest */
546
I0 := 0; // * A[1] is greatest */
551
// * test all edges of triangle 1 against the edges of triangle 2 */
552
Result := EDGE_AGAINST_TRI_EDGES(V0, V1, U0, U1, U2);
555
Result := EDGE_AGAINST_TRI_EDGES(V1, V2, U0, U1, U2);
558
Result := EDGE_AGAINST_TRI_EDGES(V2, V0, U0, U1, U2);
562
// * finally, test if tri1 is totally contained in tri2 or vice versa */
563
Result := POINT_IN_TRI(V0, U0, U1, U2);
566
Result := POINT_IN_TRI(U0, V0, V1, V2);
571
function Tri_tri_intersect(const V0, V1, V2, U0, U1,
572
U2: TAFFineFLTVector): Integer;
574
E1, E2: TAffineFLTVector;
575
N1, N2: TAffineFLTVector;
577
Du0, Du1, Du2, Dv0, Dv1, Dv2: Single;
579
Isect1: array [0 .. 1] of Single;
580
Isect2: array [0 .. 1] of Single;
581
Du0du1, Du0du2, Dv0dv1, Dv0dv2: Single;
583
Vp0, Vp1, Vp2: Single;
584
Up0, Up1, Up2: Single;
587
procedure ISECT(VV0, VV1, VV2, D0, D1, D2: Single;
588
var Isect0, Isect1: Single);
590
Isect0 := VV0 + (VV1 - VV0) * D0 / (D0 - D1);
591
Isect1 := VV0 + (VV2 - VV0) * D0 / (D0 - D2);
594
function COMPUTE_INTERVALS(VV0, VV1, VV2, D0, D1, D2, D0D1, D0D2: Single;
595
var Isect0, Isect1: Single): Integer;
599
// * here we know that D0D2<=0.0 */
600
// * that is D0, D1 are on the same side, D2 on the other or on the plane */ \
601
ISECT(VV2, VV0, VV1, D2, D0, D1, Isect0, Isect1)
602
else if (D0D2 > 0.0) then
603
// * here we know that d0d1<=0.0 */
604
ISECT(VV1, VV0, VV2, D1, D0, D2, Isect0, Isect1)
605
else if (D1 * D2 > 0.0) or (D0 <> 0.0) then
606
// * here we know that d0d1<=0.0 or that D0!=0.0 */
607
ISECT(VV0, VV1, VV2, D0, D1, D2, Isect0, Isect1)
608
else if (D1 <> 0.0) then
609
ISECT(VV1, VV0, VV2, D1, D0, D2, Isect0, Isect1)
610
else if (D2 <> 0.0) then
611
ISECT(VV2, VV0, VV1, D2, D0, D1, Isect0, Isect1)
613
// * triangles are coplanar */
614
Result := Coplanar_tri_tri(N1, V0, V1, V2, U0, U1, U2);
617
// * sort so that a<=b */
618
procedure SORT(var A: Single; var B: Single);
631
// * compute plane equation of triangle(V0,V1,V2) */
632
E1 := VectorSubtract(V1, V0);
633
E2 := VectorSubtract(V2, V0);
634
N1 := VectorCrossProduct(E1, E2);
635
D1 := -VectorDotProduct(N1, V0);
636
// * plane equation 1: N1.X+d1=0 */
638
// * put U0,U1,U2 into plane equation 1 to compute signed distances to the plane*/
639
Du0 := VectorDotProduct(N1, U0) + D1;
640
Du1 := VectorDotProduct(N1, U1) + D1;
641
Du2 := VectorDotProduct(N1, U2) + D1;
643
// * coplanarity robustness check */
644
if USE_EPSILON_TEST = TRUE then
646
if (Abs(Du0) < EPSILON) then
648
if (Abs(Du1) < EPSILON) then
650
if (Abs(Du2) < EPSILON) then
656
if (Du0du1 > 0.0) and (Du0du2 > 0.0) then
657
begin // * same sign on all of them + not equal 0 ? */
658
Result := 0; // * no intersection occurs */
662
// * compute plane of triangle (U0,U1,U2) */
663
E1 := VectorSubtract(U1, U0);
664
E2 := VectorSubtract(U2, U0);
665
N2 := VectorCrossProduct(E1, E2);
666
D2 := -VectorDotProduct(N2, U0);
667
// * plane equation 2: N2.X+d2=0 */
669
// * put V0,V1,V2 into plane equation 2 */
670
Dv0 := VectorDotProduct(N2, V0) + D2;
671
Dv1 := VectorDotProduct(N2, V1) + D2;
672
Dv2 := VectorDotProduct(N2, V2) + D2;
674
if USE_EPSILON_TEST = TRUE then
676
if (Abs(Dv0) < EPSILON) then
678
if (Abs(Dv1) < EPSILON) then
680
if (Abs(Dv2) < EPSILON) then
687
if (Dv0dv1 > 0.0) and (Dv0dv2 > 0.0) then
688
begin // * same sign on all of them + not equal 0 ? */
689
Result := 0; // * no intersection occurs */
693
// * compute direction of intersection line */
694
D := VectorCrossProduct(N1, N2);
696
// * compute and index to the largest component of D */
711
// * this is the simplified projection onto L*/
720
// * compute interval for triangle 1 */
721
COMPUTE_INTERVALS(Vp0, Vp1, Vp2, Dv0, Dv1, Dv2, Dv0dv1, Dv0dv2, Isect1[0],
724
// * compute interval for triangle 2 */
725
COMPUTE_INTERVALS(Up0, Up1, Up2, Du0, Du1, Du2, Du0du1, Du0du2, Isect2[0],
728
SORT(Isect1[0], Isect1[1]);
729
SORT(Isect2[0], Isect2[1]);
731
if (Isect1[1] < Isect2[0]) or (Isect2[1] < Isect1[0]) then
738
// ------------------ TOctree ------------------
786
// Theory on FlagMax and FlagMin:
787
// When a node is subdivided, each of the 8 children assumes 1/8th ownership of its
788
// parent's bounding box (defined by parent extents). Calculating a child's min/max
789
// extent only requires 3 values: the parent's min extent, the parent's max extent
790
// and the midpoint of the parent's extents (since the cube is divided in half twice).
791
// The following arrays assume that the children are numbered from 0 to 7, named Upper
792
// and Lower (Upper = top 4 cubes on Y axis, Bottom = lower 4 cubes), Left and Right, and
793
// Fore and Back (Fore facing furthest away from you the viewer).
794
// Each node can use its corresponding element in the array to flag the operation needed
795
// to find its new min/max extent. Note that min, mid and max refer to an array of
796
// 3 coordinates (x,y,z); each of which are flagged separately. Also note that these
797
// flags are based on the Y vector being the up vector.
799
FlagMax: array [0 .. 7] of array [0 .. 2] of Byte = ((MID, MAX, MAX),
801
(MAX, MAX, MAX), // Upper Fore Right
802
(MID, MAX, MID), // Upper Back Left
803
(MAX, MAX, MID), // Upper Back Right
805
(MID, MID, MAX), // Lower Fore Left (similar to above except height/2)
806
(MAX, MID, MAX), // Lower Fore Right
807
(MID, MID, MID), // Lower Back Left
808
(MAX, MID, MID) // Lower Back Right
811
FlagMin: array [0 .. 7] of array [0 .. 2] of Byte = ((MIN, MID, MID),
813
(MID, MID, MID), // Upper Fore Right
814
(MIN, MID, MIN), // Upper Back Left
815
(MID, MID, MIN), // Upper Back Right
817
(MIN, MIN, MID), // Lower Fore Left (similar to above except height/2)
818
(MID, MIN, MID), // Lower Fore Right
819
(MIN, MIN, MIN), // Lower Back Left
820
(MID, MIN, MIN) // Lower Back Right
823
// Design of the AABB faces, using similar method to above.. Note than normals are not
824
// correct, but not needed for current tri-intersection test.
825
// Could be removed if the tri-plane collision is replaced with a tri-box routine.
826
FlagFaces: array [0 .. 23] of array [0 .. 2] of Byte = (
828
(MIN, MAX, MAX), // Upper left corner
829
(MAX, MAX, MAX), // Upper right corner
830
(MAX, MIN, MAX), // Bottom right corner
834
(MIN, MAX, MIN), // Upper left corner
835
(MAX, MAX, MIN), // Upper right corner
836
(MAX, MIN, MIN), // Bottom right corner
840
(MIN, MAX, MAX), // Upper left corner
841
(MAX, MAX, MAX), // Upper right corner
842
(MAX, MAX, MIN), // Bottom right corner
846
(MIN, MIN, MAX), // Upper left corner
847
(MAX, MIN, MAX), // Upper right corner
848
(MAX, MIN, MIN), // Bottom right corner
852
(MIN, MAX, MAX), // Upper left corner
853
(MIN, MIN, MAX), // Upper right corner
854
(MIN, MIN, MIN), // Bottom right corner
858
(MAX, MIN, MAX), // Upper left corner
859
(MAX, MAX, MAX), // Upper right corner
860
(MAX, MAX, MIN), // Bottom right corner
865
destructor TOctree.Destroy;
873
procedure TOctree.DisposeTree;
875
procedure WalkDispose(var Node: POctreeNode);
879
if Assigned(Node) then
882
WalkDispose(Node^.ChildArray[I]);
888
WalkDispose(RootNode);
891
TriangleFiler := nil;
896
procedure TOctree.CreateTree(Depth: Integer);
898
MaxOlevel := Depth; // initialize max depth.
902
// "cuts" all the triangles in the mesh into the octree.
903
procedure TOctree.CutMesh;
905
procedure AddTriangleToNodes(N: Integer);
910
for X := 0 to High(ResultArray) do
912
P := Resultarray[X]; // Pointer to a node.
914
K := Length(P^.TriArray);
915
SetLength(P^.TriArray, K + 1); // Increase array by 1.
916
P^.TriArray[K] := N; // Store triangle # reference.
925
N: Integer; // n = triangle # in mesh
927
TriCountMesh := TriangleFiler.Count div 3;
929
while N < TriangleFiler.Count do
931
WalkTriToLeaf(RootNode, TriangleFiler.List^[N], TriangleFiler.List^[N + 1],
932
TriangleFiler.List^[N + 2]);
933
if ResultArray <> NIL then
935
AddTriangleToNodes(N);
936
Inc(TriCountOctree, 1);
942
function TOctree.GetMidPoint(Min, Max: Single): Single;
944
Result := Max / 2 + Min / 2;
945
// This formula is non-quadrant specific; ie: good.
948
function TOctree.GetExtent(const Flags: array of Byte; ParentNode: POctreeNode)
951
Emin, Emax: TAffineFLTVector;
954
Emin := ParentNode^.MinExtent; // Some easy to work with variables.
955
Emax := ParentNode^.MaxExtent;
961
Result.V[N] := Emin.V[N];
963
Result.V[N] := GetMidPoint(Emin.V[N], Emax.V[N]);
965
Result.V[N] := Emax.V[N];
972
procedure TOctree.InitializeTree(const AWorldMinExtent, AWorldMaxExtent
973
: TAffineVector; const ATriangles: TAffineVectorList;
974
const ATreeDepth: Integer);
977
Newnode: POctreeNode;
979
Self.WorldMinExtent := AWorldMinExtent;
980
Self.WorldMaxExtent := AWorldMaxExtent;
982
// set up the filer data for this mesh
983
if TriangleFiler = nil then
984
TriangleFiler := TAffineVectorList.Create;
985
TriangleFiler.Assign(ATriangles);
988
Newnode^.MinExtent := AWorldMinExtent;
989
Newnode^.MaxExtent := AWorldMaxExtent;
990
Newnode^.TriArray := NIL;
992
Newnode^.ChildArray[N] := NIL;
994
// Initialize work variables for new tree.
995
Rootnode := Newnode; // rootnode always points to root.
996
NodeCount := 0; // initialize node count
998
CreateTree(ATreeDepth);
1004
procedure TOctree.Refine(ParentNode: POctreeNode; Level: Integer);
1007
Pwork: array [0 .. 7] of POctreeNode;
1008
// Stores addresses of newly created children.
1009
Newnode: POctreeNode;
1011
if Level < MaxOlevel then
1014
begin // Create 8 new children under this parent.
1017
Pwork[N] := Newnode; // Create work pointers for the next for loop.
1019
// Generate new extents based on parent's extents
1020
Newnode^.MinExtent := GetExtent(FlagMin[N], ParentNode);
1021
Newnode^.MaxExtent := GetExtent(FlagMax[N], ParentNode);
1023
Newnode^.TriArray := nil; // Initialize.
1026
Newnode^.ChildArray[Z] := nil; // initialize all child pointers to NIL
1028
ParentNode^.ChildArray[N] := Newnode;
1029
// initialize parent's child pointer to this node
1031
for X := 0 to 7 do // Now recursively Refine each child we just made
1032
Refine(Pwork[X], Level + 1);
1038
procedure TOctree.ConvertR4(ONode: POctreeNode; const Scale: TAffineFLTVector);
1042
ScaleVector(Onode^.MinExtent, Scale);
1043
ScaleVector(Onode^.MaxExtent, Scale);
1044
if ONode^.ChildArray[0] <> NIL then
1045
begin // ie: if not a leaf then loop through children.
1048
ConvertR4(Onode^.ChildArray[N], Scale);
1055
function TOctree.PointInNode(const Min, Max, APoint: TAffineFLTVector): BOOLEAN;
1057
Result := (APoint.V[0] >= Min.V[0]) and
1058
(APoint.V[1] >= Min.V[1]) and (APoint.V[2] >= Min.V[2]) and
1059
(APoint.V[0] <= Max.V[0]) and (APoint.V[1] <= Max.V[1]) and
1060
(APoint.V[2] <= Max.V[2]);
1065
procedure TOctree.WalkPointToLeaf(Onode: POctreeNode; const P: TAffineVector);
1067
Finalize(Resultarray);
1068
WalkPointToLeafx(Onode, P);
1073
procedure TOctree.WalkPointToLeafx(Onode: POctreeNode; const P: TAffineVector);
1077
if PointInNode(Onode^.MinExtent, Onode^.MaxExtent, P) then
1079
if Assigned(Onode^.ChildArray[0]) then
1081
WalkPointToLeafx(Onode^.ChildArray[N], P)
1084
SetLength(Resultarray, Length(Resultarray) + 1);
1085
Resultarray[High(Resultarray)] := Onode;
1092
function TOctree.SphereInNode(const MinExtent, MaxExtent: TAffineVector;
1093
const C: TVector; Radius: Single): Boolean;
1094
// Sphere-AABB intersection by Miguel Gomez
1099
// find the square of the distance
1100
// from the sphere to the box
1104
if (C.V[I] < MinExtent.V[I]) then
1106
S := C.V[I] - MinExtent.V[I];
1109
else if (C.V[I] > MaxExtent.V[I]) then
1111
S := C.V[I] - MaxExtent.V[I];
1116
if D <= Radius * Radius then
1124
procedure TOctree.WalkSphereToLeaf(Onode: POctreeNode; const P: TVector;
1127
Finalize(Resultarray);
1128
WalkSphereToLeafx(Onode, P, Radius);
1133
procedure TOctree.WalkSphereToLeafx(Onode: POctreeNode; const P: TVector;
1138
if SphereInNode(Onode^.MinExtent, Onode^.MaxExtent, P, Radius) then
1140
if Assigned(Onode^.ChildArray[0]) then
1142
WalkSphereToLeafx(Onode^.ChildArray[N], P, Radius)
1145
SetLength(Resultarray, Length(Resultarray) + 1);
1146
Resultarray[High(Resultarray)] := Onode;
1151
// Cast a ray (point p, vector v) into the Octree (ie: ray-box intersect).
1152
procedure TOctree.WalkRayToLeaf(Onode: POctreeNode; const P, V: TVector);
1154
Finalize(Resultarray);
1156
WalkRayToLeafx(Onode, P, V);
1161
procedure TOctree.WalkRayToLeafx(Onode: POctreeNode; const P, V: TVector);
1166
if HitBoundingBox(Onode^.MinExtent, Onode^.MaxExtent, P, V, Coord) then
1168
if Assigned(Onode^.ChildArray[0]) then
1170
WalkRayToLeafx(Onode^.ChildArray[N], P, V)
1173
SetLength(Resultarray, Length(Resultarray) + 1);
1174
Resultarray[High(Resultarray)] := Onode;
1179
// Check triangle intersect with any of the node's faces.
1180
// Could be replaced with a tri-box check.
1181
function TOctree.TriIntersectNode(const MinExtent, MaxExtent, V1, V2,
1182
V3: TAffineVector): Boolean;
1184
F0, F1, F2, F3: TAffineFLTVector;
1186
AFace: array [0 .. 3] of TAffineFLTVector; // A face's 4 corners.
1189
begin // Do all 6 faces.
1191
begin // Do all 4 vertices for the face.
1193
begin // Do x,y,z for each vertex.
1194
if FlagFaces[O + N * 4, P] = MIN then
1195
AFace[O].V[P] := MinExtent.V[P]
1197
AFace[O].V[P] := MaxExtent.V[P];
1205
// Now check the two triangles in the face against the mesh triangle.
1206
if Tri_tri_intersect(V1, V2, V3, F0, F1, F2) = 1 then
1208
else if Tri_tri_intersect(V1, V2, V3, F2, F3, F0) = 1 then
1220
procedure TOctree.WalkTriToLeaf(Onode: POctreeNode;
1221
const V1, V2, V3: TAffineFLTVector);
1223
Finalize(Resultarray);
1224
WalkTriToLeafx(Onode, V1, V2, V3);
1229
procedure TOctree.WalkTriToLeafx(Onode: POctreeNode;
1230
const V1, V2, V3: TAffineFLTVector);
1234
if TriIntersectNode(Onode^.MinExtent, Onode^.MaxExtent, V1, V2, V3) or
1235
PointInNode(Onode^.MinExtent, Onode^.MaxExtent, V1) or
1236
PointInNode(Onode^.MinExtent, Onode^.MaxExtent, V2) or
1237
PointInNode(Onode^.MinExtent, Onode^.MaxExtent, V3) then
1239
if Onode^.ChildArray[0] <> NIL then
1241
WalkTriToLeafx(Onode^.ChildArray[M], V1, V2, V3)
1244
SetLength(Resultarray, Length(Resultarray) + 1);
1245
Resultarray[High(Resultarray)] := Onode;
1250
// RayCastIntersectAABB
1252
function TOctree.RayCastIntersect(const RayStart, RayVector: TVector;
1253
IntersectPoint: PVector = nil; IntersectNormal: PVector = nil;
1254
TriangleInfo: POctreeTriangleInfo = nil): Boolean;
1256
CInitialMinD: Single = 1E40;
1261
IPoint, INormal: TVector;
1263
WalkRayToLeaf(RootNode, RayStart, RayVector);
1265
if Resultarray = nil then
1271
MinD := CInitialMinD;
1272
for I := 0 to High(Resultarray) do
1274
P := ResultArray[I];
1275
for T := 0 to High(P^.TriArray) do
1277
K := P^.Triarray[T];
1278
if RayCastTriangleIntersect(RayStart, RayVector, TriangleFiler.List^[K],
1279
TriangleFiler.List^[K + 1], TriangleFiler.List^[K + 2], @IPoint,
1282
D := VectorDistance2(RayStart, IPoint);
1286
if IntersectPoint <> nil then
1287
IntersectPoint^ := IPoint;
1288
if IntersectNormal <> nil then
1289
IntersectNormal^ := INormal;
1290
if TriangleInfo <> nil then
1292
TriangleInfo^.Index := K;
1293
TriangleInfo^.Vertex[0] := TriangleFiler.List^[K];
1294
TriangleInfo^.Vertex[1] := TriangleFiler.List^[K + 1];
1295
TriangleInfo^.Vertex[2] := TriangleFiler.List^[K + 2];
1301
Result := (MinD <> CInitialMinD);
1304
// SphereIntersectAABB -- Walk a sphere through an Octree, given a velocity, and return the nearest polygon
1305
// intersection point on that sphere, and its plane normal.
1307
// **** This function is the result of a LOT of work and experimentation with both Paul Nettle's method
1308
// (www.fluidstudios.com) and Telemachos' method (www.peroxide.dk) over a period of about 2 months. If
1309
// you find ways to optimize the general structure, please let me know at rmch@cadvision.com. ****
1311
// TO DO: R4 conversion (routine already exists for this) for ellipsoid space.
1313
// Method for caclulating CD vs polygon: (Robert Hayes method)
1314
// ...for each triangle:
1315
// 1. Cast a ray from sphere origin to triangle's plane along plane's negative normal (set to length
1316
// of sphere radius). If no hit, skip this triangle. Otherwise this is the default plane
1317
// intersection point.
1318
// 2. If the distance is =< the sphere radius, the plane is embedded in the sphere. Go to step 6 with
1319
// plane intersection point from above step.
1320
// 3. If the distance > sphere radius, calculate the sphere intersection point to this plane by
1321
// subtracting the plane's normal from the sphere's origin.
1322
// 4. Cast a new ray from the sphere intersection point to the plane; this is the new plane
1323
// intersection point.
1324
// 5. Cast a ray from the sphere intersection point to the triangle. If it is a direct hit, then
1325
// this point is the polygon intersection point.
1326
// 6. Else, find the point on the triangle that is closest to the plane intersection point; this becomes
1327
// the polygon intersection point (ie: hit an edge or corner)
1328
// 7. Cast a ray from the polygon intersection point along the negative velocity vector of the sphere
1329
// (note: for accuracy reasons - SphereIntersect replaced with PointInSphere)
1330
// 8. If there is no intersection, the sphere cannot possibly collide with this triangle
1331
// 9. Else, save the distance from step 8 if, and only if, it is the shortest collision distance so far.
1333
// Return the polygon intersection point and the closest triangle's normal if hit.
1335
function TOctree.SphereSweepIntersect(const RayStart, RayVector: TVector;
1336
const Velocity, Radius: Single; IntersectPoint: PVector = nil;
1337
IntersectNormal: PVector = nil): Boolean;
1339
CInitialMinD2: Single = 1E40;
1340
CEpsilon: Single = 0.05;
1344
MinD2, Sd2, Radius2: Single;
1345
DistanceToTravel, DistanceToTravelMinusRadius2: Single;
1347
PNormal: TAffineVector;
1349
NEGpNormal4: TVector;
1350
SIPoint, SIPoint2: TVector; // sphere intersection point
1351
PIPoint: TVector; // plane intersection point
1352
PolyIPoint: TVector; // polygon intersection point
1353
NEGVelocity: TVector; // sphere's forward velocity
1356
P1, P2, P3: PAffineVector;
1358
// SphereSweepAABB:TAABB;
1360
// response identifiers (for future use)
1361
// destinationPoint, newdestinationPoint: TVector;
1362
// slidePlaneOrigin, slidePlaneNormal: TVector;
1363
// newvelocityVector: TVector;
1367
// Note: Current implementation only accurate for FreeForm:Radius at 1:1 ratio.
1369
Result := False; // default: no collision
1371
// quit if no movement
1372
if (Velocity = 0) or (not(VectorNorm(RayVector) > 0)) then
1374
// How far ahead to check for collisions.
1375
DistanceToTravel := Velocity + Radius + CEpsilon;
1376
DistanceToTravelMinusRadius2 := Sqr(Velocity + CEpsilon);
1377
Radius2 := Sqr(Radius);
1379
// Determine all the octree nodes this sphere intersects with.
1380
// NOTE: would be more efficient to find the bounding box that includes the
1381
// startpoint and endpoint (+sphere diameter)... especially with a large sphere
1382
// and/or a large velocity.
1383
WalkSphereToLeaf(RootNode, RayStart, DistanceToTravel);
1385
// This method may be more effective if sphere sweeps from one point to another and stops.
1386
// NOTE: If it recursively slides of planes, then WalkSphereToLeaf would probably be better, as
1387
// it will determine all possible triangles that could intersect over the whole motion
1388
// AABBFromSweep(SphereSweepAABB,rayStart,destinationPoint,Radius+cEpsilon);
1389
// GetTrianglesFromNodesIntersectingAABB(SphereSweepAABB);
1391
if not Assigned(Resultarray) then
1394
// Negative velocity vector for use with ray-sphere check.
1395
VectorScale(RayVector, -Velocity / VectorLength(RayVector), NEGVelocity);
1397
MinD2 := CInitialMinD2;
1398
for I := 0 to High(Resultarray) do
1400
P := ResultArray[I];
1401
for T := 0 to High(P^.TriArray) do
1403
K := P^.Triarray[T];
1404
// These are the vertices of the triangle to check
1405
P1 := @Trianglefiler.List[K];
1406
P2 := @Trianglefiler.List[K + 1];
1407
P3 := @Trianglefiler.List[K + 2];
1409
// Create the normal for this triangle
1410
PNormal := CalcPlaneNormal(P1^, P2^, P3^);
1412
// Ignore backfacing planes
1413
if VectorDotProduct(PNormal, PAffineVector(@RayVector)^) > 0.0 then
1416
// Set the normal to the radius of the sphere
1417
ScaleVector(PNormal, Radius);
1418
// Make some TVectors
1419
MakeVector(PNormal4, PNormal);
1420
NEGpNormal4 := VectorNegate(PNormal4);
1422
// Calculate the plane intersection point (sphere origin to plane)
1423
if RayCastPlaneIntersect(RayStart, NEGpNormal4, VectorMake(P1^), PNormal4,
1427
Sd2 := VectorDistance2(RayStart, PIPoint);
1429
// If sd <= radius, fall through to "not direct hit" code below with pIPoint
1430
// as the plane intersection point. Sphere is embedded.
1432
if Sd2 > Radius2 then
1434
// Calculate sphere intersection point (ie: point on sphere that hits plane)
1435
SetVector(SIPoint, VectorSubtract(RayStart, PNormal4));
1436
// Get new plane intersection point (in case not a direct hit)
1437
RayCastPlaneIntersect(SIPoint, RayVector, VectorMake(P1^), PNormal4,
1439
// Does the velocity vector directly hit the triangle? If yes then that is the
1440
// polygon intersection point.
1441
if RayCastTriangleIntersect(SIPoint, RayVector, P1^, P2^, P3^,
1442
@PolyIPoint, @PNormal4) then
1444
Sd2 := VectorDistance2(SIPoint, PolyIPoint);
1449
// If not a direct hit then check if sphere "nicks" the triangle's edge or corner...
1450
// If it does then that is the polygon intersection point.
1451
if not DirectHit then
1453
if not CheckPointInTriangle(AffineVectorMake(PiPoint), P1^, P2^, P3^)
1455
SetVector(PolyIPoint, ClosestPointOnTriangleEdge(P1^, P2^, P3^,
1456
PAffineVector(@PIPoint)^), 1)
1458
PolyIPoint := PiPoint;
1460
// See if this point + negative velocity vector lies within the sphere.
1461
// (This implementation seems more accurate than RayCastSphereIntersect)
1462
{ if not CheckPointInSphere(VectorAdd(polyIPoint, NEGVelocity), raystart, radius) then1464
// sd2:=0; //stops the test too soon (noticable on triangle edges in flat planes)
1465
sd2:=sqr(VectorDistance(raystart, polyIPoint)-Radius);
1467
// see if this point + negative velocity vector intersect the sphere.
1468
// (PointInSphere doesn't work for fast motion)
1469
if VectorDistance2(PolyIPoint, RayStart) > Radius2 then
1471
if RayCastSphereIntersect(PolyIPoint, VectorNormalize(NEGVelocity),
1472
RayStart, Radius, SIPoint, SIPoint2) = 0 then
1474
Sd2 := VectorDistance2(SIPoint, PolyIPoint);
1480
// Allow sphere to get close to triangle (less epsilon which is built into distanceToTravel)
1481
if Sd2 <= DistanceToTravelMinusRadius2 then
1483
Result := True; // flag a collision
1487
if IntersectPoint <> nil then
1488
IntersectPoint^ := PolyIPoint;
1489
if IntersectNormal <> nil then
1490
IntersectNormal^ := PNormal4;
1496
end; // end for t triangles
1497
end; // end for i nodes
1500
function TOctree.TriangleIntersect(const V1, V2, V3: TAffineVector): Boolean;
1504
P1, P2, P3: PAffineVector;
1506
Result := False; // default: no collision
1507
WalkTriToLeaf(RootNode, V1, V2, V3);
1508
if not Assigned(Resultarray) then
1511
for I := 0 to High(Resultarray) do
1513
P := ResultArray[I];
1514
for T := 0 to High(P^.TriArray) do
1516
K := P^.Triarray[T];
1517
// These are the vertices of the triangle to check
1518
P1 := @Trianglefiler.List[K];
1519
P2 := @Trianglefiler.List[K + 1];
1520
P3 := @Trianglefiler.List[K + 2];
1521
if Tri_tri_intersect(V1, V2, V3, P1^, P2^, P3^) <> 0 then
1526
end; // end for t triangles
1527
end; // end for i nodes
1532
function TOctree.AABBIntersect(const AABB: TAABB; M1to2, M2to1: TMatrix;
1533
Triangles: TAffineVectorList = nil): Boolean;
1535
TriList: TAffineVectorList;
1538
// get triangles in nodes intersected by the aabb
1539
TriList := GetTrianglesFromNodesIntersectingCube(Aabb, M1to2, M2to1);
1542
if Trilist.Count > 0 then
1544
Trilist.TransformAsPoints(M2to1);
1546
// run all faces and check if they're inside the aabb
1547
// uncoment the * and comment the {//} to check only vertices1548
{ // } while I < TriList.Count - 1 do1550
// *for i:= 0 to triList.count -1 do begin
1551
// * v:=VectorMake(TriList.Items[i]);
1552
// * if pointInAABB(AffinevectorMake(v), aabb) then
1553
{ // } if TriangleIntersectAABB(Aabb, TriList[I], TriList[I + 1],1554
Trilist[I + 2]) then
1557
if not Assigned(Triangles) then
1560
Triangles.Add(TriList[I], TriList[I + 1], Trilist[I + 2]);
1569
// GetTrianglesFromNodesIntersectingAABB
1571
function TOctree.GetTrianglesFromNodesIntersectingAABB(const ObjAABB: TAABB)
1572
: TAffineVectorList;
1576
procedure HandleNode(Onode: POctreeNode);
1581
AABB2.Min := Onode^.MinExtent;
1582
AABB2.Max := Onode^.MaxExtent;
1584
if IntersectAABBsAbsolute(AABB1, AABB2) then
1586
if Assigned(Onode^.ChildArray[0]) then
1589
HandleNode(Onode^.ChildArray[I])
1593
SetLength(ResultArray, Length(ResultArray) + 1);
1594
ResultArray[High(ResultArray)] := Onode;
1602
TriangleIndices: TIntegerList;
1608
SetLength(ResultArray, 0);
1609
if Assigned(RootNode) then
1610
HandleNode(RootNode);
1612
Result := TAffineVectorList.Create;
1613
TriangleIndices := TIntegerList.Create;
1615
// fill the triangles from all nodes in the resultarray to AL
1616
for I := 0 to High(ResultArray) do
1618
P := ResultArray[I];
1619
TriangleIndices.AddIntegers(P^.TriArray);
1621
TriangleIndices.SortAndRemoveDuplicates;
1622
Result.Capacity := TriangleIndices.Count * 3;
1623
for I := 0 to TriangleIndices.Count - 1 do
1625
K := TriangleIndices[I];
1626
Result.Add(TriangleFiler.List^[K], TriangleFiler.List^[K + 1],
1627
TriangleFiler.List^[K + 2]);
1630
TriangleIndices.Free;
1635
// GetTrianglesFromNodesIntersectingCube
1637
function TOctree.GetTrianglesFromNodesIntersectingCube(const ObjAABB: TAABB;
1638
const ObjToSelf, SelfToObj: TMatrix): TAffineVectorList;
1641
M1To2, M2To1: TMatrix;
1643
procedure HandleNode(Onode: POctreeNode);
1648
AABB2.Min := Onode^.MinExtent;
1649
AABB2.Max := Onode^.MaxExtent;
1651
if IntersectAABBs(AABB1, AABB2, M1To2, M2To1) then
1653
if Assigned(Onode^.ChildArray[0]) then
1656
HandleNode(Onode^.ChildArray[I])
1660
SetLength(ResultArray, Length(ResultArray) + 1);
1661
ResultArray[High(ResultArray)] := Onode;
1669
TriangleIndices: TIntegerList;
1673
// Calc Conversion Matrixes
1677
SetLength(ResultArray, 0);
1678
if Assigned(RootNode) then
1679
HandleNode(RootNode);
1681
Result := TAffineVectorList.Create;
1682
TriangleIndices := TIntegerList.Create;
1684
// fill the triangles from all nodes in the resultarray to AL
1685
for I := 0 to High(ResultArray) do
1687
P := ResultArray[I];
1688
TriangleIndices.AddIntegers(P^.TriArray);
1690
TriangleIndices.SortAndRemoveDuplicates;
1691
Result.Capacity := TriangleIndices.Count * 3;
1692
for I := 0 to TriangleIndices.Count - 1 do
1694
K := TriangleIndices[I];
1695
Result.Add(TriangleFiler.List^[K], TriangleFiler.List^[K + 1],
1696
TriangleFiler.List^[K + 2]);
1699
TriangleIndices.Free;