ISRN OrthopedicsVolume 2012 (2012), Article ID 840426, 10 pagesdoi:10.5402/2012/840426

Research Article

Geometric Structure of 3D Spinal Curves: Plane Regions and Connecting Zones

E. Berthonnaud,1,2 R. Hilmi,¹ and J. Dimnet³

¹Hôpital Nord Ouest site de Villefranche, BP436, 69655 Villefranche Saone, France²Laboratoire de Physiologie de l'Exercice, Université de Lyon, 42023 Saint Etienne, France³Group of Applied Research in Orthopaedic, 69005 Lyon, France

Received 3 October 2011; Accepted 31 October 2011

Academic Editors: M. Kawakami and E. Pola

Copyright © 2012 E. Berthonnaud et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

This paper presents a new study of the geometric structure of 3D spinal curves. The spine is considered as an heterogeneous beam, compound of vertebrae and intervertebral discs. The spine is modeled as a deformable wire along which vertebrae are beads rotating about the wire. 3D spinal curves are compound of plane regions connected together by zones of transition. The 3D spinal curve is uniquely flexed along the plane regions. The angular offsets between adjacent regions are concentrated at level of the middle zones of transition, so illustrating the heterogeneity of the spinal geometric structure. The plane regions along the 3D spinal curve must satisfy two criteria: (i) a criterion of minimum distance between the curve and the regional plane and (ii) a criterion controlling that the curve is continuously plane at the level of the region. The geometric structure of each 3D spinal curve is characterized by the sizes and orientations of regional planes, by the parameters representing flexed regions and by the sizes and functions of zones of transition. Spinal curves of asymptomatic subjects show three plane regions corresponding to spinal curvatures: lumbar, thoracic and cervical curvatures. In some scoliotic spines, four plane regions may be detected.

1. Introduction

Sagittal balance of asymptomatic subjects and patients has been studied from sagittal radiographic images [1–3]. Spinal curves show generally three regions with opposite curvatures: lumbar, thoracic, and cervical curvatures. Adjacent regional curves are bounded either by specific vertebrae in case of clinical applications or by points of inflexion when the sagittal curve is represented by a mathematical function [4, 5]. Angles of regional curvatures (Cobb angles) may be directly measured upon radiographic images or calculated from the sagittal curve equation.

The mathematical equation of the sagittal spinal projection has been obtained from using a new technique [6]. The projected spinal structure is delimited by two border curves. Each of them is represented by a B-spline curve passing a restricted number of points recorded along the border. The sagittal spinal curve representing the beam projection is defined as the median curve drawn between the two border lines. It is modeled by a B-spline curve. In sagittal projections, each spinal region showing homogeneous curvature is described by parameters defining its shape and orientation [7]. The technique defined for sagittal radiographic projections has been extended to frontal images of scoliotic spines [8], frontal and sagittal radiographic images being shot independently.

The photogrammetric technique reconstructs points in space from their two images in projection plane. The photogrammetry applied to radiographic images has been described by Suh [9]. The photogrammetric reconstruction of spinal curves from simultaneous biplanar radiography has been first presented by Brown et al. [10].

The biplanar radiography involves the setting of specific devices in case of simultaneous exposures [11–13]. Successive exposures have been used for clinical applications. They need a standard radiographic system with one X-ray source and one plate. A rotating plateau is interposed between X-ray source and plate. The patient stands motionless upon the plateau. Patient and plateau are rotated for successive frontal and sagittal radiographic images [14]. On frontal and sagittal projections median spinal curves are drawn. Points of frontal and sagittal curves are related together for the photogrammetric reconstruction of the 3D spinal curve. The relation is based on the use of epipolar planes [15]. 3D spinal curves may be projected in fixed plane. For clinical applications, three projections are available: upper view and frontal and sagittal projections for validating the technical procedure through the comparison between spinal curve projections and radiographic images.

The goal of this paper is to study the geometric structure of spinal curves defined from a series of points reconstructed from photogrammetric techniques. Stagnara et al. [16] demonstrated that scoliotic spines show highest curvatures when shot under specific radiographic incidences. A plane of maximum curvatures of scoliotic spines has been determined from points of spinal curve [14]. This technique was based on the determination of one plane located so that linear distances between 3D spinal curve and plane would be minimum. The maximum distance value governs the relevance of the geometric model representing the spinal curve as a unique plane. If the maximum distance value is too high (up to 10 mm) the unique plane model is maladapted.

The principle of local planes fitted to regions of the 3D rough spinal curve has been retained. Two criteria must be satisfied locating these planes: (i) the maximum distance between the 3D curve facing the regional plane and the unknown plane must be low (2 mm) and (ii) the regional curve must be continuously plane without local anomalies. A great number of spinal curves have an heterogeneous structure. They show plane regions, where the spine curve is purely flexed, and zones of transition between adjacent plane regions. The plane regions often correspond to lumbar, thoracic, and cervical curvatures. In asymptomatic subjects the regional planes are oriented closely. They may be strongly spreaded in patients with spinal deformities. Examples are shown for illustrating the method. Results are presented and discussed.

2. Materials

Frontal and sagittal images of subjects and patients have been obtained from two different radiographic systems.(i)New set-up shooting simultaneously from frontal and sagittal incidences. 125 members of paramedical staffs of a local hospital accepted a biplanar radiographic examination in the frame of a clinical research studying relations between spinal features and chronic low back pain.(ii)A standard radiographic system completed by a rotating plateau has been set in an hospital of Lyon specialized in back pathologies. All patients having back deformities are submitted to a biplanar radiographic examination with successive X-rays with the objective of following the 3D spinal feature evolutions versus time. The access to the geometric structure of spinal curves at different steps of their evolution brings new information for clinicians.

3. Methods

3.1. The Detection of Plane Regions along the Spinal Curve

The detection is based upon two criteria: the criterion of proximity and the criterion of continuous regional planeity.(i)The criterion of proximity to a plane is satisfied when all points of a local section of the rough spinal curve are close to the unknown regional plane. The maximum distance between rough curve and the facing plane is very low (2 mm). Accuracy reconstructing points using biplanar radiography coupled with photogrammetry is 1 mm.(ii)The criterion of continuous regional planeity controls that the plane region is continuously plane. Series of three points, located close together along the spinal curve are successively considered. The normal vector to the local plane passing through the three points is calculated. The bundle of normal vectors to the regional curve is defined. This bundle must be strongly concentrated. Aberrant points along the spinal curve introduce aberrant normal vector directions. These points are eliminated. The direction of the median normal vector is calculated from the bundle of normal vectors. Limit points bounding the plane region are then determined.

Limit points bounding each plane region are located along a diagram representing the unrolled spine from the L₅ lower plate center to the C₃ upper plate center. Each vertebra and each intervertebral disc are located along the unrolled spine using normalized values defined statistically [8]. A correlation is defined between the position of limit points along the unrolled spine and the level of vertebrae or discs.

3.2. Parametring the Orientation of the Normal Vectors to Regional Planes with Respect to the Fixed Set of Axis

The axes of the fixed frame are: vertical axis, posteroanterior axis and transverse axis. The normal vector to the regional plane is defined by its components , , in frame . The direction might be referred to using the standard technique of Euler angles. This technique involves elementary rotations about moving axes. A sequence of two rotations about fixed axes has been preferred. The normal vector is moved from by two successive rotations: an axial rotation about and a lateral flexion θ about θ and are defined by the following relations:

This technique introduces rotations about fixed axes. The direct comparison between orientations of the different plane regions does not depend on the position of the moving axes involved by the Euler calculation.

3.3. Parameters Characterizing the Geometric Pattern of the 3D Spinal Curve

Spinal curves show generally three (or more) plane regions. The planes are defined, with respect to the fixed frame , by the position of their origin (starting point of the plane region) and by the two angles and locating the vector , normal to . All points of the spinal curve , facing , are projected in on . The plane curve passing all points is then obtained.

The geometric structure of the spinal curve is heterogeneous. It has several plane regions , , of normal vectors , , , where the 3D curve is very close to the projected curves , , . Each plane curve is purely flexed. Zones of transitions between sacrum and, between and , and between and are characterized by angular offsets between normal vectors and , and , and and . These angular offsets include axial rotations and abductions needed for moving the sagittal plane to , to , to .

The plane curves are characterized by values of parameters describing their shape [14].(i)Positions of limit points starting points and arrival bounding the plane region .(ii)Length and orientation of the linear segment .(iii)Position of point of maximum regional offset referred to the segment . The tangent in to is parallel to .(iv)Angles and included between segment and tangents to in and respectively.(v)Global flexion angle .

All parameters characterizing the shape of each flexed region are calculated. But main parameters representing the geometric structure of spine, pelvis in upstanding patients, are presented.

The pelvis in erected patient is represented by two parameters.(i)One parameter represents its morphology: the pelvic incidence.(ii)One parameter describes its orientation in erected posture: the pelvic tilt.

The spinal structure is first characterized by the number of its plane flexed regions. Each regional plane is located using the two angles ( axial rotation about and θ lateral flexion about ) locating its normal vector. The limits of each region are identified by their vertebral level. The global flexion describes the regional curvature of the plane curve . All presented results show frontal, sagittal and horizontal projections of the 3D spinal curve with its flexed regions and zones of transitions. Tables of main parameter values are associated to the design of the spinal structure.

4. Results

The goal of this paper is to describe the geometric structure of 3D spinal curves assumed to be previously calculated. Figure 1 illustrates the technique used for accessing spinal curve from biplanar radiography, coupled with photogrammetric reconstruction.

Figure 1: The techniques involved in the calculation of the 3D spinal curve from biplanar radiography coupled with photogrammetric reconstruction.

The results of two different treatments are then compared. The first one assumes that the global spinal curve entirely belongs to a unique plane: the plane of maximum curvature. The second one detects several strictly plane regional sections along the spinal curve. Local geometric discontinuities are then exhibited. In asymptomatic subjects, the spinal structure is compound of three plane regions corresponding to lumbar, thoracic, and cervical curvatures. Short zones of transition are interposed between adjacent regions oriented differently. The modeling technique has been applied to patients suffering from low back pain. Local structural discontinuities of the spinal curve are observed. Scoliosis modifies significantly the spinal structure. An example of evolutive scoliosis is presented. Effects of implant on the spinal structure are shown.

4.1. Illustration of the New Technique

25 asymptomatic adults accepted a biplanar global examination (sagittal and frontal radiographic incidences). Photogrammetric reconstructions allowed to represent 3D models of pelvis and spinal curve. The pelvis is represented as a triangle based on femoral head centers and sacral plate center. 3D spinal curve and pelvic triangle are projected on sagittal plane, frontal plane, and horizontal plane. The spinal curve is compound of three plane regions and three zones of transition.

A case (male 35 years) is presented in Figure 2 the spinal curve is:(i)zone 2 from to : lumbar curvature,(ii)zone 4 from to : thoracic curvature,(iii)zone 6 from to : cervical curvature.

Figure 2: The technique allowing the access to the geometric structure of the spinal curve is illustrated by the results obtained with an asymptomatic subject (male 35 years). Are presented: the sagittal and frontal radiographic images, the sagittal, frontal, and horizontal global projections of the 3D spinal curve and the pelvic modeled as a triangle, the projections of the three plane regions, and the short zones of transition. A table gives numerical values: pelvic incidence PI and pelvic tilt PT, main parameters of the regional plane zones 2, 4, 6: components of the normal vector to the regional plane (ψ axial rotation about , θ lateral flexion about ), vertebral levels of regional limits, magnitude of the flexion angle bending the regional curve, main parameters of the zones of transition 1, 3, 5: vertebral levels of the zone, offsets of angular displacements (Δψ and Δθ).

Three zones of transition are interposed between adjacent flexed regions:(i)zone 1: disc,(ii)zone 3: disc,(iii)zone 5: disc.

The plane regions are characterized by the orientation of the normal vector to the regional plane (axial rotation about and lateral flexion θ about ) and by the angle flexing the regional section of the spinal curve.

Angular offsets between normal vectors to adjacent regional planes: Δ (axial rotation offset) and Δθ (lateral flexion offset) move adjacent regional planes at level of zones of transition. When offsets are small, the zone of transition is restricted to a unique intervertral disc.

4.2. Geometric Spinal Structure Described Using Two Different Techniques

These techniques locate the unique plane of maximum curvature and the heterogeneous structure of spine compound of flexed regions and zones of transition.

Two examples of asymptomatic subjects are presented (Figure 3) for illustrating the information brought by each of the two techniques:(i)the plane of greatest curvature (PGC) is characterized by the orientation of the global normal vector to the plane (angles about and θ about ). The linear fitting coefficient of the rough global curve to the PGC is represented by the maximum normal distance δ eg from 3D curve to PGC,(ii)the regional planes are represented by the orientations of their normal vector ( and θ) and by regional fitting linear δ lr and angular δ ar coefficients. They represent linear and angular maximum offsets between 3D curve and regional planes.

Figure 3: Cases of two asymptomatic adults (40 and 30 years). The comparison between modeling the spinal curve structure using a unique plane of greatest curvature PGC and the model involving regional planes and zones of transition. For each asymptomatic subject, are drawn: the PGC and the three regional planes and the zones of transition. The two PGCs have similar positions, but the regional structures are strongly different.

The presented asymptomatic subjects have different spinal curve features: their PGC is not significantly rotated versus the sagittal plane . In each case, the three regional planes are oriented differently, and their limits are different. The angular contribution of zones of transition may explain that these zones are lengthened.

4.3. Spinal Feature of a Patient Suffering from a Low Back Pain

90 patients suffering from low back pain accepted a biplanar radiographic examination. Figure 4 presents radiographic images of an adult male (60 years old) suffering from low back pain after a long practice of different sports (rugby, golf). The flexed region corresponding to the lumbar curvature is strongly rotated with respect to the anatomical sagittal plane. High level of axial rotation offsets between lumbar curvature versus sacrum and lumbar curvature versus thoracic one are distributed along zones of transition ( and , resp.) which include several intervertebral discs. The local heterogeneity of the spinal structure may be related to the patient pathology.

Figure 4: Patient suffering from a low back pain. The geometric structure is characterized by a strong rotation of the spinal lumbar plane plane region and by a correlative increasement of the sizes of zones of transition.

4.4. Scoliotic Patients

4.4.1. The Spinal Structure Deformations due to the Pathology (Figure 5)

Figure 5: The scoliosis alters significantly the geometric structure of the spine.

The radiographic file of an adult patient is presented for illustrating different effects of the pathology on the structure of the spinal curve. The whole spinal structure is strongly rotated versus the patient sagittal plane, main part of the rotation concerning the lumbar plane region. The thoracic curvature (from to in asymptomatic subject) is shared in two regions of opposite curvatures: – and –. This last flexed region includes the cervical curvature. Zones of transition including several intervertebral discs are related to the high values of axial rotation of regional planes.

4.4.2. The Structural Analysis of Evolutive Scoliosis (Figure 6)

Figure 6: The different steps of scoliotic evolutions. The new technique allows medical doctors to quantify local effects of deforming pathologies.

A scoliotic patient accepted four biplanar radiographic examinations after orthopaedic treatment along a period of eight years (before, 2 years, 5 years, 8 years). The spinal structure keeps comparable global features versus time, but spinal deformations are related to pelvic tilting. After eight years, the set pelvis spine recovers its initial feature. It can be observed that regional planes slightly move between two successive radiographic examinations but corresponding flexion angles keep roughly constant. The regional curvatures are submitted to small rotations versus time, but flexed shapes are kept.

4.4.3. Effect of a Surgical Operation Upon the Scoliotic Global Structure (Figure 7)

Figure 7: Illustration of the effects of a surgical operation upon a deformed spine.

Three radiographic files describe the structural changes brought by a surgical operation: initial state, 6 months (before the surgical operation), and one year (after the operation).

The two first radiographic files show the structural evolution of the spinal curve before surgical operation: the pelvic tilt is the same, the lumbar regional plane slightly increases its axial rotation versus the patient sagittal plane. The angular offset between lumbar and thoracic regions is increased.

The surgical operation consisted in implanting a rod between and , the zones of rod implantation corresponding to the limits of lumbar regional plane. The size of the lumbar regional curvature is reduced (from –, it becomes –). The zones of transition (where abduction and axial rotation are concentrated) are moved: before operation, they are , and , after operation, new values are , , .

The surgical operation reduces significantly the axial rotation value of the lumbar region (28.4° versus 68.9°). Axial rotation offsets between adjacent plane regions strongly decrease: from −24.1° to −3.9° for the lumbar/thoracic transition zone and from −24.5° to −6.9° for the zone of transition between thoracic and cervical regions.

The surgical operation acts directly upon the orientation of the lumbar region, however the lumbar flexion angle is roughly constant. This operation has indirect effects upon the global spinal structure: initial discontinuities between plane regions are strongly decreased.

5. Discussion

In perfectly balanced subjects, spinal curves belong entirely to the anatomical sagittal plane. Lumbar, thoracic, and cervical curvatures are bounded by points of inflexion of the plane spinal curve. Our objective was to describe the spinal geometric structure of spines weakly or strongly deformed.

The mathematical representation of spinal curves is a prerequisite for accessing this goal. In the present study, it involved the access to radiographic systems giving images of the motionless patient shot under different incidences. The radiographic technique is coupled with photogrammetric reconstruction of points from their two images. Spinal curves are defined from series of isolated points.

Plane regions are detected along the rough spinal curve. This detection requires two successive steps, the first one locates its position, the second one eliminates aberrant points and restrict the regional extent so that two criteria would be satisfied. The first one checks that linear distances between points of the spinal curve and facing plane are lower than a maximum value. The second criterion verifies that planes passing through three adjacent points along the curve are very closed to the regional plane. The regional normal vector to the regional plane is calculated. The rough spinal curve facing each regional plane is projected on this plane. The regional plane projection of the 3D rough curve is uniquely flexed. The flexion angle bending this region is then obtained.

Zones of transition between two adjacent plane regions must transmit the angular offsets between the two plane orientations.

Directions of normal vectors to regional planes are related to anatomical axes (fixed axes) using two angles. The normal vector direction is obtained from two successive rotations moving the transverse axis to its actual location: one axial rotation about the fixed vertical axis followed by one rotation about the fixed posteroanterior axis . This solution has been preferred to the Euler sequence, because it involves rotations about moving axes. The rotation values about fixed axes may be licitly used when comparing the different orientations of regional planes for each examination so as their evolution versus time.

Geometric structures of spinal curves have been first estimated using different techniques. Experimental testing demonstrated the existence of an optimal radiographic incidence: the spinal curve shows its greatest curvatures when projected on the plane perpendicular to this specific incidence [17]. The plane of greatest curvature (PGC) has been then defined mathematically [14]. The presented technique best fits planes to local sections of the global spinal curve. The fitting calculation is based on the same principle as used for locating the plane of greatest curvature PGC: it minimizes the normal distance value between rough curve and fitted plane. A criterion of continuous planeity completes, in regional plane technique, the effects of the linear distance minimization. It was impossible to apply this criterion of continuous planeity to the global spinal curve. At level of points of inflexion, between opposite curvatures, local planes passing through three next points of the spinal curve are inaccurately located. This phenomenon is eluded with the technique of regional planes because changes in regional plane orientations are located at level of zones of transition.

The results show that asymptomatic subjects have closely oriented regional planes. The drawing of a unique plane of greatest curvature PGC may be licitly proposed. In other asymptomatic subjects, the proposal of a unique plane of greatest curvature is not obvious. The 3D location of unknown plane PGC is always possible, but linear distances between spinal curve and PGC may reach high values, not consistent with an accurate geometric modeling. The proposed geometric model of the spinal structure displays its heterogeneity. Even in asymptomatic subject, some regional sections of the spinal curve belong to planes oriented differently. In these regions the spinal curve is purely flexed. Axial rotation and abduction components due to the different orientations of adjacent regional planes are concentrated at level of zones of transition. These zones are restricted to a unique intervertebral disc when the adjacent regional planes are close together. They are longer in patients showing high regional discontinuities.

The proposed model has been applied to asymptomatic subjects. Local unbalances are displayed. It has also been tested with patients suffering from back pathologies, chronic low back pain and scoliosis. These pathologies entail spinal structure deformations. Geometric description of spinal curve features may allow clinical people to propose well-adapted correcting techniques.

A large diffusion of this technique involves the access to the mathematical representation of the spinal curve. Local solutions have been found from using biplanar radiography with successive or simultaneous exposures coupled with photogrammetric reconstructions. But this modeling technique may be adapted to any imaging system giving the accurate representation of spinal curves.

Conflict of Interest Statement

E. Berthonnaud, R. Hilmi, and J. Dimnet disclose any financial and personal relationships with other people or organisations that could inappropriately influence their work

References

R. P. Jackson, M. D. Peterson, A. C. McManus, and C. Hales, “Compensatory spinopelvic balance over the hip axis and better reliability in measuring lordosis to the pelvic radius on standing lateral radiographs of adult volunteers and patients,” Spine, vol. 23, no. 16, pp. 1750–1767, 1998.

M. Bernhardt and K. H. Bridwell, “Segmental analysis of the sagittal plane alignment of the normal thoracic and lumbar spines and thoracolumbar junction,” Spine, vol. 14, no. 7, pp. 717–721, 1989.

R. Vedantam, L. G. Lenke, J. A. Keeney, and K. H. Bridwell, “Comparison of standing sagittal spinal alignment in asymptomatic adolescents and adults,” Spine, vol. 23, no. 2, pp. 211–215, 1998.

M. S. Goldberg, B. Poitras, N. E. Mayo, H. Labelle, R. Bourassa, and R. Cloutier, “Observer variation in assessing spinal curvature and skeletal development in adolescent idiopathic scoliosis,” Spine, vol. 13, no. 12, pp. 1371–1377, 1988.

C. Boulay, C. Tardieu, J. Hecquet et al., “Sagittal alignment of spine and pelvis regulated by pelvic incidence: standard values and prediction of lordosis,” European Spine Journal, vol. 15, no. 4, pp. 415–422, 2006.

E. Berthonnaud and J. Dimnet, “Analysis of structural features of deformed spines in frontal and sagittal projections,” Computerized Medical Imaging and Graphics, vol. 31, no. 1, pp. 9–16, 2007.

E. Berthonnaud, J. Dimnet, P. Roussouly, and H. Labelle, “Analysis of the sagittal balance of the spine and pelvis using shape and orientation parameters,” Journal of Spinal Disorders and Techniques, vol. 18, no. 1, pp. 40–47, 2005.

E. Berthonnaud and J. Dimnet, “Fast calculation ofparameters ofscolioses infrontal view forclinical applications,” ITBM-RBM, vol. 27, no. 2, pp. 56–66, 2006.

C. H. Suh, “The fundamentals of computer aided X ray analysis of the spine,” Journal of Biomechanics, vol. 7, no. 2, pp. 161–169, 1974.

R. H. Brown, A. H. Burstein, C. L. Nash, and C. C. Schock, “Spinal analysis using a three dimensional radiographic technique,” Journal of Biomechanics, vol. 9, no. 6, pp. 355–365, 1976.

A. Le Bras, S. Laporte, D. Mitton, J. A. De Guise, and W. Skalli, “3D detailed reconstruction of vertebrae with low dose digital stereoradiography,” Studies in Health Technology and Informatics, vol. 91, pp. 286–290, 2002.

J. Dubousset, G. Charpak, I. Dorion et al., “A new 2D and 3D imaging approach to musculo-skeletal physiology and pathology with low-dose radiation and the standing position: the EOS system,” Bulletin de l'Academie Nationale de Medecine, vol. 189, no. 2, pp. 287–300, 2005.

B. Andre, J. Dansereau, and H. Labelle, “Optimized vertical stereo base radiographic setup for the clinical three-dimensional reconstruction of the human spine,” Journal of Biomechanics, vol. 27, no. 8, pp. 1023–1035, 1994.

E. Berthonnaud, J. Dimnet, and R. Hilmi, “Classification of pelvic and spinal postural patterns in upright position. Specific cases of scoliotic patients,” Computerized Medical Imaging and Graphics, vol. 33, no. 8, pp. 634–643, 2009.

E. Berthonnaud, B. Moyen, and J. Dimnet, “Stereoradiographic techniques applied to three-dimensional clinical measurements,” Automedica, vol. 18, pp. 321–363, 2000.

P. Stagnara, J. C. De Mauroy, and G. Dran, “Reciprocal angulation of vertebral bodies in a sagittal plane: approach to references for the evaluation of kyphosis and lordosis,” Spine, vol. 7, no. 4, pp. 335–342, 1982.

J. Peloux, R. Fauchet, B. Faucon, and P. Stagnara, “Le plan d'élection pour l'examen radiologique des cypho-scolioses,” Revue de Chirurgie Orthopedique et Reparatrice de l'Appareil Moteur, vol. 51, no. 6, pp. 517–524, 1965.