Due to the turbulent atmosphere the resolution of conventional astrophotography is limited to ∼1 sec of arc. However, the speckle-masking method can yield diffraction-limited resolution, i.e., 0.03 sec of arc with a 3.6-m telescope. Speckle masking yields true images of general astronomical objects. No point source is required in the isoplanatic field of the object. We present the theory of speckle masking; it makes use of triple correlations and their Fourier counterparts, the bispectra. We show algorithms for the recovery of the object from genuine astronomical bispectra data.
© 1983 Optical Society of America
The resolution of conventional astrophotography with large telescopes is not limited by the diffraction on the finite aperture. The influence of the turbulent atmosphere is worse than aperture diffraction. Typically, the size of a turbulence bubble is 10 cm. Diffraction on these bubbles causes an angular blur of λ/d = 500 nm/10 cm ≈ 1 sec of arc). Therefore the resolution of conventional astrophotography is of the order of 1 sec of arc, which is much lower than the diffraction-limited resolution of large telescopes. For example, the diffraction-limited resolution of a 3.6-m telescope is 0.03 sec of arc for λ = 500 nm. This high angular resolution can be achieved by various interferometric methods– in spite of image degradation by the atmosphere. One of these methods is the speckle-masking method which has recently been applied to astronomical data. The raw data for speckle masking and also for most of the other speckle methods are speckle interferograms. Speckle interferograms are short-exposure photographs recorded in the image plane of large telescopes. The exposure time has to be shorter than ∼0.1 sec in order to freeze the atmosphere. In this case the short-exposure photographs consist of an interferometric fine structure called speckles. This interferometric fine structure contains the high-resolution object information.
In this paper we shall describe a detailed theory and astronomical applications of the speckle-masking method. Speckle masking has the following properties: (a) speckle masking yields true images instead of autocorrelations, (b) diffraction-limited resolution is obtained, (c) image degradation by the atmosphere and by telescope aberrations are overcome, (d) general objects can be measured, (e) no unresolvable star is required in the isoplanatic neighborhood of the object.
The features (b)–(e) are the same as for Labeyrie's speckle interferometry, which stands at the beginning of all recent activities in high-resolution astrophotography. The new aspect of speckle masking is feature (a). Speckle masking goes beyond speckle holography– by way of feature (e).
Originally, speckle masking was developed for high-resolution imaging with large ground-based telescopes. However, speckle masking is also useful for high-resolution imaging with large telescopes in space. In the case of large space telescopes the speckle point spread functions are produced by permanent telescope aberrations., The combination of speckle masking with the recently developed roll deconvolution method, or the aberration plate method is useful for diffraction-limited imaging with large, aberrated space telescopes. The application of aberration plates was also investigated by Dainty. The combination of speckle masking with speckle spectroscopy is useful for reconstructing objective prism spectra with diffraction-limited spatial resolution. Speckle masking can also be applied to speckle interferograms of faint objects. The speckle interferograms of such objects consist of individual photon events. In this case it may be advantageous to convolve each photon-counting speckle interferogram with an Airy disk before the speckle-masking process is applied.
In Sec. II we shall introduce the concept of bispectra of speckle interferograms in order to explain the theory of speckle masking. It will be shown that a generalized transfer function of speckle masking can be defined. We have investigated this transfer function theoretically and experimentally and found that it is greater than zero for all frequencies up to the telescope cutoff frequency.
In Sec. III we shall describe two different methods for reconstructing the diffraction-limited image of the object: (1) a method which is performed in the correlation plane, (2) a recursive method which is performed in the Fourier domain. We shall show applications of both methods to astronomical data recorded with the European Southern Observatory's (ESO) 3.6-m telescope.
II. Theory of Speckle Masking based on the Concept of Bispectra of Speckle Interferograms
In this section we shall give a theoretical explanation of the astronomical speckle-masking method. The theory will be performed mainly in the Fourier domain using the concept of bispectra of speckle interferograms. In a preceding article an intuitive description of speckle masking is given, which is based on triple correlations. Triple correlations and bispectra are Fourier transform pairs (see Appendix A).
In speckle masking a large number of speckle interferograms are processed to obtain a diffraction-limited image of the object. The 2-D intensity distribution In(x) of the nth. speckle interferogram can be described by the following incoherent, quasi-monochromatic, and space-invariant imaging equation:
For the explanation of speckle masking it may be useful to compare it with speckle interferometry since most readers are probably familiar with speckle interferometry. The first step of the speckle interferometry process is the determination of the ensemble average power spectrum 〈|Ĩn(u)|2〉 of many speckle interferograms In(x):Appendix A). The brackets 〈…〉 denote ensemble average over many statistically independent speckle interferograms. The function is called the speckle interferometry transfer function, which is positive and nonzero for all frequencies up to the diffraction limit. This speckle interferometry transfer function can be obtained experimentally by observing a single star, or it can be computed based on a reasonable model of the atmosphere. In any case, the speckle interferometry transfer function is smooth, practically constant over a wide range of spatial frequencies. Therefore from Eq. (2) the object power spectrum |Õ(u)|2 can be reconstructed. However, the phase of the object Fourier transform Õ(u) is lost in speckle interferometry. Fourier transforming the object power spectrum yields the object autocorrelation. Therefore we can also say that in speckle interferometry the diffraction-limited object autocorrelation is obtained from the average autocorrelation of many speckle interferograms, i.e.,
In speckle masking the speckle interferograms are evaluated in a different way in order to preserve the phase of the object Fourier transform. In speckle masking the average triple correlation is calculated, i.e., we evaluate the quantityFig. 1). The meaning of Eq. (3) is (i) the speckle interferogram In(x) is multiplied with the same but shifted speckle interferograms In(x + x′) yielding a product mask, (ii) the product mask In(x) · In(x + x′) is cross correlated with In(x), (iii) the ensemble average of the resulting triple correlations is calculated.
The subsequent processing steps in the speckle-masking method are such that the diffraction-limited object triple correlation O(3)(x,x′) is obtained from the average triple correlation of many speckle interferograms, i.e.,
Fortunately the object triple correlation O(3)(x,x′) of speckle masking contains complete information about the object intensity distribution as discussed below and illustrated in Fig. 2.
In Fig. 2 the usefulness of the object triple correlation is illustrated by a simple example. The first step in the illustration is the multiplication of the object O(x) with the shifted object O(x + s). The shift vector x′ = s is called the masking vector. This vector s is selected such that the product mask O(x) · O(x + s) is approximately a δ function. In this case the cross correlation of the product mask with the object yields
For a detailed explanation it seems to be advantageous to express the processing steps by bispectra rather than by triple correlations . Both descriptions are equivalent since triple correlations and bispectra are Fourier transform pairs (see Appendix A).
In Appendix A it is shown that the bispectrum can be calculated as the triple productEq. (5) and Appendix A the bispectrum is separable in u, v, and (−u − v) functions.
In speckle masking we determine the ensemble average bispectrum which is due to In = O * Pn and due to the convolution theorem (see Appendix A) given by
Equation (6) shows that a generalized speckle-masking transfer function can be defined. We have investigated this transfer function theoretically (see Appendix B) and experimentally [see Fig. 3(a)] and found that is real and larger than zero up to the telescope cutoff frequency. Due to the reality of the phase of the complex bispectrum of the object is identical to the phase of the average bispectrum of the object speckle interferograms, i.e.,Appendix B.
In the remaining part of this section we will show some experimental results (Fig. 3), which illustrate Eq. (6) and especially the generalized speckle-masking transfer function . In this experiment we have evaluated speckle interferograms recorded with the ESO 3.6-m telescope. To avoid a 4-D transfer function, we have produced 1-D projections of the 2-D speckle interferograms. Then the generalized speckle-masking transfer function is a 2-D plane of the 4-D transfer function , since a projection in the correlation domain corresponds to a subplane in the Fourier inverse domain, as well known from computed tomography.
Figure 3(a) shows a transfer function evaluated from 1000 projected speckle interferogramss of a point source. The characteristic shape of demonstrates the inherent symmetry of bispectra (see Appendix A). Figure 3(b) shows the average bispectrum of 1-D projections of 300 speckle interferograms of the close spectroscopic double-star Omega Leonis. Figure 3(c) is the bispectrum of Omega Leonis calculated by , i.e., by compensating the transfer function . For comparison, Fig. 4(b) shows the bispectrum Õ(3)(u,υ) of a computer simulated double star [Fig. 4(a)].
III. Two Image-Reconstruction Methods
Up to now we have discussed the basic principle of speckle masking. In this section we shall describe the reconstruction of a true image O(x) of the object from the average triple correlation of many speckle interferograms In(x). In Sec. III.A we present a method which is performed in the space domain and in Sec. III.B a method which is performed in the frequency domain.
A. Image Reconstruction in the Space Domain
To obtain the diffraction-limited image O(x) of the object, we have to reconstruct (1) the object triple correlation O(3)(x,x′) from the average triple correlation and (2) the image O(x) from its triple correlation O(3)(x,x′).
Step (1) can be achieved in two different ways:
- (b) By evaluating uncorrelated speckle interferograms of the object, as described in the text below (Gaussian speckle-masking method).
The latter method has the advantage that no astronomical point source has to be measured and that it can be performed in the correlation domain as well as in the Fourier domain. In Appendix B [Eq. (B11)] it is shown that for a Gaussian model of the atmosphere we obtainAppendices A and B.
Equation (9) can be used to calculate the triple correlation O(3)(x,x′ = s) for any masking vector s. If s is selected suitably, a true image of the object is reconstructed. In the case of complicated objects it is useful to choose a set of many different masking vectors in order to improve the signal-to-noise ratio. The information about suitable masking vectors is obtained from the object autocorrelation.
In Fig. 5 an application of the above described Gaussian speckle-masking method (Sec. III.A(b)] to speckle interferograms of the spectroscopic double-star Omega Leonis is shown. Figure 5(a) is one of 300 speckle interferograms which were evaluated. In Fig. 5(b) the diffraction-limited true image of Omega Leonis is shown. The raw data for the experiment were recorded with the ESO 3.6-m telescope in La Silla, Chile. In our institute we digitized the speckle interferograms with 256 × 256 pixels/frame and recentered them with respect to their center of gravity. Recentering improves the compensation of the speckle-masking transfer function. We also found that high-pass filtering of the speckle interferograms is useful in order to compensate the afterglow effect of the image intensifier.
To derive Eq. (9) (see Appendix B) it was assumed that the complex light amplitude coming from an astronomical point source has Gaussian statistics in the pupil plane of the telescope. We found that the Gaussian assumption is a good approximation, although it is well known that a lognormal model is a better description of atmospheric turbulence. We found experimentally that there must be a weighting constant C, which is almost 1, in front of the first term of Eq. (9), if we have nearly Gaussian statistics but not exactly. If we have exactly Gaussian statistics, C is equal to 1. To determine C(≈1) we calculate Eq. (9) using a masking vector s′ which is selected such that O(x) and O(x + s′) do not overlap. For such an s′ we obtain [O(x) · O(x + s′)] ⊗ O(x) = 0 and therefore C can be determined. The processing with the wrong masking vector s′ is a good test of whether the Gaussian assumption is fulfilled. If Gaussian statistics is not a good approximation, the processing will yield an intensity distribution consisting of two peaks, separated by the vector s′. If there is no constant C(≈1) such that O(3)(x,x′ = s′) becomes equal to zero or at least equal to a constant value, Eq. (9) may not be applied to compensate the speckle-masking transfer function. In this case the compensation of the speckle-masking transfer function has to be performed as described in Sec. III.A(a).
B. Recursive Image Reconstruction in the Frequency Domain
In this section we shall discuss an image-reconstruction method which unravels the complex object spectrum Õ(u) from the object bispectrum Õ(3)(u,v). How Õ(3)(u,v) can be obtained has been discussed above. Let us assume that a sampled version of the object bispectrum is available,
In the following the modulus and the phase of the complex object spectrum Õp will be calculated separately. Therefore Õp is split into
In the recursive image-reconstruction method the modulus of the object Fourier transform is obtained as in speckle interferometry. Speckle interferometry data are produced by setting p = 0 (or equivalently by q = 0 or p = − p) in , i.e.,
The next processing step is the reconstruction of the phase of the complex spectrum of the object. To this end Eq. (11) is inserted into Eq. (10). The resulting equation is split into two equations, one concerning the moduli and one concerning phase factors. The phase factors fulfill the equationEq. (13) the phase factors exp[iϕr] can be calculated recursively. It is sufficient to calculate ϕr for positive r since Õr is Hermitian and therefore ϕr = − ϕ−r. Setting q = 1 in Eq. (13) yields Eq. (14) with r = 1. Both ϕ0 and β0,1 are zero due to the reality of O(x) and O(3)(x,x′). At this point we find that ϕ1 remains undeterminable, a fact that will be discussed after writing down the next iteration steps (r = 2,r = 3…):
The recursion given in Eq. (14) uses only the phase information contained in a single line (q = 1) of the object bispectrum . Additional phase information can be obtained by setting q = 2…N in Eq. (13). Therefore each phase ϕr has (r − 1)/2 independent representations if r is odd and r/2 representations if r is even. These different representations of the ϕr can be averaged. We have used Eq. (13) to find the following recursion formula which is very insensitive to noise because of the summation:Appendix A and illustrated in Fig. 6.
The final step in the recursive image-reconstruction method is to combine the phase factors found by Eq. (15) with the modulus of the object Fourier transform (Eq. (12)]. Inverse Fourier transformation of the complex object spectrum yields the true image of the object.
In Fig. 7 we show the result of an astronomical application of the recursive image-reconstruction method. Figure 7 is a 1-D image of the spectroscopic double-star Omega Leonis. It was reconstructed from the bispectrum of Omega Leonis shown in Fig. 3(b). Since 1-D projections of the speckle interferograms were evaluated (see Sec. II), the reconstructed image is also a projection of the object intensity distribution as discussed below. This is no limitation since a double star is a quasi-1-D object. For example, all information about the object is contained in two perpendicular projections.
In the case of general 2-D objects the recursive image-reconstruction method can be extended to two dimensions. The total algorithm consists of the following steps:
- (i) Collect the raw data In(x).
- (ii) Reconstruct the modulus of the object Fourier transform Õ(u) as in speckle interferometry, this procedure includes the compensation of the speckle interferometry transfer function.
- (iii) Reconstruct the phase of Õ(u) from the phase of the 4-D bispectrum , in this processing step no turbulence defects (transfer function) have to be compensated.
- (iv) Combine the modulus and the phase information in order to reconstruct the object.
Only if the calculation of the 4-D bispectrum exceeds the storage capability of the digital equipment, we suggest the following tomography procedure. First, projections Tn(x,αi) of the raw data In(x) are derived. The projection angles αi are spread out uniformly over the range 0 to π. These projections are processed as described above, yielding projections of the object. Finally an inverse projection procedure is applied in order to reconstruct the 2-D object O(x), as known from computed tomography.
We have described:
- (a) A detailed triple correlation or bispectrum theory of image reconstruction by speckle masking.
- (b) A method for reconstructing a true diffraction-limited image from the object bispectrum.
- (c) A calculation of the generalized speckle-masking transfer function. It is shown that, fortunately, the transfer function is greater than zero up to the cutoff frequency of the telescope and that the transfer function is in essence independent of telescope aberrations.
- (d) Applications of the speckle-masking method to astronomical data recorded with the ESO 3.6-m telescope.
Appendix A: Properties of Triple Correlations and Bispectra
In this Appendix we shall list the definitions and some properties of triple correlations and bispectra. Triple correlations and bispectra have already been used in digital and optical signal processing,,– but a list of their properties is not available and may be useful.
The triple cross correlation of three intensity distributions In(x), Im(x), and Ik(x) is defined by
The Fourier transform of the triple cross correlation is defined byEqs. (A1) and (A2) yields
The bispectrum of In(x) is defined byEq. (A3) several symmetry relations for bispectra can be deduced
Now we want to derive a convolution theorem for triple correlations. IfEq. (A4′) we obtain Equations (A6) and (A7) refer to autobispectra and autotriple correlations. A generalization from auto to cross is straightforward.
1. Calculation of the Generalized Speckle-Masking Transfer Function
In this Appendix we shall calculate the generalized transfer function of speckle masking. We use a 1-D description for the sake of simplicity without sacrificing generality.
The transfer function of a space-invariant incoherent image forming system is known since Duffieux to be the autocorrelation of the pupil function :
In astronomical imaging the light from an astronomical point source propagates through the turbulent atmosphere. Therefore the complex amplitude H(u) may be split into a product of two functions, one representing the effect of the random medium [A(u)] and one representing the pupil function of the telescope [H0(u)]:Eq. (B2) into Eq. (B1) yields
From Eq. (B3) the transfer function of ordinary long-exposure imaging can be found by taking the ensemble average . In the case of large telescopes falls off to zero for frequencies much lower than the telescope cutoff frequency, as shown in [Ref. 28], for example.
Image processing in the speckle-masking method is governed by the transfer functionEq. (B3) into Eq. (B4) yields
To calculate the sixth-order moment in Eq. (B5) we assume that the real and the imaginary parts of A(w) have zero-mean Gaussian distribution. For a complex Gaussian process the momentum theorem can be used to show that
Since the complex amplitude A(u) is a stationary random variable its autocorrelation function may be defined by
Since for usual seeing conditions the complex amplitude A(u) in the pupil plane has a very fine structure compared to the telescope aperture, A(u) can be assumed to be δ correlated, i.e.,
This approximation is not necessary for the calculation, but it makes the result easier to discuss. At the end of this section we shall discuss the general case that A(u) is not δ correlated; we shall find very similar results.
Note that the term in Eq. (B8) only depends on the modulus square of the pupil function of the telescope. Therefore it is independent of both, of image degradation due to atmospheric turbulence and telescope aberrations. The term is positive and nonzero for all frequencies up to the diffraction-limit of the telescope. Terms , , , and  in Eq. (B8) contain the long-exposure transfer function as a factor and therefore they are nonzero only close to the axes u = 0, υ = 0, u = − υ and to the point u = υ = 0, respectively. Terms – cause the characteristic structure of the generalized transfer function of speckle masking, as shown in Figs. 3(a) and 8.
Finally we note that the transfer function is a real function,. i.e., it has zero phase. This can be verified by calculating terms – in Eq. (B8). Using the momentum theorem, Eqs. (B3) and (B7′), one can show that terms – can be expressed in terms of the modulus square |H0(u)|2 of the pupil function H0(u).
In the calculation given above we assumed that the complex amplitude A(u) in the pupil plane is δ correlated. If we do not make this assumption we obtain for the transfer function a result which differs from Eq. (B8) only in the term . We find that in Eq. (B8) the function must be substituted by a function which is given by
2. Compensation of the Generalized Speckle-Masking Transfer Function
In this section we shall discuss how the transfer function can be compensated by evaluating uncorrelated speckle interferograms. The compensation can be performed in the correlation domain or in the Fourier domain.Eq. (B9)]. Therefore is an aberration-free diffraction-limited transfer function. Certainly, can be compensated since it is known [see Eq. (B9)] and it is greater than zero up to the telescope cutoff frequency. After this compensation the exact object bispectrum Õ(3)(u,υ) is obtained.
Inverse Fourier transformation of Eq. (B10) yields a compensation of the generalized speckle-masking transfer function which can be performed in the correlation domainFig. 9. For an imaging system with a finite aperture of size a, the width of the function B(3)(x,x′) is approximately λ · f/a, where λ and f denote wavelength and focal length.
1. Triple Correlations and Additive Noise
In this Appendix it is briefly discussed how triple correlations are influenced by additive noise. Let a signal I(x) be contaminated by signal-independent additive noise N(x), i.e.,
2. Additive Noise in the Gaussian Speckle-Masking Method
It is shown in this section that the Gaussian speckle-masking method (Sec. III.A(b)) essentially compensates the influence of additive noise in the speckle interferograms.
If noisy speckle interferograms Jn(x) = In(x) + Nn(x) are processed in the Gaussian speckle-masking method, Eq. (9) yieldsEq. (C3) is equal to zero. We shall consider three special cases:
Finally we consider the case of additive, nonzero-mean noise, which can be described by , where is assumed to be zero mean. In this case we obtain
We thank the German Science Foundation (DFG) for financial support, the European Southern Observatory for observing time, and K.-D. Förster for very helpful discussions.
The work in this paper is based on data collected at the European Southern Observatory, La Silla, Chile.
1. A. Labeyrie, Astron. Astrophys. 6, 85 (1970).
2. C. Y. C. Liu and A. W. Lohmann, Opt. Commun. 8, 372 (1973). [CrossRef]
3. R. H. T. Bates, P. T. Gough, and P. J. Napier, Astron. Astrophys. 22, 319 (1973).
5. C. R. Lynds, S. P. Worden, and J. W. Harvey, Astrophys. J. 207, 174 (1976). [CrossRef]
6. K. T. Knox and B. J. Thompson, Astrophys. J. Lett. 193, L45 (1974). [CrossRef]
7. D. C. Ehn and P. Nisenson, J. Opt. Soc. Am. 65, 1196 (1975).
8. G. P. Weigelt, Opt. Commun. 21, 55 (1977). [CrossRef]
9. R. H. T. Bates, M. O. Milner, G. I. Lund, and A. D. Seager, Opt. Commun. 26, 22 (1978). [CrossRef]
10. K. von der Heide, Astron. Astrophys. 70, 777 (1978).
12. B. J. Brames and J. C. Dainty, J. Opt. Soc. Am. 71, 1542 (1981). [CrossRef]
15. W. J. Cocke, Proc. Soc. Photo-Opt. Instrum. Eng. 231, 11 (1980).
16. C. Roddier, F. Roddier, and J. Vernin, in Proceedings, ESO Conference on Scientific Importance of High Angular Resolution at IR and Optical Wavelengths, M. H. Ulrich and K. Kjär, Eds. (ESO, Garching, Germany, 1981), p. 165.
17. A. M. J. Huiser, Opt. Commun. 42, 226 (1982). [CrossRef]
18. A. H. Greenaway, Opt. Commun. 42, 157 (1982). [CrossRef]
19. R. H. T. Bates and F. M. Cady, Opt. Commun. 32, 365 (1980). [CrossRef]
21. A. W. Lohmann and G. Weigelt, in Proceedings, ESA/ESO Workshop on Astronomical Uses of the Space Telescope,F. Macchetto, F. Pacini, and M. Tarenghi, Eds. (ESO, Geneva, 1979), p. 353.
22. G. Weigelt, Proc. Soc. Photo-Opt. Instrum. Eng. 332, 284 (1982).
23. G. Weigelt, in Proceedings, Conference on Image Processing in Astronomy,G. Sedmak, M. Capaccioli, and R. J. Allen, Eds. (Osservatorio Astronomico di Trieste, Trieste, Italy, 1979), p. 422.
24. J. C. Dainty, Opt. Commun. 7, 129 (1973). [CrossRef]
25. G. Weigelt, in Proceedings, ESO Conference on Scientific Importance of High Angular Resolution at IR and Optical Wavelengths,M. H. Ulrich and K. Kjar, Eds. (ESO, Garching, Germany, 1981), p. 95.
26. T. R. Bader, Proc. Soc. Photo-Opt. Instrum. Eng. 185, 140 (1979).
27. G. Baier and G. Weigelt, Astron. Astrophys. 121, 137 (1983).
28. J. C. Dainty, in Laser Speckle and Related Phenomena,J. C. Dainty, Ed. (Springer, Berlin, 1975).
29. H. H. Barrett and W. Swindell, Radiological Imaging: Theory and Image Formation, Detection and Processing (Academic, New York, 1981), Vols. 1 and 2.
30. K. Sasaki, T. Sato, and Y. Namamura, IEEE Trans. Sonics Ultrasón. SU-24, 193 (1977). [CrossRef]
31. H. Gamo, J. Appl. Phys. 34, 875 (1963). [CrossRef]
33. J. Cohen, Proc. Soc. Photo-Opt. Instrum. Eng. 180, 134 (1979).
34. P. Kellman, Opt. Eng. 19, 370 (1980). [CrossRef]
35. J. W. Goodman, Aust. J. Electr. Electron. Eng. 2, 140 (1982).
36. I. S. Reed, IRE Trans. Inf. Theory IT-8, 194 (1964).