OSA's Digital Library

Optics Express

Optics Express

  • Editor: C. Martijn de Sterke
  • Vol. 18, Iss. 17 — Aug. 16, 2010
  • pp: 17819–17833
« Show journal navigation

Tensor factorization for model-free space-variant blind deconvolution of the single- and multi-frame multi-spectral image

Ivica Kopriva  »View Author Affiliations


Optics Express, Vol. 18, Issue 17, pp. 17819-17833 (2010)
http://dx.doi.org/10.1364/OE.18.017819


View Full Text Article

Acrobat PDF (5423 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

The higher order orthogonal iteration (HOOI) is used for a single-frame and multi-frame space-variant blind deconvolution (BD) performed by factorization of the tensor of blurred multi-spectral image (MSI). This is achieved by conversion of BD into blind source separation (BSS), whereupon sources represent the original image and its spatial derivatives. The HOOI-based factorization enables an essentially unique solution of the related BSS problem with orthogonality constraints imposed on factors and the core tensor of the Tucker3 model of the image tensor. In contrast, the matrix factorization-based unique solution of the same BSS problem demands sources to be statistically independent or sparse which is not true. The consequence of such an approach to BD is that it virtually does not require a priori information about the possibly space-variant point spread function (PSF): neither its model nor size of its support. For the space-variant BD problem, MSI is divided into blocks whereupon the PSF is assumed to be a space-invariant within the blocks. The success of proposed concept is demonstrated in experimentally degraded images: defocused single-frame gray scale and red-green-blue (RGB) images, single-frame gray scale and RGB images blurred by atmospheric turbulence, and a single-frame RGB image blurred by a grating (photon sieve). A comparable or better performance is demonstrated in relation to the blind Richardson-Lucy algorithm which, however, requires a priori information about parametric model of the blur.

© 2010 OSA

1. Introduction

The rest of this paper is organized as follows: in Section 2 tensor notation, the Tucker3 tensor model and HOOI decomposition are introduced. The multichannel linear mixture model of blurred image tensor is introduced in Section 3. The proposed approach to various types of blind image deconvolution problems is exemplified in Section 4 on an: (i) experimental de-focused single-frame gray scale and red-green-blue (RGB) images, (ii) experimental multi-frame gray scale and RGB images degraded by atmospheric turbulence and, (iii) experimental single-frame RGB image degraded by a grating. Conclusions are given in Section 5.

2. Basics of tensor notation, Tucker 3 model and HOOI decomposition

A tensor is a multi-way array of data: G¯I1×I2...×IN, where N denotes the number of modes or order of the tensor G¯ and {In+}n=1Ndenotes dimensions of each mode, where +denotes a set of positive integers. A particular element of the tensor G¯ is given with gi1i2...iNwhere i 1 = 1,..., I 1, i 2 = 1,..., I 2, and iN = 1, ..., IN. This is the standard notation adopted for use in multi-way analysis [29

29. H. A. L. Kiers, “Towards a standardized notation and terminology in multiway analysis,” J. Chemometr. 14(3), 105–122 (2000). [CrossRef]

]. The Tucker3 model represents tensor G¯ as an n-mode product between the core tensor R¯J1×J2...×JNand array factors {A(n)In×Jn}n=1N [27

27. A. Cichocki, R. Zdunek, A. H. Phan, and S. I. Amari, Nonnegative Matrix and Tensor Factorization (John Wiley & Sons, 2009).

,30

30. L. R. Tucker, “Some mathematical notes on three-mode factor analysis,” Psychometrika 31(3), 279–311 (1966). [CrossRef] [PubMed]

,31

31. A. Cichocki, and A. H. Phan. Fast Local Algorithms for Large Scale Nonnegative Matrix and Tensor Factorizations. IEICE Trans Fundamentals 2009; E92-A(3): 708–721.

]:
G¯R¯×1A(1)×2A(2)...×NA(N)
(1)
or on the component level

gi1i2...iNj1=1J1j2=1J2...jN=1JNrj1j2...jNai1j1(1)ai2j2(2)...aiNjN(N)
(2)

D[G¯G¯^]=G¯R¯×1A(1)×2A(2)...×NA(N)22
(7)

Due to orthogonality constraints the core tensor R¯can also be represented as
R¯G¯×1A(1)T×2A(2)T...×NA(N)T
implying that instead of minimizing Eq. (7) the cost function
D(A(1),A(2),...,A(N))=G¯×1A(1)T×2A(2)T...×NA(N)T22
(8)
can be maximized where only the array factors {A(n)}n=1Nare unknown. This cost function is maximized in alternating least square fashion where in each iteration only one array factor is optimized while keeping the others fixed. Results reported in experimental section were obtained by using an implementation of the HOOI algorithm provided through tucker_als function made available as a part of MATLAB Tensor Toolbox [32

32. B. W. Bader, and T. G. Kolda, MATLAB Tensor Toolbox version 2.2.http://csmr.ca.sandia.gov/~tkolda/TensorToolbox.

]. Note that although the original image is nonnegative, this constraint was not applied to F¯^ in Eq. (6). The reason is that orthogonality constraints, imposed by HOOI algorithm on the core tensor and array factors in tensor model [Eq. (1)], and nonnegativity constraints cannot be satisfied simultaneously.

3. Multi-dimensional linear mixture model of degraded image

Linear mixture model of a single-frame gray scale image degraded by space-invariant blur is based on an image forming equation: a convolution of space-invariant PSF H with an original source image F:
G(i1,i2)=s=MMt=MMH(s,t)F(i1s,i2t)
(9)
where M denotes half of the PSF support size. In Eq. (9) presence of the additive noise is ignored in order to focus on an essential issue: model-free blind image deconvolution. It is however clear that some form of image de-noising is performed, if necessary, prior to executing deconvolution. It is also assumed that the unknown original image F is n th order smooth implying that it is n-times differentiable at the origin (i 1,i 2), whereupon n represents the order of spatial derivatives (terms) in the Taylor expansion of F(i1s,i2t)around the origin (i 1,i 2). An implicit Taylor expansion of the original image F(i 1 -s,i 2 -t) around (i 1,i 2), is used to convert blind image deconvolution into blind source separation [20

20. S. Umeyama, “Blind deconvolution of blurred images by use of ICA,” Electron Commun. Jpn 84(12), 1–9 (2001).

,5

5. I. Kopriva, “Single-frame multichannel blind deconvolution by nonnegative matrix factorization with sparseness constraints,” Opt. Lett. 30(23), 3135–3137 (2005). [CrossRef] [PubMed]

,6

6. I. Kopriva, “3D tensor factorization approach to single-frame model-free blind-image deconvolution,” Opt. Lett. 34(18), 2835–2837 (2009). [CrossRef] [PubMed]

,17

17. I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]

], yielding:
F(i1s,i2t)=F(i1,i2)sFi1(i1,i2)tFi2(i1,i2)                                                             +12s2Fi1i1(i1,i2)+12t2Fi2i2(i1,i2)+12stFi1i2(i1,i2)...
(10)
where Fi1, Fi2 Fi1i1, Fi2i2 and Fi1i2are first- and second-order spatial derivatives in i 1 and i 2 directions respectively. By using Eq. (10), Eq. (9) can be re-written as:

G(i1,i2)=a1F(i1,i2)a2Fi1(i1,i2)a3Fi2(i1,i2)+                                       a4Fi1i1(i1,i2)+a5Fi2i2(i1,i2)+a6Fi1i2(i1,i2)...
(11)

The unknown weighting coefficients in Eq. (11) are straightforward to derive and are given in [5

5. I. Kopriva, “Single-frame multichannel blind deconvolution by nonnegative matrix factorization with sparseness constraints,” Opt. Lett. 30(23), 3135–3137 (2005). [CrossRef] [PubMed]

,6

6. I. Kopriva, “3D tensor factorization approach to single-frame model-free blind-image deconvolution,” Opt. Lett. 34(18), 2835–2837 (2009). [CrossRef] [PubMed]

,17

17. I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]

]. It is of great importance to note that blurred image [Eq. (11)] represents a linear combination of the original image and its spatial derivatives. The unknown weighting coefficients a 1 to a 6 absorbed into themselves the coefficients of the PSF: {H(s,t)}s,t=MM, including the support size parameter: M. Hence, blind image deconvolution could be converted to blind source separation provided that a multi-channel version of the blurred image [Eq. (11)] is available. As in [20

20. S. Umeyama, “Blind deconvolution of blurred images by use of ICA,” Electron Commun. Jpn 84(12), 1–9 (2001).

,5

5. I. Kopriva, “Single-frame multichannel blind deconvolution by nonnegative matrix factorization with sparseness constraints,” Opt. Lett. 30(23), 3135–3137 (2005). [CrossRef] [PubMed]

,6

6. I. Kopriva, “3D tensor factorization approach to single-frame model-free blind-image deconvolution,” Opt. Lett. 34(18), 2835–2837 (2009). [CrossRef] [PubMed]

,17

17. I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]

], this is achieved by applying a bank of 2D Gabor filters to the blurred image [Eq. (11)]. The 2D Gabor filter bank consists of filters with two spatial frequencies and four orientations, where the real and imaginary parts of the filter are applied separately. A more detailed description of the 2D Gabor filter bank can be found in [20

20. S. Umeyama, “Blind deconvolution of blurred images by use of ICA,” Electron Commun. Jpn 84(12), 1–9 (2001).

,21

21. J. G. Daugman, ““Complete Discrete 2-D Gabor Transforms by Neural Networks for Image Analysis and Compression,” IEEE Trans. Acoust. Speech Signal Process. 36(7), 1169–1179 (1988). [CrossRef]

,17

17. I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]

]. Applying it to degraded image [Eq. (11)] yields a new set of observed images:
Gi3(i1,i2)=ai31F(i1,i2)   ai32Fi1(i1,i2)ai33Fi2(i1,i2)+                                                                   ai34Fi1i1(i1,i2)   +ai35Fi2i2(i1,i2)   +ai36Fi1i2(i1,i2)...
(12)
where expressions for weighting coefficients are straightforward to derive and can also be found in [5

5. I. Kopriva, “Single-frame multichannel blind deconvolution by nonnegative matrix factorization with sparseness constraints,” Opt. Lett. 30(23), 3135–3137 (2005). [CrossRef] [PubMed]

,6

6. I. Kopriva, “3D tensor factorization approach to single-frame model-free blind-image deconvolution,” Opt. Lett. 34(18), 2835–2837 (2009). [CrossRef] [PubMed]

,17

17. I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]

]. Since the Gabor filter bank is composed of filters with two spatial frequencies and four orientations, and since real and imaginary parts are used as separate filters, there will be 16 new images after filtering a blurred image with the Gabor filter bank. Since the multi-channel version of a blurred image also contains a blurred image itself, in overall it will contain I 3 = 17 images. Thus, by together merging a blurred image [Eq. (11)] with filtered images [Eq. (12)], the multichannel image forms a 3D tensorG¯      0+I1×I2×I3. In relation to the N-order Tucker3 tensor model [Eq. (1)], this scenario is characterized by N = 3. Based on Eq. (6), it follows that an estimate of the tensor of source images is obtained by means of the HOOI-based tensor factorization as F^¯=G¯×3(A(3)), where an estimate of the original image is obtained by setting the index that corresponds with the third dimension of F^¯ to j = 1 i.e. the estimate of the original image corresponds to the first source image.

The linear mixture model [Eqs. (11) and (12)] represents a single-frame gray scale image blurred by a space-invariant blur. Its generalization to more complex scenarios elaborated in Table 1 is obtained by expanding the 2D tensor [Eq. (11)] by adding new indices that correspond to the newly introduced modalities. For example, the space-variant deconvolution of a single-frame gray scale image assumes that the blurred image is divided into I 3 blocks, where PSF is assumed to be a space-invariant within each block. Then, each of the I 3 image blocks is filtered by means of a 2D Gabor filter bank yielding in overall a 4D multi-channel image tensor G¯      0+I1×I2×I3×I4. The estimate of the tensor of source images is obtained as F^¯=G¯×4(A(4)), where an estimate of the 3D tensor of the blocks of the original image is obtained by setting an index that corresponds with the fourth dimension of F^¯ to j = 1. The original image itself is then obtained by a proper rearranging of the blocks into a matrix.

Following the same line of reasoning space-invariant deconvolution of a single-frame multi-spectral image adds a third dimension to [Eq. (11)], where I 3 represents the number of spectral bands. Filtering each of the I 3 images by means of a 2D Gabor filter bank yields a 4D multi-channel image tensor G¯      0+I1×I2×I3×I4. The estimate of the tensor of source images is obtained as F^¯=G¯×4(A(4)), where an estimate of the 3D tensor of original multi-spectral image is obtained by setting an index that corresponds with the fourth dimension of F^¯ to j = 1. The physical interpretation of tensor modes associated with blind deconvolution scenarios considered in this paper is provided in the Table 1.

4. Experimental results for images blurred by de-focus, atmospheric turbulence and grating

In this section, the success of the tensor-factorization approach in solving various types of blind image deconvolution problems, according to the classification given in Table 1, is demonstrated in three experimental scenarios. As discussed previously, the proposed approach should be considered rather complementary than competing to physically constrained iterative blind image deconvolution algorithms for situations when physically based constraints are difficult or impossible to be defined. Therefore, a comparison of the proposed method against application-specific blind deconvolution methods (for example, against methods developed for deconvolution of turbulence-degraded image in astronomy) is not fair, because these methods will almost surely perform better when a priori information is provided to them.

Experimental scenarios include: (i) de-focused RGB image shown in Fig. 1
Fig. 1 RGB experimental image obtained by digital camera in manually defocused mode.
, as well as its gray scale version, with dimensions of 384 × 512 pixels. Solutions of blind deconvolution problems classified in Table 1 as 2, 3 and 4 have been demonstrated on this image. The image has been recorded by digital camera in a manually de-focused mode and has been used previously in [5

5. I. Kopriva, “Single-frame multichannel blind deconvolution by nonnegative matrix factorization with sparseness constraints,” Opt. Lett. 30(23), 3135–3137 (2005). [CrossRef] [PubMed]

,6

6. I. Kopriva, “3D tensor factorization approach to single-frame model-free blind-image deconvolution,” Opt. Lett. 34(18), 2835–2837 (2009). [CrossRef] [PubMed]

,17

17. I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]

].; (ii) RGB and gray scale multi-frame images of the Washington monument blurred by atmospheric turbulence. Four randomly chosen frames with dimensions of 160 × 80 pixels are shown in Fig. 6a
Fig. 6 a) Four frames of experimental multi-frame RGB image degraded by atmospheric turbulence. b) left: average of the four frames from 6a; right: edges extracted from gray scale version of averaged image with Canny's algorithm and threshold set to 0.21.
. The gray scale version of this image sequence has been used previously in [33

33. I. Kopriva, Q. Du, H. Szu, and W. Wasylkiwskyj, “Independent Component Analysis Approach to Image Sharpening in the Presence of Atmospheric Turbulence,” Opt. Commun. 233(1-3), 7–14 (2004). [CrossRef]

] in formulation of an independent component analysis-based approach to sharpening an image blurred by atmospheric turbulence. Solutions of blind deconvolution problems classified in Table 1 as 5 and 7 have been demonstrated on these images.; (iii) single-frame RGB image blurred by a grating. The image with dimensions of 301 × 351 pixels is shown in the leftmost picture of Fig. 10
Fig. 10 From left to right: RGB image degraded by a grating (photon sieve); image restored by 5D tensor factorization space-variant blind deconvolution; image restored by 4D tensor factorization space-invariant blind deconvolution; image restored by blind R-L algorithm after 10 iterations and radius of the circular blur equal to 2 pixels.
. The grating-caused blur does not have a real physical equivalent such as turbulence, de-focus or motion. Therefore, it is virtually impossible to define a proper model or other type of a priori information required by physically constrained iterative blind deconvolution algorithms such as exemplified by [7

7. F. Li, X. Jia, and D. Fraser, “Superresolution Reconstruction of Multispectral Data for Improved Image Classification,” IEEE Geosci. Remote Sens. Lett. 6(4), 689–693 (2009). [CrossRef]

14

14. M. Šorel and J. Flusser, “Space-variant restoration of images degraded by camera motion blur,” IEEE Trans. Image Process. 17(2), 105–116 (2008). [CrossRef] [PubMed]

]. Thus, this grating-caused blur represents a good example for demonstrating general-purpose character of a proposed model-free approach to blind image deconvolution. Solutions of blind deconvolution problems classified in Table 1 as 3 and 4 have been demonstrated on this image.

All results obtained on the experimental images by the tensor factorization approach are compared against results obtained by a blind R-L algorithm [25

25. D. A. Fish, A. M. Brinicombe, E. R. Pike, and J. G. Walker, “Blind deconvolution by means of the Richardson-Lucy algorithm,” J. Opt. Soc. Am. A 12(1), 58–65 (1995). [CrossRef]

,26

26. D. S. C. Biggs and M. Andrews, “Acceleration of iterative image restoration algorithms,” Appl. Opt. 36(8), 1766–1775 (1997). [CrossRef] [PubMed]

], that is implemented by a MATLAB function deconvblind. The blind R-L algorithm is model-based and requires a priori information about the blur in the form of a parametric model. As other model-based deconvolution methods, it is sensitive to model misspecifications and a miss-estimation of the values of model parameters [1

1. R. L. Lagendijk, and J. Biemond, Iterative Identification and Restoration of Images (KAP, 1991).

]. As can be seen in Figs. 4
Fig. 4 Experimental RGB image restored by blind R-L algorithm after 5 iterations and radius of the circular blur equal to 2 pixels (left) and 3 pixels (right) .
, 9
Fig. 9 Images restored from averaged gray scale version of four blurred frames shown in Fig. 6a by blind R-L algorithm and Gaussian PSF with kernel width of 18 pixels and a) sigma = 1.3 pixels; b) sigma = 1.9 pixels. Shown edge maps were obtained with Canny's algorithm and threshold set to 0.21.
and 10, the blind R-L algorithm yielded results of equal or lower quality than those obtained by the proposed tensor factorization model-free approach to blind deconvolution, despite the fact that the optimal choice of a blur model and model parameters were supplied to the algorithm. Figures 4 and 9 also illustrate sensitivity of the blind R-L method to a slight over-estimation of model parameters. In Fig. 10 the best result for blind R-L is shown after many attempts to choose the blur model and model parameters.

4.1 Space-variant blind deconvolution of de-focused single-frame gray scale image

According to Table 1, this is a blind deconvolution problem 2. This problem results in 4D tensor factorization. The tensor model [Eq. (1)] of the gray scale version of a de-focused RGB image shown in Fig. 1 is characterized with: I 1 = 48, I 2 = 64, I 3 = 64, and I 4 = 17. The image has been divided into 64 blocks with the size of 48 × 64 pixels and PSF is assumed to be the space-invariant within each block. Figure 2
Fig. 2 Experimental gray scale version of RGB image restored by left: space-variant single-frame blind deconvolution 4D tensor factorization algorithm; right: single-frame blind deconvolution algorithm where each block is restored through one 3D tensor factorization.
-left shows a result of the 4D tensor factorization approach. Figure 2-right shows a corresponding result, with poorer resolution, obtained when each block has been deconvolved separately through 3D tensor factorization. This result is explained by the fact that there is only one scaling indeterminacy associated with the solution of the 4D tensor factorization problem, while there are I 3 = 64 scaling indeterminacies associated with the solutions of 3D tensor factorization problems. White artifact lines and a sharp difference between the few blocks are visible in restored image and are the consequence of an adopted block-wise approach to space-variant blind deconvolution. However, the super-resolution effect of the deconvolution method is demonstrated. I consider the block-wise induced artifacts as an open problem that can possibly be reduced through some continuity (smoothness) constraints either directly through factorization or in a form of post-processing. Choosing the optimal size of the image blocks is also an open issue to be explored in the future.

4.2 Space-invariant blind deconvolution of de-focused single-frame multi-spectral image

According to Table 1, this is a blind deconvolution problem 3. It results in 4D tensor factorization. Tensor model [Eq. (1)] of multi-channel de-focused RGB image is characterized with: I 1 = 384, I 2 = 512, I 3 = 3 and I 4 = 17. Figure 3
Fig. 3 Experimental RGB image restored by left: space-invariant single-frame blind deconvolution 4D tensor factorization algorithm; right: single-frame blind deconvolution where each spectral channel is restored through one 3D tensor factorization problem.
-left shows the result of a 4D tensor factorization approach to model-free blind deconvolution. A corresponding result, with poorer resolution, is presented in Fig. 3-right, and also in [7

7. F. Li, X. Jia, and D. Fraser, “Superresolution Reconstruction of Multispectral Data for Improved Image Classification,” IEEE Geosci. Remote Sens. Lett. 6(4), 689–693 (2009). [CrossRef]

], where each spectral image has been deconvolved separately through 3D tensor factorization. This result is explained by the fact that there is only one scaling indeterminacy associated with the solution of the 4D tensor factorization problem, while there are I 3 = 3 scaling indeterminacies associated with the solutions of 3D tensor factorization problems. The result obtained by model-free blind deconvolution has served as a reference to optimize parameters for a blind R-L algorithm, the results of which are shown in Fig. 4. Knowing that the RGB image was de-focused, a disk with the proper radius was supplied to the algorithm as a blur model. The optimal result is shown in Fig. 4-left. As shown in Fig. 4-right, over-estimating the disk radius for 1 pixel only deteriorated the performance of the blind R-L method significantly.

4.3 Space-variant blind deconvolutin of a de-focused single-frame multi-spectral image

According to Table 1, this is a blind deconvolution problem 4. To demonstrate this case, the RGB image shown in Fig. 1 has been divided into 64 blocks with the size of 48 × 64 pixels. PSF is assumed to be space-invariant within each block. The 5D tensor of multi-channel blurred image [Eq. (1)] is characterized with: I 1 = 48, I 2 = 64, I 3 = 64, I 4 = 3 and I 5 = 17. Figure 5
Fig. 5 Experimental RGB image restored by left: 5D tensor factorization space-variant single-frame blind deconvolution; right: 4D tensor factorization single-frame blind deconvolution, where each block at each spectral channel is restored through a 3D tensor factorization.
-left shows a result obtained with a proposed tensor factorization approach. Note, in direct comparison with Fig. 4-left and Fig. 3-left, that the book title and names of the authors are more readable (details are better resolved) in Fig. 5-left. Figure 5-right shows a corresponding result, with poorer resolution, obtained by deconvolving each of I 3 = 64 blocks in each of I 4 = 3 spectral bands separately through 3D tensor factorizations. This result is explained by the fact that there is only one scaling indeterminacy associated with the solution of the 5D tensor factorization problem, while there are I 3 × I 4 = 192 scaling indeterminacies associated with the solutions of 3D tensor factorization problems.

4.4 Space-invariant blind deconvolution of a multi-frame multi-spectral image blurred by atmospheric turbulence

According to Table 1, this is a blind deconvolution problem 7. To demonstrate this case, the multi-frame RGB image of the Washington monument is used. Four frames chosen randomly are shown in Fig. 6a. The 5D tensor [Eq. (1)] of a blurred multi-channel multi-frame and multi-spectral image is characterized with: I 1 = 160, I 2 = 80, I 3 = 3, I 4 = 4 and I 5 = 17. Figure 6b-left shows an average of the four frames shown in Fig. 6a, while Fig. 6b-right shows an edge map extracted from Fig. 6b-left by Canny's algorithm with the threshold set to 0.21. Figure 7a
Fig. 7 a) Time evolution of RGB source image extracted by 5D tensor factorization space-invariant multi-frame blind deconvolution algorithm from experimental dynamic RGB shown in Fig. 6a. b) left: average of the four source images from 7a; right: edges extracted from the gray scale version of averaged image with Canny's algorithm and threshold set to 0.21.
shows four source image frames restored with a proposed tensor factorization approach, while Fig. 7b-left shows the average of these four source image frames. Figure 7b-right shows an edge map extracted from Fig. 7b-left by Canny's algorithm with the threshold set again to 0.21. Details like windows and the top of the monument are reconstructed. Some of these details (windows) are missed on edges extracted from turbulence degraded image. They can be recovered by reducing the threshold (increasing the sensitivity) of Canny's algorithm. However, in this case, turbulence-induced artifacts would be picked up as well.

4.5 Space-invariant blind deconvolution of multi-frame gray scale image blurred by atmospheric turbulence

According to Table 1, this is a blind deconvolution problem 5. To demonstrate this case, a gray scale version of the RGB image sequence from section 4.4 is used. 4D tensor [Eq. (1)] of a blurred image is characterized with: I 1 = 160, I 2 = 80, I 3 = 4 and I 4 = 17. Figure 8a
Fig. 8 a) Time evolution of gray scale source image restored by 4D tensor factorization space-invariant multi-frame blind deconvolution algorithm from gray scale version of the RGB shown in Fig. 6a. b) left: average of the four source images from 8a; right: edges extracted from averaged image with Canny's algorithm and threshold set to 0.21.
shows four source image frames obtained with the proposed tensor factorization approach, while Fig. 8b-left shows an average of the four source image frames shown in Fig. 8a. Figure 8b-right shows an edge map extracted from Fig. 8b-left by Canny's algorithm with the threshold set again to 0.21. Again, details like windows and the top of the monument are reconstructed. Figure 9a-left shows the image restored from a gray scale version of the average of the four blurred frames, shown in Fig. 6a, by means of a blind R-L algorithm. The corresponding edge map extracted by Canny's method with a threshold set again to 0.21 is shown in Fig. 9a-right. Knowledge of the type of degradation has been used to specify the correct model of the blur that is required by blind R-L algorithm: 2D Gaussian with a support size of 18 pixels and a standard deviation of 1.3 pixels has been used to obtain the shown result. The result obtained by the model-free approach, shown in Fig. 8b, has been used as a reference to obtain optimal values of the parameters of the 2D Gaussian. Yet, the blind R-L algorithm with parameters optimized for a blur model yielded a result that is not better than the one obtained by a model-free blind deconvolution approach. Figure 9b shows the result obtained by a blind R-L algorithm when the standard deviation is over-estimated to 1.9 pixels. This, again, illustrates sensitivity of the blind R-L algorithm to miss-estimation of model parameters. No such problems arise with the proposed model-free approach to blind deconvolution.

4.6 Space-(in)variant blind deconvolutin of single-frame multi-spectral image blurred by a grating

According to Table 1, these are blind deconvolution problems 3 and 4. Space-variant blind deconvolution problem 4, results in 5D tensor factorization. The tensor (1) of the blurred image is characterized with: I 1 = 150, I 2 = 117, I 3 = 6, I 4 = 3 and I 5 = 17. A grating blurred RGB image is shown in the leftmost picture in Fig. 10. The true image is composed of a palette of color pens, and a painting on the white board. Since there is no “real life” physical analogy to this grating-caused blur, it is virtually impossible to select a specialized method to perform blind deconvolution of this grating-blurred image. The picture second from the left in Fig. 10 shows the result obtained by a solution of space-variant blind deconvolution problem 4. The picture second from the right in Fig. 10 shows the result obtained by a space-invariant blind deconvolution problem 3 that results in 4D tensor factorization. The tensor (1) of the blurred image is characterized with: I 1 = 301, I 2 = 351, I 3 = 3 and I 4 = 17. The best result obtained by blind R-L algorithm for a space-invariant problem is shown in the rightmost picture of Fig. 10. It has been obtained after many attempts to find both the proper model and parameters of the blur. Still, obtained results is not better than the one show in the picture second from right and obtained by a solution of space-invariant blind deconvolution problem 3. It is, however, visible that the space-variant version of blind deconvolution really improved the spatial resolution of the restored image significantly in relation to images restored by space-invariant blind deconvolution methods. This result has been obtained truly without any a priori information provided to the blind deconvolution algorithm.

5. Conclusion and future work

The HOOI-based tensor factorization approach has been proposed for the space-(in)variant model-free blind deconvolution of a single- and multi-frame multi-spectral image. This is achieved by converting blind deconvolution into blind source separation using the implicit Taylor expansion of the original image in the convolution image-forming equation. Herein, the sources represent the original image and its spatial derivatives. In addition to generality, the two major contributions of the proposed approach to blind image deconvolution are: (i) the HOOI-based factorization of the tensor of the blurred image is essentially unique with no hard constraints imposed on source images compared to matrix factorization based methods; (ii) neither model nor size of the support of the point spread function is required to be a priori known or estimated. However, use of the implicit Taylor expansion implies a certain level of smoothness of the original image. This might limit the performance of the proposed approach to blind image deconvolution when the blurring process is strong or the original image contains sharp boundaries. Nevertheless, the proposed method is expected to be useful in scenarios when a priori information required by physically constrained iterative blind deconvolution methods are difficult or impossible to define. The two fundamental issues are considered to be important exploring in the future work: optimal selection of the size of the image blocks and neutralization of block-wise induced artifacts associated with space-variant deconvolution; sequence partitioning-based methods as a possible replacement of Gabor filter bank approach to single-channel blind source separation.

Appendix: Mode-n unfolding

The mode-n unfolding of a tensor G¯I1×I2...×INis denoted by G(n)and arranges the mode-n fibers of the tensor G¯into a matrix G(n), see also definition 1.4 in [27

27. A. Cichocki, R. Zdunek, A. H. Phan, and S. I. Amari, Nonnegative Matrix and Tensor Factorization (John Wiley & Sons, 2009).

], as well as [31

31. A. Cichocki, and A. H. Phan. Fast Local Algorithms for Large Scale Nonnegative Matrix and Tensor Factorizations. IEICE Trans Fundamentals 2009; E92-A(3): 708–721.

]. A tensor element (i 1, i 2,...,i N) maps onto a matrix element (i n, j), where

j=1+pn(ip1)Jp,            with      Jp={1,                if    p=1      or    p=2      and    n=1mnp1Im,              otherwise.

For example, a third order tensor G¯I1×I2×I3, with elements gi1i2i3and indices (i 1,i 2,i 3) has a corresponding position (i n, j) in the mode-n unfolded matrix G(n) (n=1,2,3) as follows: mode-1: j = i 2 + (i 3 - 1)I 2; mode-2: j = i 1 + (i 3 - 1)I 1; mode-3: j = i 1 + (i 2 - 1)I 1.

Acknowledgments

This work was supported by the Ministry of Science, Education and Sports, Republic of Croatia, grant 098-0982903-2558. I also thank the anonymous reviewers for their very constructive comments. Professor Wasyl Wasylkiwskyj's and Dr. Jordi Sancho-Parramon's help in proofreading the manuscript is also gratefully acknowledged.

References and links

1.

R. L. Lagendijk, and J. Biemond, Iterative Identification and Restoration of Images (KAP, 1991).

2.

M. R. Banham and A. K. Katsaggelos, “Digital Image Restoration,” IEEE Signal Process. Mag. 14(2), 24–41 (1997). [CrossRef]

3.

D. Kundur and D. Hatzinakos, “Blind Image Deconvolution,” IEEE Signal Process. Mag. 13(3), 43–64 (1996). [CrossRef]

4.

P. Campisi, and K. Egiazarian, eds., Blind Image Deconvolution (CRC Press, Boca Raton, 2007).

5.

I. Kopriva, “Single-frame multichannel blind deconvolution by nonnegative matrix factorization with sparseness constraints,” Opt. Lett. 30(23), 3135–3137 (2005). [CrossRef] [PubMed]

6.

I. Kopriva, “3D tensor factorization approach to single-frame model-free blind-image deconvolution,” Opt. Lett. 34(18), 2835–2837 (2009). [CrossRef] [PubMed]

7.

F. Li, X. Jia, and D. Fraser, “Superresolution Reconstruction of Multispectral Data for Improved Image Classification,” IEEE Geosci. Remote Sens. Lett. 6(4), 689–693 (2009). [CrossRef]

8.

H. Ji and C. Fermüller, “Robust wavelet-based super-resolution reconstruction: theory and algorithm,” IEEE Trans. Pattern Anal. Mach. Intell. 31(4), 649–660 (2009). [CrossRef] [PubMed]

9.

T. J. Schulz, “Multiframe blind deconvolution of astronomical images,” J. Opt. Soc. Am. A 10(5), 1064–1073 (1993). [CrossRef]

10.

R. Fergus, B. Singh, A. Hertzmann, S. T. Roweis, and W. T. Freeman, “Removing Camera Shake from a Single Photograph,” ACM Trans. Graph. 25(3), 787–794 (2006). [CrossRef]

11.

Q. Shan, J. Jia, and A. Agarwala, “High-quality Motion Deblurring from a Single Image,” ACM Trans. Graph. 27(3), 1 (2008). [CrossRef]

12.

S. Cho and S. Lee, “Fast Motion Deblurring,” ACM Trans. Graph. 28(5), 1 (2009). [CrossRef]

13.

J. Miskin, and D. J. C. MacKay, “Ensemble Learning for Blind Image Separation and Deconvolution,” in: Advances in Independent Component Analysis (Springer, London, 2000).

14.

M. Šorel and J. Flusser, “Space-variant restoration of images degraded by camera motion blur,” IEEE Trans. Image Process. 17(2), 105–116 (2008). [CrossRef] [PubMed]

15.

J. Bardsley, S. Jefferies, J. Nagy, and R. Plemmons, “A computational method for the restoration of images with an unknown, spatially-varying blur,” Opt. Express 14(5), 1767–1782 (2006). [CrossRef] [PubMed]

16.

E. F. Y. Hom, F. Marchis, T. K. Lee, S. Haase, D. A. Agard, and J. W. Sedat, “AIDA: an adaptive image deconvolution algorithm with application to multi-frame and three-dimensional data,” J. Opt. Soc. Am. A 24(6), 1580–1600 (2007). [CrossRef]

17.

I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]

18.

L. De Lathauwer, B. De Moor, and J. Vandewalle, “A multilinear singular value decomposition,” SIAM J. Matrix Anal. Appl. 21(4), 1253–1278 (2000). [CrossRef]

19.

L. De Lathauwer, B. De Moor, and J. Vandewalle, “On the best rank-1 and rank-(R1,R2,…,RN) approximation of higher-order tensors,” SIAM J. Matrix Anal. Appl. 21(4), 1324–1342 (2000). [CrossRef]

20.

S. Umeyama, “Blind deconvolution of blurred images by use of ICA,” Electron Commun. Jpn 84(12), 1–9 (2001).

21.

J. G. Daugman, ““Complete Discrete 2-D Gabor Transforms by Neural Networks for Image Analysis and Compression,” IEEE Trans. Acoust. Speech Signal Process. 36(7), 1169–1179 (1988). [CrossRef]

22.

J. Lin and A. Zhang, “Fault feature separation using wavelet-ICA filter,” NDT Int. 38(6), 421–427 (2005). [CrossRef]

23.

M. E. Davies and C. J. James, “Source separation using single channel ICA,” Signal Process. 87(8), 1819–1832 (2007). [CrossRef]

24.

H. G. Ma, Q. B. Jiang, Z. Q. Liu, G. Liu, and Z. Y. Ma, “A novel blind source separation method for single-channel signal,” Signal Process. 90(12), 3232–3241 (2010). [CrossRef]

25.

D. A. Fish, A. M. Brinicombe, E. R. Pike, and J. G. Walker, “Blind deconvolution by means of the Richardson-Lucy algorithm,” J. Opt. Soc. Am. A 12(1), 58–65 (1995). [CrossRef]

26.

D. S. C. Biggs and M. Andrews, “Acceleration of iterative image restoration algorithms,” Appl. Opt. 36(8), 1766–1775 (1997). [CrossRef] [PubMed]

27.

A. Cichocki, R. Zdunek, A. H. Phan, and S. I. Amari, Nonnegative Matrix and Tensor Factorization (John Wiley & Sons, 2009).

28.

A. Cichocki, and S. Amari, Adaptive Blind Signal and Image Processing (John Wiley, New York, 2002).

29.

H. A. L. Kiers, “Towards a standardized notation and terminology in multiway analysis,” J. Chemometr. 14(3), 105–122 (2000). [CrossRef]

30.

L. R. Tucker, “Some mathematical notes on three-mode factor analysis,” Psychometrika 31(3), 279–311 (1966). [CrossRef] [PubMed]

31.

A. Cichocki, and A. H. Phan. Fast Local Algorithms for Large Scale Nonnegative Matrix and Tensor Factorizations. IEICE Trans Fundamentals 2009; E92-A(3): 708–721.

32.

B. W. Bader, and T. G. Kolda, MATLAB Tensor Toolbox version 2.2.http://csmr.ca.sandia.gov/~tkolda/TensorToolbox.

33.

I. Kopriva, Q. Du, H. Szu, and W. Wasylkiwskyj, “Independent Component Analysis Approach to Image Sharpening in the Presence of Atmospheric Turbulence,” Opt. Commun. 233(1-3), 7–14 (2004). [CrossRef]

OCIS Codes
(100.1830) Image processing : Deconvolution
(100.3010) Image processing : Image reconstruction techniques
(100.3190) Image processing : Inverse problems
(100.6640) Image processing : Superresolution
(100.6890) Image processing : Three-dimensional image processing

ToC Category:
Image Processing

History
Original Manuscript: June 21, 2010
Revised Manuscript: July 23, 2010
Manuscript Accepted: July 24, 2010
Published: August 3, 2010

Citation
Ivica Kopriva, "Tensor factorization for model-free space-variant blind deconvolution of the single- and multi-frame multi-spectral image," Opt. Express 18, 17819-17833 (2010)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-18-17-17819


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. R. L. Lagendijk, and J. Biemond, Iterative Identification and Restoration of Images (KAP, 1991).
  2. M. R. Banham and A. K. Katsaggelos, “Digital Image Restoration,” IEEE Signal Process. Mag. 14(2), 24–41 (1997). [CrossRef]
  3. D. Kundur and D. Hatzinakos, “Blind Image Deconvolution,” IEEE Signal Process. Mag. 13(3), 43–64 (1996). [CrossRef]
  4. P. Campisi and K. Egiazarian, eds., Blind Image Deconvolution (CRC Press, Boca Raton, 2007).
  5. I. Kopriva, “Single-frame multichannel blind deconvolution by nonnegative matrix factorization with sparseness constraints,” Opt. Lett. 30(23), 3135–3137 (2005). [CrossRef] [PubMed]
  6. I. Kopriva, “3D tensor factorization approach to single-frame model-free blind-image deconvolution,” Opt. Lett. 34(18), 2835–2837 (2009). [CrossRef] [PubMed]
  7. F. Li, X. Jia, and D. Fraser, “Superresolution Reconstruction of Multispectral Data for Improved Image Classification,” IEEE Geosci. Remote Sens. Lett. 6(4), 689–693 (2009). [CrossRef]
  8. H. Ji and C. Fermüller, “Robust wavelet-based super-resolution reconstruction: theory and algorithm,” IEEE Trans. Pattern Anal. Mach. Intell. 31(4), 649–660 (2009). [CrossRef] [PubMed]
  9. T. J. Schulz, “Multiframe blind deconvolution of astronomical images,” J. Opt. Soc. Am. A 10(5), 1064–1073 (1993). [CrossRef]
  10. R. Fergus, B. Singh, A. Hertzmann, S. T. Roweis, and W. T. Freeman, “Removing Camera Shake from a Single Photograph,” ACM Trans. Graph. 25(3), 787–794 (2006). [CrossRef]
  11. Q. Shan, J. Jia, and A. Agarwala, “High-quality Motion Deblurring from a Single Image,” ACM Trans. Graph. 27(3), 1 (2008). [CrossRef]
  12. S. Cho and S. Lee, “Fast Motion Deblurring,” ACM Trans. Graph. 28(5), 1 (2009). [CrossRef]
  13. J. Miskin and D. J. C. MacKay, “Ensemble Learning for Blind Image Separation and Deconvolution,” in: Advances in Independent Component Analysis (Springer, London, 2000).
  14. M. Šorel and J. Flusser, “Space-variant restoration of images degraded by camera motion blur,” IEEE Trans. Image Process. 17(2), 105–116 (2008). [CrossRef] [PubMed]
  15. J. Bardsley, S. Jefferies, J. Nagy, and R. Plemmons, “A computational method for the restoration of images with an unknown, spatially-varying blur,” Opt. Express 14(5), 1767–1782 (2006). [CrossRef] [PubMed]
  16. E. F. Y. Hom, F. Marchis, T. K. Lee, S. Haase, D. A. Agard, and J. W. Sedat, “AIDA: an adaptive image deconvolution algorithm with application to multi-frame and three-dimensional data,” J. Opt. Soc. Am. A 24(6), 1580–1600 (2007). [CrossRef]
  17. I. Kopriva, “Approach to blind image deconvolution by multiscale subband decomposition and independent component analysis,” J. Opt. Soc. Am. A 24(4), 973–983 (2007). [CrossRef]
  18. L. De Lathauwer, B. De Moor, and J. Vandewalle, “A multilinear singular value decomposition,” SIAM J. Matrix Anal. Appl. 21(4), 1253–1278 (2000). [CrossRef]
  19. L. De Lathauwer, B. De Moor, and J. Vandewalle, “On the best rank-1 and rank-(R1,R2,…,RN) approximation of higher-order tensors,” SIAM J. Matrix Anal. Appl. 21(4), 1324–1342 (2000). [CrossRef]
  20. S. Umeyama, “Blind deconvolution of blurred images by use of ICA,” Electron Commun. Jpn 84(12), 1–9 (2001).
  21. J. G. Daugman, ““Complete Discrete 2-D Gabor Transforms by Neural Networks for Image Analysis and Compression,” IEEE Trans. Acoust. Speech Signal Process. 36(7), 1169–1179 (1988). [CrossRef]
  22. J. Lin and A. Zhang, “Fault feature separation using wavelet-ICA filter,” NDT Int. 38(6), 421–427 (2005). [CrossRef]
  23. M. E. Davies and C. J. James, “Source separation using single channel ICA,” Signal Process. 87(8), 1819–1832 (2007). [CrossRef]
  24. H. G. Ma, Q. B. Jiang, Z. Q. Liu, G. Liu, and Z. Y. Ma, “A novel blind source separation method for single-channel signal,” Signal Process. 90(12), 3232–3241 (2010). [CrossRef]
  25. D. A. Fish, A. M. Brinicombe, E. R. Pike, and J. G. Walker, “Blind deconvolution by means of the Richardson-Lucy algorithm,” J. Opt. Soc. Am. A 12(1), 58–65 (1995). [CrossRef]
  26. D. S. C. Biggs and M. Andrews, “Acceleration of iterative image restoration algorithms,” Appl. Opt. 36(8), 1766–1775 (1997). [CrossRef] [PubMed]
  27. A. Cichocki, R. Zdunek, A. H. Phan, and S. I. Amari, Nonnegative Matrix and Tensor Factorization (John Wiley & Sons, 2009).
  28. A. Cichocki and S. Amari, Adaptive Blind Signal and Image Processing (John Wiley, New York, 2002).
  29. H. A. L. Kiers, “Towards a standardized notation and terminology in multiway analysis,” J. Chemometr. 14(3), 105–122 (2000). [CrossRef]
  30. L. R. Tucker, “Some mathematical notes on three-mode factor analysis,” Psychometrika 31(3), 279–311 (1966). [CrossRef] [PubMed]
  31. A. Cichocki, and A. H. Phan. Fast Local Algorithms for Large Scale Nonnegative Matrix and Tensor Factorizations. IEICE Trans Fundamentals 2009; E92-A(3): 708–721.
  32. B. W. Bader and T. G. Kolda, MATLAB Tensor Toolbox version 2.2. http://csmr.ca.sandia.gov/~tkolda/TensorToolbox .
  33. I. Kopriva, Q. Du, H. Szu, and W. Wasylkiwskyj, “Independent Component Analysis Approach to Image Sharpening in the Presence of Atmospheric Turbulence,” Opt. Commun. 233(1-3), 7–14 (2004). [CrossRef]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.


« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited