## Primary and secondary superresolution by data inversion

Optics Express, Vol. 14, Issue 2, pp. 456-473 (2006)

http://dx.doi.org/10.1364/OPEX.14.000456

Acrobat PDF (353 KB)

### Abstract

Superresolution by data inversion is the extrapolation of measured Fourier data to regions outside the measurement bandwidth using postprocessing techniques. Here we characterize superresolution by data inversion for objects with finite support using the twin concepts of primary and secondary superresolution, where primary superresolution is the essentially unbiased portion of the superresolved spectra and secondary superresolution is the remainder. We show that this partition of superresolution into primary and secondary components can be used to explain why some researchers believe that meaningful superresolution is achievable with realistic signal-to-noise ratios, and other researchers do not.

© 2006 Optical Society of America

## 1. Introduction

## 2. Forward and inverse models

_{o},x

_{o}], the system PSF h(x), and the additive noise n(x), is given by

_{o}(in cycles/unit length) and n(x) is stationary zero-mean Gaussian noise with variance σ [6

6. C. K. Rushforth and R. W. Harris, “Restoration, resolution, and noise,” J. Opt. Soc. Am. **58**, 539–545 (1968). [CrossRef]

_{o}f

_{o}[8]. These eigenfunctions and their associated eigenvalues can be used to rewrite Eq. (1) in a form that permits solving for o(x) in terms of the measured data i(x). Although o(x) can be estimated from ‘raw’ values of i(x), lower-noise estimates of o(x) can be obtained if only the portion of i(x) associated with data inside the measurement bandwidth is included in the estimation process [11

11. M. Bertero and E. R. Pike, “Resolution in diffraction-limited imaging, a singular value analysis I. The case of coherent illumination,” Opt. Acta **29**, 727–746 (1982). [CrossRef]

_{e}(f) is the estimated (and superresolved) Fourier spectrum. By restricting the integration limits in Eq. (2) to ±f

_{o}, only the signal-bearing portion of I(f) is included in the estimate of O(f). For this reason, the estimator in Eq. (2) has the same noise properties as the image-domain estimator proposed in [11

11. M. Bertero and E. R. Pike, “Resolution in diffraction-limited imaging, a singular value analysis I. The case of coherent illumination,” Opt. Acta **29**, 727–746 (1982). [CrossRef]

_{m}(x)} of Eq. (2) are scaled versions of the eigenfunctions of the forward model and their associated eigenvalues {λ

_{m}} are their energies inside the measurement bandwidth [9]. Although, strictly speaking, the PSWFs in Eq. (2) are defined only on [-f

_{o},f

_{o}], they are analytic functions and can be extended to (-∞,∞). In the following analysis, this extension is assumed to have occurred. Because ϕ

_{m}(f) is non-zero almost everywhere for each value of m, M in Eq. (2) must be set to infinity to obtain an unbiased estimate of O

_{e}(f) at any frequency, in general. Although any practical application can include only a finite number of terms, analyzing the infinite-summation case is important because it provides insight into the SNR properties of unbiased estimates. From Eq. (2), it is easily seen that the expected value of O

_{e}(f) for M=∞ is just the true object spectrum. To obtain the SNRs of O

_{e}(f) at any frequency, an expression for the variances of O

_{e}(f) are needed and can be obtained in a straightforward manner from Eq. (2), yielding

_{e}(f) depend on the properties of the PSWFs and their eigenvalues. For these reasons, these properties are explored next.

*superresolving PSWF*for this reason. To gain insight into additional properties of PSWFs, it is useful to look at a representative set of PSWFs as a function of frequency and index value. Using the method of Latham and Tilton [13

13. W. P. Latham and M. L. Tilton, “Calculation of prolate functions for optical analysis,” Appl. Opt. **26**, 2653–2658 (1987). [CrossRef] [PubMed]

_{m}≈10

^{-30}.

_{m}in frequency space monotonically increases as a function of m. This implies that the energies of the superresolving PSWFs become increasingly concentrated at higher and higher frequencies. Recall, however, that the energy of any PSWF inside the measurement bandwidth is equal to its eigenvalue; thus, from Fig. 2, it is clear that every PSWF has non-zero energy inside the measurement bandwidth. As a result, noise-free knowledge of the true object’s Fourier spectrum inside the measurement bandwidth permits superresolving this knowledge to all frequencies.

_{e}(f) for all frequencies when M = ∞, recall that the expected value of O

_{e}(f) is just the true object spectrum. Therefore, for the SNRs to remain finite, var{O

_{e}(f)} must remain finite. From. Eq. (3), it can be seen that the question of whether or not var{O

_{e}(f)} is finite for any frequency is determined by the index dependence of the ratio ϕ

^{2}

_{m}(f)/λ

_{m}. From Fig. 1, it can be seen that ϕ

^{2}

_{m}(f)→0 as m → ∞ for any frequency f, and from Fig. 2 that λ

_{m}→ 0 as m → ∞. As a result, their relative rates of convergence to zero determine whether or not var{O

_{e}(f)} remains finite. Because there does not exist an analytical expression for ϕ

_{m}(f), Eq. (3) cannot be used to definitively answer this question. However, we can use the PSWFs and eigenvalues calculated for Figs. 1 and 2 to calculate var{O

_{e}(f)} for values of λ

_{m}down to ≈10

^{-30}to see if it appears to be converging to a finite number for any frequency. A movie of var{O

_{e}(f)} as a function of the upper limit of the summation in Eq.(3) is shown in Fig. 3. A key property seen in this movie is that var{O

_{e}(f)} does not seem to be converging to a finite value for any frequency, even those inside the measurement bandwidth, although the rate of increase of var{O

_{e}(f)} inside the measurement bandwidth is significantly less than outside. Because of this fact, it appears that the estimator given by Eq.(2) can provide an unbiased estimate of the object spectrum at any frequency only at the cost of an infinite variance.

## 3. Cramér-Rao bounds

_{o},x

_{o}] for some arbitrary number of locations. Let

**α**be a vector containing these locations and let

**y**,

**θ**and

**η**be vectors that contain the values of i(

**α**), o(

**α**), and n(

**α**), respectively. In addition, let

**H**be the matrix associated with h(

**α**) [15]. This permits rewriting Eq. (1) as a matrix-vector equation given by

**y**conditioned on

**θ**is given by

_{η}(

**η**-

**Hθ**) is the PDF of

**η**but with mean

**Hθ**.

**F**(

**θ**) corresponding to

**θ**. The element of

**F**(

**θ**) in the p

^{th}row and the q

^{th}column is given by

**α**rather than in terms of the elements of

**H**to enhance the reader’s understanding of how we calculate

**F**(

**θ**). Notice that each element of

**F**(

**θ**) is the inverse of the noise variance multiplying a discrete approximation to the integral of the product of two shifted versions of the PSF. This discrete integration approximation must be accurate to greater than one part in 10

^{30}to generate results comparable in accuracy to the noise variance calculations in Section 2. One approach to obtaining this accuracy is to densely sample the support of the object; however, that leads to large dimensions for

**F**(

**θ**) and prohibitively long calculations. A second approach, and the one we implemented, is to carry out the full integration implied by the summation in Eq. (7). We used a quad-precision Romberg integration routine based upon the Numerical Recipes’ qromb routine [16].

**F**(

**θ**)

^{-1}. Prior to the calculation, the eigenvalues of

**F**(

**θ**) were obtained to determine the stability of the inverse. A set of eigenvalues for the same set of parameters as the results in Section 2 (SBP of 40, vector size of 1024 elements, and an object support size of 161 elements) is plotted in Fig. 4 as a function of the FIM eigenvalue index. From Fig. 4, it can be seen that the FIM eigenvalue behavior as a function of the FIM eigenvalue index is quite similar to the PSWF eigenvalue behavior as a function of the PSWF index; specifically, the eigenvalues are close to one for values of the FIM eigenvalue index ≳ SBP and decay at approximately the same rate as do the PSWF eigenvalues for larger values of the index. This similarity can be explained by looking again at Eq. (7). It can be seen that this equation is a discrete approximation to the convolution of a sinc function with itself, which results in a sinc function. Since each row of the FIM is the previous row circularly-shifted by one pixel, the FIM is just a discrete version of the convolution operator in Eq. (1). Thus the eigenvectors and eigenvalues of the FIM are just PSWFs and their eigenvalues. This implies that, even though the FIM is invertible for any (finite) discretization of the interval [-x

_{o},x

_{o}], the condition number of the FIM approaches infinity as the number of elements in

**α**approaches infinity. Thus unbiased estimates of o(x) with finite variances do not exist. However, this does not necessarily imply that unbiased estimators of O(f) do not exist, at least for some set of frequencies. The existence of unbiased estimators of O(f) depends upon the properties of the transformation of

**F**(

**θ**)

^{-1}into the inverse of the FIM corresponding to O(f). This is discussed next.

**F**(

**θ**)

^{-1}into the inverse of the FIM corresponding to O(f) for any desired set of frequencies and extract the CRBs that are along its diagonal. This transformation is achieved by premultiplying

**F**(

**θ**)

^{-1}by

**G**, the Jacobian of the Fourier transform operation for the desired set of frequencies, and postmultiplying by

**G**

^{T}. Because the eigenvalues of

**F**(

**θ**) converge to zero, and because computer calculations have finite precisions, the matrix product

**GF**

_{n}(

**θ**)

^{†}

**G**

^{T}must be computed in a way to clearly illustrate at what frequencies, if any, unbiased estimates of O(f) can be obtained with finite variances. Our approach is to calculate the sequence of matrix products {

**GF**

_{n}(

**θ**)

^{†}

**G**

^{T}} as a function of n=1,2,…, where

**F**

_{n}(

**θ**)

^{†}is the pseudoinverse of

**F**(

**θ**) associated with the n largest eigenvalues of

**F**(

**θ**), and to examine the associated CRBs as a function of n. We believe that, if the sequence of CRBs at a given frequency are a non-converging and increasing function of n, it is plausible to assume that they approach infinity as n approaches infinity. Because of the finite precision of computer calculations, our results are limited to approximately thirty orders of magnitude; however, for all practical purposes, the CRBs are infinite at this point anyway.

**F**

_{n}(

**θ**)

^{†}are displayed in the movie in Fig. 5 for the same set of parameters as used for generating Fig. 4. Each frame of the movie corresponds to adding one more eigenvalue to the calculation of

**F**

_{n}(

**θ**)

^{†}than in the previous frame. Two important observations can be made from this movie. The first is that none of the CRBs appear to be converging to a finite value. Clearly, the CRBs outside the measurement bandwidth are dramatically increasing in value as a function of n. Even the CRBs inside the measurement bandwidth do not appear to be converging. This slow increase in the CRBs inside the measurement bandwidth is undoubtedly due to the “leakage” of the large CRB values outside the measurement bandwidth into the measurement bandwidth brought about by the Fourier-domain correlations enforced by the support constraint [17

17. C. L. Matson, “Variance reduction in Fourier spectra and their corresponding images with the use of support constraints,” J. Opt. Soc. Am. A **11**, 97–106 (1994). [CrossRef]

## 4. Primary superresolution

_{e}(f) at that frequency. In addition, the movie of the PSWFs in Fig. 1 shows that the energies of superresolving PSWFs are negligible in a connected region outside of but adjacent to the measurement bandwidth, and that the size of this region is an increasing function of the PSWF index. The combination of these two facts implies that the bias at a given superresolved frequency is a function of the energies at that frequency contained in PSWFs not included in the summation in Eq. (2) and that this bias is a decreasing function of the summation upper limit. At a given superresolved frequency, if the energies of the PSWFs not included in the estimator are negligible, the superresolved spectrum at this frequency can be said to be essentially unbiased.

19. P. J. Sementilli, B. R. Hunt, and M. S. Nadar, “Analysis of the limit to superresolution in coherent imaging,” J. Opt. Soc. Am. A **10**, 2265–2276 (1993). [CrossRef]

20. C. L. Matson, “Fourier spectrum extrapolation and enhancement using support constraints,” IEEE Trans. Sig. Process. **42**, 156–163 (1994). [CrossRef]

*added*to the reconstructed object. If, however, superresolution is viewed in terms of the

*percentage increase*in the number of degrees of freedom, as is done in [1], the increase in resolution for poorly-resolved objects is much greater than for well-resolved objects because poorly-resolved objects start out with very little information content. We argue that the best way to view primary superresolution is in terms of adding degrees of freedom, rather than in terms of percentage increase in the number of degrees of freedom, because primary superresolution is an inherently additive phenomenon. In the next section, where secondary superresolution is explored, it will be shown that an additive degrees-of-freedom approach is also an appropriate way to view superresolution in the interior of the support region; however, the definition of degree of freedom changes. At the edges of the support region, secondary superresolution will be shown to be a function of the SBP and thus not additive.

## 5. Secondary superresolution

_{o}f

_{s}, where f

_{s}is the largest frequency included in the image reconstruction. Because f

_{s}=∞ for any reconstruction where superresolving PSWFs are included, the classical degree-of-freedom approach states that there also are an infinite number of degrees of freedom in the reconstruction. Clearly this is in error. Thus the classical degree-of-freedom model is not valid to characterize the resolution properties of secondary superresolution. An alternate model must be used.

^{th}PSWF can be characterized by the average distance between these zeros [1], [6

6. C. K. Rushforth and R. W. Harris, “Restoration, resolution, and noise,” J. Opt. Soc. Am. **58**, 539–545 (1968). [CrossRef]

^{-4}(that is, if the data SNR ≈ 100), the percent increase in the resolution throughout the entire support region due to primary superresolution would be 2.5% because primary superresolution only provides a single additional degree of freedom. The resolution increase including both primary superresolution and secondary superresolution in the center region of the support would be 5%, a factor of two increase over what primary superresolution provides, but still small. Finally, at the edges of the support region, the resolution increase including both primary and secondary superresolution would be 30%, a large and noticeable increase. To demonstrate visually the impact of these numbers, three noise-free reconstructions of the triple star whose Fourier amplitude spectrum is plotted in Fig. 8 are displayed in Fig. 13. The first reconstruction uses just the Fourier data inside the measurement bandwidth, the second reconstruction uses both the measured Fourier data and the primary component of the superresolved spectrum, and the third reconstruction uses all the superresolved spectrum along with the measured Fourier data. The two rightmost components of the triple star are contained inside the center 75% of the support; as a result, the amount of superresolution possible is predicted to be ≤ 5% for all three reconstructions and thus essentially negligible. The reconstructions in Fig. 13 confirm this prediction. However, the leftmost component of the triple star is near the edge of the support region and benefits from the increased resolution brought about by secondary superresolution. As predicted by the plots in Fig. 7, the primary superresolution reconstruction of this component differs very little from the no-superresolution reconstruction. However, as predicted by the plots in Fig. 12, the reconstruction that uses all the superresolved spectrum has an increased resolution of ~30%. The increase in resolution is greater on the side of the component near the edge of the support as compared to the side nearer the center of the support because of the decreasing zero spacings of PSWFs as one moves toward the support boundaries. Although the increased resolution is encouraging, it can be seen that the reconstruction using all the superresolved spectrum produces less accurate estimates of the relative magnitudes of the triple star components.

## 6. Discussion and future work

7. J. J. Green and B. R. Hunt, “Improved restoration of space object imagery,” J. Opt. Soc. Am. A **16**, 2859–2865 (1999). [CrossRef]

## Acknowledgments

## References and links

1. | M. Bertero and C. De Mol, “Super-resolution by data inversion,” in |

2. | S. Bhattacharjee and M. K. Sundareshan, “Mathematical extrapolation of image spectrum for constraint-set design and set-theoretic superresolution,” J. Opt. Soc. Am. A |

3. | B. R. Hunt, “Super-resolution of images: algorithms, principles, and performance,” Int. J. Imaging Syst. Technol. |

4. | H. Liu, Y. Yan, Q. Tan, and G. Jin, “Theories for the design of diffractive superresolution elements and limits of optical superresolution,” J. Opt. Soc. Am. A |

5. | V. F. Canales, D. M. de Juana, and M. P. Cagigal, “Superresolution in compensated telescopes,” Opt. Lett. |

6. | C. K. Rushforth and R. W. Harris, “Restoration, resolution, and noise,” J. Opt. Soc. Am. |

7. | J. J. Green and B. R. Hunt, “Improved restoration of space object imagery,” J. Opt. Soc. Am. A |

8. | B. R. Frieden, “Evaluation, design, and extrapolation methods for optical signals based on the use of prolate functions,” in |

9. | D. Slepian and H. O. Pollak, “Prolate spheroidal wave functions, Fourier analysis and uncertainty – I,” Bell Syst. Tech. J. |

10. | H. J. Landau and H. O. Pollak, “Prolate spheroidal wave functions, Fourier analysis and uncertainty – II,” Bell Syst. Tech. J. |

11. | M. Bertero and E. R. Pike, “Resolution in diffraction-limited imaging, a singular value analysis I. The case of coherent illumination,” Opt. Acta |

12. | M. C. Roggemann and B. Welsh, |

13. | W. P. Latham and M. L. Tilton, “Calculation of prolate functions for optical analysis,” Appl. Opt. |

14. | B. Porat, |

15. | R. C. Gonzalez and R. E. Woods, |

16. | W. H. Press, B. P. Flannery, S. A. Teukolsky, and W. T. Vetterling, |

17. | C. L. Matson, “Variance reduction in Fourier spectra and their corresponding images with the use of support constraints,” J. Opt. Soc. Am. A |

18. | H. J. Landau and H. O. Pollack, “Prolate spheroidal wave functions, Fourier analysis and uncertainty – III: the dimension of the space of essentially time-and band-limited signals,” Bell Syst. Tech. J. |

19. | P. J. Sementilli, B. R. Hunt, and M. S. Nadar, “Analysis of the limit to superresolution in coherent imaging,” J. Opt. Soc. Am. A |

20. | C. L. Matson, “Fourier spectrum extrapolation and enhancement using support constraints,” IEEE Trans. Sig. Process. |

21. | Y. L. Kosarev, “On the superresolution limit in signal reconstruction,” Sov. J. Commun. Technol. Electron. |

**OCIS Codes**

(100.2980) Image processing : Image enhancement

(100.3190) Image processing : Inverse problems

(100.6640) Image processing : Superresolution

**ToC Category:**

Focus issue: Signal recovery and synthesis

**Citation**

Charles L. Matson and David W. Tyler, "Primary and secondary superresolution by data inversion," Opt. Express **14**, 456-473 (2006)

http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-14-2-456

Sort: Journal | Reset

### References

- M. Bertero and C. De Mol, "Superresolution by data inversion," in Progress in Optics XXXVI, E. Wolf, ed. (Elsevier, Amsterdam, 1996), 129-178
- S. Bhattacharjee and M. K. Sundareshan, "Mathematical extrapolation of image spectrum for constraint-set design and set-theoretic superresolution," J. Opt. Soc. Am. A 20, 1516-1527 (2003). [CrossRef]
- B. R. Hunt, "Super-resolution of images: algorithms, principles, and performance," Int. J. Imaging Syst. Technol. 6, 297-304 (1995). [CrossRef]
- H. Liu, Y. Yan, Q. Tan, and G. Jin, "Theories for the design of diffractive superresolution elements and limits of optical superresolution," J. Opt. Soc. Am. A 19, 2185-2193 (2002). [CrossRef]
- V. F. Canales, D. M. de Juana, and M. P. Cagigal, "Superresolution in compensated telescopes," Opt. Lett. 29, 935-937 (2004). [CrossRef] [PubMed]
- C. K. Rushforth and R. W. Harris, "Restoration, resolution, and noise," J. Opt. Soc. Am. 58, 539-545 (1968). [CrossRef]
- J. J. Green and B. R. Hunt, "Improved restoration of space object imagery," J. Opt. Soc. Am. A 16, 2859-2865 (1999). [CrossRef]
- B. R. Frieden, "Evaluation, design, and extrapolation methods for optical signals based on the use of prolate functions," in Progress in Optics IX, E. Wolf, ed. (North-Holland, Amsterdam, 1971), 313-407
- D. Slepian and H. O. Pollak, "Prolate spheroidal wave functions, Fourier analysis and uncertainty - I," Bell Syst. Tech. J. 40, 43-63 (1961).
- H. J. Landau and H. O. Pollak, "Prolate spheroidal wave functions, Fourier analysis and uncertainty - II," Bell Syst. Tech. J. 40, 65-84 (1961).
- M. Bertero and E. R. Pike, "Resolution in diffraction-limited imaging, a singular value analysis I. The case of coherent illumination," Opt. Acta 29, 727-746 (1982). [CrossRef]
- M. C. Roggemann and B. Welsh, Imaging Through Turbulence (CRC Press, Boca Raton, 1996), 49.
- W. P. Latham and M. L. Tilton, "Calculation of prolate functions for optical analysis," Appl. Opt. 26, 2653-2658 (1987). [CrossRef] [PubMed]
- B. Porat, Digital Processing of Random Signals, Theory and Methods (Prentice-Hall, Englewood Cliffs, 1994), 65-67.
- R. C. Gonzalez and R. E. Woods, Digital image processing (Addison-Wesley, Reading, 1992), chap. 5.
- W. H. Press, B. P. Flannery, S. A. Teukolsky, and W. T. Vetterling, Numerical Recipes in FORTRAN, 2nd ed.,(Cambridge Press, Cambridge, 1996), 134-135.
- C. L. Matson, "Variance reduction in Fourier spectra and their corresponding images with the use of support constraints," J. Opt. Soc. Am. A 11, 97-106 (1994). [CrossRef]
- P. J. Sementilli, B. R. Hunt, and M. S. Nadar, "Analysis of the limit to superresolution in coherent imaging," J. Opt. Soc. Am. A 10, 2265-2276 (1993). [CrossRef]
- C. L. Matson, "Fourier spectrum extrapolation and enhancement using support constraints," IEEE Trans. Signal Process. 42, 156-163 (1994). [CrossRef]
- Y. L. Kosarev, "On the superresolution limit in signal reconstruction," Sov. J. Commun. Technol. Electron. 35, 90-108 (1990).

## 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.

### Figures

Fig. 1. |
Fig. 2. |
Fig. 3. |

Fig. 4. |
Fig. 5. |
Fig. 6. |

Fig. 7. |
Fig. 8. |
Fig. 9. |

Fig. 10. |
Fig. 11. |
Fig. 12. |

Fig. 13. |
||

### Supplementary Material

» Media 1: MOV (191 KB)

» Media 2: MOV (87 KB)

» Media 3: MOV (89 KB)

» Media 4: MOV (37 KB)

« Previous Article | Next Article »

OSA is a member of CrossRef.