This article introduces a new adaptive method for image interpolation. In order to obtain a high resolution (HR) image from its low resolution (LR) counterpart (original image), an interpolator function (array) is used, and the main focus of this manuscript is to formulate and define this function. By applying this interpolator function to each row and column of a LR image, it is possible to construct its HR counterpart. One of the main challenges of image interpolation algorithms is to maintain the edge structures while developing an HR image from the LR replica. The proposed approach overcomes this challenge and exhibits remarkable results at the image edges. The peak signal to noise ratio and structural similarity criteria by using this innovative technique are notably better than those achieved by alternative schemes. Also, in terms of implementation speed, this method displays a clear advantage and outperforms the high performance algorithms in the ability to decrease the artifact results of image enlargement such as blurring and zigzagging.
Keywords:Image interpolation; Linear least square; Low resolution; High resolution
Image interpolation in this article refers to algorithms that transform a low resolution (LR) image to high resolution (HR) one and has many applications. In medical applications, it is highly preferred for images to have high resolution while medical equipments cannot produce images with resolutions higher than a specific standard. In another application, with the development of digital monitors, we need to transform Standard Definition TeleVision (SDTV) video frames to High Definition TeleVision (HDTV) equivalents. Another application of this technology is to find missing pixels or blocks in an image.
Traditional methods such as bilinear, bicubic, and cubic convolution [1,2] are the common methods in image enlargement. Although the merits of these techniques are simplicity and fast implementation, they suffer from visual degradations such as jagged edges, blurring, and ringing around the image edges.
Human visual system is sensitive to the edges of the objects in an image. As a result, there are several mechanisms which improve the edges of an image. Undoubtedly, as the quality of image enhances, the computational complexity increases. Hence, in all of these methods, there is a trade off between performance and affordable complexity. In the last few decades, several authors have carried out researches on this topic. Along the same line, we briefly introduce [3-10] as follows.
In NEDI , the aim is to maintain the edge structure in HR images. This is a nonlinear algorithm whose basic idea is to estimate the covariance of HR image from the covariance of LR counterpart. The estimation is based on the geometric duality between low and high resolution. The algorithm presented in  which is based on directional filtering and data fusion to find the missing pixel, considers two orthogonal sets which produce an estimation of the pixel value. In SAI , using a moving window in the input LR image, one can find the model parameters. The pixel structure obtained in this way can be led to a block of both available and estimated pixels using a soft-decision estimation process. In WZP-CS , wavelet transformation is applied to LR image. Using sub-bands of the resultant image, one can estimate sub-bands of HR image. Afterwards, the HR image can be obtained by applying the inverse wavelet transform. In Optimized spline , by taking advantages of calculus of variations, the optimization problem is simplified from a nonlinear infinite dimensional case to a linear finite dimensional one in order to design compact-support interpolation kernels for a class of signals. The method presented in  tries to create up-scaled images that are suitable for real time applications and appear natural to human observer. In this algorithm a two step grid filling is used. Then, an iterative correction is obtained by minimizing an objective function depending on the second order directional derivatives of the image intensity.
In , in order to enlarge the image, it is down-sampled by using the bilinear method and then, the missing pixels are estimated by using a combination of the other pixels, while the desired coefficients are optimized by using a least mean square technique. In , the authors have extended the algorithm presented in  to display products such as video restoration.
In our approach, we work out the parameters that are used to enlarge an image from its LR counterpart. These parameters depend on the original image and the factor which we want to enlarge the image by. In order to enlarge an image or find the missing pixels in it, an interpolation array which is constructed by the mentioned parameters is convolved by the matrix that contains the values of the image pixels.
The rest of the article is organized as follows: The following section describes the proposed interpolation method for one-dimensional vectors. The two-dimensional application is proposed in section “Two-dimensional algorithm”. Simulation results and comparisons are given in section “Simulation results” , and section “Conclusion” concludes the paper.
In order to describe the proposed algorithm for image interpolation, we first apply it to one-dimensional vectors. Let XLR be a 1×N vector. Suppose we wish to enlarge XLR to a kN sized vector XHRusing the following steps: We put k−1 zeros between any two successive entries of vector XLR. In order to achieve these kN−Nnew entries of XHR, we choose 2M of N entries of XLR. The new value will be a linear combination of these 2M entries.
In order to find the 2M indeterminate multipliers of this linear combination, we should follow the procedure stated below, bearing in mind that for each fixed vector of any size the procedure and the multipliers applied to enlarge the vector are fixed. This means that, the multipliers used for obtaining a 1×N vector from its counterpart are the same as the ones used to calculate a vector from its counterpart.Step(1-1): As presented in (1) and (2), we down-sample the original vector XLR by a factor of k to find .Step(1-2): Again we down-sample by a factor of k to obtain vector.Step(1-3): Now, by constructing vector from vector of size , we can find the multipliers used to obtain the enlarged image. To approach this, we must zero-pad by a factor of k to build vector . This means that, k−1 zeros are inserted between every two successive entries of .
As it was stated earlier, by using these multipliers, we can find a new vector XHR from .
In order to evaluate the enlarged image, we take advantages of the criteria peak signal to noise ratio (PSNR) and Structural SIMilarity (SSIM)  by down-sampling the original vector by a factor of k and interpolating it using the proposed algorithm.
down-sample by a factor of k↓
down-sample by a factor of k↓
Zero padding ↓
Step(1-4): In order to compute vector from , we must insert a linear combination of ’s; instead of each zero in . For this purpose, we can convolve with an interpolator vector, A (which is determined by (5)).
In this equation, A is a 1×(1 + 2(M−1) + 2M(k−1)) vector. The 1 + (M−1) + M(k−1)st entry, which is located at the middle of the vector, is 1. Except for the mentioned entry, the entries with the indices that are multiples of k are 0, and there are (k−1)ai″s and (k−1)ais between any two successive zeros at the left and right side of the entry 1, respectively which can be calculated by using the proposed algorithm.
In (5), the parameter M can take an optional value which is set to 2 in our experiments. However, by setting M higher than this value, implementation speed decreases while image quality increases.
In (6), and are known and A is unknown, and there are equations and 2M(k−1) unknown variables. This statement means that, if (6) is to have a solution, must be greater than or equal to 2M(k−1).Step(1-5): To solve (6), we have the following matrix equation:
where , and A denotes the interpolator matrix. Step(1-6): In order to find matrix B, We add (Mk−1) zeros to the beginning and the end of vector to obtain vector .
where is the ith entry of vector .Step(1-7): By using , we can find matrix B as follows:
where Bi∗is the ithrow of matrix B.
Now, following Step(1-5) and then Step(1-4), we can obtain the multipliers used for image enlargement.
To understand the above algorithm, we can illustrate the above relation with an example in which N=36,K=3,M=2.Step(1-1):
A typical image is a matrix of N rows and N columns which can be enlarged by a factor of k to kN×kNby applying the one-dimensional algorithm discussed in the previous section to each row and column separately. Enlargement requires that we first down-sample the given matrix by factors of k1 and k2, then up-sample the modified matrix by factors of k1and k2 to obtain a new N×Nmatrix. In order to illustrate the mentioned statement, we assume a 16×16 matrix F. We also define a matrix with entries that are the odd entries of F. Down-sampling by a factor of k=2, we can obtain , and up-sampling by a factor of k=2 leads us to .
In this case, we know all the entries of F. Equation (7) can be applied to each row and column of . As an illustration, let us consider the following statement for the first row: C is the first row of matrix , B is obtained using (9), and matrix A is unknown. As a result, we should use (7) eight times, four times for the rows, and four times for the columns of matrix F. In as much as (7), itself, contains n equations, where n is the size of vector C, we have 4×8=32 equations (where four is the number of equations: ). We use matrices B and C that are obtained from (7) as the blocks of new matrices named and , respectively. Then, these 32 equations are solved using linear least square method; moreover, A can be calculated as follows:
Now, using A, we zero pad matrix by factors of k1and k1, through the following steps:
Figure 1 shows an image that is zero-padded by factors of k1=k2=2.
Figure 1. A zero-padded image by factors of k1 = k2 = 2.
The gray blocks in Figure 1 contain the pixels of the main N×Nimage and the remaining blocks contain zero. Now, in order to apply the proposed method to this figure, we need to take the following steps:
We first inspect the odd rows and then apply the rule mentioned before to each one to determine the amounts which must be replaced with zeros in the white blocks in Figure 1, using the amounts of adjacent blocks.
Next, we apply step (1) to the odd columns.
As can be seen from Figure 1, even rows and even columns contain zeros; therefore by using steps (1) and (2), the values that some of these blocks contain can be calculated. In Figure 2, blocks that contain circles are obtained by applying step (1), and those containing diamonds are obtained by applying step (2) to Figure 1.
Consequently, only the blocks in even rows and columns hold zeros (white blocks in Figure 2), and the goal is to determine the values of these blocks. To overcome this challenge, we can use one of the following approaches:
After applying step (1), the blocks that contain circles in Figure 2 are identified. Now, by applying step (2) to the same columns, we can obtain the remaining unknown blocks.
After applying step (2), the blocks that contain diamonds in Figure 2 are spotted. Then, by applying step (1) to the same rows, the remaining unknown blocks are extracted.
We can average the amounts obtained through methods (1) and (2) to each block to work out a new value.
In this way, the proposed method can provide us a high resolution image from a low resolution one. We can find the interpolator multipliers by employing the downsampled version of the original image. Then, the HR image can be estimated by using the original LR image and the computed multipliers.
The proposed image interpolation algorithm was implemented and its performance was compared with some existing methods. The results of the proposed method simulation are benchmarked against some classic methods such as bilinear, bicubic, and some recent methods such as SAI  and WZP-CS  that are based on wavelet transform, and the results are tabulated in Tables 1 and 2. Table 1 shows the PSNR of implementing the mentioned methods on nine 256×256 images, in order to enlarge them by a factor of 2, and the SSIM of them is indicated in Table 2. The simulations are implemented using Matlab 188.8.131.529(R2010a) software on a VAIO SR 490 DCB laptop with an Intel Core2 Duo P8800 2.66 GHz CPU and a 4 GB RAM.
As can be seen from these two tables, the results achieved by our proposed approach outperform those from other solutions. The bold values given in each row of these two tables indicate that the proposed method outperforms the other compared methods when applied to the images of the mentioned rows. The average PSNR of implementing the proposed method, is 29.0030 dB, which is higher than the other values of the compared solutions. On the other hand, the average SSIM of implementing the proposed method on the mentioned images is 0.8229 dB which is higher than the average SSIM for the other compared algorithms. In order to compare the methods objectively, Figure 3 shows that, applying the bilinear method yields blurring and also has less implementation speed than the proposed method. In addition, the scheme based on WZP-CS causes damages at the edges of the image. In the Cameraman image (Figure 4), applying the proposed method produces a better visual quality compared to other known models.
Figure 3. Comparison of different methods for the Lena image: (a) the original image, (b) bilinear interpolation, (c) bicubic interpolation, (d) WZP cycle-spinning, (e) SAI, (f) NEDI, (g) optimized spline, and (h) the proposed method.
Also, since we do not use blocks in this presentation, unlike the SAI method, there is no blocking result in the obtained image.
Moreover, one of the most important criteria that should be taken into account for comparison is computational complexity of the algorithms, which can be computed by considering the number of operations in the algorithms. In the proposed method, the number of multiplications and summations in order to calculate the matrix A from expression (15) for an N×N image to be enlarged by a factor of K are and , respectively; where 2M is the number of entries of the LR image that are used for interpolation of each pixel. By assuming M=K=2 in our experiments, there will be (28N + 5201) multiplications and (28N + 2499) summations for calculating the matrix A. On the other hand, after calculating the matrix A, the number of multiplications and summations to estimate the value of each pixel are and , respectively. By assuming M=K=2 in our experiments, there will be 3N2multiplications and summations for estimating the value of each pixel. The number of operations to estimate the value of each pixel is approximately on the same order in the compared methods, and the difference is in the number of operations in order to calculate the multipliers of the linear combination, hence, in order to compare the proposed method with the other mentioned ones, we consider the number of operations in order to calculate the matrix A, where, this number of operations is acceptable for our purpose.
In addition to this, since all algorithms are simulated in the same situation, running time may have a better sense of complexity. As it is indicated in Table 3, the proposed method outperforms the high performance algorithms such as SAI and optimized spline in terms of the average time that is needed to simulate the algorithm, and features a higher implementation speed.
Table 3. Average calculation time for enlarging a 256 * 256 image by a factor of 2
The linear least square study unveiled in this article also offers an advantage when transmitting an image where we have memory limitation; i.e. instead of transmitting the enlarged version of the original image by a factor of k, we can transmit the original image and the interpolator matrix A, and then use the interpolator matrix to enlarge the original image by a factor of k at the receiving end. Accordingly, less memory is used in transmitting the image. In the images that are simulated, the interpolator matrix entries are decimal numbers with three digits to the right of the decimal point; hence, ten bits would be sufficient in order to send the interpolator matrix as side information.
In this article, we presented a novel image interpolation method. For each image to be interpolated, we found an interpolator function which was applied to each row and column of the image separately by convolving this interpolator function with the vector considered for interpolation. PSNR and SSIM criteria were then used to compare our results versus those obtained by other algorithms. Experimental results demonstrated that the solution revealed in this article achieves better results both objectively and subjectively.
Another key feature of the proposed method is that it allows enlargement of the image by any factor of k=2,3,4,5,… by simply setting this parameter to the desired value.
The authors declare that they have no competing interests.
The authors would like to thank Seyyed Amirhossein Hosseini for editing this article and his helpful comments.
R Keys, Cubic convolution interpolation for digital image processing. IEEE Trans. Acoust. Speech Signal Process 29(6), 1153–1160 (1981). Publisher Full Text
HS Hou, HC Andrews, Cubic splines for image interpolation and digital filtering. IEEE Trans. Acoust. Speech Signal Process 26(6), 508–517 (1978). Publisher Full Text
L Zhang, X Wu, An edge guided image interpolation algorithm via directional filtering and data fusion. IEEE Trans. Image Process 15(8), 2226–2238 (2006). PubMed Abstract
A Temizel, T Vlachos, W Visioprime, Wavelet domain image resolution enhancement using cycle-spinning. Electron. Lett 41(3), 119–121 (2005). Publisher Full Text
L Shao, H Zhang, G de Haan, An overview and performance evaluation of classification based least squares trained filters. IEEE Trans. Image Process 17(10), 1772–1782 (2008). PubMed Abstract | Publisher Full Text
Z Wang, AC Bovik, HR Sheikh, EP Simoncelli, Image quality assessment: from error visibility to structural similarity. IEEE Trans. Image Process 13(4), 600–612 (2004). PubMed Abstract | Publisher Full Text