Michael Lahanas and Kostas Karouzakis

In
medical applications three-dimensional (3D) anatomical models are
usually produced from contours from cross-sectional images.
Conventional algorithms are described which are used to generate
points inside these objects which are required for example for
volumetric purposes or for the calculation of the dose distributions
in radiotherapy or brachytherapy [LB00],
see Fig. 1. These algorithms use triangles obtained from a
triangulation of the contours. Using a preprocessing method we
generate from the contours so called 2D *virtual polygons*. This
enables us to reduce the 3D point in a polyhedron test into a 2D
point in a virtual polygon test. Additional we use optimized bounding
boxes [LK00] which further improves the efficiency of this method.

**Fig.
1** 100000 sampling Points in a part of the body which includes a
rib tumor with 50000 sampling points used for volumetric and
dosimetric purposes.

We
consider a 2D contour as a closed polygon made up of *N*
vertices (*x _{i},y_{i}*) i, Î
1,...,

**The
Intersection number test.**

In
this method we count the number of intersections of a ray with all
edges of the polygon. If the number of edges of the polygon
intersected by the ray is even then the point **P** is outside the
polygon, and if the number of intersections is odd then the point **P**
lies inside the polygon. This is a consequence of the Jordan curve
theorem in topology and a proof is given by Courant *et al*.
[CR41] There
are some special cases when a ray passes exactly through a vertex.
These situations can be considered using tie breaking rules or by a
method such as described by Edelsbrunner *et al*.
[ED90]

**Fig 2**
Schematic diagram illustrating a generalized point **P**
in relation to a 2D cross-section through a polygon (a) with and (b)
without a critical structure inside the cross-section. The number of
intersections is shown for each point.

**The winding
number test in 2D.**

An
alternative test for whether a point is inside or outside a 2D
contour can be achieved using the property that if a point **P**
is inside a polygon such as that described by the contour, then the
sum of all angles at **P** which are made between the adjacent
lines from each vertex to the point **P** equals 2p
otherwise it is 0. This method
is algorithmic more complex and therefore more time consuming than
the intersection number method.

We also assume that the volume is described as a polyhedron. As for the 2D case there exist a similar intersection and winding number test.

**The
Intersection number test with the Polyhedron surface.**

We have used a variant of the intersection algorithm described by O’Rourke:

*Algorithm:
Point in Polyhedron*

*Compute
bounding radius R*

*loop
forever*

*r0
= random ray of length R*

*r
= q + r0.*

*crossings
= 0.*

*for
each face F of polyhedron P do*

*SegFInt(F,q,r)*

*if
degenerate intersection *

*then
Go back to loop*

*else
increment crossings appropriately.*

*if
crossings odd*

*then
q is inside P. *

*else
q is outside P.*

*Exit.*

O’Rourke Point in Polyhedron Algorithm, see [OR98].

This algorithm avoids degenerate cases such as if the ray crosses an edge or a polyhedron vertex. If such a case occurs a new ray is generated and the test is repeated. Usually it is common to triangulate the facets, so that a line triangle intersection routine is used which is simpler than a general polygon line intersection test.

**B.
The winding number test.**

The winding number test can be extended in three dimensions[CA95]:

Project each face of the polyhedron onto a unit sphere of center P and compute the signed area of the spherical polygon thus determined. The sign is positive if the spherical polygon has counterclockwise orientation and negative otherwise.

The sum S of the signed areas of the projections of all facets of the simple polyhedron P onto the unit sphere p is necessarily 0, 4p or 4p.

**Fig. 3**
Spherical Triangle from the projection from a point P of
a triangular polyhedron facet on a unit sphere.

If the facets of P have clockwise or counterclockwise orientation and the point P is inside the polyhedron then S = -4p or 4p respectively. If the point is outside the polyhedron then regardless of the orientation S = 0.

The area of a spherical triangle can be obtained by Girard s formula:

S = A + B + C = 2p

Where A,B,C are the spherical angles at each vertex. The spherical angle a at a given vertex A is the angle determined by the tangents to the two great circles corresponding to the sides that cross at A. a is also the angle formed by the planes containing those two great circles. Thus a can be determined from the normal vectors to each of these planes, which can be computer without projecting the vertices onto the sphere.

**Generation
of the intermediate polygons.**

For
the generation of the polygons we have to consider that the triangles
are not necessary stored in the appropriate order. We take a random
triangle and calculate the intersection point P_{1} of this
triangle with the plane between the two consecutive slices. Then we
take the next side of the triangle and the next intersection point
P2. Then we search the triangle which has the same edge A2A3. This is
the triangle number 2 with the intersection point P3. We continue
this procedure until we arrive again to the point P1. The
intersection points P1,P2,P3 define the polygon.

**Fig.
4** Generation of a polygon P1P2 from the
triangulated surface between two slices (dotted line).

**Generation
of the virtual polygons.**

We have developed a method which avoids the generation of polygons from the intersection of the triangulated surface and the plane through the point to be tested which is computation complex. We use a preprocessing method where so called virtual polygons are generated and stored.

For
each consecutive pair of actual contours a virtual polygon is
generated. The vertices of this polygon are functions of z, i.e. the
i^{th} vertex is given as P_{i}(x,y) = (f_{i}(z),
g_{i}(z)). A polygon is therefore a function of z and we call
therefore these polygons *virtual polygons*. We use the same
method as shown in Fig. 4. but this time we do not store the actual
x,y values for each polygon vertex but the function (f_{i}(z),
g_{i}(z)).

**Fig.
5** Example of an object defined by three
contours and the two corresponding virtual polygons shown in red
color generated for some z coordinates.

**Fig. 6**
Extension of point in 3D Polyhedron considering objects
inside the polyhedron.

**VI. REFEReNCES**

[LB00]^{
}M. Lahanas, D. Baltas, N. Milickovic, S. Giannouli, and N.
Zamboglou,Generation of uniformly distributed dose points for
anatomy-based-three-dimensional dose optimization in brachytherapy,
Med. Phys. 27, 1034-1046, 2000. **Abstract
**

[LK00] M. Lahanas, T. Kemmerer, N. Milickovic, D. Baltas, N. Zamboglou," Optimized bounding boxes for three-dimensional treatment planning in brachytherapy," accepted for publication in Med. Phys. 2000.

[CR41] R. Courant and H. Robbins, "What is Mathematics?," Oxford University Press, (1941).

[ED90]
H. Edelsbrunner and E. P. Mücke, "Simulation of Simplicity:
A technique to cope with degenerate cases in geometric algorithms,"
ACM Transactions on Graphics **9**, 66-104 (1990).

[CA95]^{
}P. C. P. Carvaiho, P. R. Cavalcanti, Point in Polyhedron
Testing Using Spherical Polygons, Graphics Gems V edited by A. W.
Paeth, Academic Press 1995

[OR98]J. O’Rourke, Computational Geometry in C , Cambridge, Second Edition, Academic Press 1998