The use of sparse CT datasets for auto-generating accurate FE models of the femur and pelvis
Article Outline
Abstract
The finite element (FE) method when coupled with computed tomography (CT) is a powerful tool in orthopaedic biomechanics. However, substantial data is required for patient-specific modelling. Here we present a new method for generating a FE model with a minimum amount of patient data. Our method uses high order cubic Hermite basis functions for mesh generation and least-square fits the mesh to the dataset. We have tested our method on seven patient data sets obtained from CT assisted osteodensitometry of the proximal femur. Using only 12 CT slices we generated smooth and accurate meshes of the proximal femur with a geometric root mean square (RMS) error of less than 1
mm and peak errors less than 8mm. To model the complex geometry of the pelvis we developed a hybrid method which supplements sparse patient data with data from the visible human data set. We tested this method on three patient data sets, generating FE meshes of the pelvis using only 10 CT slices with an overall RMS error less than 3mm. Although we have peak errors about 12mm in these meshes, they occur relatively far from the region of interest (the acetabulum) and will have minimal effects on the performance of the model. Considering that linear meshes usually require about 70–100 pelvic CT slices (in axial mode) to generate FE models, our method has brought a significant data reduction to the automatic mesh generation step. The method, that is fully automated except for a semi-automatic bone/tissue boundary extraction part, will bring the benefits of FE methods to the clinical environment with much reduced radiation risks and data requirement.
Keywords: Finite element modelling, Automatic mesh generation, CT, Radiation dosage
1. Introduction
Since its introduction in 1972 (Brekelmans et al., 1972) the finite element (FE) method has been widely used in orthopaedic biomechanics from stress strain calculation to micro-structural modelling of bones (Huiskes and Chao, 1983; Huiskes and Hollister, 1993). However, it has not yet reached widespread clinical use mainly due to the intensive labour needed in generating accurate FE meshes. There have been a number of attempts to develop an automated method for generating FE models. Keyak and coworkers developed an efficient and robust automatic mesh generation (AMG) method (Keyak et al., 1990) that has proved to be useful in many research areas (Keyak et al., 1998, Keyak et al., 2001; Cody et al., 1999). Viceconti and coworkers developed an automatic method that generates patient-specific FE models of various different types of bones other than the femur (Viceconti et al., 2004).
All the aforementioned methods use computed tomography (CT) to obtain patient specific geometry and material properties. The data requirements are quite intensive, typically requiring very dense FE meshes with thousands of small elements. This will inevitably lead to increased radiation exposure to subjects/patients.
There have been a number of concerns about the possible adverse effects of CT radiation exposure (Rehani et al., 2000; Pierce and Preston, 2000; Schmidt et al., 2002). As a result, significant efforts have been made to reduce radiation dose mainly by controlling scanning parameters such as tube current, voltage, and scanning modes (Kalra et al., 2003). There have also been some efforts to optimize scanning protocols to reduce radiation exposure (Lattanzi et al., 2004; Zannoni et al., 1998). However, to the authors’ knowledge, there is no published study that develops an AMG technique with sparse CT data.
Moreover, it is common for CT datasets to be limited to the part of a bone that is a region of clinical interest such as the acetabulum of the pelvis. Data describing a structure such as the pelvic blade may be very limited or missing. However, to have a full FE model of a bone for mechanical simulations, it is vital to have structural information for the whole bone.
In this paper we provide a new mesh generation method for models of the femur and pelvis from sparse CT datasets. We also provide a hybrid method for generating pelvic FE models, which supplements sparse patient pelvic data with geometric data from the visible human dataset (Ackerman, 1998). The models use high order cubic Hermite elements which have the advantage of capturing complex geometry using few elements (Bradley et al., 1997).
2. Materials and methods
Our AMG method uses serial cross sections to create a mesh composed of eight-noded hexahedrons and six-noded pentahedrons in the following three steps: (1) boundary extraction from CT slices; (2) fitting a predefined grid to boundaries in the CT slices, where the grid will be defined by a user as numbers of rows and columns in the grid; (3) connecting fitted grids to form a three-dimensional mesh.
The CT slices are segmented using Deriche's edge detection algorithm (Deriche, 1990) to find bone/soft tissue boundaries. Rectangular and elliptical grids (Fig. 1) are used in generating meshes of the pelvis and femur, respectively. Since we use high order cubic Hermite basis functions for the FE meshes local element coordinate 
directions need to be kept consistent.
Fig. 1
Two types of grids used in our mesh generation. The grid on the left was used for the pelvis and the one on the right was used for the femur. Arrows show that
directions are kept consistent.
The chosen grid is fitted to the image. Firstly, the size of the object is determined by sweeping the image with a straight line. While sweeping the image, we also search for some special features in the image such as a hole or multiple objects and save their locations. We then divide the grid into a given number of evenly distributed rows. If the image has special features such as holes or multiple objects, we place rows at their boundaries rather than evenly distributing them. When the number of rows given is not enough to describe the special features in the image we add additional nodes as shown in Fig. 2. After distributing rows, the width of the object that each row covers is evenly divided into columns and nodes are placed on them. We store the location of each node but if it is placed on an empty space such as a hole, we store it as an empty node. Additional nodes are also introduced when the given number of columns are not enough to describe the special features in the image (Fig. 2).
Fig. 2
Fitting of a lattice to a given boundary. The flowchart on the left shows the steps in our mesh generation method and the pictures on the right explain the steps. (1) Shows cross sections and their initial boundaries. Images are from the illium of the pelvis, the shaft of the femur and the bottom condyle of the femur. Note that different types of grids were used for each case and for the hole boundary. The rows have been moved to coincide with the start and end of the hole. (2) Shows how grids are fitted and nodes are distributed. For simple cross-sections like the illium in 2.1, nodes are evenly distributed. However, for hole and bifurcation, rows are placed at their boundaries and additional nodes are introduced to accommodate geometry change. In 2.2, the number of columns given is two and the hole inside the object cannot be described with two columns. Hence additional nodes are introduced. In 2.3, the number of rows given is 6 and the gap between the two objects cannot be described accurately with 6 rows. Hence additional nodes are introduced. (3) Shows the 3D elements formed by connecting two neighbouring grids. Since the initial grid over the whole area is covered by object and only non-empty nodes are used in forming an element holes and bifurcations are created. (a) Mesh generation flow chart. (b) Examples of each step in mesh generation.
Once all the nodes are placed, the grids are connected to their neighbouring grids to form a three-dimensional mesh. We only use non-empty nodes to form an element and if more than six neighbouring nodes are not empty, then we connect them to form an element. Holes and bifurcations are created in this way as shown in Fig. 2. The number of nodes available to form an element determine the type of the element as either a hexahedron or a pentahedron. After the initial mesh is generated, we apply a least-squares fitting algorithm using all data points from the serial cross sections (Fernandez et al., 2003). The resulting meshes give an accurate and smooth depiction of the bone surface geometry as seen in Fig. 3 where we applied our mesh generation method on the femur and pelvis from the Visible Human male data set.
Fig. 3
Initial linear mesh and final fitted mesh. The left column shows linear meshes generated from our automatic mesh generation method. The right column shows the final fitted meshes resulting from the least-squares fitting. (a) Linear mesh before fitting. (b) Cubic mesh after fitting.
When detailed patient-specific data is not available and to accommodate the complexities of pelvic geometry, we supplement the patient CT slices (concentrated in the region of interest, the acetabulum) using CT data from the Visible Human. In this hybrid method, Visible Human CT slices are transformed to match the size and orientation of the patient CT sets. The transformation is done using the landmark-based method of Challis (1995). We extended this method so that the corresponding landmark points are automatically selected from CT slices. Our method extracts the skeleton of the object in a CT slice and then uses it to select landmark points (Fig. 4).
Fig. 4
Steps in skeletonization to extract landmark points from a boundary. (a) The boundary of a section in the acetabulum; (b) the skeleton of the boundary in dotted line; (c) the boundary, skeleton and landmark points found. The solid lines drawn from the skeleton to the boundary are rays. The points where rays meet with the boundary are landmark points.
To get the accurate transformation from the Visible Human data to our patient specific one, we divide the pelvis into four regions according to their common geometric characteristics as shown in Fig. 5. Each region is described by a transformation matrix calculated from the patient CT slice at the start/end of the region and its matching Visible Human slice, which is obtained by the procedure explained in Fig. 6.
Fig. 5
Pelvic regions and the way that extra slices were placed. Extra slices are shown as lines. (a) Four principal pelvis regions assumed in our study. Those four regions are the blade region comprising the wide section of the illium just above the greater sciatic notch, the isthmus region from the greater sciatic notch to just above the acetabulum, the acetabulum region, and the foramen obturatum region surrounded by the pubis and the ischium. (b) A case where the foramen obturatum region was divided into two. Note that two extra slices were used; one at the start of the blade and the other at the middle of the foramen obturatum region. The isthmus region was described by the last slice of the acetabulum. (c) A case where the foramen obturatum and the blade regions were divided into two. Note that three extra slices were used. (d) A case where the foramen obturatum region is divided into two and the blade was divided into three. Four extra slices were used.
Fig. 6
Overall process of generating a hybrid dataset composed of Visible Human CT slices as well as patient CT slices. The flow chart on the left shows the overall steps in the method and the one on the right shows how the matching Visible Human slices were selected.
Once all of the matching Visible Human CT slices are found, we determine which left-over Visible Human CT slices should fill the areas that are not covered by patient CT slices. These slices are transformed using the transformation matrix of the region and form a complete hybrid data set along with the patient slices. Then we use our AMG method described previously to create meshes of the pelvis. Each region will have its own transformation matrix and this can result in discontinuities between regions that are removed by introducing a Sobolev smoothing term on our least-squares fitting algorithm (Young et al., 1989). Our method requires two inputs from the user: (1) how many extra slices should be used; (2) how the empty regions in the patient dataset should be filled.
The AMG method was tested using femoral CT datasets from seven patients with unilateral hip joint osteoarthritis, who had an uncemented total hip replacement inserted. The CT scans were collected 10 days after the operation in axial mode and the slice thickness 2mm. The table feed was 5mm proximal to the lesser trochanter and 10mm distal to it. Five different meshes were developed from each patient data set by using different numbers of CT slices in order to investigate the minimum number of CT slices needed to generated accurate FE meshes (Table 1).
Table 1. Locations and numbers of CT slices used for the five different mesh types tested
| Below lesser trochanter | Intertrochanteric region | Greater trochanter | Average number of slices used | |
|---|---|---|---|---|
| Mesh 1 | 4cm gap | 2cm gap | 2cm gap | 7 |
| Mesh 2 | 2cm gap | 2cm gap | 2cm gap | 9 |
| Mesh 3 | 2cm gap | 1cm gap | 1cm gap | 12 |
| Mesh 4 | 1cm gap | 1cm gap | 1cm gap | 17 |
| Mesh 5 | 1cm gap | 1cm gap | 0.5cm gap | 22 |
The hybrid method for generating pelvic meshes was tested using three patient CT (Somatom Plus 4, Siemens) data sets of the abdomen and pelvis (2 male, 1 female) from a routine abdominal scanning. The slide thickness was 2mm and the table feed was 8mm. We generated four different types of mesh by utilizing only a limited number of CT slices among all the available slices as shown in Fig. 5. Details on how many slices were used are shown in Table 2.
Table 2. Four different types of pelvis mesh generated
| No. of slices below the acetabulum | No. of slices in the acetabulum | No. of slices above the acetabulum | |
|---|---|---|---|
| Mesh 1 | 1 | All of the slices (4–6) | 1 |
| Mesh 2 | 1 | All of the slices (4–6) | 2 |
| Mesh 3 | 1 | All of the slices (4–6) | 3 |
| Mesh 4 | 1 | All of the slices (4–6) | 4 |
The accuracy was tested by calculating the RMS error, which was the distance between the data points collected from all of the CT slices (11,142 points on average for the femur and 10,454 points on average for the pelvis) and the fitted surface. The maximum error and its location were also investigated.
3. Results
When all the available slices were used in generating the femoral meshes, the average RMS error after fitting was less than 0.5mm. Meshes generated with only a limited number of CT slices had larger errors, and as the number of slices used decreased, the error increased. The overall trend in the RMS error is shown in Fig. 8 and the final fitted meshes are shown in Fig. 7.
Fig. 7
Final fitted meshes generated. (a) Mesh 1 (7 slices); (b) mesh 2 (9 slices); (c) mesh 3 (11 slices); (d) mesh 4 (17 slices); (e) mesh 5 (22 slices). The bottom row shows the error plot where the maximum error was depicted both in colour and with the length of the error bar. Maximum error is indicated by red lines and the colour scale is on the bottom right corner.
As can be seen from Fig. 8, the RMS error dropped significantly (more than 50%) when 10 or more slices were used. The error plot from Fig. 7 also shows a large maximum error in the femoral neck region reduced dramatically as we used 10 or more slices. Detailed results of the error in terms of the RMS and maximum errors are given in Table 3 for all of the meshes generated. All the meshes have similar RMS errors but maximum errors vary among patients especially when less than 10 slices were used in mesh generation. However, when more than 10 CT slices were used, the localised large maximum error dramatically reduced as can be seen from Fig. 7.
Fig. 8
Error vs number of slices used in mesh generation.
Table 3. Error analysis for all the femurs used
| Mesh 1 | Mesh 2 | Mesh 3 | Mesh 4 | Mesh 5 | ||||||
|---|---|---|---|---|---|---|---|---|---|---|
| RMS | Max | RMS | Max | RMS | Max | RMS | Max | RMS | Max | |
| Femur 1 | 1.86 | 15.00 | 1.75 | 7.14 | 0.69 | 6.46 | 0.67 | 6.55 | 0.33 | 4.63 |
| Femur 2 | 2.41 | 15.03 | 2.53 | 15.00 | 0.93 | 8.15 | 0.83 | 8.16 | 0.38 | 4.46 |
| Femur 3 | 1.75 | 10.05 | 1.67 | 10.00 | 0.59 | 5.06 | 0.52 | 5.06 | 0.37 | 3.64 |
| Femur 4 | 1.54 | 10.04 | 1.39 | 9.91 | 0.70 | 5.66 | 0.70 | 5.33 | 0.39 | 3.67 |
| Femur 5 | 2.15 | 10.45 | 1.95 | 10.42 | 0.78 | 7.99 | 0.46 | 4.73 | 0.30 | 3.54 |
| Femur 6 | 1.99 | 12.89 | 1.99 | 12.98 | 0.88 | 6.04 | 0.72 | 5.17 | 0.45 | 5.09 |
| Femur 7 | 1.89 | 10.64 | 1.68 | 10.59 | 1.05 | 6.07 | 0.61 | 5.14 | 0.31 | 2.49 |
The RMS errors from the pelvic side, where we used the hybrid method were larger than the femoral side. When compared with the whole patient data set, including those that were not used in mesh generation, the error increased from 2.7mm to 4.85mm as we decreased the number of extra slices above and below the acetabulum. When we used 10 or more slices (four or more extra slices) the average RMS error became less than 3mm. Visual comparison of the meshes generated from the full patient dataset and from the limited dataset is given in Fig. 9 and the trend of average RMS errors with respect to the number of slices used is shown in Fig. 11. Final meshes generated from all three patient datasets are shown in Fig. 10, which also shows the maximum errors and their locations. Maximum errors are localised in different spots for each mesh and they occur close to the area where the geometry changes rapidly such as the greater sciatic notch and the pubic crest (Figs. 10, 11).
Fig. 9
Comparison of resulting meshes from the hybrid dataset as described in Section 2.3 to the one generated with the full patient dataset. (a) Mesh generated with 11 CT slices; (b) mesh generated with all of the available CT slices (27 slices).
Fig. 10
Pelvic meshes generated from three patient datasets and the locations of the maximum errors. Note that the maximum error is denoted as red on the mesh and the scale bars on the bottom right corner shows the colour scales.
Fig. 11
Error vs number of slices used in mesh generation.
Detailed results of the maximum errors and the RMS errors are shown in Table 4 for all of the pelvic meshes generated. As with the femoral cases, the RMS errors are quite consistent across patients while the maximum errors varied, especially when less than 10 CT scans were used in mesh generation. However, when four extra slices were used, the maximum errors were similar across three patients and were always less than 13mm.
Table 4. Error analysis for all the hips used
| Mesh 1 Error | Mesh 2 Error | Mesh 3 Error | Mesh 4 Error | |||||
|---|---|---|---|---|---|---|---|---|
| RMS | Max | RMS | Max | RMS | Max | RMS | Max | |
| Pelvis 1 | 5.01 | 19.08 | 3.48 | 13.67 | 2.85 | 13.23 | 2.77 | 11.11 |
| Pelvis 2 | 4.695 | 31.72 | 3.876 | 25.50 | 2.87 | 18.34 | 2.77 | 11.20 |
| Pelvis 3 | 3.828 | 12.61 | 3.99 | 18.76 | 3.187 | 18.49 | 2.48 | 12.97 |
4. Discussion
We have developed a new method for automatic FE mesh generation that is computationally efficient and accurate and that reduces radiation exposure from CT. There are a number of practical applications associated with our method such as large scale clinical studies and preoperative planning and computer assisted surgery.
The reduction in data requirement was possible due in part to the computational efficiency of the cubic Hermite basis functions. These basis functions provide continuity of first derivatives as well as nodal values, which means that fewer elements are needed to describe complex geometries (Nielsen et al., 1991). As cubic Hermite elements have higher degrees of freedom (DOF) per node (24 per node) than linear meshes (3 per node), we can achieve better accuracy than linear meshes with much fewer elements.
In terms of radiation reduction, the direct comparison of our method to those of other researchers is impossible because of different CT scanning parameters. However, since the total irradiated dose is proportional to the number of CT slices, it follows that by reducing the number of slices required, we will directly reduce the dose, assuming that slice thickness stays the same.
A more practical application of our work will be to use it in generating a patient-specific FE mesh using an incomplete patient dataset. Since patient CT scans are usually limited to only a part of a bone we have shown here how limited datasets can be used in generating FE models by supplementing them with a complete dataset such as the Visible Human. Moreover, the automatic nature of our method makes it an attractive candidate for processing a large number of patient datasets to examine the behaviour and effects of a hip prosthesis after surgery.
We also acknowledge the limitations of our study. Firstly the underlying assumption of our hybrid method is that the angular orientation of the patient pelvis is the same as that of the Visible Human pelvis, which was scanned while lying supine on the table. This is due to the fact that we compare the patient dataset with the Visible Human on a slice by slice basis, hence the slices from both datasets should run almost parallel to each other. However, we believe that this will not restrict our method because lying supine on the table is the most common positioning method in abdominal CT scanning. Moreover, one can resample the Visible Human dataset at another angle if necessary to match the degree of tilting in the patient. Secondly, our method cannot be applied if a severe deformity is present as we merge the patient data with that from a normal person (VH dataset). However, it would be unwise to reduce the amount of CT data if deformities are present and if one wants to visualise and model the deformities. Therefore, we do not suggest to use limited or incomplete CT datasets in such circumstances. Thirdly the generality of our method will not be as broad as the methods of Viceconti et al. (2004) or Keyak et al. (1990) which can even be applied to micro CT data. However, our method is general enough to accommodate the patient specificity of a known bone geometry (pelvis and femur), which is most frequently required in clinical settings. Moreover, it can be easily expanded to other bones such as the tibia or fibula. Finally, relatively high peak errors occur in some regions of the pelvis (Fig. 10). These peak errors get worse as the mesh becomes coarser. It is unknown whether these peak errors may compromise the overall accuracy of the mesh at this stage. However, it should also be noted that locations of high peak errors are relatively far away from the region of interest, which is the acetabulum. Therefore, we are optimistic that effects of these errors will not dominate the overall performance of the mesh in predicting stress/strain in the acetabulum.
Bone tissue material properties can be obtained from CT scans by using an empirical relationship between bone density and modulus such as the one by Carter and Hayes (1977) These can be assigned to our mesh using a spatially varying field where each Gauss point in elements is assigned with a different material property value. Thus, each Gauss point (not element) will have a different modulus value.
Gauss points work as material points where the constitutive relationship is evaluated. Hence, when using the Gaussian quadrature scheme to evaluate integrals, we use the modulus value that we saved at each Gauss point rather than using one modulus for all the Gauss points in the elements. Therefore, our elements can have as many material property values as Gauss points in the element. This effectively captures the variation in material property that occurs within a single element when its size is large. This method and implementation has been successfully applied in various biomechanical situations (Stevens et al., 2003; Fernandez and Hunter, 2005; Shim et al., 2004) using a FE software called CMISS1. Moreover, it has also been shown that meshes with this type of element give an equivalent result as having a global mesh refined to the same extent (Edidin and Taylor, 1991).
The material properties themselves can be obtained in a similar manner to the geometry—the region of interest from the patient dataset and the rest from the Visible Human dataset—VH CT slices are scaled for use in assigning material properties to the areas where the patient data is not available. The result of convergence studies to validate the element quality will be reported in a follow-up paper.
In conclusion, it can be said that our method offers a new and efficient way of generating patient specific FE models that can be applied to clinical studies for preoperative planning and computer-assisted surgery.
References
- . The visible human project. Proceedings of the IEEE. 1998;86(3):504–511
- . Geometric modeling of the human torso using cubic hermite elements. Annals of Biomedical Engineering. 1997;25:96–111
- . A new method to analyse the mechanical behaviour of skeletal parts. Acta Orthopaedica Scandinavica. 1972;43(5):301–317
- . The compressive behavior of bone as a two-phase porous structure. Journal of Bone and Joint Surgery. 1977;59:954–962
- . A procedure for determining rigid body transformation parameters. Journal of Biomechanics. 1995;28(6):733–737
- . Femoral strength is better predicted by finite element models than qct and dxa. Journal of Biomechanics. 1999;32(10):1013–1020
- . Fast algorithm for low-level vision. IEEE Transactions on Pattern Analysis and Machine Intelligence. 1990;12(1):78–87
- . Use of variable stiffness elements is equivalent to global mesh refinement in a bone-alone proximal femoral finite element model. Transactions of the annual meeting of the Orthopaedic Research Society. 1991;16:309
- Fernandez, J., Hunter, P., 2005. An anatomically based patient-specific finite element model of patella articulation: towards a diagnostic tool. Biomechanics and Modelling in Mechanobiology, 2(3), 139–155.
- . Anatomically based geometric modeling of the musculo-skeletal system and other organs. Biomechanics and Modelling in Mechanobiology. 2003;4(1):20–38
- . A survey of finite element analysis in orthopedic biomechanics—the first decade. Journal of Biomechanics. 1983;16(6):385–409
- . From structure to process, from organ to cell: recent developments of fe-analysis in orthopaedic biomechanics. Journal of Biomechanical Engineering. 1993;115:520–527
- Kalra, M.K., Maher, M., Saini, S., 2003. Ct radiation exposure: rationale for concern and strategies for dose reduction. In: Proceedings from the 26th Annual Course of the Society of Computed Body Tomography and Magnetic Resonance.
- . Automated three-dimensional finite element modelling of bone: a new method. Journal of Biomedical Engineering. 1990;12(5):389–397
- . Prediction of femoral fracture load using automated finite element modeling. Journal of Biomechanics. 1998;31(2):125–133
- . Prediction of fracture location in the proximal femur using finite element models. Medical Engineering & Physics. 2001;23(9):657–664
- . Specialised ct scan protocols for 3-d pre-operative planning of total hip replacement. Medical Engineering & Physics. 2004;26:237–245
- . Mathematical model of geometry and fibrous structure of the heart. American Journal of Physiology. 1991;260:H1365–H1378
- . Radiation-related cancer risks at low doses among atomic bomb survivors. Radiation Research. 2000;154:178–186
- . Managing X-ray dose in computed tomography. ICRP special task force report. Annals of the ICRP. 2000;30:7–45
- . A computed tomography assessment of femoral and acetabular bone changes after total hip arthroplasty. International Orthopaedics (SICOT). 2002;26:299–302
- Shim, V., Pitto, R., Streicher, R., Hunter, P., Anderson, I., 2004. Generation of finite element models of patients with total hip replacement surgery: mechanical simulations with cell based model. In: Medical Science Congress, 28th November–1st December 2004, Queenstown, New Zealand.
- . Ventricular mechanics in diastole: material parameter sensitivity. Journal of Biomechanics. 2003;36:737–748
- . Automatic generation of accurate subject-specific bone finite element models to be used in clinical studies. Journal of Biomechanics. 2004;37(10):
- . Epicardial surface estimation from coronary angiograms. Computer Vision, Graphics, and Image Processing. 1989;47:111–127
- . Optimal ct scanning plan for long-bone 3-d reconstruction. IEEE Transactions on Medical Imaging. 1998;17:663–666
- 1 An interactive computer program developed by the Bioengineering Institute for Continum Mechanics, Image analysis, Signal processing and System identification. Freely available for academic use.
PII: S0021-9290(05)00535-X
doi:10.1016/j.jbiomech.2005.11.018
© 2006 Elsevier Ltd. All rights reserved.
