PIECEWISE
MATHEMATICAL REPRESENTATION ARTICULAR SURFACES* P. K.
SCHERRER+
OF
and B. M. HILLIWKRI.:
School of Mechanical Engineering, Purdue University. West Lafayette, IN 47907. U.S.A. Abstract -Spatial linkages have been used to measure the relative position of anatomical bodies in vim If an accurate method is used for mathematically representing the articular surfaces, the relative position of which is measured by a linkage, in ciao areas of contact in anatomical joints can be determined. Articular surfaces can be represented mathematically in a piecewise continuous fashion using parametric surface patches. This paper describes the use of Coons’ bicubic surface patches, with null twist vectors, to obtain the piecewise mathematical representation of the articular surfaces of the proximal humerus and scapula from the shoulder of a German Shepherd Dog. The accuracy of the representation is examined qualitatively and quantitatively.
PARAMETRIC SURFACE PATCHES
INTRODUCTION In order to gain some insight into the dynamics of the canine shoulder, Kinzel (1972. 1973) developed an instrumented spatial linkage for measuring the relative position of two rigid bodies in space. With this linkage, Kinzel (1973, 1976) measured the relative position of the scapula and humerus in the canine shoulder, it! cico. for various activities (i.e. walking and trotting over a flat surface, a step, a hurdle and a ramp). Kinzel then plotted the periphery and two orthogonal curves from the articular surface of the scapula and the articular surface of the proximal humerus in their relative position, based on the data from the instrumented linkage. These plots provided valuable insight into the motion between the two articular surfaces. If the articular surfaces were represented mathematically, the points on one articular surface which are in contact with the mating surface could be determined. Unfortunately, anatomical surfaces are not the classic geometric surfaces readily described mathematically, such as a sphere or cylinder. Therefore, a method was sought which would provide the best approximation to the articular surfaces. Due to the shape of articular surfaces, the method used to represent them would have to be capable of representing three-dimensional, doubly-curved surfaces. A surface of arbitrary shape, that cannot be represented by a single analytical expression, can be represented in a piecewise continuous fashion using surface patches. If parametric surface patches are used, the mathematical representation will be independent of the location of the coordinate axes, and the surface can
be single or multivalued.
In this paper, the mathematical representation of articular surfaces using parametric surface patches is presented. * Receiced 7 Murch 1978. t Former Graduate Research Assistant. : Professor.
Literature dealing with parametric surface patches and their application began to appear in the mid 1960’s. Two papers in particular, one by Ferguson (1964) and the other by Coons (1967), seem to have kindled interest in parametric patches. Coons’ method is a more general method than Ferguson’s, and contains Ferguson’s method as a special case. Bezier (1968, 1972) developed a parametric method which he used for the design of automobile bodies. Sadeghi et al. (1974) compared Bezier’s and Ferguson’s methods. They pointed out that although the mathematical basis is similar for both methods, Ferguson’s method is basically for surface fitting, whereas Bezier’s method is for designing with patches. MacCallum (1970) and Lee (1969) have developed parametric representations specifically for interactive computer-aided design systems rather than surface fitting. Forrest (in Bezier, 1972) pointed out that a Coon’s bicubic surface patch with null twist vectors is equivalent to Ferguson’s F-surface. It is also equivalent to a Coons’ F-surface with cubic blending functions and cubic boundary curves. This form of surface description is widely used, Johnson (1966), Hogue (1966), Sabin (1971), Hart (1971), Almond (1972) and South (1965). It is this type of surface patch, Coons’ bicubic surface patch with null twist vectors, which was used in this study to obtain the mathematical representation of the articular surfaces.
COOS!? BICUBIC SURFACE PATCH
Each Coons’ bicubic surface patch is bounded by four parametric cubic curves and, as pointed out by Almond (1972): “A cubic is the lowest order polynomial that can symbolize a curve that twists through space and, therefore, form patch boundaries that are nonplanar”. Another reason for choosing this form of Coons’ surface patches was due to the fact that the boundary conditions needed to define patches con301
302
P. K.
SCHERRER
and 8. M. HILLBERRY
where the superscript t is used to indicate the transpose of a matrix, and the base vectors U and Ware the following cubic vectors :
(a)
lJ = (u’u2u1)
(2)
w = (w3wW).
(3)
The constant matrix M (called the “magic matrix” by Lee, 1969) is as follows : 2
-2
M/-3
(b)
tinuous in position and slope are easily obtained from a mesh of points describing the surface. The mesh of points shown in Fig. l(a) would be used to form the multipatch surface of Fig l(b). A point on a surface has two degrees of freedom, thus enabling one to represent the equation of the surface as a function of two variables. Let the two variables (parametric coordinates) for a patch be u and w. To make it easy to bound a surface patch, Coons makes the stipulation that the surface patch in parametric space is a unit square, that is 0 -< II 5 1 and 0 I w I 1. Coons uses the single variable uw to represent a vector containing the X, Y and Z coordinates of a point on a surface patch at the parametric coordinates (u, w). Therefore, a typical surface patch would have the four boundary cures I&, ul, Ow and lw, which intersect at the points 00, 10,Ol and 11 as shown in Fig. 2. The term u0 represents a vector containing the X, Y and Z coordinates of the boundary curve generated when w is zero and u is allowed to vary between 0 and 1. The equation for a point on a bicubic surface patch at the parametric coordinates (a, w) can be expressed in the following compact form, Coons (1967): uw = UMBM’W,
Fig. 2. A Coons’ surface patch.
B =
1
0
0
1
0’
1
0
0
01
The B matrix, called the”boundary has the following form : Fig. 1. A multipatch surface obtained from a mesh of points. (a) Mesh of surface points. (b) Multipatch surface formed from the mesh of points shown in (a).
1 3-2-l
condition matrix”,
00
01
00,
OlJ
10 ---
11
10, ----.
11.
00” 10, L
01” 11,
ww
01,.
rI
(4)
(5)
10,“. 1l,,
The subscripts in the B matrix indicate differentiation with respect to the subscribed variable. Thus 00, is equal to d(uw)/dw evaluated at u=O and w =O. Similarly, Ol,, is equal to 8(uw)/du dw evaluated at u=O and w = 1. The B matrix contains the four items that must be specified at each of the four comer points of a Coons’ bicubic surface patch. The four elements in the upper left hand portion of the partitioned B matrix contain the coordinates of the four corner points. The upper right and lower left four items in the partitioned matrix contain the comer slopes (tangents) for the patch as shown in Fig. 3. The four elements in the lower right hand side of the partitioned B matrix contain the cross derivatives (which Coons calls twist vectors and which are sometimes referred to as twists) at the corner points. Since the elements of the B matrix are vectors and there are sixteen elements in B (each with three components), forty-eight numbers are needed to describe a bicubic patch. The ith component of the vectors in B are used to obtain the ith component of the Cartesian coordinates at the parametric coordinates (u, w).
(1)
Fig. 3. Comer point tangent vectors.
Piecewise mathematical representation
of articular surfaces
303
slopes across the boundary and the twists at the end points of the boundary curve are the same for both patches. With Coon’s method, three-sided (singular) patches can be used. This type of patch is obtained by letting one of the four boundary curves be a point. However, this causes small discontinuities in the slope between a three sided patch and its neighbors across the boundary curves of the singular patch having the degenerate curve as an end point. The procedure to be used in obtaining the tangent Fig. 4. Two adjacent surface patches vectors for the surface patches from a mesh of surface points is described in the next section. Due to the difficulty encountered in determining the twist vectors The mathematical representation of a complex for Coons’ surface patches, some researchers, Hart surface is obtained from a composite of surface (1971). Almond (1972), have set the twist vectors equal patches. In forming the surface with patches defined by to zero. They pointed out that in doing this, some equation 1, it is easy to have adjacent patches conflattening of the surface patches occurs. To decrease tinuous in position and slope along their common the effect of setting the twist vectors equal to zero, the edge. Consider the two surface patches A and B, shown size of the patches should be decreased. Almond (1972) in Fig. 1, for which it is desired to have the two patches stated that they “deliberately chose to let the twist continuous in position along their common edge. In vectors be zero, realizing that first-order slope conorder for the lw boundary curve of patch A to be tinuity at patch boundaries alone could leave a residue continuous with the Ow curve of patch B, 1~” = Ov?‘, of second-order quasi-flattening at patch corners. In the following conditions must hold: our experiment, the decision was justified in that propellers machined.. exhibit no measurable flatten1. 10” = oo* ing.” For the work described in this paper, the twist 2. 11’ = 018 vectors will also be assumed to be null vectors. 3. 10;. = 00; 4
11: = Olf,.
Therefore, two adjacent patches will have zero order continuity along their common edge, if the end points for the curves and the tangents along the boundary curves at the end points are the same for both patches. Now consider what is necessary to insure that two adjacent patches are continuous in slope across their common boundary. In order for the slope across the 1~. boundary curve of patch A to be equal to the slope across the OH.curve of patch B. 1~‘: = Of, the following conditions must hold : 1. 10; = 00; 7 _.
11: = 01;
3.
lo& = 00;”
3.
11;“. = 01;*.
Therefore, two adjacent patches will havecontinuity of slope across their common boundary if both the
DETERMINING THE TANGENT VECTORS FOR COONS’ BICLBIC SURFACE PATCHES
In the previous section, the equation for a Coons’ bicubic surface patch was presented. One element in the matrix representation was the boundary condition matrix, B, which contains the coordinates, tangents and twists for a patch. Assuming that the twists are zero, only the comer point coordinates and tangents would have to be specified to define a patch. In representing an existing surface, the coordinates of the corner points can be measured directly, but the tangent vectors must be determined. If a smooth curve is formed through the corner points that lie on a mesh line, as shown in Fig. 5, the derivative of the curve at the corner points can be used as an estimate for the corner point tangent vectors. Since the boundary curves for the patches are parametric cubin, it would seem logical to use a piecewise parametric cubic spline
Segment
Y
Segment
n-IA
j
Fig. 5. Points along a mesh line of a muitipatch surface which are used in computing the tangent vectors.
P. K.
303
SCHERRER and B. M. HILLBERRY
to obtain the smooth curve from which the tangent vectors are obtained. The points used to form the spline would be corner points for the patches as shown in Fig. 5. With a piecewise parametric cubic spline, a parametric cubic is formed between adjacent points, subject to the condition that the segments are continuous in position, as well as first and second derivatives. South and Kelly (1965). Almond (1972) and Adams (1974) used the same method for fitting a spline. This method requires the solution of a system of linear equations and permits a number of different end conditions for open curves. The procedure that was used in this work for fitting a piecewise cubic spline to a set of points is based on the method of South and Kelly (1965) as reported by Almond (1972). With this method, the parametric cubic for segment i of the spline will be of the following form: $(y) = dj c B,y + C,$ + b+/‘, (6) for y in the interval (O,R,), where Rj is the distance between points pj and Pj,, which is referred to as the chord length of segment j. Since the segments are to have continuous second derivatives at the points where the segments join,$_ ,(q) must be equal to s;‘(q) at Pi. Therefore, equating $_ , (R ;_ 1) and $(O.O), one obtains cj= cj-t + 3di_,Ri_,. (7) Using thecontinuity conditions for position-and slope at the points where the segments join (Sj(O) = Pi, sj(Ri) = pj+l, s;(O) = i‘j and $(Rj) = Tj+,), the coefficients for equations 6 and 7 can be obtained in terms of the end points and end point tangents. Replacing thecoefficients in equation (7), the following equation is obtained :
RjTj-1 +Z(R;+ R,_,)?, + Rj_,T,+, = 3(Rj_,(P,+,
- Pj) + Rj(P, - P,_ ,)),‘R,R,_ ,.
For an open spline containing :V points, equation (8) holds for i = 2.3 , . , N - 1 since thereare N - 2 points at which the segments are joined. Thus there are .V - 2 linear equations and X unknowns (the tangent xctors). The two additional equations which are needed can be obtained by imposing boundary conditions on the spline, such as requiring the first and last seaments to contain known points or specifying I‘, and-T,. (A table of the different boundary conditions that can be used in obtaining the two additional equations is given by Adams (1974) and Scherrer (197i).) If a smooth closed curve is desired, and points 1 and N are the same, there will be N - 1 linear equations oftheform of equation 8 and there will be N - 1 unknown tangent vectors. The tangent vectors determined using the method described above are approximately unit vectors, Almond (1972) and they are only directly applicable for a spline where the parametric coordinate for segment i varies from 0 to Rj. The parametric coordinates along a boundary curve of a bicubic surface patch, however, vary from 0 to 1. It is therefore necessary to modify the magnitude of the tangent vectors before they can be used for the surface patches. In practice, the tangent vector at a point is modified by the shortest chord length for the segments containing that point, before being used for the patches. Modifying the magnitude of the tangent vectors in this way will produce boundary curves free from undesirable loops and bulges, although some Rattening may result, and which lie in the triangle formed by the tangent vectors and chord, South and Kelly (1965).
The XY plane iS the viewing plane
Fig. 6. The projection
of a surface
(8)
point onto the viewing
plane,
White (1974).
Piece*ise
mathematical
representation
of articular
surfaces
Surface to be viewed is defined X, ,q ,Z, coordinate system
305
in the
Fig. 7. The viewing coordinate system for a surface is different from the coordinate system in which rhe surface is defined.
DETERMINING THE ACCURACY Oc‘ THE MULTIPATCH REPRESENTATION
After all the tangent vectors at the corner points of the patches have been estimated, a perspective view of the composite patchescan be plotted in order to obtain a qualitative check on the mathematical representation. A perspective view of a point is obtained by finding the point of intersection of a plane (the viewing plane) and a line from the eye of the viewer to the point on the surface being viewed. The location of the projected points on the plane are dependent on the distance the points on the objects are from the viewer, and the location of the viewing plane. Let it be assumed that the eye of the viewer lies DlST units along the positive Z axis, the object to be viewed is on the negative side of the X Y plane, and that the XY plane is the viewing plane. Referring to Fig. 6, it can be seen that based on similar triangles the following equations hold: XJDIST
= Xj(DlST
the original coordinate system to the viewing coordinate system can be determined. The surface patches can be transformed into the viewing coordinate system, and then points on the patches can be generated and projected onto the viewing plane using equation (9). (The details of the transformation are given by Scherrer (1977).) A major problem with a perspective view is that all the surface points generated are plotted. Thus points on the object which are hidden from the viewer are plotted. This probiem is evident in Fig. S(a) which shows a perspective view of a multipatch representation of the articular surface of the humeral head
- 2,)
YJDIST = YJ(DIST - Z,) ’
(9)
where X, and Y, are the coordinates of the point to be plotted for a perspective view of the point on the surface with coordinates (X,, Y,, Z,), White (1974). Equation (9) cannot be used directly if the viewer does not lie on the positive Z axis, or if the object is not on the negative side of the X Y plane (viewing plane). To overcome this problem and still make use of equation (9), it is necessary to set up a viewing coordinate system as shown in Fig. 7, White (1974). The location of the viewing coordinate system is defined by the points (X,, Y,,Z,) and (X,, Y,,Z,), where (X,, Y,, Z,) is the point from which the surface is to be viewed and (X,, Y,, Z,) is a point on the desired viewing plane. The viewing plane will be assumed to be perpendicular to the line from (X,, Y,,Z,) to (X,, Y,, Z,) and should be located between the object being viewed and the viewer. Based on the coordinates of the points (X,, Y,,Z,) and (X,, Y,,Z,), the coordinate transformation used to transform a point from
Fig. 8. Perspective of a multipatch surface. (a) All the patches plotted. (b) Patches on the front of the surface pLotted
306
P. K. SC~~ERRER and B. M. HILLBERRY
from a canine shoulder. To overcome the problem of hidden points, a procedure was sought for determining which points are on the back side of an object so they would not be plotted. Since the perspective view is obtained with the viewer located on the positive Z axis, a vector perpendicular to the surface at a point on the back side of the surface will have a negative Z component. Therefore, if vectors perpendicular to the surface at the comers of a patch are determined and they all have a positive Z component, the patch will be on the front of the surface and is considered visible, whereas the patch is considered hidden if all the Z components are negative. If any, but not all, of the corner points of a patch are found to lie on the backside of the surface, a vector perpendicular to the surface patch is generated at each point to be plotted to determine if it is hidden. If the Z component for the normal vector at a point is positive, the point is plotted. Figure 8(b) shows the same surface as the one in Fig. 8(a) except the procedure described above was used to remove the points which lie on the back of the surface from the plot. A Coons’ surface patch contains the four corner points for the patch; therefore, there is no error in the representation at these points. At any other location on a patch, a quantitative measure of the error in the patch can be obtained by determining the perpendicular distance between the patch and an actual surface point that the patch would contain if it were an exact representation of the surface. Let ps be a point on the wtual surface, and P, be the point on a surface patch at which a vector perpendicular to the patch intersects ps. The distance from P,v to ps can be used as a measure of the error in the patch representation. The point P, is located at the point on the surface patch where the vector from P,v to ps is perpendicular to the tangent vectors in the u and w parametric directions as illustrated in Fig. 9. Since the scalar product of two perpendicular vectors is zero, the parametric coordinates at which P, is located can be determined by solving the following two simultaneous nonlinear equations in u and w, Scherrer (1977) :
Fig. 9. Point on a Patch t?om which a normal vector intersects a known point.
* Dow-Coming, Midland, Michigan, U.S.A. t Ransom and Randolph Company, Toledo. Ohio, U.S.A.
E&u, w) PsP(f4, w) . ---&= 0
(10)
a&u, w) psp(u, w) ‘~ZO, dW
where PsP(u,w) = &u, w) - lis and P(u,w) is the coordinates of a point on the surface patch at the parametric coordinates (u, w). Equation (10) can be easily solved, using Newton’s method, to find the parametric coordinates at which P,v is located. If P,v does not lie on the patch for which the B matrix was formed, one or both of the parametric coordinates at which P,v is located will exceed the range O-l. MOLDING
THE ARTICULAR
SURFACES
In order to obtain the mathematical representation of a surface using surface patches, it is necessary to know thecoordinates ofthecorner points for the patches. To draw patches on the articular surfaces of the actual bones would have been difficult since the surfaces must be kept moist, and they could not be damaged since they were to be examined histochemically. Another problem with working with the actual articular surface is that it would have been necessary to keep the surface moist while the coordinates of surface points were being measured. Therefore, a silicone rubber mold and then a plaster casting of the articular surfaces were made. Patches weredrawn on thecastings obtained for the articular ends of the scapula and the proximal humerus, and the coordinates of the corner points were measured. Figures 10 and 11 show the patches drawn on the castings for the articular surfaces ofthe shoulder of a 30 kg, 2 year old female German Shepherd Dog. (Test subject GS071576.) In order to avoid loss in accuracy, care had to be taken in choosing the materials to be used in making the mold and casting. The molds were made of Dow Coming Silastic 382 Medical Elastomer* (silicone rubber) and the castings were made of Ransom and Randolpht dental plaster. The dental plaster produced castings that were both accurate and durable. The setting expansion (linear only) for the plaster was between 0.05 and 0.09%. The silicone rubber used was a two-component (catalyst and base), room temperature vulcanizing (RTV) type, The amount of time needed for the silicone rubber to dry depended on the amount of catalyst used. In the molds that were made of the articular surfaces of the canine shoulder, the silicone rubber dried in less than five minutes. Medical grade Silastic was chosen to avoid severe damage to the cartilage. (It appears that the Silastic did affect the c+on the sufface of the cartilage, Kincaid (1977)) Silicone rubber has the desired flexibility needed to remove the articular ends of the bones from the molds, but the shrinkage characteristics over an extended period of time are not very good. According to Braley (1970), at 77°F the shrinkage (linear) is 0.4% after 3 days, 0.6% after 6 days and 0.7% after 14 days. It was
Piecewise
Fig. 10. The surface Fip. Il. Thz surface
mathematical
patches patches
drawn drawn
representation
~2 artiiular
on the casting
of rhe left humerus.
on the casting
surfaces
of rhe left scapul:l.
Piecewise mathematical representation
of articular surfaces
found, however, that shrinkage effects would be neghgible if the castings were made immediately after the bones were removed from the molds. This was determined by making a mold of a U-shaped block consisting of two 25.4 mm square uprights, 50.8 mm high with 1 in. between them. The maximum difference between the dimensions of the block and its casting was found to be 0.05 mm. Thus the molding procedure appears quite accurate. MAPPING THE SURFACES
After the patches were drawn on the surface of the casting (see Figs. 10 and ll), the coordinates of the corner points for the patches were measured using a Brown and Sharpe Validator SO+.The Validator has a granite table, air bearings on all axes, and digital readout in inches or millimeters. ft has a resolution of 0.002 mm (0.0001 in.), and the overall accuracy is k0.012mm (+O.OOOSin.), with an accuracy in a 300 mm (12 in.) square of +_0.008 mm (+_0.0003 in.). The repeatability is given as f 0.005 mm (+O.OOOt in.). In mounting the casting of the articular surface of the humerus. it was necessary to use three set-ups so that the probe carrier was approximately perpendicular to the surface for all measurements. Only one set up was required for the scapula casting. After the coordinates of the surface points were measured, the casting and clamping blocks were removed from the table and then remounted to check the repeatability of the coordinate measurements. The coordinate measurements were then checked by setting the X and Y coordinates of the probe on the Validator to the X and Y coordinates previously measured at a point and comparing the Z coordinate. The maximum difference in the readings for the Z coordinates from the first and second set-up was 0.012 mm (0.0005 in.).
Cramal
Idol
Fig. 13. Dorsal view of the multipatch representation of the proximal articular surface of the humerus.
THE MATHEMATICAL REPRESENTATION OF THE ARTICULAR SURFACES OF THE CANINE SHOULDER
Once the corner points for the patches were measured, the corner point tangent vectors for the surface patches were generated using a piecewise parametric cubic spline (Program TNGNTS was used to generate the tangent vectors for the patches, Scherrer (1977).) Since Coons’ bicubic surface patches with null twist vectors were used to describe the surfaces, the mathematical representations were complete once the tangent vectors were obtained. The coordinates of the corner points and the components of the tangent vectors for the patches were used in a computer program (program PSP, Scherrer, 1977) to obtain perspective views of the multipatch repr-esentation of the surfaces. This provided a qualitative check on the mathematical representations. Perspective views of the articular surfaces of the left shoulder of a German Shepherd Dog (test subject GS071576) are shown in Figs. 12-15. The solid lines in these figures represent the patch boundaries. Figures 12 and 13 show the proximal humeral articular surface, and Figs. 14 and 15 show the articular surface of the
Croniol
Caudal
Fig. 12. Caudal-lateral view of the multipatch representation of the proximal articular surface of the humerus. l
Brown and Sharpe, North Kingstown, RI., U.S.A.
Media\
Fig. 14. End-on view, from a mediil aspect, of the multipatch representation of the articular surface of the scapula.
P. K. SCHERRERand B. .%I.HILLBERRY
310 Lateral
Medial
Fig. 15. Direct end-on view of the multipatch representation of the articular surface of the scapula.
scapula. Due to a concavity in the articular surface of the proximal humerus near the medial osteochondral border, a portion of the medial-caucial articular surface was not mathematically represented. The missing portion of the articular surface is not shown in Figs. 12 and 13. Figures 12-15 can bz compared to Figs. 10 and 11 which show the castings of the articular surfaces of the proximal humerus and scapula with the patches drawn on them. (The sides of the patches on both the humerus and scapula were approximately 3.2 mm (0.125 in) long.) The pronounced curve in the boundary curves of some of the patches of the multipatch surfaces shown in Figs. 12-15 ace not a result of coding errors in any computer routine. For end conditions on the splines used to generate the tangent vectors, the splines were required to contain known points. Thecurves would have been smoother if greater care had been taken in choosing the surface points used for the boundary conditions. The accuracy of the representations were examined quantitatively by determining the perpendicular distance between points on the actual surfaces and the multipatch representations. (This was done using program ERR, Scherrer, 1977.) For the scapula, the coordinates of one point at the center of each of the 48 patches was measured, as well as the coordinates of 116 other points located along the boundary curves of the patches.* It was found that, for the points at the center of the patches of the scapula, the error in the representation ranged from 0.19 mm (0.0075 in.) to -0.21 mm (-0.0084 in.), with the average error being 0.063 mm (0.0025 in.). The average error is defined here as the sum of the absolute value of the error at the l The coordinates of a number of points between the ends of the boundary curves for the patches on the scapular surface were originally measured for use in an experimental computer routine. This routine contained a procedure for varying the magnitude of the tangent vectors at the knots of a piecewise cubic splint to minimize the error in the spline. The error was defined as the perpendicular distance between the spline and points located between the spline and points located between the knots. Since these coordinates had been measured, they were used in examining the error along the patch boundaries. It should also be noted that the routine for minimizing the, error did not appear to make an improvement in the fit.
points, divided by the number of points. In descrrbing the range of the error, a positive distance indicates that the actual surface point was above the representation. There was one point for the scapula at which the computed error was 0.70 mm (0.027-1in.). It appears that this large error was a result of an error in recording the coordinates of the point, and therefore is not representative of the error in the mathematical representation. However. the 0.70 mm (0.0274 in.) error was used in determining the average error of 0.06 mm (0.0025 in.) for the scapula. For the 116 points located along the boundary curves of the scapula, the error ranged from 0.26 mm (0.0103 in.) to -0.22 mm (-0.0088 in.), with an average error of 0.04mm (0.0016 in.). Based on theerror in the 164 points on the scapula that were not corner points. the average error was 0.05 mm (0.0019 in.). The coordinates ofone point at the center of each of the 65 patches on the casting of the humerus were also measured. The maximum error was 0.21 mm (0.0083 in.), the minimum error was -0.152 mm (-0.0060 in.) and the average error was 0.04 mm (0.0017 in.). After determining the error, increased accuracy can be obtained by using smaller patches. After the error was determined at the points mentioned above, it was found that the error at other points should have been examined. According to Sabin (1978), the effect of using zero twist vectors is expected to be largest around the parametric coordinates (0.2, 0.2), (0.2, 0.8). (0.8, 0.2) and (0.8, 0.8). Thus the error should have been examined at these parametric coordinates.
SU,MMARY
Dental plaster castings of the articular surfaces from the shoulder of a German Shepherd Dog were made directly from the shoulder components. Surface patches were drawn on the castings of the articular surfaces, and the coordinates of the corner points on the patches were measured. (As noted previously, patches with sides approximately 3.2 mm (0.125 in.) long were used on both the scapular and humeral castings.) Using the coordinates of the corner points, the tangent vectors used to describe the surface with surface patches were obtained. Perspective views were generated to verify the multipatch representations of the surfaces. The error in the representation was examined by determining the perpendicular distance between the representation and actual surface points. Based on 164 points, the error in the representation of the articular surface of the scapula ranged from 0.26 mm (0.0103 in.) to -0.21 mm (-0.0088 in.), with an average error of 0.04mm (0.0016 in.). For the representation of the articular surface of the proximal humerus, the maximum error was 0.21 mm (0.0083 in.) and the minimum error was -0.15 mm (-0.0060 in.), with the averageerror being0.04 mm (Oo.0017in.). This was based on 65 points. Increased accuracy can b-e obtained by using smaller patches.
Piecewise mathematical representation Ackno%ledgemenrs - This project was sponsored by the National Institute of Arthritis, Metabolism and Digestive Diseases. HEW, under Grant Number 5 ROI 17482. The authors wish to thank Dr. D. C. Van Sickle, DVM, Ph.D. and Dr. G. J. Pijnowski, DVM, School of Veterinary Science and Medicine and Professor J. R. Rice, Department of Computer Science, Purdue University for their assistance on this project.
REFERENCES
Adams, J. A. (1974) Cubic spline curve fitting with controlled end conditions. Comp. Aided Desiyn 6, 2-9. Almond, D. B. (1972) Numerical control for machining complex surfaces. fBM Sysr. J. 150-168. Bezier. P. (1968) How Renault uses numerical control for car body design and tooling. SAE Paper 680010. Bezier, P. (1972) Numerical control in automobiledesign and manufacture of curved surfaces. pp. 44-48, Proc. Curved Surfaces in Engny Conf Churchill College, Cambridge, England. Bezier, P. (1972) Numerical control. Mathematics and Applicarions. Wiley, New York. Braley, S. (1970) The chemistry and properties of the medicalgrade silicons. J. Macromol. Sci.-Chem. 529-544. Coons, S. A. (1967) Surface ofcomputer-aided design of space forms. MIT Project MAC-TR-41. Ferguson, J. (1964) Multivalued curve interpolation. J. ACM 1I, 221-228. Hart, W. B. (1971) Computer graphics aid fuselage design. Enyny Mat. Des. 205-209. Hope, W. M. (1966) Computer-aided design in body engineering. SAE Paper 660153.
B.Y. 12,4-E
of articular surfaces
311
Johnson, W. L.. Sanders, J. W. and South, N. E. (1966) Analytic surfaces for computer-aided design SAE Paper 660152. Kincaid, S. A. (1977) Personal communication. Purdue University. Kinzel, G. L., Hall, A. S.. Jr. and Hillberry, B. ?+t. (1972) Measurement of the total motion between two body segments, Part I - analytical development. J. Biumech. 3, 93- 105. Kinzel, G. L. (1973) On the design of instrumented linkages for the measurement of relative motion between two rigid bodies. Ph.D. Thesis, Purdue University. Kinzel, G. L.. Van Sickle, D. C., Hillberry, 6. M. and Hall. A. S., Jr, (1976) A preliminary study of the in viuo motion in the canine shoulder, Am. J. cer. Res. 37, 1505-15 IO. . Lee, T. M. P. (1969) Three-dimensionof Curves and Surfaces for Ropid Computer Display. Management Information Service, Detroit, Michigan. MacCullum, K. J. (1970) Surfaces for interactive graphical design. Camp. J. 13, 352-358. Sabin, M. A. (1971) An existing system in the aircraft industry. The British Aircraft Corporation numerical master geometry system. Proc. R. Sot. (A)321. 197-205. Sabin, M. A. (1978) Personal communications. Sadehgi, M. M. and Gould, S. S. (1974) A comparison of two parametric surface patch methods. Comp. Aided Design 6. South, N. E. and Kelly, J. P. (1965) Analytic Surface Methods. Numerical Conrrol Decelopment L’nit. Ford \lotor Company, Dearborn, Michigan. White, P. (1974) Notes from ME 573 course. Purdue University. Scherrer. P. K. (1977) Determining the distance between the surface of two rigid bodies using parametric surface _ patches and an instrumented spatial linkage. Ph.D. Thesis, Purdue University.
MATHEMATICAL REPRESENTATION ARTICULAR SURFACES* P. K.
SCHERRER+
OF
and B. M. HILLIWKRI.:
School of Mechanical Engineering, Purdue University. West Lafayette, IN 47907. U.S.A. Abstract -Spatial linkages have been used to measure the relative position of anatomical bodies in vim If an accurate method is used for mathematically representing the articular surfaces, the relative position of which is measured by a linkage, in ciao areas of contact in anatomical joints can be determined. Articular surfaces can be represented mathematically in a piecewise continuous fashion using parametric surface patches. This paper describes the use of Coons’ bicubic surface patches, with null twist vectors, to obtain the piecewise mathematical representation of the articular surfaces of the proximal humerus and scapula from the shoulder of a German Shepherd Dog. The accuracy of the representation is examined qualitatively and quantitatively.
PARAMETRIC SURFACE PATCHES
INTRODUCTION In order to gain some insight into the dynamics of the canine shoulder, Kinzel (1972. 1973) developed an instrumented spatial linkage for measuring the relative position of two rigid bodies in space. With this linkage, Kinzel (1973, 1976) measured the relative position of the scapula and humerus in the canine shoulder, it! cico. for various activities (i.e. walking and trotting over a flat surface, a step, a hurdle and a ramp). Kinzel then plotted the periphery and two orthogonal curves from the articular surface of the scapula and the articular surface of the proximal humerus in their relative position, based on the data from the instrumented linkage. These plots provided valuable insight into the motion between the two articular surfaces. If the articular surfaces were represented mathematically, the points on one articular surface which are in contact with the mating surface could be determined. Unfortunately, anatomical surfaces are not the classic geometric surfaces readily described mathematically, such as a sphere or cylinder. Therefore, a method was sought which would provide the best approximation to the articular surfaces. Due to the shape of articular surfaces, the method used to represent them would have to be capable of representing three-dimensional, doubly-curved surfaces. A surface of arbitrary shape, that cannot be represented by a single analytical expression, can be represented in a piecewise continuous fashion using surface patches. If parametric surface patches are used, the mathematical representation will be independent of the location of the coordinate axes, and the surface can
be single or multivalued.
In this paper, the mathematical representation of articular surfaces using parametric surface patches is presented. * Receiced 7 Murch 1978. t Former Graduate Research Assistant. : Professor.
Literature dealing with parametric surface patches and their application began to appear in the mid 1960’s. Two papers in particular, one by Ferguson (1964) and the other by Coons (1967), seem to have kindled interest in parametric patches. Coons’ method is a more general method than Ferguson’s, and contains Ferguson’s method as a special case. Bezier (1968, 1972) developed a parametric method which he used for the design of automobile bodies. Sadeghi et al. (1974) compared Bezier’s and Ferguson’s methods. They pointed out that although the mathematical basis is similar for both methods, Ferguson’s method is basically for surface fitting, whereas Bezier’s method is for designing with patches. MacCallum (1970) and Lee (1969) have developed parametric representations specifically for interactive computer-aided design systems rather than surface fitting. Forrest (in Bezier, 1972) pointed out that a Coon’s bicubic surface patch with null twist vectors is equivalent to Ferguson’s F-surface. It is also equivalent to a Coons’ F-surface with cubic blending functions and cubic boundary curves. This form of surface description is widely used, Johnson (1966), Hogue (1966), Sabin (1971), Hart (1971), Almond (1972) and South (1965). It is this type of surface patch, Coons’ bicubic surface patch with null twist vectors, which was used in this study to obtain the mathematical representation of the articular surfaces.
COOS!? BICUBIC SURFACE PATCH
Each Coons’ bicubic surface patch is bounded by four parametric cubic curves and, as pointed out by Almond (1972): “A cubic is the lowest order polynomial that can symbolize a curve that twists through space and, therefore, form patch boundaries that are nonplanar”. Another reason for choosing this form of Coons’ surface patches was due to the fact that the boundary conditions needed to define patches con301
302
P. K.
SCHERRER
and 8. M. HILLBERRY
where the superscript t is used to indicate the transpose of a matrix, and the base vectors U and Ware the following cubic vectors :
(a)
lJ = (u’u2u1)
(2)
w = (w3wW).
(3)
The constant matrix M (called the “magic matrix” by Lee, 1969) is as follows : 2
-2
M/-3
(b)
tinuous in position and slope are easily obtained from a mesh of points describing the surface. The mesh of points shown in Fig. l(a) would be used to form the multipatch surface of Fig l(b). A point on a surface has two degrees of freedom, thus enabling one to represent the equation of the surface as a function of two variables. Let the two variables (parametric coordinates) for a patch be u and w. To make it easy to bound a surface patch, Coons makes the stipulation that the surface patch in parametric space is a unit square, that is 0 -< II 5 1 and 0 I w I 1. Coons uses the single variable uw to represent a vector containing the X, Y and Z coordinates of a point on a surface patch at the parametric coordinates (u, w). Therefore, a typical surface patch would have the four boundary cures I&, ul, Ow and lw, which intersect at the points 00, 10,Ol and 11 as shown in Fig. 2. The term u0 represents a vector containing the X, Y and Z coordinates of the boundary curve generated when w is zero and u is allowed to vary between 0 and 1. The equation for a point on a bicubic surface patch at the parametric coordinates (a, w) can be expressed in the following compact form, Coons (1967): uw = UMBM’W,
Fig. 2. A Coons’ surface patch.
B =
1
0
0
1
0’
1
0
0
01
The B matrix, called the”boundary has the following form : Fig. 1. A multipatch surface obtained from a mesh of points. (a) Mesh of surface points. (b) Multipatch surface formed from the mesh of points shown in (a).
1 3-2-l
condition matrix”,
00
01
00,
OlJ
10 ---
11
10, ----.
11.
00” 10, L
01” 11,
ww
01,.
rI
(4)
(5)
10,“. 1l,,
The subscripts in the B matrix indicate differentiation with respect to the subscribed variable. Thus 00, is equal to d(uw)/dw evaluated at u=O and w =O. Similarly, Ol,, is equal to 8(uw)/du dw evaluated at u=O and w = 1. The B matrix contains the four items that must be specified at each of the four comer points of a Coons’ bicubic surface patch. The four elements in the upper left hand portion of the partitioned B matrix contain the coordinates of the four corner points. The upper right and lower left four items in the partitioned matrix contain the comer slopes (tangents) for the patch as shown in Fig. 3. The four elements in the lower right hand side of the partitioned B matrix contain the cross derivatives (which Coons calls twist vectors and which are sometimes referred to as twists) at the corner points. Since the elements of the B matrix are vectors and there are sixteen elements in B (each with three components), forty-eight numbers are needed to describe a bicubic patch. The ith component of the vectors in B are used to obtain the ith component of the Cartesian coordinates at the parametric coordinates (u, w).
(1)
Fig. 3. Comer point tangent vectors.
Piecewise mathematical representation
of articular surfaces
303
slopes across the boundary and the twists at the end points of the boundary curve are the same for both patches. With Coon’s method, three-sided (singular) patches can be used. This type of patch is obtained by letting one of the four boundary curves be a point. However, this causes small discontinuities in the slope between a three sided patch and its neighbors across the boundary curves of the singular patch having the degenerate curve as an end point. The procedure to be used in obtaining the tangent Fig. 4. Two adjacent surface patches vectors for the surface patches from a mesh of surface points is described in the next section. Due to the difficulty encountered in determining the twist vectors The mathematical representation of a complex for Coons’ surface patches, some researchers, Hart surface is obtained from a composite of surface (1971). Almond (1972), have set the twist vectors equal patches. In forming the surface with patches defined by to zero. They pointed out that in doing this, some equation 1, it is easy to have adjacent patches conflattening of the surface patches occurs. To decrease tinuous in position and slope along their common the effect of setting the twist vectors equal to zero, the edge. Consider the two surface patches A and B, shown size of the patches should be decreased. Almond (1972) in Fig. 1, for which it is desired to have the two patches stated that they “deliberately chose to let the twist continuous in position along their common edge. In vectors be zero, realizing that first-order slope conorder for the lw boundary curve of patch A to be tinuity at patch boundaries alone could leave a residue continuous with the Ow curve of patch B, 1~” = Ov?‘, of second-order quasi-flattening at patch corners. In the following conditions must hold: our experiment, the decision was justified in that propellers machined.. exhibit no measurable flatten1. 10” = oo* ing.” For the work described in this paper, the twist 2. 11’ = 018 vectors will also be assumed to be null vectors. 3. 10;. = 00; 4
11: = Olf,.
Therefore, two adjacent patches will have zero order continuity along their common edge, if the end points for the curves and the tangents along the boundary curves at the end points are the same for both patches. Now consider what is necessary to insure that two adjacent patches are continuous in slope across their common boundary. In order for the slope across the 1~. boundary curve of patch A to be equal to the slope across the OH.curve of patch B. 1~‘: = Of, the following conditions must hold : 1. 10; = 00; 7 _.
11: = 01;
3.
lo& = 00;”
3.
11;“. = 01;*.
Therefore, two adjacent patches will havecontinuity of slope across their common boundary if both the
DETERMINING THE TANGENT VECTORS FOR COONS’ BICLBIC SURFACE PATCHES
In the previous section, the equation for a Coons’ bicubic surface patch was presented. One element in the matrix representation was the boundary condition matrix, B, which contains the coordinates, tangents and twists for a patch. Assuming that the twists are zero, only the comer point coordinates and tangents would have to be specified to define a patch. In representing an existing surface, the coordinates of the corner points can be measured directly, but the tangent vectors must be determined. If a smooth curve is formed through the corner points that lie on a mesh line, as shown in Fig. 5, the derivative of the curve at the corner points can be used as an estimate for the corner point tangent vectors. Since the boundary curves for the patches are parametric cubin, it would seem logical to use a piecewise parametric cubic spline
Segment
Y
Segment
n-IA
j
Fig. 5. Points along a mesh line of a muitipatch surface which are used in computing the tangent vectors.
P. K.
303
SCHERRER and B. M. HILLBERRY
to obtain the smooth curve from which the tangent vectors are obtained. The points used to form the spline would be corner points for the patches as shown in Fig. 5. With a piecewise parametric cubic spline, a parametric cubic is formed between adjacent points, subject to the condition that the segments are continuous in position, as well as first and second derivatives. South and Kelly (1965). Almond (1972) and Adams (1974) used the same method for fitting a spline. This method requires the solution of a system of linear equations and permits a number of different end conditions for open curves. The procedure that was used in this work for fitting a piecewise cubic spline to a set of points is based on the method of South and Kelly (1965) as reported by Almond (1972). With this method, the parametric cubic for segment i of the spline will be of the following form: $(y) = dj c B,y + C,$ + b+/‘, (6) for y in the interval (O,R,), where Rj is the distance between points pj and Pj,, which is referred to as the chord length of segment j. Since the segments are to have continuous second derivatives at the points where the segments join,$_ ,(q) must be equal to s;‘(q) at Pi. Therefore, equating $_ , (R ;_ 1) and $(O.O), one obtains cj= cj-t + 3di_,Ri_,. (7) Using thecontinuity conditions for position-and slope at the points where the segments join (Sj(O) = Pi, sj(Ri) = pj+l, s;(O) = i‘j and $(Rj) = Tj+,), the coefficients for equations 6 and 7 can be obtained in terms of the end points and end point tangents. Replacing thecoefficients in equation (7), the following equation is obtained :
RjTj-1 +Z(R;+ R,_,)?, + Rj_,T,+, = 3(Rj_,(P,+,
- Pj) + Rj(P, - P,_ ,)),‘R,R,_ ,.
For an open spline containing :V points, equation (8) holds for i = 2.3 , . , N - 1 since thereare N - 2 points at which the segments are joined. Thus there are .V - 2 linear equations and X unknowns (the tangent xctors). The two additional equations which are needed can be obtained by imposing boundary conditions on the spline, such as requiring the first and last seaments to contain known points or specifying I‘, and-T,. (A table of the different boundary conditions that can be used in obtaining the two additional equations is given by Adams (1974) and Scherrer (197i).) If a smooth closed curve is desired, and points 1 and N are the same, there will be N - 1 linear equations oftheform of equation 8 and there will be N - 1 unknown tangent vectors. The tangent vectors determined using the method described above are approximately unit vectors, Almond (1972) and they are only directly applicable for a spline where the parametric coordinate for segment i varies from 0 to Rj. The parametric coordinates along a boundary curve of a bicubic surface patch, however, vary from 0 to 1. It is therefore necessary to modify the magnitude of the tangent vectors before they can be used for the surface patches. In practice, the tangent vector at a point is modified by the shortest chord length for the segments containing that point, before being used for the patches. Modifying the magnitude of the tangent vectors in this way will produce boundary curves free from undesirable loops and bulges, although some Rattening may result, and which lie in the triangle formed by the tangent vectors and chord, South and Kelly (1965).
The XY plane iS the viewing plane
Fig. 6. The projection
of a surface
(8)
point onto the viewing
plane,
White (1974).
Piece*ise
mathematical
representation
of articular
surfaces
Surface to be viewed is defined X, ,q ,Z, coordinate system
305
in the
Fig. 7. The viewing coordinate system for a surface is different from the coordinate system in which rhe surface is defined.
DETERMINING THE ACCURACY Oc‘ THE MULTIPATCH REPRESENTATION
After all the tangent vectors at the corner points of the patches have been estimated, a perspective view of the composite patchescan be plotted in order to obtain a qualitative check on the mathematical representation. A perspective view of a point is obtained by finding the point of intersection of a plane (the viewing plane) and a line from the eye of the viewer to the point on the surface being viewed. The location of the projected points on the plane are dependent on the distance the points on the objects are from the viewer, and the location of the viewing plane. Let it be assumed that the eye of the viewer lies DlST units along the positive Z axis, the object to be viewed is on the negative side of the X Y plane, and that the XY plane is the viewing plane. Referring to Fig. 6, it can be seen that based on similar triangles the following equations hold: XJDIST
= Xj(DlST
the original coordinate system to the viewing coordinate system can be determined. The surface patches can be transformed into the viewing coordinate system, and then points on the patches can be generated and projected onto the viewing plane using equation (9). (The details of the transformation are given by Scherrer (1977).) A major problem with a perspective view is that all the surface points generated are plotted. Thus points on the object which are hidden from the viewer are plotted. This probiem is evident in Fig. S(a) which shows a perspective view of a multipatch representation of the articular surface of the humeral head
- 2,)
YJDIST = YJ(DIST - Z,) ’
(9)
where X, and Y, are the coordinates of the point to be plotted for a perspective view of the point on the surface with coordinates (X,, Y,, Z,), White (1974). Equation (9) cannot be used directly if the viewer does not lie on the positive Z axis, or if the object is not on the negative side of the X Y plane (viewing plane). To overcome this problem and still make use of equation (9), it is necessary to set up a viewing coordinate system as shown in Fig. 7, White (1974). The location of the viewing coordinate system is defined by the points (X,, Y,,Z,) and (X,, Y,,Z,), where (X,, Y,, Z,) is the point from which the surface is to be viewed and (X,, Y,, Z,) is a point on the desired viewing plane. The viewing plane will be assumed to be perpendicular to the line from (X,, Y,,Z,) to (X,, Y,, Z,) and should be located between the object being viewed and the viewer. Based on the coordinates of the points (X,, Y,,Z,) and (X,, Y,,Z,), the coordinate transformation used to transform a point from
Fig. 8. Perspective of a multipatch surface. (a) All the patches plotted. (b) Patches on the front of the surface pLotted
306
P. K. SC~~ERRER and B. M. HILLBERRY
from a canine shoulder. To overcome the problem of hidden points, a procedure was sought for determining which points are on the back side of an object so they would not be plotted. Since the perspective view is obtained with the viewer located on the positive Z axis, a vector perpendicular to the surface at a point on the back side of the surface will have a negative Z component. Therefore, if vectors perpendicular to the surface at the comers of a patch are determined and they all have a positive Z component, the patch will be on the front of the surface and is considered visible, whereas the patch is considered hidden if all the Z components are negative. If any, but not all, of the corner points of a patch are found to lie on the backside of the surface, a vector perpendicular to the surface patch is generated at each point to be plotted to determine if it is hidden. If the Z component for the normal vector at a point is positive, the point is plotted. Figure 8(b) shows the same surface as the one in Fig. 8(a) except the procedure described above was used to remove the points which lie on the back of the surface from the plot. A Coons’ surface patch contains the four corner points for the patch; therefore, there is no error in the representation at these points. At any other location on a patch, a quantitative measure of the error in the patch can be obtained by determining the perpendicular distance between the patch and an actual surface point that the patch would contain if it were an exact representation of the surface. Let ps be a point on the wtual surface, and P, be the point on a surface patch at which a vector perpendicular to the patch intersects ps. The distance from P,v to ps can be used as a measure of the error in the patch representation. The point P, is located at the point on the surface patch where the vector from P,v to ps is perpendicular to the tangent vectors in the u and w parametric directions as illustrated in Fig. 9. Since the scalar product of two perpendicular vectors is zero, the parametric coordinates at which P, is located can be determined by solving the following two simultaneous nonlinear equations in u and w, Scherrer (1977) :
Fig. 9. Point on a Patch t?om which a normal vector intersects a known point.
* Dow-Coming, Midland, Michigan, U.S.A. t Ransom and Randolph Company, Toledo. Ohio, U.S.A.
E&u, w) PsP(f4, w) . ---&= 0
(10)
a&u, w) psp(u, w) ‘~ZO, dW
where PsP(u,w) = &u, w) - lis and P(u,w) is the coordinates of a point on the surface patch at the parametric coordinates (u, w). Equation (10) can be easily solved, using Newton’s method, to find the parametric coordinates at which P,v is located. If P,v does not lie on the patch for which the B matrix was formed, one or both of the parametric coordinates at which P,v is located will exceed the range O-l. MOLDING
THE ARTICULAR
SURFACES
In order to obtain the mathematical representation of a surface using surface patches, it is necessary to know thecoordinates ofthecorner points for the patches. To draw patches on the articular surfaces of the actual bones would have been difficult since the surfaces must be kept moist, and they could not be damaged since they were to be examined histochemically. Another problem with working with the actual articular surface is that it would have been necessary to keep the surface moist while the coordinates of surface points were being measured. Therefore, a silicone rubber mold and then a plaster casting of the articular surfaces were made. Patches weredrawn on thecastings obtained for the articular ends of the scapula and the proximal humerus, and the coordinates of the corner points were measured. Figures 10 and 11 show the patches drawn on the castings for the articular surfaces ofthe shoulder of a 30 kg, 2 year old female German Shepherd Dog. (Test subject GS071576.) In order to avoid loss in accuracy, care had to be taken in choosing the materials to be used in making the mold and casting. The molds were made of Dow Coming Silastic 382 Medical Elastomer* (silicone rubber) and the castings were made of Ransom and Randolpht dental plaster. The dental plaster produced castings that were both accurate and durable. The setting expansion (linear only) for the plaster was between 0.05 and 0.09%. The silicone rubber used was a two-component (catalyst and base), room temperature vulcanizing (RTV) type, The amount of time needed for the silicone rubber to dry depended on the amount of catalyst used. In the molds that were made of the articular surfaces of the canine shoulder, the silicone rubber dried in less than five minutes. Medical grade Silastic was chosen to avoid severe damage to the cartilage. (It appears that the Silastic did affect the c+on the sufface of the cartilage, Kincaid (1977)) Silicone rubber has the desired flexibility needed to remove the articular ends of the bones from the molds, but the shrinkage characteristics over an extended period of time are not very good. According to Braley (1970), at 77°F the shrinkage (linear) is 0.4% after 3 days, 0.6% after 6 days and 0.7% after 14 days. It was
Piecewise
Fig. 10. The surface Fip. Il. Thz surface
mathematical
patches patches
drawn drawn
representation
~2 artiiular
on the casting
of rhe left humerus.
on the casting
surfaces
of rhe left scapul:l.
Piecewise mathematical representation
of articular surfaces
found, however, that shrinkage effects would be neghgible if the castings were made immediately after the bones were removed from the molds. This was determined by making a mold of a U-shaped block consisting of two 25.4 mm square uprights, 50.8 mm high with 1 in. between them. The maximum difference between the dimensions of the block and its casting was found to be 0.05 mm. Thus the molding procedure appears quite accurate. MAPPING THE SURFACES
After the patches were drawn on the surface of the casting (see Figs. 10 and ll), the coordinates of the corner points for the patches were measured using a Brown and Sharpe Validator SO+.The Validator has a granite table, air bearings on all axes, and digital readout in inches or millimeters. ft has a resolution of 0.002 mm (0.0001 in.), and the overall accuracy is k0.012mm (+O.OOOSin.), with an accuracy in a 300 mm (12 in.) square of +_0.008 mm (+_0.0003 in.). The repeatability is given as f 0.005 mm (+O.OOOt in.). In mounting the casting of the articular surface of the humerus. it was necessary to use three set-ups so that the probe carrier was approximately perpendicular to the surface for all measurements. Only one set up was required for the scapula casting. After the coordinates of the surface points were measured, the casting and clamping blocks were removed from the table and then remounted to check the repeatability of the coordinate measurements. The coordinate measurements were then checked by setting the X and Y coordinates of the probe on the Validator to the X and Y coordinates previously measured at a point and comparing the Z coordinate. The maximum difference in the readings for the Z coordinates from the first and second set-up was 0.012 mm (0.0005 in.).
Cramal
Idol
Fig. 13. Dorsal view of the multipatch representation of the proximal articular surface of the humerus.
THE MATHEMATICAL REPRESENTATION OF THE ARTICULAR SURFACES OF THE CANINE SHOULDER
Once the corner points for the patches were measured, the corner point tangent vectors for the surface patches were generated using a piecewise parametric cubic spline (Program TNGNTS was used to generate the tangent vectors for the patches, Scherrer (1977).) Since Coons’ bicubic surface patches with null twist vectors were used to describe the surfaces, the mathematical representations were complete once the tangent vectors were obtained. The coordinates of the corner points and the components of the tangent vectors for the patches were used in a computer program (program PSP, Scherrer, 1977) to obtain perspective views of the multipatch repr-esentation of the surfaces. This provided a qualitative check on the mathematical representations. Perspective views of the articular surfaces of the left shoulder of a German Shepherd Dog (test subject GS071576) are shown in Figs. 12-15. The solid lines in these figures represent the patch boundaries. Figures 12 and 13 show the proximal humeral articular surface, and Figs. 14 and 15 show the articular surface of the
Croniol
Caudal
Fig. 12. Caudal-lateral view of the multipatch representation of the proximal articular surface of the humerus. l
Brown and Sharpe, North Kingstown, RI., U.S.A.
Media\
Fig. 14. End-on view, from a mediil aspect, of the multipatch representation of the articular surface of the scapula.
P. K. SCHERRERand B. .%I.HILLBERRY
310 Lateral
Medial
Fig. 15. Direct end-on view of the multipatch representation of the articular surface of the scapula.
scapula. Due to a concavity in the articular surface of the proximal humerus near the medial osteochondral border, a portion of the medial-caucial articular surface was not mathematically represented. The missing portion of the articular surface is not shown in Figs. 12 and 13. Figures 12-15 can bz compared to Figs. 10 and 11 which show the castings of the articular surfaces of the proximal humerus and scapula with the patches drawn on them. (The sides of the patches on both the humerus and scapula were approximately 3.2 mm (0.125 in) long.) The pronounced curve in the boundary curves of some of the patches of the multipatch surfaces shown in Figs. 12-15 ace not a result of coding errors in any computer routine. For end conditions on the splines used to generate the tangent vectors, the splines were required to contain known points. Thecurves would have been smoother if greater care had been taken in choosing the surface points used for the boundary conditions. The accuracy of the representations were examined quantitatively by determining the perpendicular distance between points on the actual surfaces and the multipatch representations. (This was done using program ERR, Scherrer, 1977.) For the scapula, the coordinates of one point at the center of each of the 48 patches was measured, as well as the coordinates of 116 other points located along the boundary curves of the patches.* It was found that, for the points at the center of the patches of the scapula, the error in the representation ranged from 0.19 mm (0.0075 in.) to -0.21 mm (-0.0084 in.), with the average error being 0.063 mm (0.0025 in.). The average error is defined here as the sum of the absolute value of the error at the l The coordinates of a number of points between the ends of the boundary curves for the patches on the scapular surface were originally measured for use in an experimental computer routine. This routine contained a procedure for varying the magnitude of the tangent vectors at the knots of a piecewise cubic splint to minimize the error in the spline. The error was defined as the perpendicular distance between the spline and points located between the spline and points located between the knots. Since these coordinates had been measured, they were used in examining the error along the patch boundaries. It should also be noted that the routine for minimizing the, error did not appear to make an improvement in the fit.
points, divided by the number of points. In descrrbing the range of the error, a positive distance indicates that the actual surface point was above the representation. There was one point for the scapula at which the computed error was 0.70 mm (0.027-1in.). It appears that this large error was a result of an error in recording the coordinates of the point, and therefore is not representative of the error in the mathematical representation. However. the 0.70 mm (0.0274 in.) error was used in determining the average error of 0.06 mm (0.0025 in.) for the scapula. For the 116 points located along the boundary curves of the scapula, the error ranged from 0.26 mm (0.0103 in.) to -0.22 mm (-0.0088 in.), with an average error of 0.04mm (0.0016 in.). Based on theerror in the 164 points on the scapula that were not corner points. the average error was 0.05 mm (0.0019 in.). The coordinates ofone point at the center of each of the 65 patches on the casting of the humerus were also measured. The maximum error was 0.21 mm (0.0083 in.), the minimum error was -0.152 mm (-0.0060 in.) and the average error was 0.04 mm (0.0017 in.). After determining the error, increased accuracy can be obtained by using smaller patches. After the error was determined at the points mentioned above, it was found that the error at other points should have been examined. According to Sabin (1978), the effect of using zero twist vectors is expected to be largest around the parametric coordinates (0.2, 0.2), (0.2, 0.8). (0.8, 0.2) and (0.8, 0.8). Thus the error should have been examined at these parametric coordinates.
SU,MMARY
Dental plaster castings of the articular surfaces from the shoulder of a German Shepherd Dog were made directly from the shoulder components. Surface patches were drawn on the castings of the articular surfaces, and the coordinates of the corner points on the patches were measured. (As noted previously, patches with sides approximately 3.2 mm (0.125 in.) long were used on both the scapular and humeral castings.) Using the coordinates of the corner points, the tangent vectors used to describe the surface with surface patches were obtained. Perspective views were generated to verify the multipatch representations of the surfaces. The error in the representation was examined by determining the perpendicular distance between the representation and actual surface points. Based on 164 points, the error in the representation of the articular surface of the scapula ranged from 0.26 mm (0.0103 in.) to -0.21 mm (-0.0088 in.), with an average error of 0.04mm (0.0016 in.). For the representation of the articular surface of the proximal humerus, the maximum error was 0.21 mm (0.0083 in.) and the minimum error was -0.15 mm (-0.0060 in.), with the averageerror being0.04 mm (Oo.0017in.). This was based on 65 points. Increased accuracy can b-e obtained by using smaller patches.
Piecewise mathematical representation Ackno%ledgemenrs - This project was sponsored by the National Institute of Arthritis, Metabolism and Digestive Diseases. HEW, under Grant Number 5 ROI 17482. The authors wish to thank Dr. D. C. Van Sickle, DVM, Ph.D. and Dr. G. J. Pijnowski, DVM, School of Veterinary Science and Medicine and Professor J. R. Rice, Department of Computer Science, Purdue University for their assistance on this project.
REFERENCES
Adams, J. A. (1974) Cubic spline curve fitting with controlled end conditions. Comp. Aided Desiyn 6, 2-9. Almond, D. B. (1972) Numerical control for machining complex surfaces. fBM Sysr. J. 150-168. Bezier. P. (1968) How Renault uses numerical control for car body design and tooling. SAE Paper 680010. Bezier, P. (1972) Numerical control in automobiledesign and manufacture of curved surfaces. pp. 44-48, Proc. Curved Surfaces in Engny Conf Churchill College, Cambridge, England. Bezier, P. (1972) Numerical control. Mathematics and Applicarions. Wiley, New York. Braley, S. (1970) The chemistry and properties of the medicalgrade silicons. J. Macromol. Sci.-Chem. 529-544. Coons, S. A. (1967) Surface ofcomputer-aided design of space forms. MIT Project MAC-TR-41. Ferguson, J. (1964) Multivalued curve interpolation. J. ACM 1I, 221-228. Hart, W. B. (1971) Computer graphics aid fuselage design. Enyny Mat. Des. 205-209. Hope, W. M. (1966) Computer-aided design in body engineering. SAE Paper 660153.
B.Y. 12,4-E
of articular surfaces
311
Johnson, W. L.. Sanders, J. W. and South, N. E. (1966) Analytic surfaces for computer-aided design SAE Paper 660152. Kincaid, S. A. (1977) Personal communication. Purdue University. Kinzel, G. L., Hall, A. S.. Jr. and Hillberry, B. ?+t. (1972) Measurement of the total motion between two body segments, Part I - analytical development. J. Biumech. 3, 93- 105. Kinzel, G. L. (1973) On the design of instrumented linkages for the measurement of relative motion between two rigid bodies. Ph.D. Thesis, Purdue University. Kinzel, G. L.. Van Sickle, D. C., Hillberry, 6. M. and Hall. A. S., Jr, (1976) A preliminary study of the in viuo motion in the canine shoulder, Am. J. cer. Res. 37, 1505-15 IO. . Lee, T. M. P. (1969) Three-dimensionof Curves and Surfaces for Ropid Computer Display. Management Information Service, Detroit, Michigan. MacCullum, K. J. (1970) Surfaces for interactive graphical design. Camp. J. 13, 352-358. Sabin, M. A. (1971) An existing system in the aircraft industry. The British Aircraft Corporation numerical master geometry system. Proc. R. Sot. (A)321. 197-205. Sabin, M. A. (1978) Personal communications. Sadehgi, M. M. and Gould, S. S. (1974) A comparison of two parametric surface patch methods. Comp. Aided Design 6. South, N. E. and Kelly, J. P. (1965) Analytic Surface Methods. Numerical Conrrol Decelopment L’nit. Ford \lotor Company, Dearborn, Michigan. White, P. (1974) Notes from ME 573 course. Purdue University. Scherrer. P. K. (1977) Determining the distance between the surface of two rigid bodies using parametric surface _ patches and an instrumented spatial linkage. Ph.D. Thesis, Purdue University.
