Craniofacial Reconstruction Using Rational Cubic Ball Curves

This paper proposes the reconstruction of craniofacial fracture using rational cubic Ball curve. The idea of choosing Ball curve is based on its robustness of computing efficiency over Bezier curve. The main steps are conversion of Digital Imaging and Communications in Medicine (Dicom) images to binary images, boundary extraction and corner point detection, Ball curve fitting with genetic algorithm and final solution conversion to Dicom format. The last section illustrates a real case of craniofacial reconstruction using the proposed method which clearly indicates the applicability of this method. A Graphical User Interface (GUI) has also been developed for practical application.


Introduction
Craniofacial refers to the study and treatment of certain congenital malformations or facial injuries. Craniofacial fractures are common with the principal causes by sports-related injuries, gunshot wounds, and motor vehicle accidents. Frequently encountered craniofacial fractures hold certain distinct patterns. Sometimes the fracture is small and at times it makes a large cavity in the skull. A direct surgical reconstruction is complex, time consuming and expensive, where the surgeon relies on personal experience to reconstruct the structure back to its normal contour.
Ball curve was first introduced by Alan Ball at British Aircraft Corporation (BAC) which was featured in BAC Computer Aided Design (CAD) system [11]. It has a number of similar features as the world renowned Bezier curves such as symmetry, satisfies convex hull property, variation diminishing property, coordinate system independence, invariant under affine transformation and interpolates end points.
The distinct advantages of Ball basis functions over Bezier basis functions are twofold. Firstly, a robust algorithm has been developed to evaluate the Ball curve which suits interactive design environment [12][13][14][15]. Secondly, generalized Ball basis suits much better in degree elevation and reduction which eases data portability and curve approximation in CAD systems [16,17]. We have chosen rational Ball curve for the reconstruction of craniofacial due to these reasons. Other advantages include the flexibility, easy shape tweaks by changing control points and extra free parameters which can be used to satisfy design constraints.
This paper proposes the construction of missing portion of traumatised parts of the craniofacial region using rational cubic Ball curves and the final solution is represented in Dicom format for direct application. The next section discusses on the representation of rational cubic Ball basis functions. It is followed by an explanation on boundary extraction and corner detection, and then on genetic algorithm (GA) which is used as an optimization engine to obtain the best fit for derived corner points. An algorithm on craniofacial reconstruction using rational Ball curve is proposed after the explanation of GA configuration. A real case of craniofacial fracture has been chosen to show the applicability of proposed algorithm. These data are obtained from Universiti Sains Malaysia (USM) Health Campus. Matlab has been used to develop a Graphics User Interface (GUI) for craniofacial reconstruction using the proposed algorithm, which can be used by surgeons who may employ this tool without detailed understanding of its mathematical knowledge.

The Rational Cubic Ball Interpolant
The cubic Ball polynomial basis was first proposed by Ball [11] for CAD systems application.  S 1 ðyÞ ¼ 2yð1 À yÞ 2 ; S 2 ðyÞ ¼ 2y 2 ð1 À yÞ; where X 3 i¼0 S i ðyÞ ¼ 1: These functions were further generalized by Said to arbitrary degree [12]. He further proposed a recursive formula to compute Ball curve which were proven to be more efficient than Bezier curves. A rational cubic Ball curve for craniofacial reconstruction is given as follows: where and 0 ¼ y 0 < y 1 < y 2 < y 3 . . . : < y n ¼ 1; Eq (1) satisfies the following conditions

Boundary Extraction and Corner Detection
The first step on craniofacial reconstruction is to identify the boundary of missing part from a scanned image which is in Dicom format. The boundary of the original image is obtained using mathematical morphology defined by β(A) = A− (AΘB) where A is the set of all black pixel, B is the 3×3 structuring element and β(A) is the boundary of set A. Θ and − represent the erosion and difference operator, respectively. The corner points are used to divide the boundaries into smaller segments. This paper employs Sarfraz et al.'s method [18] to identify the corner points.

Parameterization
Chord length parameterization is used to evaluate the values of θ i related to points f i = D i , where D i represent the data points of segments, It can be observed that θ i is a normalized form and varies from 0 to 1 with

Tangent Vectors
In general, the tangent vectors d i at f i are defined as follows: For open curve, We use both open and closed curves. To reconstruct the complete skull we use closed curves and to reconstruct curves of the missing parts for each slice, we use open curves.

Normalized Mean Squares Error
A normalized mean squares error is used as cost function which is employed to find the errors of free parameters of rational Ball interpolant. The normalized mean squares error is computed by where D i represents the pixel values of the segments and θ is the parameterized chord length.
The goal is to find the values of unknown parameters a i and b i in Eq (1) in order to minimize the normalized mean squares error. GA is used to optimize the value of a i and b i for craniofacial reconstruction.

Genetic Algorithm
Once the boundaries and corners are identified, the next step is to fit rational cubic Ball curves to reconstruct the missing craniofacial parts. For this we have to minimize Eq (4) for best fitting rational Ball curve. Two free parameters a i and b i are involved in Eq (4). We use a continuous genetic algorithm (GA) as an engine to optimize free parameters a i and b i to obtain a suitable fit from given boundaries and corner points.
The GA employed in this study is similar to binary GA where the primary difference is variables are no longer represented by bits of zeros and ones. We start the GA process on E 2 by defining a chromosome as an array of variables to be optimized. Since there are two variables for our interpolant, the chromosome is represented as an array of 1×2 elements i.e chromosome = [a i ,b i ], i = 1, ...n.
The GA is an optimization engine which must be constrained to explore a reasonable region of variable space. Sometimes this is done by imposing a constraint on the interpolant. For example, if an initial search region is unknown, then we define a suitable search space to start the GA process. It is then followed by stating an initial population of chromosome; a normal practise is the initial population is varied from 2k to 4k, where k is the number of chromosomes. In our case, initial population size is set to be 8. The chromosomes are now ready for evaluation and selection of suitable interpolant.
In the next step, we select chromosomes in the initial population which fits enough to survive and possibly reproduce a better offspring in the next generation. For this, the cost and its associated chromosomes are ranked from lowest to highest. Only the top 50% of population with lowest cost will be selected for crossover which are represented as new set of parents. This process of natural selection occurs at each iteration of the algorithm.
The simplest method for crossover is to choose one or more points in the selected chromosomes to produce various solutions. The variables at these points are merely swapped between the two parents. Mutation is the final step where we change some chromosome values randomly within the search space; usually the mutation rate is set to 20%. This process continues iteratively for a given number of generations or until a result obtained is less than a defined value.

Curves to Dicom Format
After fitting the missing part with rational Ball curves of each Dicom image, the next step is to convert the missing part into a Dicom format. For this, we take the corner points of both curves and then join the initial and final points by using the following linear form: where A i and B i are initial and final points of curves, and t varies from 0 to 1. Our next step is to define a black image (all zero values) of the same size as the original Dicom image slices. The fitted rational Ball curve is then converted into white form by equaling this data to 1. Finally, we use Matlab command imfill to fill the area bounded by the curves. To note, ezDicom is a software which is used to read the constructed Dicom image.

Graphical user interface(GUI)
A graphical user interface (GUI) is graphical display which is used to carry out interactive tasks and consists on one or more windows having controls. These are called components. This interface eases out user from creating and executing the commands in order to complete tasks. In contrast to coding programs GUI user does not have to go into details to understand the procedure of execution of commands. This interface include menu, start and stop buttons, boxes etc.
Each control, and the GUI itself, has one or more user written routines (executable MATLAB code) known as callbacks, named for the fact that they call back to MATLAB to ask it to do things. The execution of each callback is triggered by a particular user action such as pressing a screen button, clicking a mouse button, selecting a menu item, typing a string or a numeric value, or passing the cursor over a component. The GUI then responds to these events. In this article GUI is used to control the curve by using free parameters. The input is the end points of missing part depending on the users. The free parameters allow the user to control the curve.

Proposed Algorithm
This section simplifies the algorithm for craniofacial fracture reconstruction. 8. Finally, convert the area enclosed by the curves into Dicom format using Matlab, and ezDicom software is used to visualize the constructed Dicom image 9. GUI is used to construct the missing part with out understanding the mathematics behind  image is the combination of first two images, which is the required Dicom solution. This process is repeated for various CT scan slices as shown in Figs 7-12. The given Dicom data for skull are in sequence and consistent for different frames or contours. The constructed curves of the missing parts using rational cubic Ball interpolant will be consistent between different frames since these curves obey convex hull property, variation diminishing property (VDP), independent of coordinate system and satisfy the end point conditions. By convex hull property all constructed curves will lie within the convex hull of control polygon. Fig 13 represents the graph of normalized mean squares error of slice 176. Fig 14 is the display of Graphical User Interface (GUI).

Conclusion
This paper proposes craniofacial reconstruction using cubic rational Ball interpolant. A patient with head injury as a case study is illustrated to show its applicability. The computed errors and final outputs indicate that the proposed interpolant is suitable for solving craniofacial reconstruction problem. In addition, a user friendly GUI is developed which surgeons may use for practical applications.