# MULTIFOCAL IMAGES FUSION

Ondrej KOVÁČ*, Iveta GLADIŠOVÁ**

*Department of Technologies in Electronics, Faculty of Electrical Engineering and Informatics, Technical University of Košice, Letná 9, 042 00 Košice, Slovak Republic, Tel.: +421 55 602 3013, E-mail: ondrej.kovac @tuke.sk

**Department of Electronics and Multimedia Telecommunications, Faculty of Electrical Engineering and Informatics, Technical University of Košice, Letná 9, 042 00 Košice, Slovak Republic, Tel.: +421 55 602 2940, E-mail: iveta.gladisova @tuke.sk

## ABSTRACT

The paper deals with the issue of multifocal image fusion. Two techniques of the hierarchical image representation are analyzed herein. The first one is hierarchical decomposition based on the Laplacian pyramid and the second one is based on integer discrete wavelet transformation. After the analysis of the above-mentioned decompositions techniques the algorithm of multifocal image fusion is described. Attention is given to the detection of the focused parts of approximated images. From the detected parts, the approximation of an image is obtained. The analyzed algorithm is experimentally verified and the accuracy with regard to the decomposition method and level of decomposition is evaluated. The fusion is done with real microscopic images.

**Keywords:** DWT, Laplacian pyramid, image fusion, multifocus

## 1. INTRODUCTION

Capturing focused images of a scene which contains multiple objects located at different distances from the camera is a common issue in the field of image processing. This task cannot be done by one focal length. This issue is even more significant in microscopic imaging. Therefore, it is necessary to capture multiple images by using different focal lengths of the camera. These images are called multifocal images [1]. The multifocal images, focused on different objects of the scene, have to be merged into one in the whole area of the focused image. This image fusion can be done by several techniques. In general, they can be divided into spatial based and transform based techniques. The spatial based techniques include pixel level fusion, feature level fusion and decision level fusion [2, 3]. Most widely used transform techniques [4] are Discrete Wavelet Transformation (DWT), Gabor transformation as well as Contourlet transformation, which is able to capture intrinsic geometric structures [5]. This paper deals with two techniques of multifocal image fusion. They are described in the chapter two and three. The first one is based on decomposition of images into the Laplacian Pyramid (LP) [6, 7] and the second one is based on integer DWT [8]. Both mentioned methods of decomposition are from the group of transforms of hierarchical representation [9]. This image representation allows to divide an image into several sublevels which have different resolution. The algorithm of fusion itself is described in chapter four. This algorithm is identical for decomposition of both images. The results of the experiment are discussed in chapter five. The experiment was performed on the images which were known but in the defined regions were purposely blurred. The goal of this experiment was to compare the effectiveness of the described methods. Also experiments on the real multifocal microscopic images of the printed circuit board (PCB) were done.

## 2. DISCRETE WAVELET TRANSFORM

Nowadays, DWT is implemented in many areas of signal and image processing. Mostly it is used for its compress and decorrelation properties. DWT is also used for the analysis, denoising and enhancement of images as well as in digital photomontages, contrast, white balancing etc. The presented technique of multifocal image fusion is one of the image enhancement techniques. The benefit of DWT is the fact that DWT allows representation of images in a hierarchical structure of the scaled subimages. The principle of DWT consists of image decomposition into the approximation and wavelets (details) subimages. These subimages are created by splitting an image into two branches where in the first one the image is filtered by a low-pass filter and in the second one by a high-pass filter. After that the filtered images are decimated by factor two. The need of convolution filtration is the main drawback of the classical DWT implementation [10]. The results of DWT are decomposition coefficients which are in the domain of real numbers which is another drawback.

### 2.1. Common notes Lifting implementation of DWT

The above-mentioned drawbacks can be eliminated using Integer Discrete Wavelet Transform (IDWT). The difference between IDWT and DWT is that IDWT does not use convolution between the input sequence and pulse response of the filter. IDWT is also known as a lifting implementation (L) DWT; it uses only operations of prediction, correction, summation and rounding. One decomposition stage of LDWT for biorthogonal bank of filters BF (5, 3) [11] is shown in Fig. 1.

The decomposition coefficients are obtained from the input sequence as follows. The input sequence of samples $c_{j+1}(n)$ is divided into odd $a_{j+1}(k)$ and even $b_{j+1}(k)$ samples into two branches. On the basis of the odd sequence, samples of the even sequence are predicted by the predictor P. These prediction samples are subtracted from the even sequence. Thus, the sequence of detail
coefficients $d_j(k)$ is created. From $d_j(k)$ the correction sequence is obtained by the corrector C. These correction samples are summed up with the odd sequence $a_{j+1}(k)$ and the approximation sequence $c_j(k)$ is obtained.

<figure><img src="image_1.png" /><figcaption>Fig. 1 One decomposition stage of LDWT</figcaption></figure>

From Fig. 1 it is clear that values of the correction and prediction are rounded. By rounding R, the integer output of the LDWT is achieved. It is obvious that the recovered sequence $c_{j+1}(n)$ is not influenced by the rounding on the side of the synthesis. The sequences $c_j(k)$ and $d_j(k)$ are given by the following equations (1, 2).

$$d_j(k) = b_{j+1}(k) - \left\lfloor \frac{1}{2} \left( a_{j+1}(k) + a_{j+1}(k+1) \right) \right\rfloor \quad (1)$$

$$c_j(k) = a_{j+1}(k) - \left\lfloor \frac{1}{4} \left( d_j(k-1) + d_j(k) \right) \right\rfloor \quad (2)$$

The backward transform is simply performed using reordering operations of the prediction, correction and summation. By using properties of the separated transformation core, a two dimensional (2D) transform is achieved by 1D LDWT which is at first applied on rows of the input image [12]. This leads to creation of two subimages. In the next step, these subimages are decomposed by 1D LDWT applied on columns and the final decomposition consisting of four subimages is created. Block diagram of 2D LDWT is shown in Fig. 2.

<figure><img src="image_2.png" /><figcaption>Fig. 2 The block diagram of 2D LDWT</figcaption></figure>

### 3. LAPLACIAN PYRAMID

A pyramid representation of images is also based on the hierarchical image representation. Fig. 3 shows the block diagram of the hierarchical image representation by pyramid analysis and synthesis. As can be seen from Fig. 3, the pyramid analysis and synthesis is based on the decomposition of an input image into the Gaussian Pyramid (GP) which represents low-pass (approximation) part. The Gaussian images in Fig. 3 are presented as $G_0, G_1, G_2, \dots, G_n$, where n is the number of decomposition

levels. In the hierarchical image representation, the image is given in several levels with different resolution. In this case, the subimage $G_0$ represents the level with the highest and $G_n$ represents the level with the lowest resolution.

<figure><img src="image_3.png" /><figcaption>Fig. 3 The block diagram of four level pyramid analysis and synthesis</figcaption></figure>

Fig. 3 shows that the resolution of the input image is reduced by quadratic value of 2. The resolution reduction and filtration of the high frequencies are performed by decimation filters (DF) with the decimation factor two. This could be written as (3).

$$G_i = Df_2(G_{i-1}) \quad (3)$$

where $Df_2(G_{i-1})$ represents operation of the decimation by factor 2 and filtration. In this hierarchical representation, LP represents a high-pass (detailed) part of the input image. The Laplacian images are presented as $L_0, L_1, \dots, L_n$. The i-level of the Laplacian pyramid $L_i$ is given by subduction of the Gaussian image $G_i$ and interpolated Gaussian image of the next level $G_{i+1}$ that can be written as (4).

$$L_i = G_i - If_2(G_{i+1}) \quad (4)$$

where $If_2(G_{i+1})$ represents operation of the interpolation by factor 2 and the filtration. The interpolation (IF) provides equal resolution of the subimages. The subimages at the highest level of the analysis for GP and LP are identical. Hence, $If_2(G_{i+1})$. The reconstruction or synthesis is performed by creation of the Gaussian pyramid from the Laplacian subimages. The image with the highest resolution is the final reconstructed image. In general, it can be written as (5).

$$\tilde{G}_i = L_i + If_2(\tilde{G}_{i+1}) \quad (5)$$

It has to be mentioned that in the process of interpolation the raster is extended by inserting zero rows and columns. It causes quadruple energy reduction with regard to the original Gaussian image. For this reason, the interpolated and filtrated image has to be multiplied by the constant 4.
## 4. IMAGE FUSION

Thanks to the hierarchical representation of the images regardless of the decomposition method, images can be divided into the approximation and detail subimages. The approximation part represents low-pass decomposition coefficients and details, in images edges, presented by high-pass decomposition coefficients. In the process of multifocal image fusion attention is focused on the detail decomposition coefficients. In the next subchapter, the algorithm of the multifocal images fusion in LDWT and LP domain is analysed.

### 4.1. The algorithm of synthesis

The multifocal image fusion is based on the synthesis of a focused image from IDWT detail coefficients or coefficients of the Laplacian pyramid. These coefficients may be processed by the algorithm given by (6). At first, the input subimages are decomposed into IDWT or Laplacian subimages. The subimages which represent details are joined as follows:

$$d_{i,j} = \begin{cases} d_{i,j}^A, & \text{if } d_{i,j}^A > d_{i,j}^B \\ d_{i,j}^B, & \text{if } d_{i,j}^A \le d_{i,j}^B \end{cases} \qquad (6)$$

where $d_{i,j}$ represents a decomposition coefficient of the given decomposition level of a joined image, $d_{i,j}^A$ and $d_{i,j}^B$ are decomposition coefficients of the input multifocal images A and B. The detail coefficients for IDWT are from subimages LH, HL, and HH and detail coefficients for LP are subimages $L_0, L_1, ..., L_{n-1}$, where n is the number of levels. The question arises as to the approximation subimage of which input image is to be used in the process of synthesis. In general, the approximation image of arbitrary input images can be used. But in the case of images which are significantly blurred, the low-pass part could disturb the final focused image. This may be partly eliminated by averaging of the approximation coefficients of both input images. The averaging is given as (7).

$$a_{i,j} = \frac{a_{i,j}^A + a_{i,j}^B}{2} \qquad (7)$$

But the best fusion result is achieved only if the approximation part is also from the focused areas of the input images. This could be achieved by segmentation of an image into smaller blocs of which the ones are chosen that are focused to subsequently construct the approximation subimage. Whether the image A in a particular block is more focused than image B can be determined by the edge detector. For this purpose, the simplest detector is sufficient – the first order difference. The focused image has more edges than the blurred image. Hence, the energy of the differentiated version of the focused image must be higher than the energy of the blurred one. By simple comparison of the blocks energy of the differentiated images, the approximation image used for the image fusion can be obtained. The input microscopic multifocal images of PCB are shown in Fig. 4a), b). The images show Integrated Circuit (IC) with high package and SMD components. The image where the

microscope is focused on the surface of SMDs and PCB is shown in Fig. 4a), while in Fig. 4b), on the contrary, the microscope is focused on the surface of IC. Fig. 4c), d) shows the results of edge detection for the first level of GP. Fig. 5 represents a block diagram of the complete algorithm of the multifocal image fusion. Fig. 5 shows that input images (A and B) are decomposed into subimage of IDWT or LP. Subimages of details are directly fused by application of (6). Thus, the combined subimages are obtained. The approximation part is processed separately. Approximation subimages are segmented into orthogonal blocks and subsequently edge detection is performed. Based on the comparison of energy of the detected edges, the final approximation image is assembled. The final focused image is obtained by synthesis of the combined details and approximation subimages.

<figure><img src="image_4.png" /><figcaption>Fig. 4 Image of PCB focused on the a) SMD components, b) IC's surface, and c), d) edge detections</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 5 The block diagram of described algorithm</figcaption></figure>
## 5. EXPERIMENTAL RESULTS

The three images shown in Fig. 6 were modified for the purpose of evaluating the effectiveness and accuracy of the method. The modifications were done by application of a low-pass filter in the selected areas of the images. These modified images are shown in Fig. 7. Subsequently the modified images were fused into one image by the above-mentioned algorithm. The fused images were compared with the original images. The difference between them was expressed by the signal to noise ratio (SNR). The fusions were performed using IDWT and LP. The influence of the decomposition level on SNR was also examined. The achieved results are listed in Table 1.

<figure><img src="image_6.png" /><figcaption>Fig. 6 The experimental images a) Apples, b) Baboons, c) F1</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 7 Locally blurred experimental images: a), b) Apples, c), d) Baboons, e), f) Formula 1</figcaption></figure>

According to Table 1, for all of the experimental images better results were achieved using LP decomposition, and for LP one level decomposition was sufficient. In the case of DWT, the best results for all the tested images were achieved by decomposition into two levels. With regard to the informative value of Table 1,

<table><caption>Table 1 The comparison of SNRs for IDWT and LP method of decomposition with regard to the level of decomposition</caption><thead><tr><th>Image</th><th>Method</th><th>Level</th><th>SNR (dB)</th></tr></thead><tbody><tr><td rowspan="6">Apples</td><td rowspan="3">IDWT</td><td>1</td><td>35,938</td></tr><tr><td>2</td><td><b>36,390</b></td></tr><tr><td>3</td><td>36,228</td></tr><tr><td rowspan="3">Laplacian</td><td>1</td><td>39,415</td></tr><tr><td>2</td><td>42,076</td></tr><tr><td>3</td><td><b>42,555</b></td></tr><tr><td rowspan="6">Baboons</td><td rowspan="3">IDWT</td><td>1</td><td>28,929</td></tr><tr><td>2</td><td><b>33,387</b></td></tr><tr><td>3</td><td>32,121</td></tr><tr><td rowspan="3">Laplacian</td><td>1</td><td><b>42,334</b></td></tr><tr><td>2</td><td>40,706</td></tr><tr><td>3</td><td>39,423</td></tr><tr><td rowspan="6">Formula 1</td><td rowspan="3">IDWT</td><td>1</td><td>17,635</td></tr><tr><td>2</td><td><b>24,438</b></td></tr><tr><td>3</td><td>22,809</td></tr><tr><td rowspan="3">Laplacian</td><td>1</td><td><b>37,468</b></td></tr><tr><td>2</td><td>33,676</td></tr><tr><td>3</td><td>32,723</td></tr></tbody></table>

the fused images are not given. The next experiment was based on fusion of real multifocal images obtained by microscopic imaging of PCB. The microscope was at first

<figure><img src="image_8.png" /><figcaption>Fig. 8 The input multifocal microscopic images of DPSs</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 9 Fused images for a), b) IDWT and c), d) LP method</figcaption></figure>
focused on the elements with small packages and then on the SMD elements and the surface of PCB. The multifocal images are shown in Fig. 8a) - d). Images focused by image fusion for IDWT (level 2) are given in Fig. 9a), b) and for LP (level 1) in Fig. 9c), d). Because the experimental images were taken by camera with different focus, they weren't created from the known image as were images in Fig. 6, only a subjective quality level can be evaluated. The SNR cannot be computed.

## 6. CONCLUSIONS

The techniques of image decomposition by DWT and LP were analyzed in the paper. DWT was performed using its integer lifting implementation. Subsequently, the issue of the multifocal image fusion using the above-mentioned transform was analyzed. The well-known algorithm of multifocal image fusion was enhanced by automatic detection of approximation parts of the images which contain more details. From the detected focused parts of these subimages, the final approximation subimages were created. The fusion was applied to images which were blurred in the selected parts. The results were evaluated by SNR. From the achieved SNRs it was obvious that the multifocal fusion based on LP achieved better results than fusion based on IDWT. Additionally, sufficient results, in the case of LP, were achieved already at the first level of decomposition. On the other hand, implementation with IDWT needed two levels of decomposition. Besides fusion of blurred images, the multifocal images acquired by a microscope were fused as well. The subjective quality of the fused images was better in the case of LP than IDWT. This observation was in accordance with the results of the quantitative quality evaluation.

## REFERENCES

[1] ZHANG, Z.: Camera Calibration, Chapter 2, in G. Medioni and S. B. Kang, eds., Emerging Topics in Computer Vision, Prentice Hall Professional Technical Reference, 2004, pp. 4–44.

[2] KHARE, A. – SRIVASTAVA, R. – SINGH, S.: Edge Preserving Image Fusion Based on Contourlet Transform, in Image and Signal Processing, Springer Berlin Heidelberg, 2012, pp. 93–102.

[3] SIMONE, G. – FARINA, A. – MORABITO, F. C. – SERPICO, S. B. – BRUZZONE, L.: Image fusion techniques for remote sensing applications, Information Fusion, Vol. 3, No. 1, 2002, pp. 3–15.

[4] EHLERS, M. – KLONUS, M. – ASTRAND, P. J. – ROSSO, P.: Multi-sensor image fusion for pansharpening in remote sensing, International Journal of Image and Data Fusion, Vol. 1, No. 1, 2010, pp. 25–45.

[5] DO, M. N. – VETTERLI, M.: The contourlet transform: an efficient directional multiresolution image representation, IEEE Transactions on Image Processing, Vol. 14, No. 12, 2005, pp. 2091–2106.

[6] JAYARAMAN, R.: Digital Image Processing, Tata McGraw-Hill Education, 2011.

[7] PARIS, S. – HASINOFF, S. W. – KAUTZ, J.: Local Laplacian filters: edge-aware image processing with a Laplacian pyramid, Communications of the ACM, Vol. 58, No. 3, 2015, pp. 81–91.

[8] SRUTHY, S. – PARAMESWARAN, L. – SASI, A. P.: Image Fusion Technique using DT-CWT, Automation, Computing, Communication, Control and Compressed Sensing, 2013 International Multi-Conference on (iMac4s), March 2013, pp.160–164.

[9] PARVATIKAR, M. V. – PHADKE, G. S.: Comparative study of different image fusion techniques, International Journal of Scientific Engineering and Technology, Vol. 3, No. 4, 2014, pp. 375–379.

[10] PUTTARAJU, M – ASWATHA, A. R.: FPGA Implementation of 5/3 Integer DWT for Image Compression, International Journal of Advanced Computer Science and Applications, Vol. 3, No. 10, 2012.

[11] MIHALÍK, J. – ZAVACKÝ, J.: Filter Banks (Banky filtrov), Košice: LČSOV FEI TU Košice, 2012.

[12] KHEMIRI, R. – SAYADI, F. E. – ATRI, M.: MatLab-GPU-based 2D-DWT Acceleration for JPEG2000 with Single and Double-Precision", Indian Journal of Science and Technology, Vol. 9, No. 1, 2016.

Received January 12, 2017, accepted July 3, 2017

## BIOGRAPHIES

**Ondrej Kováč** was born on 1986. In 2011 he received the Ing. degree in Multimedia telecommunications from the Faculty of Electrical Engineering, Technical University of Košice. His PhD. thesis topic was focused to texture generating, 3D modeling and coding of human head. Since June 2015, he has been working as assistant professor at the Department of Technologies in Electronics, Technical University of Košice.

**Iveta Gladišová** graduated from the Technical University of Košice in 1984. She is employed at same university of Košice, where received her PhD. degree in Radioelectronics in 1997. She is Assistant Professor at the Department of Electronics and Multimedia Telecommunications, Faculty of Electrical Engineering and Informatics, Technical University of Košice. Her research interests include entropy coding, segmentation, vector quantization and digital image processing.