Material-solid modeling of human body: A heterogeneous B-spline based approach

Computer-Aided Design 41 (2009) 586–597

Contents lists available at ScienceDirect

Computer-Aided Design journal homepage: www.elsevier.com/locate/cad

Material-solid modeling of human body: A heterogeneous B-spline based approach Ravi M. Warkhedkar, Amba D. Bhatt ∗ Mechanical Engineering Department, M.N. National Institute of Technology, Allahabad, UP 211004, India

article

info

Article history: Received 13 July 2007 Accepted 20 October 2008 Keywords: Heterogeneous solid model Reverse engineering Material representation Human model B-spline

abstract Human body is a natural heterogeneous object optimized in its structure and composition by natural processes over years of evolution. In this paper, B-spline solid representation method is extended to represent material composition to develop a heterogeneous model of the human body. Two different approaches, namely surface fairing and surface fit, are used to create a slice by slice model from CT scan data. Both the approaches are compared for different regression parameters. This methodology has the potential to represent internal body details accurately with fewer digital input data,with applications in FEM analysis, freeform-fabrication and tissue engineering. © 2008 Elsevier Ltd. All rights reserved.

1. Introduction Recent advances in information technology and bio-medicine have prompted Computer Aided Design (CAD) to find many novel applications in biomedical engineering. CAD with the help of medical imaging, and Rapid Prototyping (RP) technologies has the capacity to create anatomic models which have diagnostic, therapeutic and rehabilitatory medical applications [1–4]. These models also have non-medical applications in the field of passenger safety design and crash analysis [5–7]. Generally, human body models are prepared by reverse engineering process using non-invasive imaging techniques like CT/MRI. These models are used for creating biomedical implants and tissue scaffolds or as a surgical/diagnostic aid [2,3,8–12]. Although non-invasive medical imaging is able to produce three dimensional descriptions of the internal body details, the capture of internal detail of the body is essentially a voxel based representation. This representation requires huge storage memory and lacks geometric topological relations. Hence they are unsuitable for biomedical analysis and simulations activity which requires a vector based CAD solid modeling environment [2,3]. Effective methods for the conversion of CT data in a CAD solid model is still in a developmental stage. Another issue is regarding the internal material representation of the body part. The available CAD models represent the human body as a solid, with a homogeneous material distribution. Heterogeneous objects, also called as Functionally Graded Materials

∗ Corresponding address: Mechanical Engineering, National Institute of Technology Hamirpur, 177 005 Hamirpur, Himachal Pradesh, India. Tel.: +91 1972 254748. E-mail address: [email protected] (A.D. Bhatt). 0010-4485/$ – see front matter © 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.cad.2008.10.016

(FGM) are objects made of more than one constituent material. The human body is a heterogeneous object which has been optimized in its structure and composition by natural processes over years of evolution. Qian and Dutta [13] have proposed B-spline tensor solid representation for heterogeneous material representation. The method was applied for designing and the analysis of a turbine blade, by using heterogeneous B-spline lofting [14,15]. Material distribution in the human body is quite complex compared to man made materials. Observations on a human-CT scan data reveal that usually there is a smooth variation of material density across the slice (Fig. 1). The B-spline curves and surfaces are known to represent such freeform objects geometry closely. The Bspline based method has excellent model coverage due to the large number of control points. Therefore, in this paper, B-spline tensor product methodology is extended to represent geometry as well as material distribution inside the human body by reverse engineering CT scan data. The resultant model being a analytical geometry-material model, it is expected to facilitate easy manufacture and analysis of the bio-objects. In earlier work [16], an application of heterogeneous B-spline tensor product modeling, namely heterogeneous B-spline fairing and heterogeneous B-spline fit, was developed, and heterogeneous surface models were created for the proximal femur, skull and mandible slices using this approach. The present work is an attempt to extend a B-spline tensor product modeling approach to create heterogeneous solid model of the human body. Salient points of this modeling approach are as follows:

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

587

Fig. 1. Typical femur CT slice and variation of radio intensity along the profile line. There is linear relationship between radio intensity and material density.

Fig. 3. Image cropping (femur).

Fig. 2. Flowchart for the B-spline heterogeneous human body modeling procedure.

• Iso-material contours and iso-material surfaces can easily be extracted from the model to manufacture a complete object by the Solid Freeform Fabrication (SFF) machine. • Material gradient information can easily be extracted from the model. The Finite Element (FE) meshing of this model is also expected to be easier, especially based on the material gradient variation. • Data compression up to 30% of original CT scan data can be achieved by this work, with acceptable error on reconstruction for visualization purpose. As complete inside details of a human body can be captured with this model, a better analysis of the human-body, especially related to crash–impact analysis, can be done. This can be used for safer design of automobiles and other high-speed vehicles. • This model can also be used for designing tissue scaffolds by locally modifying the model, as local editing of B-spline based curves or surfaces is possible.

of the B-spline scheme to represent a heterogeneous object is presented in Section 3.1 and representation of heterogeneous properties with this scheme is given in Section 3.2. The Section 4 describes modeling methodology steps. It includes Section 4.1 on data preprocessing and Section 4.2 on the model development. The results are presented in Section 5. After presenting potential applications of the present modeling approach in Section 6, the paper is concluded in Section 7.

Flowchart for the B-spline heterogeneous human body modeling procedure is shown in Fig. 2. Organization of the remaining paper is as follows. Section 2 is dedicated to literature review. In Section 3, B-spline is introduced as the basic representation scheme in this work. The extension

Heterogeneous model development is carried out in two steps by using the captured data from CT images and then converting it to a heterogeneous model. The literature pertaining to Reverse Engineering (RE) of a human body and heterogeneous modeling is reviewed in brief below.

Fig. 4. Image organization (Femur).

2. Literature review

588

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

Fig. 5. Typical CT slice after thresholding (A) Femur (B) Mandible (C) Skull.

Fig. 6. B-spline faired curve.

2.1. Human body models by Reverse Engineering RE processes use non-invasive medical imaging for preparing human body models. Serially sectioned CT image data of an entire human body has been made available to researchers through various research projects like Visible Human Project (VHP), Chinese Visible Human (CVH) and Visible Korean Human [17]. Human body models created from such data have applications in creating artificial tissues and bio-medical implants [2,3,8–10]. Anatomical models were also reported for post operative implant geometry reconstruction [18], 3D interactive virtual implant manipulation [11] and shape optimization of prosthesis [12]. 2.2. Heterogeneous human body modeling Various attempts were made to define material heterogeneity on geometric models, especially for non-medical applications. Recently Kou and Tan [19] have presented a review of representation methodologies for heterogeneous objects. Depending on the material specification modalities, such schemes can be classified into four categories as voxel based, explicit function based, implicit function based and control point based modeling methods. In the voxel based method, the entire model space is discretized into small elemental volume blocks called voxels [20, 21]. Geometrical coordinates, as well as material composition, are specified at the voxel. Jackson et al. [18] proposed a volumemeshed tetrahedron model where material information inside the tetrahedron is represented as an interpolation function. Although voxel-based representation has good model coverage, it has disadvantages in terms of a large storage memory requirement. The local editing of the model is also cumbersome and inefficient. In an explicit function based method, analytical functions are used to define material heterogeneity [22,23]. Siu and Tan [24] introduced the concept of ‘grading source’ and material composition function to represent material variation. Grading source is defined as the ‘origin of material variation’. Material composition functions are linear or non linear real-valued mathematical functions. They are defined in term of the distance from the ‘grading source’. As the mathematical functions like linear, exponential, parabolic and power functions are used

for material distribution with explicit functions, the complex heterogeneities are difficult to express by such explicit functions. While the explicit function based method employs boundary representation (B-rep) for geometry and explicit functions for material representations, the implicit function based method uses functional representation (F-rep) for both geometry and material representation [25]. The implicit functions (i.e. f (m) = 0) are used for defining material heterogeneity. Adzheiv [26] utilized F-rep to model heterogeneity of geological structures and adaptive finite element mesh generation. Biswas et al. [27] used R-function with distance field, where the inverse distance weighing function with transfinite interpolation were employed with a multiple material feature. Model coverage is poor for this method and, as such, it is difficult to find an implicit function for a domain enclosed by NURBS or spline curves. In the control point based method, geometry-material values are specified at control points and interpolated by shape functions. The shape functions usually used are B-spline [13–15], Bezier [28] and NURBS [29]. Due to the large number of control points, this method has excellent model coverage but a high degree of freedom. On standardization of heterogeneous object representation, Lalit et al. [30,31] suggested an information model to represent heterogeneous object for ISO 10303. While all these methods are available for heterogeneous representation, generally voxel-based methods are used to represent heterogeneous bio-modeling. In the bio-modeling area, most of the attempts are directed toward bone modeling. Attempts are made to model compact bone as a molecular composite [32]. For analysis of bone, the heterogeneity of cortical and cancellous bone is commonly accounted for by considering several sets of elastic moduli for bone material [33]. Muller et al. [34] developed the three dimensional finite element (FE) model of the human tibia bone using the p-version of the FE method. Recently Cheng et al. proposed the heterogeneous solid representation scheme based on material feature [35]. Lian et al. [36] suggested a three dimensional concentric microstructure model with structural extensionality and an optimization feature for artificial bone. Although some researchers presented the bio-object as a composite [37] and others provided theoretical directions for heterogeneous bio-modeling [35], a heterogeneous model has not been reported in practice, to our knowledge. Bio-models created as surface models, use reverse engineering of medical imaging data. Outer body contours are created by processing CT data and triangulated, B-spline or NURBS patches are fitted across the outer shape of the model. Although good for visualization, these models lack internal details. Voxel based representation suffers from problems like aliasing, loss of geometric information after voxelization, huge memory requirements and inconvenient model editing. Implicit and explicit functions, though

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

589

Fig. 7. Pseudo-code for the development of BHSOLID(fair).

easy to manipulate, are unable to address geometrical and material complexity of a bio model. The B-splines are known to represent the freeform objects closely. The control point based model is a direct extension of parametric curves, surfaces and volume where, along with geometrical information, addition material information is stored at control points. Both these pieces of information can easily be controlled by changing control points. Hence the B-spline based control point method is selected for preparing human body models in this work. Also, the method differs from other methods, in that the 3D grid nature of the CT data can directly be used for 3D parameterization for a building model. Hence the method is intuitive, easy and quick. In our previous work, B-spline based heterogeneous, analytical, 2D models of human body parts were created from the CT scan images [16]. These models, created by the fairing and fit method, are able to represent hard and soft tissues. They are suitable for significant data reduction and have potential applications in SSF and FE analysis of bio-objects.

Fig. 8. B-spline fitted curve.

3.1. B-spline representation of heterogeneous hyperpatch Extending the modeling philosophy to represent a point, Qian and Dutta [13] proposed B-spline tensor product representation of a heterogeneous object. A point p = {x, y, z }T in a parametric domain is represented by B-spline as follows.

3. Representation of B-spline heterogeneous solid object Heterogeneous object modeling is essentially geometry-material modeling. A bounded region A3 in Euclidean space E 3 is defined for denoting the geometry of an object as [24,35]; A3 = {P ∈ E 3 |Pg

inside the object|}

(1)

3

where P is a point in E . The Euclidean space E 3 can be extended to E n to represent material composition m at any point p; p = {x, y, z}T by introducing fourth vector m as follows. p = {x, y, z , m} ; T

[u, v, w] ∈ [0, 1]

nm

i=1

mi = 1.

(3)

n X m X L X

Ni,p (u)Nj,q (v)Nk,r (w)Pi,j,k

(4)

i=1 j=1 k=1

where Pi,j.k = {xi,j,k , yi,j,k , zi,j,k , mi,j,k }T are control points for the heterogeneous solid volume; p, q, r are the order of the Bspline basis functions Ni,p , Nj,q , Nk,r , in the direction of u, v, w respectively. B-spline basis function, Ni,k is defined as follows. Ni,s (u) = 1

if t1 ≤ u < ti+1

= 0 otherwise

(2)

where all the elements x, y, z, and m are functions of (u, v, w) The vector m represents the vector of primary materials within the geometric solid. E n is the n dimensional space with the dimensionality of 3 + nm − 1. nm is the number of normalized materials such that

X

p(u, v, w) =

(5)

and Ni,s (u) =

(u − ti )Ni,s−1 (u) (ti+s − u)Ni+1,s−1 (u) + ti+s−1 − ti ti+s − ti+1

(6)

where ti are the knot values satisfying the relation ti ≤ ti+1 . The B-spline heterogeneous solid can be constructed using Eq. (4) from the tridirectional net of n ∗ m ∗ l control points in the u, v and w directions, respectively.

590

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

Fig. 9. BHSOLID created from unsegmented skull data is seen in (A). Isomaterial surfaces are extracted from the same model for (A) soft tissue and (C) hard tissues.

Fig. 10. BHSOLID created from unsegmented proximal femur data shown as (A) longitudinal section (B) isomaterial surface (C) isomaterial surface with section.

Fig. 11. BHSOLID created from unsegmented mandible data shown as (A) BHSOLID (B) isomaterial surface (C) primary sections through BHSOLID.

The parametric derivatives of a B-spline solid are given by formally differentiating Eq. (4) as follows. p (u, v, w) = u

n X m X L X

(u)Nj,q (v)Nk,r (w)Pi,j,k

(7)

Ni,p (u)Njv,q (v)Nk,r (w)Pi,j,k

(8)

Ni,p (u)Nj,q (v)Nkv,r (w)Pi,j,k

(9)

Niu,p

i=1 j=1 k=1

pv (u, v, w) =

n X m X L X i=1 j=1 k=1

pw (u, v, w) =

n X m X L X i=1 j=1 k=1

where pu , pv , and pw are the parametric derivatives of B-spline solid with respect to u, v and w respectively. The present work makes use of Eq. (4) for the specific purpose of defining the geometry-material along the parametric variables u, v and w within a CT scan slices. The Eqs. (7)–(9) are used for obtaining a material gradient in parametric directions u, v and w respectively. The methodology proposed for reverse engineering of the human body is as follows.

3.2. B-spline representation of heterogeneous material property By controlling material composition, we control material properties indirectly, as can be seen from Eq. (4). Many times, it is more convenient to control heterogeneous properties directly than to control material composition. B-spline representation provides us with the convenient means to directly control material properties [13]. For example, suppose ρa , ρb are the densities at a point for two composite materials with volume fraction ma and mb respectively. Then the effective density at a point is given by

ρ = ρa ma + ρb mb .

(10)

In general E = f (Ea , Eb , ma , mb )

(11)

where E is the property of interest, Ea , Eb are the material properties for two materials [13]. For directly controlling the material property, the material vector m in Eq. (2) can be replaced by E. Relationship between density and material composition, as shown in the example is linear, but a more complex relationship may exist for other material properties.

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

591

Fig. 12. (A) X-ray image of typical proximal femur (B) BHSOLID longitudinal section (C) Material variation inside the solid is represented with color code. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 13. CT slice data input can be varied (a) by varying resolution in x and y direction along each slice (b) by varying data points number in z direction (i.e. by varying number of slices) or (c) by combination of (a) and (b).

4. Material and method The primary source of reverse engineering data for the present work is CT scan data. The CT image is a pixel map of linear Xray attenuation coefficients of the tissue. The pixel values are scaled such that the linear X-ray attenuation coefficient (h) of air equals −1024 and that of water equals zero. This scale is called the Hounsfield scale with the unit HU. Apparent density (ρ ) of tissue is the mass of mineralized tissue per total tissue volume. The relationship between apparent density ρ and h is linear and is given by the empirical relation for the human body [38]:

ρ = a + bh.

(12)

The constants a and b depend on the type of the hard tissue. For example, for typical proximal femur, a = 131, b = 1.067 and h is CT intensity in HU [38]. Linear correlation between Xray attenuation coefficient and density is also reported for the soft tissues and the fat present in the body [39–41]. In this work, for preparing a model, all negative CT intensities (h) are converted to positive values by transformation. All h are divided by maximum intensity (hmax ) and normalized intensity hn (hn = h/hmax ) is considered as composition/property indicator at the geometric point. The B-spline interpolated hn can be reverse transformed and equation (12) can be used for calculating the density. A linear relationship is also reported between CT intensity and other material properties, like the modulus of elasticity [38]. The same procedure can be adopted to calculate values of such properties. 4.1. Data pre-processing Mimics [42], the medical based image processing software, is used for preprocessing of CT images as well as for exporting point cloud data. The procedure is carried out in the following steps.

• Image cropping : CT images are imported in Mimics environment. The slices are cropped after observing 2D and 3D views to demarcate the region of interest. This is necessary to limit point cloud data and to avoid unnecessary noise (Fig. 3).

Fig. 14. Iso-material contours extracted from BHSOLID (proximal femur) (A) transverse section (B) longitudinal section.

592

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

Fig. 15. Material gradient for proximal femur (A) for a transverse section in u parametric direction (B) for a transverse section in v parametric direction (C) for a longitudinal section in u parametric direction (D) for a longitudinal section in v parametric direction (epiphyseal plates are more clearly visible).

• Image organization: CT images are organized in the image

• B-spline fairing approach: The B-spline fairing approach is

organization interface. The images are selected one at a time for further processing (Fig. 4). • Image thresholding : Thresholding is used to ensure segmented object contain the pixels with a specified value-range. (Fig. 5). • Image data export: After thresholding, bone tissue pixels are exported as a text file which contains the geometric data, as well as the image intensity data.

explained with B-spline fared curve shown in Fig. 6. In the fairing method, data points act as control points and the curve approximates the control points. The BHSOLID(fair) makes use of the tri-parametric tensor product Eq. (4) for the specific purpose of defining the geometry-material along the parametric variables u, v and w within a group of CT scan slices. The input data points are treated as control points (Pi,j,k ) which have four components, three representing geometrics dimensions x, y and z and fourth specifying material composition m. The 3D grid of such input data points forms a control polyhedron. The BHSOLID(fair) approximates geometrical and material co-ordinates specified at the control points of the control polyhedron. Pseudo code for the creation of BHSOLID(fair) is shown in Fig. 7. • B-spline fit approach: When a curve or surface actually passes through the data point, it is said to fit the data. B-splines can be used for fitting the data point. The method is explained with the help of the B-spline fitted curve shown in Fig. 8. In the fairing method, the B-spline solid is generated from the control polyhedron. Surface fit approach is the reverse problem. From a known set of data in a solid, control polyhedron, is evaluated so as to fit the surface-points. The methodology for surface fitting given by Rogers et al. [43] is extended to heterogeneous solid fitting as follows:

4.2. Model development The preprocessed slice data is arranged in a rectangular grid of n × m points in u and v parametric directions respectively. Stack of l number of such slices gives a hexahedral grid of n × m × l points in u, v and w parametric directions. This hexahedral grid acts as a control polyhedron. Gridded data points in control polyhedron act as control points where geometry is defined by Cartesian components x, y, z and material composition by normalized Hounsfield unit hn . B-spline heterogeneous solid model, or BHSOLID in short, is developed for a proximal femur from the group of 34 CT slices with 2 mm separation. Additional BHSOLIDs are prepared for the skull (46 CT slices with 2 mm separation) and mandible (38 slices with 1 mm separation) to ensure the robustness of the program and consistency of the results. BHSOLIDs are created by two methods, namely B-spline fairing and B-spline fit, and are referred to as BHSOLID(fair) and BHSOLID(fit) in the remainder of the paper. The methods are explained below.

Let D = {x, y, z }T be known data points for a ∗ b ∗ c topologically rectangular set of 3-D data, Ni,p , Ni,q and Nk,r are basis functions

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

593

Table 1 R2 plots for BHSOLID for three data variation methods. Type of the CT data

Data variation along x–y direction of a CT slice (a)

Data variation along z direction of a CT slice (b)

Data variation by combination of (a) and (b) (c)

Unsegmented data (proximal femur)

Segmented data (proximal femur)

Unsegmented data (Skull)

Segmented data (Skull)

Table 2 Regression parameter values for the data variation methods by (a) Data variation along x–y direction of a CT slice (b) Data variation along z direction of a CT slice (c) Data variation by the combination of (a) and (b). Regression parameter

Proximal femur (unsegmented data) (a)

2

Maximum R Minimum range Minimum standard deviation Minimum volumetric error (%)

Proximal femur (segmented data)

(b)

(c)

D(u, v, w) =

Fair

Fit

Fair

Fair

Fit

Fit

Fair

Fit

Fair

Fair

0.979 0.5 0.02 −0.2

1 0 0 0

0.98 0.48 0.02 0.1

1 0 0 0

0.98 0.45 0.01 −0.02

1 0 0 0

0.96 0.82 0.04 0.45

1 0 0 0

0.92 0.8 0.03 0.4

1 0 0.01 0

.96 0.9 0.032 0.5

1 0 0 0

Ni,p (u)Nj,q (v)Nk,r (w)Pi,j,k .

(13)

For each known data point, the above equation provides a linear equation with unknown control polyhedron Pi,j,k . The Eq. (11) can be put in matrix form as: (14)

where Ni,j,k = Ni,p (u).Nj,q (v), Nkr (w). For a ∗ b ∗ c rectangular set of data, [D] is a a ∗ b ∗ c × 4 matrix containing three geometric point dimensions and the remaining dimension representing material composition at that point. [N ] is a ∗ b ∗ c × n ∗ m ∗ l matrix of the products of the B-spline basis functions and [P ] is a n ∗ m ∗ l × 4 matrix of the four-dimensional coordinates of the required control polyhedron points. If the matrix [N ] is square:

[P ] = [N ]−1 [D].

(c)

Fit

i=0 j=0 k=0

[D] = [N ][P ]

(b)

Fair

and Pi,j,k = {xi,j,k , yi,j,k , zi,j,k , mi,j,k }T are control points. Eq. (4) can be written as n X m X L X

(a)

(15)

For non-square matrix [N ]:

[P ] = [[N ]T [N ]]−1 [N ]T [D].

(16)

The program is developed to find out the control polyhedron from input data points which have four components, three representing geometrics dimensions x, y and z and a fourth specifying material composition by applying the B-spline fit methodology. For the control polyhedron thus found, B-spline heterogeneous solid model, BHSOLID(fit) is generated by using the algorithm given in Fig. 7. 5. Results Initially the B-spline heterogeneous surface models (BHSURF) were created by reducing the Eq. (4) to a bi-parametric tensor product equation [16]. Results of these surface models showed that the open knot vector with third order B-spline basis functions, are most suitable for an accurate model representation. Hence the open knot vector with third order B-spline basis functions, is used for the development of BHSOLIDs. BHSOLIDs are created in the form of hyperpatch from the hexahedral control point grid formed

594

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

Table 3 Comparison between fairing and fit method for BHSOLID (proximal femur). Regression parameter

Unsegmented data

Segmented data

R2

Range

Standard deviation

Volumetric error (%)

by gridded CT slice data. The three geometrical dimensions are represented along X , Y and Z as usual, and the fourth material dimension is represented by color in the figure. Both the BHSOLID(fair) and BHSOLID(fit) are tested for the following types of data.

• Unsegmented CT scan data: Point cloud data are read without thresholding the CT image intensities. This gives a complete range of radio intensities, including the soft issues, the hard tissues and the bone marrows • Segmented CT scan data: Point cloud data are read after thresholding the CT image intensities. The model has the capability to represent both the hard and the soft tissues (Fig. 9). We can extract different isomaterial surfaces from the same BHSOLID (Fig. 9) or show them as a sectional view for the visualization (Figs. 10 and 11). When compared with the actual X-ray of a typical femur, it can be seen that the trabecular region and the cortical region in femur are well represented with our model (Fig. 12). BHSOLIDs are compared with the actual CT data for assessing the error in interpolation of the material variation. The models are tested for various regression parameters like residue (e), range, volumetric error (representative of deviation of volume enclosed by interpolated material surface from the volume enclosed by actual material surface) and the coefficient of determination (R2 ). These regression parameters are calculated by increasing CT data in steps. The CT data input was changed for a set of CT slices in three ways (Fig. 13) by

(a) varying data points simultaneously in x and y direction (b) varying data points in z direction (c) varying data points in both (a) and (b). Quality of fit is found to improve with the increase in input data points. R2 values have no dependence on the choice of data variation method, except at the very low input data (Table 1). This holds good for other regression parameters, also. The minimum values of errors for different data variation methods are shown in Table 2. As the choice of data variation method has little bearing on the interpolation, we choose the third method in which data is varied in equal proportions in x, y and z directions, for further discussion. The fairing and the fit models can be compared for their efficacy of material representation by referring to Table 3. The accuracy of the model increases with input CT data for both the fairing and fit method. The fit model was found to give perfect interpolation at 100% CT data. This is obvious as the hyperpatch passes through all geometry-material points of CT slices for a fit method at maximum data. The fairing method, being the data approximation method, gives some error even at the full input data. Except for the range parameter, the difference in error with the fairing method and fit method is marginal, especially at a higher percentage of the CT data. While the segmented data reduces the data storage requirement, unsegmented data has the ability to represent soft and hard tissues, together. The models prepared with these two types of input data are compared for accuracy of material representation (Table 4). Generally, unsegmented input data gives good fit even if the

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

595

Table 4 Comparison between unsegmented data and segmented data for BHSOLID (femur). Regression parameter

Fairing

Fit

R2

Range

Standard deviation

Volumetric error (%)

data is as low as 30%. This is because the variation in the material density is more gradual in unsegmented data than in case of segmented data. Another advantage with the unsegmented data is that both the fairing and fit results are compatible with each other, as the deviations in the regression parameters are small. R-square value is above 90% for most of the range of unsegmented CT data for both fairing and fit methods. Segmented data show higher range in both the fairing and fit model, compared to unsegmented data. However, segmented data used with the fit method can compete with the unsegmented data model if the input data is more than 70%, as the R2 is nearly 90% and volumetric error is less than 1% for this data range. 6. Applications This work has applications in the following areas.

• Solid free-form fabrication (SFF): A SFF machine with capability of handling heterogeneous materials can get the slice by slice material-geometry information from this model. Isomaterial information (Fig. 14) can easily be extracted from this model to manufacture a complete object with the help of an SFF machine.

• Data compression: As indicated, data compression up to 50% of the original CT scan data can be achieved by this work, with acceptable error on reconstruction for visualization purpose. • Analysis: Since this model is based on B-spline hyperpatch tensor product model, it is easier to get material gradient information (Fig. 15) across the slice. Material gradients can be co-related with the material property-gradients. Material gradient information may be used for visualizing structural features more clearly, which otherwise are not distinct in normal model (Fig. 15. (D)). The FE meshing of this model is also expected to be easier, especially the meshing based on material gradient variation, as material gradient information is easily extracted from this model. As complete inside details of a human body can be captured with this model, a better mechanical stress analysis of the human-body, especially related to crash–impact analysis can be done. This can be used for safer design of automobiles and other high-speed vehicles. This model can also be used for evaluating the customized position of an astronaut during take off of a space craft. This model can also be used to study the effect of implants/ prosthesis on the healthy body tissues.

596

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597

• Tissue engineering: This model can also be used for designing tissue scaffolds, by locally modifying the model, as local editing of B-spline based curves or surfaces is possible. 7. Conclusions In this paper, we presented novel application of the B-spline hyperpatch tensor product approach to represent the human body. The method is analytical, intuitive, has a wide model coverage, and has the capability of exact material representation. The model creation is also easy and has several practical applications. Two different types of modeling approaches, namely B-spline heterogeneous fairing and B-spline heterogeneous fit, used in this work, are capable of representing both the hard and the soft tissues. In general, the B-spline fit model is found to represent human body parts more accurately compared to the B-spline fairing model. However, it is computationally more extensive. A model based on unsegmented data input is suitable for data reduction up to 50% of the input data. These two approaches, when used with different proportions of CT data, give user modeling flexibility. The ultimate choice may depend on the model accuracy requirement and allowable data size for a particular bio-engineering approach. For real time medical applications and for non-medical applications like human body testing for crash analysis, complete body part data, including hard tissues and soft tissues, are needed. The unsegmented data model can be useful for such application for the purpose of analysis. Another advantage of the method is that one can directly extract iso-material information and material gradient information from the model. Iso-material information can be handy in the process planning of SSF manufacture of body implants and tissue scaffolds. As the material features are included within the model, the present approach is expected to evolve to a better analysis model. Currently, the FE analysis research based on this work is under progress. References [1] Magnenat-Thalmann N, Cordier F. Construction of a human topological model. IEEE transaction on Information Technology in Biomedicine 2000;4(2): 137–43. [2] Sun Wie, Lal Pallavi. Recent development on computer aided tissue engineering — a review. Computer Methods and Programs in Biomedicine 2002;87:85–103. [3] Sun W, Starley B, Nam J, Darling A. Bio-CAD modeling and its application in computer-aided tissue engineering. Computer-aided Design 2005;37: 1097–114. [4] Taun Ho, Saey Hutchmacher, Deitmar W. Application of micro CT and computation modeling in bone tissue engineering. Computer-Aided Design 2005;37:1151–61. [5] Johnson EAC, Young PG. On the use of patient-specific rapid-prototyped model to simulate the response of the human head to impact and comparison with analytical and finite element models. Journal of Biomechanics 2005;38:39–45. [6] Roberts JC, Merkle AC, Biermann PJ, Ward EE, Carkhuff BG, Cain RP, O’Connor JV. Computational and experimental models of the human torso for nonpenetrating ballistic impact. [7] Aare Magnus, Kleiven Svein. Evaluation of head response to ballistic helmet impacts using the finite element method. International Journal of Impact Engineering 2007;34:596–608. [8] Jaecques SVN, Oosterwyck H, Van Mararu L, Cleynenbreugel T, Van Smet E, De Wevers M, et al. Individualised, micro CT-based finite element modeling as a tool for biomechanical analysis related to tissue engineering of bone. Biomaterials 2004;25:1683–96. [9] Sarti Alessandro, Gori Roberto, Lamberti Claudio. A physically based model to simulate maxillo-facial surgery from 3D CT images. Future Generation Computer Systems 1999;15:217–21. [10] Viceconti Marco, Zannoni Cinzia, Pierotti Luisa, Casali Massimiliano. Spatial positioning of an hip stem solid model within the CT data set of the host bone. Computer Methods and Programs in Biomedicine 1999;58:219–26. [11] Seipel Stefan, Wagner Ina-Veronika, Koch Sabine, Schneide Werner. Oral implant treatment planning in a virtual reality environment. Computer Methods and Programs in Biomedicine 1998;57:95–103. [12] Pawlikowski M, Skalski K, Haraburda M. Process of Hip joint prosthesis design including bone remodeling phenomenon. Computers and Structures 2003;81: 887–93.

[13] Qian X, Dutta D. Physics-based modeling for heterogeneous objects. Journal of Mechanical Engineering 2003;125:416–27. [14] Qian X, Dutta D. Design of heterogeneous turbine blade. Computer-Aided Design 2003;35:329–39. [15] Yang Pinghai, Qian Xiaoping. A B-spline-based approach to heterogeneous objects design and analysis. Computer-Aided Design 2007;39:95–111. [16] Bhatt Amba D, Warkhedkar Ravi M. Reverse Engineering of human body: A Bspline based heterogeneous modeling approach. Computer-Aided Design and Applications 2008;5(1–4):194–208. [17] Park Jin, Seo Chung, Min Suk, Hwang Sung, Bae Lee, Sook Yong, et al. Visible Korean human: Improved serially sectioned images of the entire body. IEEE Transactions on Medical Imaging 2005;24(1):352–60. [18] Jackson TR, Patrikalakis NM, Sachs EM, Cima MJ. Modeling and designing components with locally controlled composition. Materials and Design 1999; 20:63–75. [19] Kou XY, Tan ST. Heterogeneous object modeling: A review. Computer-Aided Design 2007;39:284–301. [20] Chandru Vijay, Manohar Swami, Prakash , Edmond C. Voxel –based Modeling for Layered Manufacturing. IEEE Computer Graphics and Applications 1995; 42–7. [21] Cho JR, Ha DY. Optimal tailoring of 2D volume fraction distributions for heat resisting functionally graded materials using FDM. Computer Methods in Applied Mechanics and Engineering 2002;3195–211. [22] Shin Ki-Hoon, Dutta Debasish. Constructive Representation of Heterogeneous Objects. Journal of Computing and Information Science in Engineering 2001; 1:205–17. [23] Kumar Vinod, Dutta Debashish. An approach to modeling multi-material object. ACM conference Solid Modeling 2007;336–45. [24] Siu YK, Tan ST. ‘Source-based’ heterogeneous solid modeling. Computer-Aided Design 2002;34:41–5. [25] Rvachev VL, Sheiko TI, Shapiro V, Tsukanov I. Transfinite interpolation over implicitly defined sets. Computer Aided Geometric Design 2001;18:195–220. [26] Adzhiev Valery, Kartasheva Elena, Kunii Tosiyasu, Pasko Alexander, Schmitt Benjamin. Hybrid cellular-functional modeling of heterogeneous objects. Journal of Computing and Information Sciences in Engineering 2002; 2:312–21. [27] Biswas Arpan, Shapiro Vadim, Tsukanov Igor. Heterogeneous material modeling with distance field distance fields. Computer Aided Geometric Design 2004;21:215–42. [28] Huang Jinhua, Fadel M George, Blouin Vincent Y, Grujicic Mica. Bi-objective optimization design of functionally gradient materials. Materials and Design 2002;23:657–66. [29] Martin William, Cohen Elaine. Representation and extraction of volumetric attributes using trivariate splines: A mathematical framework. In: Proceedings of the sixth ACM symposium on Solid modeling and applications. 2001. p. 234–40. [30] Patil Lalit, Dutta Debasish, Bhatt AD, Jurrens K, Lyons K, Pratt MJ, et al. A proposed standards-based approach for representing heterogeneous objects for layered manufacturing. Rapid Prototyping Journal 2002;8(3):134–46. [31] Pratt MJ, Bhatt AD, Dutta D, Lyons KW, Patil L, Sriram RD. Progress towards an international standard for data transfer in rapid prototyping and layered manufacturing. Computer-Aided Design 2002;34:1111–21. [32] Zhang Dajun. Modeling compact bone as molecular composite. Bioengineering conference ASME 2001;50:1–2. [33] Chang Chih-Han. Effect of material inhomogeneous for femoral finite element analysis. Journal of Medical and Biological Engineering 2002;22(3):121–8. [34] Muller-Karger CM, Rank E, Cerrolaza M. P-version of the finite –element method for highly heterogeneous simulation of human bone. Finite Elements in Analysis and Design 2004;40:757–70. [35] Cheng Jie, Lin Feng. Approach of heterogeneous bio-modelling based on material features. Computer-Aided Design 2005;37:1115–26. [36] Lian Qin, Li Di-Chen, Tang Yi-Ping, Zhang Youg-Rui. Computer modeling approach for a novel internal architecture of artificial bone. Computer-Aided Design 2006;38:507–14. [37] Milkowski LM, Gervasi VR, Kumaresan S, Crockett RS. Development of a mechanically similar composite bone replica. In: Proceedings of the first joint BMES/EMBS conf Atlanta. 1999. p. 495. [38] Rho JY, Hobatho MC, Ashman RB. Relations of mechanical properties to density and CT numbers in bone. Medical Engineering and Physics 1995;17:347–55. [39] Lukaski Henry. Sarcopenia: Assessment of muscle mass. The Journal of Nutrition 1997;127:994–7. [40] Zhang Futang, Peck ChristopherC, Hannam AlanG. Mass properties of the human mandible. Journal of Biomechanics 2002;35:975–8. [41] Teo JeremyCM, Si-Hoe Kuan, Ming Keh, Justin EL, Teoh SweeHin. Relationship between CT intensity, micro-architecture and mechanical properties of porcine vertebral cancellous bone. Clinical Biomechanics 2006;21:235–44. [42] Mimics 9.1, reference guide Materialize software Belgium. [43] Roger DF, Adams JA. Mathematical elements of computer graphics. II ed New Delhi: Tata McGraw-Hill Publication; 2002.

R.M. Warkhedkar, A.D. Bhatt / Computer-Aided Design 41 (2009) 586–597 Ravi M. Warkhedkar has received an M.E. in Mechanical Engineering from Amravati University, India in 1997 and a Ph.D. in Mechanical Engineering from MNNIT Allahabad (India) in 2008. He is working as a lecturer in Government Engineering College Jalgaon, India. His interest includes CAD and Bio-mechanical Engineering.

597

Amba D. Bhatt is a Professor in Mechanical Engineering at the National Institute of Technology, Hamirpur, HP, India. He received a B.E. in Mechanical engineering and M.E. in Design from M.N. National Institute of Technology, Allahabad in the year 1979 and 1992, respectively, and a Ph.D. from I.I.T. Kanpur, India in the year 2000. He was also a guest researcher at the Manufacturing engineering Laboratory at the NIST, Gaitherburg MD. He was a Design engineer at HINDALCO, Mirzapur India and Assistant Manager at Steel Authority of India Limited (1980–1988). His research interests include Product design, Rapid Prototyping, Heterogeneous modeling and Computer Graphics.