·Home ·Table of Contents ·Aeronautics and Aerospace

Automatic mosaicing with error equalisation of non-destructive testing images for aerospace industry inspection Bogdan J. Matuszewski , Lik-Kwan Shark , Martin R. Varley Department of Engineering and Product Design, University of Central Lancashire Preston, PR1 2HE, United Kingdom John P. Smith BAE SYSTEMS Warton, Preston, PR4 1AX, United Kingdom Contact

Abstract

This paper presents a special mosaic construction method developed for merging multiple non-destructive testing (NDT)images. Described in the paper are different geometrical image transformation models, estima tion methods of their unknown parameters, and a new method for equalisation of the transformation errors. The proposed new method is shown to reduce the final alignment error to no bigger, and in most cases smaller than the error introduced by the transformation parameter estimation methods. The performance of the proposed method is validated by simulation data as well as real ultrasonic and x-ray NDT data of aircraft components.

1. Introduction

Due to the relatively large size of aircraft components, NDT systems are usually capable of inspecting a comparatively small area of the component under inspection at a time. It is therefore necessary to perform the image mosaicing operation to merge a large number of sub-images acquired during the inspection process into one complete image of the entire component under inspection. However, as there exists an uncertainty in the component position and orientation with respect to the sensor, the mosaic construction requires special methods to be applied to ensure the correct alignment of all images in the mosaic.
In the development of an automatic mosaic construction method, three key problems have to be addressed:
1. selection of a suitable geometrical image transformation model for image alignment;
2. estimation of transformation parameters; and
3. blending of overlapping sub-images.
While all the three problems are discussed in this paper, the main focus of the paper is on error equalisation in mosaic construction using the rigid body transformation model for the image displacement. The purpose of error equalisation is to overcome mosaic inconsistency, which is one of the main difficulties in the construction of a mosaic with a large number of images, i.e. the local transformation estimated for an overlapping image pair based on their image contents differs from the corresponding transformation computed based on their positions and orientations in the mosaic. By numerical minimisation of a specially designed cost function, mosaic construction with error equalisation is shown, using simulated data as well as real NDT data, to offer an average misalignment error lower than the estimated local transformation.

2 Geometrical image transformation models

To construct a mosaic, the correspondence between different images has to be established first such that matched sub-areas of the images represent the same component point (component area). The class of geometrical transformation needed to align images depends on the component geometry, sensor movement as well as the inspection technique used. If the component point with co-ordinates [X,Y,Z] is represented by I_k(x,y) in the first image and by I_l(x,y) in the second image, the geometrical transformation which moves point (x,y) in the first image to the corresponding point in the second image can be represented by a vector function f :

(1)

In some cases, there may be no analytical function to model this geometrical transformation with sufficient accuracy, or the function is too complex. In these cases, the movement of each image pixel has to be described separately (known as dense motion representation) [1]. This means that the number of parameters to be estimated can be twice the number of pixels in the overlapping area. The number of parameters to be estimated can be reduced by approximation of the function f using two-dimensional splines controlled by a smaller number of displacement parameters defined on a coarse spline control grid [2]. A further and significant reduction in the number of parameters to be estimated can be achieved by using simpler image transformation models. One of the most common models is the two-dimensional projective transformation with eight parameters for function f. Using homogeneous co-ordinates this function is defined as:

(2)

The corresponding Cartesian co-ordinates are [x/w ,y/w ] ^T , and [x/w ,y/w ] ^T . It can be shown that the two-dimensional projective transformation is an accurate model of the geometrical image transformation for arbitrary sensor movements and planar objects [3]. If the sensor motion is restricted (e.g. panning or rotation only), the number of parameters to be estimated can be further reduced. When the component surface is not planar but can be approximately represented as a piecewise planar, the two-dimensional projective transformation model can be used locally [4].

One of the simplest image displacement models is the rigid transformation given by:

(3)

This model has three parameters, namely, the translation vector co-ordinates T=[Tx,Ty]^T and rotation angle a. The meaning of these parameters is illustrated in Figure 1 which shows a rigid displacement between two images.

Fig 1: Rigid transformation model of image displacement (Ti,j represents translation vector between images `j' and `i' and ai,j rotation angle of image j in image i co-ordinate system)

Although simple, the rigid transformation is an accurate image displacement model for many NDT inspection methods, especially when the component to be inspected has a relatively simple geometry. For example, this model can be used for mosaic construction of ultrasonic, shearographic and thermographic images if the component surface is perpendicular (or almost perpendicular) to the optical axis of the camera. The mosaic construction of x-ray images is also possible using this model if the component to be tested is relatively thin and has a simple internal structure.

3 Estimation of image transformation parameters

Having selected the image displacement model, the next step is to estimate its unknown parameters. Many methods are available. One of the possibilities is to detect predefined features [5] (for example corners, points or line segments) in the image pairs and then find the correspondence between these features. Using matched features the unknown parameters of the geometrical image transformation model can be estimated and subsequently used to align the images. However, this feature based parameter estimation method can not be used for the case of the dense motion fields. In addition, it requires presence of well-defined features (geometrical structures) in the image, which can be detected and subsequently used for parameter estimation.
When there are problems with feature detection in the images (or dense motion field has to be found), the methods based on pixel values are used. The most often used methods of this type are based on the sum of squared differences given by:

(4)

where I _j (x _k , y _k), is the displaced k-th pixel of the j-th image in the co-ordinate system of the i-th image. This method can be used to estimate the parameters of the model described by equation (2) which includes the rigid transformation. The estimation procedure usually requires performing numerical minimisation in order to determine the values for the unknown parameters [4]. Variations of this method with extra constraints added and images modelled as random fields have also been used to estimate dense motion fields [6].
When image displacement is modelled by the rigid transformation, methods based on the cross-correlation measure are often used. Many different variants of this measure can be employed [7] depending on the image structure, robustness required, accuracy and computation cost. In the work described in this paper, the cross-correlation coefficient function, which measures the similarity between pixels in different images covering the same component area, was used. This function is defined by:

(5)

where I _j (x _k , y _l) and I _i ( x_k , y_l) denote two images with the respective image mean values _j and _i. Since a higher value of the correlation coefficient function indicates a higher degree of similarity between the two images, it was applied to the overlapping areas of neighbouring images to estimate the transformation parameters. The translation vector and rotation angle are estimated iteratively. Following the initial estimation of the translation vector, the rotation angle is then estimated by moving the second image to the position defined by the translation vector found and transforming the overlapping images' sub-areas to polar co-ordinates (as rotation in Cartesian co-ordinates becomes translation in polar co-ordinates). The whole process is repeated until no significant changes in the estimated parameters are detected. To reduce computation cost, the calculations are performed in a hierarchical manner using image pyramid decomposition [8] with the translation vector and the rotation angle found for a lower resolution level used as the initial estimates for a higher resolution level.
A review of other methods available for the estimation of the image displacement parameters can be found in [9]. The estimation yields a set of the parameter values which define displacement transformations (local transformations) for all possible image pairs with sufficient overlapping. These parameters are denoted by

(6)

where and are respectively the translation vector and rotation angle between images i and j, matrix M defines adjacent image pairs with M _i,j =1 and M _i,j =0 indicating respectively sufficient and insufficient overlapping between images i and j, and m is the number of images to be aligned. With the transformation for all possible image pairs with sufficient overlapping found, the next step in the image mosaic construction is to put all images into a common co-ordinate system. The next section explains the difficulties associated with this task and describes the authors' solution to this problem.

4 Error equalisation

Figure 2 shows a diagram representing an example of a `rigid' mosaic of five images. The necessary geometrical transformations needed to align the five images in the mosaic are translation and rotation operations only.

Fig 2: Diagram of `rigid' mosaic of five images

In the diagram shown in Figure 2, the current position and orientation of the images with respect to a common co- ordinate system in the mosaic are denoted by:

(7)

Using the differences in the positions and orientations between the images, the local transformations required can be computed for the current mosaic and are denoted by:

(8)

Since these computed local transformations are consistent (because they are obtained from the image positions and orientations in the mosaic), it implies the following relations to be true for the mosaic shown in Figure 2:

(9)

(10)

However, the relations (9) and (10) are unlikely to be true for the local transformations previously estimated, because they are obtained based on the similarity existing in an image pair without relating to their positions and orientations in the mosaic, and are also contaminated by some inherent estimation errors. Using the previously estimated local transformations based on image similarity, the images can still be placed into a mosaic as long as no closed loop is formed from them in the process. For example, in Figure (2), a mosaic could be built using the sequence of the estimated local transformations . Errors are accumulated as the local transformations are cascaded to build the mosaic. Consquently, the local transformation estimated using image similarity for images one and five, is usually not equal to the local transformation computed, (T_1,5a_1,5 ) using their positions and orientations of (x ₁ ,y ₁ ,q1 ) and (x₅ ,y₅ ,q₅ ) in the current mosaic. As a matter of the fact, for large mosaics (containing dozens or hundreds of images) the local transformation estimated based on image similarity can differ significantly from those computed based on image positions and orientations. This means that the accumulated errors of the local transformation can cause a serious degradation in the mosaic quality and the estimated local transformations based on image similarity should not be used directly for mosaic construction.

Let an optimum set of the local transformations to achieve a best alignment of all images by distributing the errors evenly across the whole mosaic be denoted by:

(11)

The optimal local transformations should not only satisfy the consistency conditions (similar to conditions (9) and (10)), but also be as close as possible to the local transformations estimated based on image similarity (because each one is a good approximation of the true local transformation). The latter constraint can be expressed as

(12)

The determination of the optimal local transformations is achieved via the minimisation of a cost function constructed based on the sum of the squared differences between the two local transformations with one estimated using the image similarity and the other computed using the image positions and orientations in the mosaic:

(13)

where

(14)

(15)

In the algorithm implementation, the position and orientation of the first image is used as the reference to anchor the whole mosaic (with (x₁,y₁,q₁ )= (0,0,0 )for example). Minimising the cost function of equation (13) yields a set of optimum position and orientation values in p for the images in the mosaic:

(16)

and the optimal local transformations are computed from using

(17)

These optimal local transformations satisfy the consistency conditions as they are computed in equation (17) using the image positions and orientations in the mosaic. Furthermore, the optimal local transformations also satisfy condition (12) as the cost function reaches the lowest possible value when these transformations equal exactly to the local transformations estimated based on image similarity.

To simplify the notation, the following mapping operation is employed to convert the two subscripts used in equation (13) to one subscript:

(18)

where n is the total number of the local transformations in the mosaic. Rewriting equation (13) using one subscript yields

(19)

As the minimisation defined by equation (16) is a non-linear problem, a numerical method based on the Levenberg-Marquardt (LM) algorithm [10] is used. This algorithm requires the gradient and an approximation of the Hessian matrix of the cost function E(p), they are calculated using the partial derivatives of for 1£ i £ 3n ,1£ j £ 3m and are given by the following equations:

(20)

(21)

(22)

All others entries of the Jacobian matrix are zero. The initial positions and orientations of the images in the mosaic can be obtained either from the NDT inspection system or if this is not possible, they can be computed from the estimated transformations between images in the acquired image sequence.

After the correct alignment of all images, mosaic composition based on the Voronoi diagram construction method is performed to combine the overlapping areas of the neighbouring images into one image. In this method, the values of each pixel in the mosaic is taken from one image only with the image having the image centre closest to the location of the pixel concerned.

5 Simulation results

Fig 3: Simulated mosaic

To investigate the performance of the proposed method of image mosaic construction with equalisation of the local transformation errors, tests were performed using a simulated data set shown in Figure 3. To generate the simulated data, the positions and orientations of each image in the mosaic were first selected randomly with the constraint that an image has to overlap at least with the previous image in the sequence. From the positions and orientation of these images, the exact local transformation are then computed using equation (15) and denoted as the true local transformations as

(23)

The simulation of the local transformations estimated based on image similarity is done by adding random errors to the exact local transformations . Figure 4 shows the histogram of the simulated local transformation errors obtained by adding the random errors with the uniform distribution over the intervals of ±1 pixel and ±1 degree to and . The local transformations estimated based on consecutive images, with 1 £ i £ m - 1 were used to compute the initial positions and orientations of the images in the mosaic. Figure 5 shows the histogram of the final errors between the exact local transformation, and the local transformation computed from the initial image positions and orientations,. While the maximum translation error is seen to about 12 pixels, the maximum rotation error is seen to be about 2.3 degrees. Comparing Figure 5 with Figure 4, these errors are significantly higher than the corresponding errors from

Fig 4: Histograms of estimated translation errors from (with errors in horizontal and vertical co-ordinates in blue and red respectively) and rotation error from

Fig 5: Histograms of translational errors from (with errors in horizontal and vertical coordinates in blue and red respectively) and rotation errors from computed using initial image positions

Using the algorithm of error equalisation described in the previous section, the histograms of the errors between the equalised local transformations and the estimated local transformations is shown in Figure 6, whereas the histograms of the errors between the equalised local transformations and the exact local transformations is shown in Figure 7. From these error histograms, improvements on the quality of the mosaic is obvious as the maximum translation error is about 0.6 pixel and the maximum rotation error is about 0.5 degree. Comparing Figure 7 with Figure 4, the errors from and j are seen to be smaller than the local transformation estimation errors from

Fig 6: Histograms of equalised translational errors from (with errors in horizontal and vertical co-ordinates in blue and red respectively)and rotation errors from

Fig 7: Histograms of equalised translation errors from (with errors in horizontal and vertical co-ordinates in blue and red respectively)and rotation errors from

6 Results of real data mosaic construction

The proposed method of the image mosaic construction was used to build mosaics for various NDT images of different aircraft components and acquired using different inspection techniques.

Fig 8: Three ultrasonic sub-images of an aircraft component

Fig 9: Constructed ultrasonic mosaic for images shown in Figure (8).

Fig 10: Mosaic of 132 x-ray images.


Figure 9 shows the mosaic constructed from three ultrasonic sub-images of an aircraft component shown in Figure 8. The mosaic was constructed without any initial constraints imposed on the images' position or orientation. The rigid transformation image displacement model described by equation (3) was used with the parameters estimated using cross-correlation coefficient function using equation (5). The proposed error equalisation method was used to compute consistent local transformation between images.

Figure (10) shows the mosaic of 132 x-ray images of an aircraft component captured by using a real-time radiographic system. In this case, translation only was used to model the image displacement with the parameters estimated using the cross-correlation coefficient function. The initial positions of images were acquired from the real-time radiographic machine. Although the rigid transformation is only an approximation of the exact image displacement model, the developed method was able to produce the satisfactory quality mosaic in this case also.

7. Concluding Remarks

The paper presents a new method for error equalisation in the automatic mosaic construction. The results show that the proposed method achieves an average misalignment error lower than that using the local transformations estimated. Due to the high precision ensured by this method, the NDT mosaic image can be used to examine the structural integrity and geometrical dimensional tolerance of large aircraft component, to register NDT image data with the component CAD model, to fuse multiple NDT data types and to perform data compression [5].

8 Acknowledgement

The work is supported by UK Engineering and Physical Science Research Council (EPSRC) under Innovative Manufacturing Initiative (Grant Ref: GR/L34464).

9 References

1. Stiller C. & Konrad J., Estimating Motion in Image Sequences, IEEE Signal Processing Magazine, July 1999, pp. 70-91.
2. Szeliski R. & Coughlan J., Hierarchical Spline-Based Image Registration, Proceedings of IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR'94), Seattle June 1994, pp. 194-201.
3. Szeliski R., Video Mosaic for Virtual Environment, IEEE Computer Graphics and Applications, March 1996, pp. 22-30.
4. Shark L.-K., Matuszewski B.J., Smith J.P., Varley M.R., Automatic Mosaic and Construction of Three-Dimensional Shearographic Surface Images, Proceedings of the IASTED International Conference on Signal and Image Processing (SIP'99), Nassau, Bahamas October 1999, pp. 183-187.
5. Matuszewski B.J., Shark L.-K., Smith J.P. & Varley M.R., Automatic fusion of multiple non-destructive testing images and CAD models for industrial inspection', IEE 7 th International Conference on Image Processing and Its Applications, Manchester 1999, pp. 661-665.
6. Andre R., Hendriks E., Biemond J., Correspondence Estimation in Image Pairs, IEEE Signal Processing Magazine, May 1999, pp. 29-46.
7. Burt P.J., Yen C., Xu X., Local Correlation Measures For Motion Analysis, A Comparative Study, Proceedings of Conference on Pattern Recognition & Image Processing, 1982, pp. 269-274.
8. Hansen M., Anandan P., Dana K., van der Wal G., Burt P., Real-time Scene Stabilisation and Mosaic Construction, Proceedings of IEEE Workshop an Application of Computer Vision, 1994, pp. 54-62.
9. Brown L.G., A Survey of Image Registration Techniques, ACM Computing Surveys, Vol. 24, No. 4, December 1992, pp.325-376.
10. Press W.H., et al., Numerical Recipes in C: The Art of Scientific Computing, 2 ^nd edition, Cambridge Univ. Press, 1995.

© AIPnD , created by NDT.net

|Home|    |Top|