# GENERAL PSF RETRIEVAL AND IMAGE DEBLURRING USING THE CEPSTRAL DOMAIN.

1. INTRODUCTION

The motivation behind this paper was to implement a blind deblurring method that relies only on natural assumptions: the PSF is present on the entire image (limiting the approach to uniform motion blur) and the image contains sufficient non-repetitive information (a small pattern that is present all over the image could be mistaken for the kernel).

The method is based on properties of the cepstrum domain, described in the chapters to come. Other recent advances on the problem at hand are presented in [35][37].

This paper concludes the work done during the 2 years of master studies at the faculty of Automatics and Computers from the "Politehnica" University of Bucharest by the author of [36] and provides a new method for reconstructing the alteration kernel from a corrupted image during the phase of blind image deblurring, thus continuing the work presented in [38].

2. PROPOSED METHOD

2.1. Cepstrum domain

The cepstrum was firstly defined in 1963 by Bogert et al. as:

power _ cepstrum = |[Fourier.sup.-1] [{log[|Fourier{f(t)}|.sup.2]}.sup.2]

Or in a more intuitive form:

Fourier- > abs- > log - > [Fourier.sup.-1]

The domain of application was human speech, which naturally doesn't emit clean sinusoids that can be visualized on a Fourier transform as a clear spike. The vocal chords vibrate not only at a single defined frequency but also generate harmonics (consonant higher frequencies). Inside the vocal tract, some harmonics are amplified and other diminished. If the frontier analysis is to be applied directly to thin signal one would have to choose from the multitude of spikes present at every moment at different frequencies making even the simplest of sounds hard to interpret. These echoes, or harmonics, form a wave-like (sinc) shape in the frequency domain. By applying logarithm and taking the Fourier transform again the waves end up in a single spot. Thus, all the harmonics generate a single impulse on the cepstral graph.

The idea that started research in the PSF estimation of a blurred image in the cepstral domain was that inherently the cepstrum transforms the multiplication from the frequency domain into a summation in the quefrency (cepstral; reversal of letters from 'frequency') domain. This is because applying logarithm on the multiplication gives the sum of logarithms.

[Fourier.sup.-1] {log(F (v)G(v))}(t) = [Fourier.sup.-1] {log(F(v))}(t) + [Fourier.sup.-1] {log(G(v))}(t)

2.2. Power cepstrum observation

After the initial definition, multiple types of cepstrum were developed and found application in other kinds of signal analysis:

* The power considers as input the amplitude of the first Fourier transform

* The real cepstrum receives as input the real part of the Fourier transform

* The complex cepstrum considers the entire output data of the first transform, thus it is fully reversible to the original signal

In our problem, the convolution signal is present on the entire image surface, thus it would be useless to try to find information in the phase of the Fourier transform. Having this in mind and establishing that the phase is forever lost in the computation, the signal must be reconstructed from the remaining elements: the amplitude of the frequencies. This is the first observation.

The second refers to the shapes of the kernels and some important shapes for the Fourier transform: Any motion blur kernel can be created from a filiform shape (movement of the camera) and a Gaussian (defocus). These two functions remain the same after two consecutive Fourier transforms. That is why a slightly modified cepstrum is used: instead of calculating the inverse of the Fourier in the second step, the direct Fourier transform will be applied again.

The importance of this observation is that the shape of the motion blur kernel (which split into the rows of columns is a square convolved with a Gaussian) is preserved in the Cepstrum transform. So this is the theory behind the geometrical cepstrum methods. What wasn't taken into account was the fact that the phase was lost, the recovered signal being a combination of transformed original kernels.

2.3. Baseline method for PSF reconstruction

Based on the observations made earlier, after applying the Power Cepstrum transform, the result is an image that has the shape of the general image but no phase.

The fact that the image is random, from the point of view of the phase, but the kernel is present on the entire surface, this means that the kernel's characteristics will be at most N times stronger than the latent (clear) image's characteristics, because it is present in all N pixels of the image and in the Cepstrum domain it will be over-imposed on the same pixels. This means that the shape of the kernel with the phase lost would be strongly evident.

In (Figure 3, center) the image is mapped almost randomly over the entire surface of the Cepstrum, but a stronger shape steps out of the picture, which resembles the actual PSF. At closer inspection of the kernel autocorrelation, it is evident that it is the same shape, and this is the second observation. Thus, the problems remaining are those of filtering the Cepstrum and resolving the de-autocorrelation.

2.4. De-autocorrelation

The autocorrelation represents the process of correlation of a signal with itself. It would seem a straightforward problem of finding a signal from the autocorrelation, as the autocorrelation is the multiplication with the complex conjugate in the frequency domain:

[R.sub.ff] (T) = [Fourier.sup.-1] (F(t) x [bar.F] (t))

In this situation, the phase becomes 0 so this is the same problem. There is an infinite number of functions that autocorrelated generate the same signal!

In conclusion, our problem becomes that of offering enough extra information in order to limit the number of possible functions. The multiple projection method [26] is used for achieving this objective.

2.5. The multiple projection method

The basic idea is this: the input is an image that lost its phase in the frequency domain but must satisfy some constraints in the spatial domain and the amplitude in the frequency domain is known. This is the same problem as generating a 3D image from multiple2D images.

So the base algorithm is this:

When a change is made in the spatial domain, a phase is generated and the amplitude is distorted. Going back to the frequency domain and changing the amplitude, the spatial image is distorted. Iteratively the image will get to the point of satisfying both conditions.

The problem is to find enough restrictions in order to limit the infinite number of functions that autocorrelated generate the observed signal.

The constraints used in the program include: the image must be real, positive, have values smaller than 255, centered in a square the dimension of the estimated kernel; all of these with an epsilon tolerance, directly proportional to the level of noise in the estimated autocorrelation. The tolerance is needed because with the noise the exact restrictions cannot be met.

2.6 Filtration techniques

As observed in (Figure 3, center), the deduced kernel has a very small signal to noise ratio rendering it impossible to use the already too unstable projection method.

2.6.1 Multiple patches

By instinct, we are tempted to use only the patches of the image in which the PSF is evident.

In order to find the best patches, an estimated PSF is correlated with the image, thus generating high-intensity spots where the PSF is the most evident.

In the frequency domain, the correlation is the multiplication with the conjugate, as seen earlier. The problem is that a light part of the image ends up being white on the correlation result, distorting the data.

The normalized cross-correlation removes the image's influence:

[mathematical expression not reproducible]

Where n is the number of pixels in the kernel and image, F is the sub-image under the kernel and K the kernel, [sigma] the standard deviation and [bar.F], [bar.K] the means. In other words, for every over imposing of the kernel on the image, calculate the standard deviations of the kernel and region from the image; their means; subtract the means from the regions and divide by the deviations and normalize by the number of pixels. Put the result in the current pixel.

As the correlation has great complexity, it cannot be used directly. The implementation idea comes from Industrial Light and Magic: [27]. The mean and the deviations of the image are calculated by a progressive sum and the actual correlation is done in the frequency domain.

As all colors have the same movement applied, 3 times more data can be used by combining the colors along with the best image patches.

The selection algorithm works by choosing the most intense points from the cross-correlation, in order, but doesn't choose overlapping windows in order not to influence an object in the image multiple times.

The combination of the patches is done using the difference between them. The PSF's characteristics must be present in all of the images because only the best patches were selected. So the more different an image is the less influence it will have in the combination because it means that it contains more image information than more kernel information.

2.6.2 Margins removal

Another problem that is evident is the cross in the middle. It is generated by the abrupt margins of the patches. A solution could be to blur the margins, but this would just generate another type of artifact because it would be present on every patch.

The idea used here is to apply logarithm on the image in order to get the lines uniform (because a square generates a sinc in the frequency domain. Applying logarithm to a sinc generates a more uniform function, thus the edges will have almost equal intensities from the center to the margins), and delete the margins from the entire cepstrum. Deleting is done by actually removing the values from the margins from the entire cepstrum patch. This implies that the estimation of the kernel must be 2px larger than the kernel itself.

2.6.3 Noise removal

The last step is removing the remaining noise. Because it isn't a frequency dependent noise, but rather random, it can't be done in Fourier space. That is noise removal will be done using the spatial domain.

Noise reduction starts with computing the Laplace pyramid. Noise only regions are deduced by automatically selecting the most uniform patch from the image. All the frequencies present in the Laplace pyramid that are smaller than those deduced in the uniform region will be removed.

In order to increase the influence of filiform kernels, a derivative of the input image can be applied before inserting it in the application.

Now the kernel's autocorrelation is highly evident and De-autocorrelation can begin.

Because the autocorrelation is symmetrical, the deduced kernel can possibly have reversed Orientation.

Finally, by means of inverse filtering, with the obtained kernel on the input image, the clear image is restored.

3. RESULTS

4. CONCLUSIONS

The presented method is a new approach for estimating the Point Spread Function of an image corrupted by convolution, a method that doesn't rely on opaque assumptions to generate the result. It is based on the work presented in detail in [36]. The different and lite assumptions made in this approach make it useful in combining the results with results from other existing approaches in voting-like systems [29][30][31] in order to obtain a better output than any individual one.

The method proves to be robust for strong motion blur and can also be used to correct small defocus aberrations. What was presented is only the base algorithm, with minor deduction improvements. As future work, all the correction steps mentioned in the related work chapter can be applied, along with some improvements in the existing implementation:

* Remove the waves by generating a gradient border around the image by utilizing the algorithm presented in the ringing artifact removal

* Use one of the better multiple projection algorithms, as this is the simplest and has the possibility to fall in a local minimum

* Add support for non-uniform motion deblurring by splitting the image into layers and selecting patches from the current layer

* Add a preprocessing phase for eliminating noise and superimposing it on the result (in order not to lose the look and feel of the image or details that were wrongly classified as noise)

* Add gamma correction and outliers handling

ACKNOWLEDGEMENT

This work was supported by a grant of the Romanian Ministry of Research and Innovation, CCCDI - UEFISCDI, project number PN-ffl-P1-1.2-PCCDI-2017-0689/,,Lib2Life- Revitalizarea bibliotecilor si a patrimoniului cultural prin tehnologii avansate" / "Revitalizing Libraries and Cultural Heritage through Advanced Technologies", within PNCDI III.

REFERENCES

[1] M. Ben-Ezra, S.K. Nayar, "Motion-Based Motion Deblurring, Pattern Analysis And Machine Intelligence", IEEE Transactions, Vol. 26(6), June 2004, pp. 689-698

[2] J.F. Caia, H. Jib, C. Liua, Z. Shenb, "Blind Motion Deblurring Using Multiple Images", Journal Of Computational Physics, Vol. 228(14), Aug. 2009, pp. 5057-5071

[3] L. Yuan, J. Sun, L. Quan, H.Y. Shum, "Image Deblurring With Blurred/Noisy Image Pairs", ACM Transactions On Graphics (TOG), Vol. 26(3), July 2007

[4] S. Dai, Y. Wu, "Motion from Blur", Computer Vision and Pattern Recognition, 2008. CVPR 2008, 23-28 June 2008, pp. 1-8

[5] R. Fergus, B. Singh, A. Hertzmann, S.T. Roweis, W.T. Freeman, "Removing Camera Shake From A Single Photograph", ACM Transactions On Graphics (TOG), Vol. 25(3), July 2006, pp. 787-794

[6] L. Xu, J. Jia, "Two-Phase Kernel Estimation For Robust Motion Deblurring", ECCV'10 Proceedings Of The 11th European Conference On Computer Vision: Part I, pp. 157-170

[7] Q. Shan, J. Jia, A. Agarwala, "High-Quality Motion Deblurring From A Single Image", ACM Transactions On Graphics (TOG), Vol. 27(3), Aug. 2008

[8] L. Yuan, J. Sun, L. Quan, H.Y. Shum, "Progressive Inter-Scale And Intra-Scale Non-Blind Image Deconvolution", ACM Transactions On Graphics (TOG), Vol. 27(3), August 2008

[9] L. Zouy, H. Zhouz, S. Chengx, C. Heyy, "Dual Range Deringing For Non-Blind Image Deconvolution", Image Processing (ICIP), 26-29 Sept. 2010, pp. 1701-1704

[10] S. Cho, J. Wang, S. Lee, "Handling Outliers In Non-Blind Image Deconvolution", Computer Vision (ICCV), 6-13 Nov. 2011, pp. 495-502

[11] J.H. Lee, Y.S. Ho, "Non-Blind Image Deconvolution With Adaptive Regularization", PCM'10 Proceedings Of The 11th Pacific Rim Conference On Advances In Multimedia Information Processing: Part I, pp. 719-730

[12] A. Jalobeanu, L. Blanc-Feraud, J. Zerubia, "Satellite Image Deconvolution Using Complex Wavelet Packets", Image Processing, 2000, Vol. 3, pp. 809-812

[13] Q. Shan, W. Xiong, J. Jia, "Rotational Motion Deblurring Of A Rigid Object From A Single Image", Computer Vision (ICCV) 2007, pp. 1-8

[14] A. Levin, A. Rav-Acha, D. Lischinski, "Spectral Matting, Computer Vision and Pattern Recognition", CVPR '07. IEEE Conference, June 2007, pp. 1-8

[15] S. Cho, Y. Matsushita, S. Lee, "Removing Non-Uniform Motion Blur From Images", Computer Vision (ICCV) 2007, 14-21 Oct. 2007 pp. 1-8

[16] B. Wohlberg, P. Rodriguez, "An L1-Tv Algorithm For Deconvolution With Salt And Pepper Noise", ICASSP '09 Proceedings Of The 2009 IEEE International Conference On Acoustics, Speech And Signal Processing, pp. 1257-1260

[17] A. Levin, R. Fergus, F. Durand, W.T. Freeman, "Deconvolution Using Natural Image Priors", ACM Trans. Graphics, Vol. 26, No. 3. (2007)

[18] M. Irani, S. Peleg, "Super Resolution From Image Sequences", 10th ICPR, Vol. 2, Jun 1990, pp.115-120

[19] A. Agrawal, R. Raskar, "Resolving Objects At Higher Resolution from A Single Motion-Blurred Image", Computer Vision and Pattern Recognition, 2007. CVPR '07. 17-22 June 2007, pp. 1-8

[20] A. Levin, R. Fergus, F. Durand, W.T. Freeman, "Image And Depth From A Conventional Camera With A Coded Aperture", ACM Transactions On Graphics (TOG), Vol. 26(3), July 2007

[21] Http://En.Wikipedia.Org/Wiki/Convolution, Accessed 20.10.2018

[22] Http://En.Wikipedia.Org/Wiki/Richardson%E2%80%93lucy_Deconvolution, Accessed 20.10.2018

[23] Http://Www.Owlnet.Rice.Edu/~Elec539/Projects99/Bach/Proj2/Inverse.Html, Accessed 20.10.2018

[24] Y. Oyamada, H. Asa, H. Saito, "Blind Deconvolution For A Curved Motion Based On Cepstral Analysis", IPSI Transactions On Computer Vision And Applications Vol. 3 (2011) pp. 32-43

[25] O. Whyte, J. Sivic, A. Zisserman, J. Ponce, "Non-Uniform Deblurring For Shaken Images", Computer Vision And Pattern Recognition (CVPR), 2010, pp. 491-498

[26] E. Osherovich, M. Zibulevsky, I. Yavneh, Signal Reconstruction From The Modulus Of Its Fourier Transform, 24/12/2008 Technion

[27] J. P. Lewis, Fast Normalized Cross-Correlation, Industrial Light & Magic

[28] Http://Www.Owlnet.Rice.Edu/~Elec539/Projects99/Bach/Proj2/Wiener.Html, Accessed 20.10.2018

[29] C.A. Boiangiu, P. Boglis, G. Simion, R. Ioanitescu, "Voting-Based Layout Analysis", Journal Of Information Systems & Operations Management, Vol. 8(1), 2014, pp. 39-47

[30] C.A. Boiangiu, R. Ioanitescu, "Voting-Based Image Segmentation", Journal Of Information Systems & Operations Management, Vol. 7(2), 2013, pp. 211-220

[31] C.A. Boiangiu, M. Simion, V. Lionte, Z. Mihai, "Voting-Based Image Binarization", Journal Of Information Systems & Operations Management, Vol. 8(2), 2014, pp. 343-351

[32] F. Manaila, C.A. Boiangiu, I. Bucur, "Super Resolution From Multiple Low Resolution Images", Journal Of Information Systems & Operations Management, Vol. 8(2), 2014, pp. 316-322

[33] C.A. Boiangiu, A.V. Stefanescu, "Target Validation And Image Calibration In Scanning Systems", In Proceedings Of The 1st International Conference On Image Processing And Pattern Recognition (IPPR '13), Budapest, 2013, pp. 72-78

[34] M. Pollefeys, M. Vergauwen, F. Verbiest, K. Cornelis, L. Van Gool, "From Image Sequences To 3d Models", Proc. Automatic Extraction Of Man-Made Objects From Aerial And Space Images (III), pp.403-410, 2001

[35] M. Zaharescu, C.A. Boiangiu, "An Investigation Of Image Deblurring Techniques", International Journal Of Mathematical Models And Methods In Applied Sciences, Vol. 8, 2014

[36] M. Zaharescu, "Deblurring-ul Unei Singure Imagini", Master Thesis, Unpublished Work, Bucharest, 2013.

[37] M. Zaharescu, C.A. Boiangiu, "Image Deblurring: Challenges And Solutions," In Proceedings Of The 12th International Conference On Circuits, Systems, Electronics, Control & Signal Processing (CSECS '13), Budapest, 2013, pp. 187-196.

[38] M. Zaharescu, C.A. Boiangiu, I. Bucur, "Blind Image Deblurring Using A Single-Image Source: The State Of The Art", Journal Of Information Systems & Operations Management, Vol. 11(1), 2017, pp. 17-28.

Mihai ZAHARESCU (1*)

Costin-Anton BOIANGIU (2)

(1*) corresponding author, Engineer, "Politehnica" University of Bucharest, 060042 Bucharest, Romania, mihai.zaharescu@cs.pub.ro

(2) Professor PhD Eng., "Politehnica" University of Bucharest, 060042 Bucharest, Romania, costin.boiangiu@cs.pub.ro