Large amplitude tip/tilt estimation by geometric diversity for multiple-aperture telescopes

TL;DR

The paper introduces ELASTIC, a new noniterative method for large-range tip/tilt estimation in multi-aperture telescopes using only two diversity images, demonstrated through simulations and experimental alignment.

Contribution

It presents ELASTIC, a novel, hardware-compatible, noniterative algorithm for large-range tip/tilt estimation in segmented telescopes, requiring only two images.

Findings

Residues within interferometric phasing range at 5000 photo-electrons/sub-aperture/frame

Method successfully aligned a 6 sub-aperture mirror experimentally

Performance optimized through simulations

Abstract

A novel method nicknamed ELASTIC is proposed for the alignment of multiple-aperture telescopes, in particular segmented telescopes. It only needs the acquisition of two diversity images of an unresolved source, and is based on the computation of a modified, frequency-shifted, cross-spectrum. It provides a polychromatic large range tip/tilt estimation with the existing hardware and an inexpensive noniterative unsupervised algorithm. Its performance is studied and optimized by means of simulations. They show that with 5000 photo-electrons/sub-aperture/frame and 1024x1024 pixel images, residues are within the capture range of interferometric phasing algorithms such as phase diversity. The closed-loop alignment of a 6 sub-aperture mirror provides an experimental demonstration of the effectiveness of the method. Author accepted version. Final version is Copyright 2017 Optical Society of America. One print or electronic copy may be made for personal use only. Systematic reproduction and distribution, duplication of any material in this paper for a fee or for commercial purposes, or modifications of the content of this paper are prohibited.

Tables

Table 1. Table 1: Variation of field coverage with δ 𝛿 \delta , and corresponding maximum aberration coefficient a m ​ a ​ x subscript 𝑎 𝑚 𝑎 𝑥 a_{max} in our simulated case.

Shift amplitude ($\delta$)	$1$	$2$	$3$	$4$
Image field coverage ($\%$)	$100$	$100$	$66,6$	$50$
$a_{max}\left(wave\right)$	$25,6$	$12,8$	$8,5$	$6,4$

Equations

$$\Pi(\boldsymbol{u})=\left\{\begin{array}[]{lll}1&\mbox{for}\,\,|\boldsymbol{u}|\leq 1,\\ 0&\mbox{elsewhere,}\end{array}\right.$$

$${\mathbf{p}}_{n,d}$$

$$\text{withleavevmodenobreak }{\Pi }_{n,d}$$

$${\mathbf{p}}_d=\sum _{n=1}^{{\mathcal{N}}_a}{\mathbf{p}}_{n,d}.$$

$${\mathbf{s}}_d(\mathbf{u})=({\mathbf{p}}_d\otimes {\mathbf{p}}_d)(\mathbf{u})=\int {\mathbf{p}}_d({\mathbf{u}}^{\prime }){\mathbf{p}}_d^{\ast }(\mathbf{u}+{\mathbf{u}}^{\prime })d{\mathbf{u}}^{\prime }.$$

$${\mathbf{s}}_d(\mathbf{u})=\sum _{n=1}^{{\mathcal{N}}_a}\sum _{n^{\prime }=1}^{{\mathcal{N}}_a}\left[({\mathbf{p}}_{n,d}\otimes {\mathbf{p}}_{n^{\prime },d})\star {\delta }_{({\mathbf{c}}_n-{\mathbf{c}}_{n^{\prime }})}\right](\mathbf{u}).$$

$${\mathbf{s}}_d(\mathbf{u})=$$

$$\text{withleavevmodenobreak }{\Lambda }_{n,d}=$$

$$n_s=\frac{\lambda F}{2R{p}_{pix}}.$$

$${\mathbf{s}}_{d,\mathbf{q}}=\sum _{n=1}^{{\mathcal{N}}_a}{\rho }_n^2{\Lambda }_{n,d}\left(\frac{n_s\mathbf{q}}{{\mathcal{N}}_p}\right)\exp \left[j\sum _{k=2}^{3}4\pi a_{kn}{\mathbf{Z}}_k\left(\frac{n_s\mathbf{q}}{2{\mathcal{N}}_p}\right)\right].$$

$${\mathbf{\delta }}_2$$

$${\mathbf{j}}_{m,\mathbf{q}}=$$

$${\mathbf{j}}_{m,\mathbf{q}}=\sum _{n=1}^{{\mathcal{N}}_a}$$

$$\times \exp \left[j\sum _{k=2}^{3}4\pi a_{k,n}({\mathbf{Z}}_k-{\Delta }_m\left[{\mathbf{Z}}_k\right])\left(\frac{n_s\mathbf{q}}{2{\mathcal{N}}_p}\right)\right].$$

$$({\mathbf{Z}}_k-{\Delta }_m\left[{\mathbf{Z}}_k\right])\left(\frac{n_s\mathbf{q}}{2{\mathcal{N}}_p}\right)=-\frac{n_s\delta }{{\mathcal{N}}_p}{\delta }_{k,m},$$

$${\mathbf{j}}_{m,\mathbf{q}}$$

$$\displaystyle\text{with\leavevmode\nobreak\ }\mbox{\huge$\curlywedge$}_{n,m}$$

$${\mathbf{x}}_m$$

$$\text{mbox}with{x}_{m,n}$$

$${\mathbf{j}}_m={\mathrm{C}}_{\mathrm{m}}{\mathbf{x}}_m,$$

$$\widehat{{\mathbf{x}}_m}={\mathrm{C}}_{\mathrm{m}}^{†}{\mathbf{j}}_m.$$

$$\left\{\begin{array}[]{r@{}l}\displaystyle\widehat{a}_{m,n}&=\frac{-\mathcal{N}_{p}}{4\pi\leavevmode\nobreak\ n_{s}\delta}\,\operatorname{Arg}(\widehat{\boldsymbol{x}}_{m,n}),\cr\widehat{\rho}_{m,n}&=|\widehat{\boldsymbol{x}}_{m,n}|^{1/4}.\end{array}\right.$$

$$RMSE=\sqrt{⟨\frac{1}{{\mathcal{N}}_a}\sum _{m,n}({\hat{a}}_{m,n}-a_{m,n}{)}^2⟩},$$

$$a_{field}=\frac{{\mathcal{N}}_p}{8n_s}.$$

$$a_{max}(\delta )=\frac{2a_{field}}{\delta }=\frac{{\mathcal{N}}_p}{4n_s\delta }.$$

$${\hat{a}}_{m,n}\simeq a_{m,n}+\frac{{\mathcal{N}}_p}{4\pi n_s\delta }\mathrm{Im}\left(\frac{{\mathbf{ϵ}}_n}{x_{m,n}}\right),$$

$$Z_1({\mathbf{u}}_n)$$

$$Z_2({\mathbf{u}}_n)$$

$$Z_3({\mathbf{u}}_n)$$

$$Z_4({\mathbf{u}}_n)$$

$$\mathbf{r}={\mathbf{r}}_n+R{\mathbf{c}}_n.$$

$${\mathbf{u}}_n\stackrel{△}{=}\frac{{\mathbf{r}}_n}{R},$$

$${\mathbf{u}}^{\prime }=\frac{2R}{B}({\mathbf{u}}_n+{\mathbf{c}}_n).$$

$$Z_4({\mathbf{u}}^{\prime })=$$

Full text

††journal: josaa

1 Introduction

The resolution of a telescope is ultimately limited by its aperture diameter, but the size of mirrors is bounded by current technology to about 10 m on the ground and to a few meters in space. To overcome this limitation, interferometry consists in making an array of sub-apertures interfere; the resulting instrument is called an interferometer or a multi-aperture telescope.

The sub-apertures can either be telescopes per se, as in current ground-based interferometers [1, 2, 3], or segments of a common primary mirror, such as in the Keck telescopes, the future extremely large telescopes [4, 5, 6, 7]) or large ground collectors (Cherenkov Telescope Array [8]).

So far, this technique has been operational only on ground-based telescopes, but interferometers have long been forecast for high-resolution spaceborne astronomy, for the soon-to-be-launched James Webb Space Telescope (JWST) [9] or other projects [10, 11, 12, 13], and for Earth observation [14, 15].

To reach the diffraction limit, all the sub-apertures of such a telescope must be phased to within a small fraction of the wavelength. A critical sub-system of interferometers is thus the Cophasing Sensor (CS), whose goal is to measure with this kind of precision the relative positioning errors between the sub-apertures (differential piston and tip/tilt), which are the specific low-order aberrations of an interferometer and the main sources of wave-front degradation (Fig. 1).

Focal-plane wavefront sensing is an elegant solution to measure the wavefront degradation. Since the near-focal images of any source taken by a 2D camera show distortions when the telescope is not aligned (Fig. 1), these misalignments can be retrieved by solving the associated inverse problem. The phase retrieval technique, based on the analysis of the sole focal-plane image, is generally not sufficient to retrieve piston and tip/tilt without ambiguity except in specific cases [16]. The phase diversity technique [17, 18], typically based on a focal and a slightly defocused images, removes all ambiguities and operates even on unknown extended sources.

The fact that phase diversity can be used as a CS on a segmented aperture telescope was recognized very early [19], and extensively studied for the JWST [20, 21, 22, 23], for Darwin [24], and for the European Extremely Large Telescope (E-ELT) [25] in particular. Additionally, in contrast with most pupil-plane-based devices, phase diversity enjoys three appealing characteristics: firstly, it is appropriate for a large number of sub-apertures, because the hardware complexity remains essentially independent of the number of sub-apertures. Secondly, the CS is included in the main imaging detector, simplifying the hardware and minimizing differential paths. Thirdly, it can be used on very extended objects. These properties are strong motivations for the choice of phase diversity as a CS, even when looking at an unresolved source.

The measurement and correction described above of piston-tip-tilt by means of a CS to within a small fraction of $\lambda$ is hereafter called the fine phasing mode. This mode assumes that the tip/tilts are smaller than typically $\lambda/8$ (see Eq.A36) and the differential pistons are within the coherence length (cf [16]). However, during the first alignment steps (after integration or deployment), the disturbances are much larger than a few micrometers.

A preliminary geometrical alignment mode is thus mandatory, to efficiently drive the telescope into the reduced capture range of an interferometric CS. This mode consists in an incoherent superimposition of the focal plane images from each sub-aperture. Only then can the relative piston errors between sub-apertures be measured and corrected for the fine phasing to operate. Using only one image, such as the common focal image of an unresolved source, the geometrical alignment is impossible at least when all sub-apertures are identical since their Point Spread Functions (or sub-PSFs) all have the same shape. Even if these sub-PSFs are well separated in the field with clearly identifiable positions, it is impossible to associate them with their respective sub-apertures.

The solution selected by JWST to identify the sub-PSFs is temporal modulation [26], but the measurement time scales with the number of sub-apertures. A method called Geometrical Phase Retrieval (GPR), with the ability to increase the JWST fine algorithm capture range, was developed [27], refined [28], experimentally demonstrated [29] and should be implemented in the JWST [30]. Based on geometrical optics, it uses typically four to six defocused images per segment. The JWST geometrical alignment mode, called image-stacking operation [31], was originally forecast to last around 1 week of commissioning time, but the GPR should provide substantial time savings [29].

In this paper, a novel method for the geometrical alignment of the subapertures is proposed. The alignment error is encoded by each sub-PSF’s shift in the focal plane image. As can be seen on Figure 1, a defocus induces an additional shift of the sub-PSFs in the diversity plane with respect to the focal plane. Each of these additional shifts is given by the defocus, and is different for each sub-aperture. Using these two images jointly, our method extracts the position and identifies the sub-aperture sub-PSF, thanks to a modified cross-correlation. Based on geometrical optics similarly to GPR, it thus only needs two images of an unresolved source and a simple data processing. Additionally, its capture range is only limited by the imaging sensor field of view.

In Section 2, a closed-form model of the multi-aperture Optical Transfer Function (OTF) is derived. Section 3 provides a solution for large amplitude tip/tilt measurement from a modified cross-spectrum of two diversity images. The method is then validated and optimized by simulation in Section 4. Finally, in Section 5, an experimental validation is performed on a segmented mirror.

2 Closed-form model for the multi-aperture OTF with geometric diversity

Throughout this paper, the considered multiple-aperture telescope is based either on a segmented telescope or on a telescope array, all feeding a common focal-plane image detector and operates at a central wavelength $\lambda$ with a focal length $F$. We assume the object is a point source at infinity. In this section, the telescope’s pupil, its disturbances and its OTF are successively introduced.

2.1 Model of the multiple aperture

The pupil of the telescope is assumed for simplicity to be made of a set of ${\mathcal{N}_{a}}$ circular sub-apertures (of index $n$ ranging from 1 to ${\mathcal{N}_{a}}$) with identical radius $R$. In the following we use the reduced coordinate $\boldsymbol{u}$, defined as the usual coordinate $\boldsymbol{r}$ in the pupil plane divided by $R$ ($\boldsymbol{u}=\boldsymbol{r}/R$). Similarly we define $\boldsymbol{c}_{n}$ as the reduced coordinates of the center of each sub-aperture. In each of the two so-called diversity planes (of index $d\in\{1,2\}$), each sub-aperture is characterized by its complex amplitude transmittance $\boldsymbol{p}_{n,d}$ whose shape is described by the disk function $\Pi$,

[TABLE]

and the pupil transmission is constant over each sub-aperture, equal to $\rho_{n}$. The phase of $\boldsymbol{p}_{n,d}$ is the sum of unknown (and seeked) phases plus known diversity phases $\Psi_{n,d}(\boldsymbol{u})$, and is expanded on a local orthonormal basis of Zernike polynomials [32]. A description of the first $4$ modes is given in Appendix 1. The amplitude, in wave unit, of the $k$th mode over the $n$th sub-aperture is defined as $a_{k,n}$. Since the specific aberrations of a multi-aperture telescope are piston ($k=1$), tip and tilt ($k=2$ and $k=3$), the maximal value of $k$ considered here is $3$. Each sub-aperture transmittance in the $d$th diversity plane can then be written as:

[TABLE]

where $\star\delta_{\boldsymbol{c}_{n}}$ denotes the shift operator by vector $\boldsymbol{c}_{n}$. The total pupil transmitance in the $d$th diversity plane is then :

[TABLE]

There are several ways to acquire these diversity images. A first solution is global temporal modulation: a known full-aperture aberration is applied prior to image acquisition; a usual example is the use of the focal image (no aberration) and a defocused image (by a longitudinal shift of the detector or of a mirror, or by lens insertion through a filter wheel [29]). A second solution is spatial modulation: the beam is split in two parts, routed to simultaneously form two distinct images with different aberrations (e.g. focus, astigmatism) on two different (parts of the) detector(s) [33]. A third solution is local temporal modulation (Fig. 2): a global defocus is approximated by its piston/tip/tilt projection on the sub-apertures, or even to the sole tip/tilt part: the latter will be called pseudo-defocus from now on. This pseudo-defocus can be introduced simply, either by moving the actuators intended to correct for the sub-aperture misalignment, or by an array of (achromatic) prisms in the filter wheel. This provides the sub-PSF’s sprawl induced by the global defocus, and allows to exclusively acquire focal plane images.

2.2 Closed-form model of the OTF

We define the reduced OTF $\boldsymbol{s}_{d}$ of the $d^{th}$ diversity plane as the autocorrelation ($\otimes$) of the reduced pupil transmittance $\boldsymbol{p}_{d}$:

[TABLE]

Because the support width is doubled by autocorrelation, $|\boldsymbol{u}|$ ranges between [math] and $2$ in the OTF planes. The OTF with the usual spatial frequency $\boldsymbol{f}$ in $rad^{-1}$ units in the object frequency space is thus $\boldsymbol{s}_{d}(\boldsymbol{f}\lambda/R)$. Merging Eqs. (2) and (5), the global OTF $\boldsymbol{s}_{d}$ in the $d^{th}$ diversity plane is then written as:

[TABLE]

The telescope’s PSF can be obtained by computing the Inverse Fourier Transform of the OTF. Eq. (6) shows that the OTF is the sum of ${\mathcal{N}_{a}}^{2}$ terms: ${\mathcal{N}_{a}}$ so-called photometric terms (when $n=n^{\prime}$, centered at $\boldsymbol{u}=\boldsymbol{0}$), and ${\mathcal{N}_{a}}({\mathcal{N}_{a}}-1)$ interferometric terms (when $n\neq n^{\prime}$, centered at $\boldsymbol{u}_{n,n^{\prime}}=(\boldsymbol{c}_{n}-\boldsymbol{c}_{n^{\prime}})$). The photometric terms are the autocorrelations of each sub-aperture, and the interferometric terms are the correlations of two different sub-apertures. In this paper dealing with large amplitude alignment, the interferometric terms will be discarded for 4 reasons. Firstly, when large tip/tilt errors are present, it is very likely that large piston errors are also present. When their amplitude is larger than the coherence length (a few micrometers for broadband observations), interferometric terms are dimmed by the coherence enveloppe. Secondly, interferometric terms are also dimmed by significant differential tip/tilts. Thirdly, at this stage of operation, there is no need to see the high frequency contents of the image, only the sub-PSF positions are of interest; the detector can thus be rebinned to reduce the number of pixels to process. And lastly, if all these natural filters are not sufficient, an explicit low-pass filter can be used to isolate the superimposed photometric terms (in the center of the Fourier plane) from the fringe peaks.

An important result is that in presence of tip/tilt errors, the photometric terms from each sub-aperture are affected by a phase ramp, with the same slope than their associated pupil transmittance $\boldsymbol{p}_{n,d}$ [16]. This can be simply understood by the fact that the Fourier Transform (FT) of the two gives the focal signal (in intensity or amplitude) at the same location. The Zernike polynomials were previously used over sub-apertures with a unit radius. In order to use the same modes (Appendix 1), Zernike polynomials $2$ and $3$ with a doubled support (but same slope) will be written $2\boldsymbol{Z}_{k}\left(\boldsymbol{u}/2\right)$. Keeping only incoherent terms, the $d^{th}$ diversity plane OTF writes:

[TABLE]

To illustrate the previous derivations, Fig. 3 presents the pupil, the PSF, the OTF modulus and the OTF phase of a 3 sub-aperture telescope. Line 1 is the no-aberration coherent case. All sub-PSFs are superimposed in the focal plane, and the OTF modulus is the sum of 9 terms: 6 interferometric terms (two for each baseline) with a chinese-hat peak shape (hence the symbol $\Lambda=\Pi\otimes\Pi$) and 3 superimposed photometric terms which form a stronger central peak. In line 2 the impact of tip/tilt aberrations on the sub-apertures is highlighted. The sub-PSFs are pulled apart in the focal plane. Although some weak fringes remain present in the image, the OTF interferometric peaks are strongly attenuated, as previously mentioned. Line 3 illustrates the effect of the same previous aberration in an incoherent case, cancelling all interferometric terms. The information in the photometric peak of the OTF is not modified with respect to the coherent case.

In conclusion, the sum of all superimposed photometric peaks in Eq. (7) creates an intricate pattern (Fig. 3 lines 2 and 3), with a different value in each diversity plane, from which the seeked unknown tip/tilts of each sub-aperture must be estimated. Such an estimator is derived in the next section.

3 ELASTIC, a sub-aperture tip-tilt estimator

Before the quantative mathematical derivation in subsection 3.2, the physical origin of this large amplitude tip-tilt estimator is explained in subsection 3.1.

3.1 Principle of the estimator

Figure 1 shows that the focal image is only sensitive to the tip/tilt from each sub-aperture, whereas in the diversity image the seeked tip/tilt is combined with an additional shift resulting from the diversity (here a defocus obtained by longitudinal shift of the detector). If the diversity image is telecentric, then this additional shift does not depend on the sub-aperture tip/tilt, but only on the introduced diversity. The same holds for a temporal diversity induced by the sub-aperture actuactors, such as a pseudo-defocus. If these diversity-induced shifts over all the sub-apertures are sufficiently different, then this differential information between the two diversity images enables the association of the sub-PSFs to their sub-aperture index. The goal is now to convert this qualitative visual process to an unsupervised quantitative measurement, using a simple algorithm that can be operated whatever the number of sub-apertures is or their relative positions are.

A first ingredient is to compute the correlation between the two diversity images to access this differential tip/tilt. Since each diversity image is composed by ${\mathcal{N}_{a}}$ sub-PSFs (Fig. 4 line 1), the image correlation (Fig. 4 line 2) contains ${\mathcal{N}_{a}}^{2}$ terms: ${\mathcal{N}_{a}}$ terms that will be called autospots in the following, obtained by the correlation of the sub-PSF from sub-aperture $n$ in one diversity image by the sub-PSF from the same aperture $n$ but in the other diversity image. These autospots lie in the central circled area in Fig. 4 line 2. Outside this circle are ${\mathcal{N}_{a}}({\mathcal{N}_{a}}-1)$ so-called interspots, obtained by the correlation of the sub-PSF from sub-aperture $n$ in one diversity image by the sub-PSF from sub-aperture $n^{\prime}\neq n$ from the other diversity image. Those two kinds of correlation terms are clearly distinguishable by their positions. Indeed, the interspots result from the differential tip/tilt between the two involved sub-apertures, which can be large during the first alignment steps. The autospots positions, near the origin, are not affected by the tip/tilt of their related aperture, but only by the deterministic diversity-induced tip/tilt. Fig. 4 line 2 shows that the autospots, whose amplitude scales with $\rho_{n}^{4}$ according to Eq. (7), allow the simple estimation of the intensity of each sub-PSF.

The second ingredient is to make these autospots essentially not overlapping. Firstly, using a sufficiently large diversity (Appendix 4), so that the sub-PSFs are shifted from one diversity plane to the other, it is possible to have all the autospots’ cores separated. The small overlapping via the diffraction rings creates a slight deterministic coupling, which decreases as the separation increases, but does not prevent the autospots to form a nearly orthogonal family. Secondly, if the relative tip/tilts between sub-apertures are large enough with respect to the diversity, autospots are only slightly affected by the diffraction rings of the interspots. Indeed, the interspots are outside the circle overplotted on Fig. 4 line 2 as soon as the distance (derived as a tip/tilt value in Appendix 5) between two spots is larger than the diversity-induced shift from one plane to the other. If these two conditions are fulfilled, each of these autospots can be isolated by projecting the correlation on a predefined filter matching each autospot. From now on, the configuration in which the sub-PSFs have the minimal distance between each other to satisfy the second condition is called parking position.

The third ingredient is to rather compute the FT of this correlation, also called cross-spectrum, which contains the same information (under a different form). This cross-spectrum can be computed as the simple product between the first diversity image FT ($\boldsymbol{s}_{1}$ from subsection 2.2) and the second diversity image conjugated FT ($\boldsymbol{s}_{2}^{*}$). The autospot and interspot FTs will be respectively called autopeaks and interpeaks. Although all the peaks are superimposed in the center of the Fourier plane (Fig. 4, bottom line), the auto/inter peaks inherit the (near) orthogonality properties of the autospots and interspots by the Parseval-Plancherel theorem.

The fourth ingredient of the method is to code the sought tip-tilts in the autopeaks. In the basic cross-spectrum computation (Fig. 5, left) the contribution of the seeked tip-tilt, in the autopeaks, is cancelled by the phase conjugation. The idea is to realize that (1) the sought tilt is a phase ramp in the FT of the sub-PSF, and (2) if we perform a slight shift of the FT of one sub-PSF, this phase ramp is transformed into itself plus a piston, whose amplitude is proportional to the input tilt coefficient and the shift offset (Fig. 5, right). The nice feature of this piston factor is that it is constant all over the autopeak, and thus can be factored as a global weighting coefficient, like its amplitude $\rho_{n}^{4}$. The used operator, detailed in appendix 6, is hereafter called the Frequency Shifted Cross-Spectrum (FSCS).

The last ingredient, in the context of a closed-loop control, is to bring the sub-PSFs from a random scattered state to the parking position. From this parking position, an open-loop pre-defined offset can then superimpose all the sub-PSFs in the focal plane.

This has led to the ELASTIC (Estimation of Large Amplitude Sub-aperture Tip/tilt from Image Correlation) algorithm, which we now explain with all the technical details.

3.2 Derivation of the algorithm

In practice, the inputs of the algorithm are the diversity images from a $\mathcal{N}_{p}\times\mathcal{N}_{p}$ pixel detector. The typical width (in pixels) of a sub-PSF is thus equal to the sampling factor $n_{s}$ defined from the pixel pitch $p_{pix}$ by:

[TABLE]

Note that $n_{s}=2$ if the sub-PSFs are Shannon-sampled; $n_{s}$ is greater in practice for the interferometric PSF to be at least Shannon sampled. To compute the FSCS without wrapping, images are zero-padded to a width of $2\mathcal{N}_{p}$ pixels before the computation of the discrete OTFs $\boldsymbol{s}_{d,\boldsymbol{q}}$ by a Discrete Fourier Transform (DFT). The support width of $\boldsymbol{s}_{d,\boldsymbol{q}}$ is thus $4\mathcal{N}_{p}/n_{s}$ pixels. Then, we approximate the discrete OTFs $\boldsymbol{s}_{d,\boldsymbol{q}}$ to a sampled version of their continuous OTF model $\boldsymbol{s}_{d}$ derived in Eq. (7). The link between the 2D index $\boldsymbol{q}$ of $\boldsymbol{s}_{d,\boldsymbol{q}}$ and the continuous reduced coordinates $\boldsymbol{u}$ of Eq. (7) is thus $\boldsymbol{u}=\boldsymbol{q}n_{s}/\mathcal{N}_{p}$ since the support half width of $\boldsymbol{s}_{d}$ is $|\boldsymbol{u}_{max}|=2$. Therefore:

[TABLE]

Let $\Delta_{m}$ the operator that performs a shift ${\boldsymbol{\delta}}_{m}$ ($m={2,3}$), with an amplitude of $\delta$ pixel(s):

[TABLE]

We define the FSCS vector $\boldsymbol{j}_{m,\boldsymbol{q}}$ for each value of $m$ (2 or 3) and pixel $q$ as:

[TABLE]

where $\cdot^{*}$ denotes complex conjugation. This computation is illustrated on the right part of Fig. 5. Inserting Eq. (10), and keeping only the ${\mathcal{N}_{a}}$ autopeaks as explained in section 3.1, Eq. (12) becomes:

[TABLE]

According to Eqs. (A26)-(A28) from Appendix 1:

[TABLE]

with $\delta_{k,m}$ the Kronecker delta.

Thus Eq. (13) becomes:

[TABLE]

Figs. 3 and 4 suggest that the FSCS computation can be limited to a predefined reduced set of indexes, the OTF central part, where the autopeaks values are significant. So let us concatenate the corresponding values of $j_{m,\boldsymbol{q}}$ into a vector $\boldsymbol{j}_{m}$. Besides, our unknowns of interest are the tilts $a_{m,n}$, and possibly the transmittances $\rho_{n}$. So we group them into a new vector of complex unknowns :

[TABLE]

Lastly, Eq. (15) is obviously linear in $\boldsymbol{x}_{m}$, so it can be rewritten in the final matrix form of the direct model:

[TABLE]

where $\mathrm{C_{m}}$ is a matrix with $2\times\,{\mathcal{N}_{a}}$ columns, a number of lines equal to the number of indices kept in $\boldsymbol{j}_{m}$ and the elements of $\mathrm{C_{m}}$ are made with the appropriate sampled values of \mbox{\huge\curlywedge}_{n,m}.

The generalized inverse $\mathrm{C_{m}^{\dagger}}$ of $\mathrm{C_{m}}$ can easily be computed by Singular Value Decomposition. Thus the resolution of the inverse problem yields the solution $\widehat{\boldsymbol{x}_{m}}$:

[TABLE]

The seeked tip/tilts and pupil amplitudes are then simply computed as:

[TABLE]

A self-consistency test can check that both amplitude estimations agree: $\forall\leavevmode\nobreak\ n,\rho_{2,n}\simeq\rho_{3,n}$. It must be noted that $\mathrm{C_{m}^{\dagger}}$ can be pre-computed, so that the only real-time operations are the computation the two (tip and tilt) frequency shifted cross-spectra from the two (possibly undersampled) diversity images (Eq. (15)), then the two projections (Eq. (20)) and arguments of Eq. (21). Moreover, for a given sub-aperture tip/tilt, the sub-PSF shifts (and thus the OTF slopes) are the same for all wavelengths. The ELASTIC algorithm can thus operate with broadband illumination, as long as the sub-OTFs $\Lambda_{n,d}$ considered in the direct model are the polychromatic OTF obtained by averaging over the spectral band the monochromatic sub-OTFs weighted by the source amplitude and the detector efficiency.

4 Algorithm optimization and performance

The optimization and performance evaluation performed in this section by numerical simulations are based on a compact 18 sub-aperture pupil chosen to mimick the JWST pupil (Fig. 6), with circular sub-apertures for compatibility with our software.

4.1 Definition of the simulation parameters

Quasi-monochomatic images of size $\mathcal{N}_{p}^{2}=1024\times 1024$ pixels are simulated from an unresolved object. The focal-plane detector is assumed to be sampled at the Nyquist frequency in the fine phasing mode. Since each sub-aperture diameter is $1/5$ of the full aperture, the sub-PSFs are oversampled and the sub-aperture sampling ratio, $n_{s}$, is $10$. An amplitude of $A_{4}=0.9$ wave (Eq. (A37)) is chosen for starters but will be detailed later.

Images are simulated with a total of $3\times 10^{5}$ photo-electrons, Poissonian photon noise and a $5$ electrons per pixel read-out-noise. Algorithm performance is quantified with the estimation of the Root Mean Square Error (RMSE) of $\mathcal{N}_{o}$ outcomes (with $\mathcal{N}_{o}=50$) defined as:

[TABLE]

where $\widehat{a}_{m,n}$ and $a_{m,n}$ are respectively the estimated and introduced aberration coefficients in waves. $\langle\cdot\rangle$ is the average over the $\mathcal{N}_{o}$ outcomes.

4.2 Tilt dynamic range and choice of $\delta$

The maximum tilt that can be measured is ultimately limited by the detector field. This maximum tilt, noted $a_{field}$, is given by Eq. (A36) assuming a $\pm\mathcal{N}_{p}/2$ shift from the central origin.

[TABLE]

Here, $a_{field}=1024/(8\times\,10)\simeq 13$ waves. According to Eq. (21), the phase of the $\boldsymbol{x}_{m,n}$ coefficients is proportional to the seeked aberration coefficients $a_{m,n}$ and to the shift amplitude $\delta$ chosen for the computation of the FSCS. Because the dynamic range of the $\operatorname{Arg}$ function is only $(-\pi;\pi]$, increasing the shift amplitude reduces the range since from Eqs. (21) and (23), the maximal unwrapped estimation is:

[TABLE]

Eq. (24) highlights that the algorithm does not limit the field when $\delta\leq 2$.

These limits, given in Table 1, are confirmed by simulation. Fig. 7 presents the ELASTIC estimation versus the introduced tip aberration for different values of the shift amplitude $\delta$. A tip slope is applied over one sub-aperture, the other sub-apertures being randomly scattered in the field. Each point on the graph is the average of $\mathcal{N}_{o}$ noise outcomes. A full-field linear response is obtained for a 1 and 2 pixel shifts. The curves corresponding to 3 and 4 pixel shifts confirm the previously mentioned phase wrapping, with the theoretical boundaries around $66\%$ or $50\%$ field coverage, respectively for a 3 or 4 pixel shift.

A 1 or 2 pixel shift are necessary for a maximum field coverage, hence can be chosen for initial error estimation. Then, as the field covered by the sub-PSFs decreases, the pixel shift amplitude can be increased to improve the estimation accuracy as seen in next section.

4.3 Optimization of the algorithm’s parameters

The ELASTIC algorithm has two free parameters: the introduced diversity and the frequency shift amplitude $\delta$. As mentioned in the previous subsection, it is required to start the alignment with $\delta=2$ if the sub-PSFs are spread over the whole detector field. The goal here is to optimize and quantify performance in the ultimate step, when the sub-PSFs are close to the parking position. Considering Eqs. (18-20), each $\widehat{x}_{m,n}$ value can be written as $x_{m,n}(1+\boldsymbol{\epsilon}_{n}/x_{m,n})$, with $\boldsymbol{\epsilon}_{n}$ the small noise. Eq. (21) becomes:

[TABLE]

with $\operatorname{Im}(\cdot)$ the imaginary part of a complex number. Because $\boldsymbol{\epsilon}_{n}$ is roughly independant of $\delta$, it comes that increasing $\delta$ reduces the error on the estimation (if $\mathcal{N}_{p}$ and $n_{s}$ are fixed). Therefore the parking position (depending on the pseudo-defocus amplitude) needs to be tightly packed for two reasons. Firstly, it allows to increase $\delta$ hence to improve the estimation accuracy. Secondly, it will minimize errors due to the uncertainty of segment displacements (for the open loop superimposition after ELASTIC). Fig. 8 presents the parking position chosen for the sub-PSFs. The positions of the sub-PSFs from the central corona were determined thanks to Appendix 5. Then, the sub-PSFs from the external corona were positioned so that the interspots do not overlap the autospots after image correlation (Fig. 8, line 2).

Fig. 9 presents the evolution of the ELASTIC estimation RMSE versus the pseudo-defocus amplitude $A_{4}$, for different values of the shift amplitude $\delta$. It is clearly evidenced that ELASTIC estimation RMSE decreases when the pseudo-defocus amplitude and the pixel shift amplitude are increased. The algorithm performance is better than $\lambda/8$ for a pseudo-defocus larger than $0.8$ to $1$ wave. This is in agreement with the typical value of $0.9$ wave given by Eq. (A37). For a larger diversity amplitude the RMSE estimation goes down to less than $\lambda/500$ for a $2$ or $3$ pixel shift. Indeed, when the diversity amplitude increases, the interspots are shifted further from the autospots. Hence their influence on the autospots decreases then the estimation’s accuracy is better. The limit to these performances is that the diversity amplitude must be chosen so that the sub-PSFs do not slip out of the unwrapped field (settled by $\delta$). This limitation is not fulfilled when the pseudo defocus amplitude is larger than 2.3 waves and 3.3 waves for respectively a 4 and 3 pixel shift amplitude: the RMSE increases. There is then a trade off between the choice of $\delta$ and the pseudo-defocus amplitude. The evolution of the condition number of $\mathrm{C_{m}}$ is also plotted. We see that it is indeed a good indicator of the matrix inversion’s sensivity to noise. As the ratio between the highest and the smallest singular values, it can help one qualify the ability to separate autospots as a function of the pseudo-defocus amplitude. The condition number decreases when the pseudo-defocus amplitude increases, and reaches an asymptote at 1. Indeed a large diversity improves the autospots separation. Therefore, the singular modes are well distinguishable for the matrix inversion. The condition number is the same for each shift amplitude case. Indeed the autospot separation in the correlation only depends on the diversity, not on the introduced frequency shift. Hence, increasing the latter makes the estimation more accurate.

4.4 Robustness to noise

In order to have the best accuracy with the largest frequency shift amplitude and a relatively small defocus amplitude, considering the sub-PSFs in parking position, a $2.2$ waves pseudo-defocus and a $4$ pixel shift amplitude are chosen for the evaluation of noise propagation. $\mathcal{N}_{o}$ noise outcomes are simulated for different image brightnesses, quantified by $\mathcal{N}_{ph}$, the total number of photo-electrons. Fig. 10 presents the estimation’s RMSE as a function of the flux.

Several regimes are to be noticed on Fig. 10. Below $10^{4}$ photo-electrons (low illumination), the estimated RMSE is constant. In this case, the distribution of $\operatorname{Arg}(\widehat{\boldsymbol{x}}_{m,n})$ behaves as a uniform noise between $-\pi$ and $\pi$. Such a distribution has a theoretical standard deviation of $\pi/\sqrt{3}$. Using Eq. (21), we deduce that the theoretical saturation of the RMSE is $\mathcal{N}_{p}/(4\leavevmode\nobreak\ n_{s}\delta\sqrt{3})=a_{max}/\sqrt{3}$. The numerical value of the saturation for our setup is around $4$ waves, and matches with the saturation observed in the simulation.

A second regime is for a total image flux greater than $10^{5}$ photo-electrons: the estimation RMSE is saturated. This last limitation is the algorithm’s bias. The bias can be explained by the small influence of the interpeaks on the autopeaks. However, the RMSE is down to less than $\lambda/100$, which is already much better than needed to enter the fine phasing mode.

These results confirm the capacity of ELASTIC to bring a misaligned multi-aperture telescope to a configuration state where the wavefront errors are $\lambda/100$ allowing the fine phasing algorithms to operate.

4.5 Robustness to phase diversity error

An interesting feature of the ELASTIC algorithm is its small sensitivity to the diversity error. This can be understood in the focal domain: if the diversity used in the numerical model does not match exactly the actual optical diversity, then the autospots in the frequency-shifted cross-correlation (cf appendix 6) do not perfectly overlap the associated data in the FT of the $\mathrm{C_{m}^{\dagger}}$ vectors. The estimated intensity will be affected by the partial overlap. But since the piston information induced by the seeked tip/tilt is constant and not degraded over the autospot, its extraction is not much affected (neglecting the overlapping of the spots associated to different sub-apertures). This is confirmed by simulations: the estimation error is less than $\lambda/8$ as long as the diversity-induced tip or tilt error is less than $\lambda/2.5$ on each sub-aperture [41].

4.6 Adaptation to atmospheric turbulence

In the case of a ground-based telescope, the atmospheric turbulence induces (among others) tip-tilt disturbances on the sub-PSFs. A solution to average the random atmospheric tip/tilt is to record long-exposure images, leading to sub-PSFs of size $\lambda/r_{0}$ (instead of $\lambda/D$) with $r_{0}$ the Fried parameter [38]. Thus, the ELASTIC algorithm can be adapted by changing the photometric peaks $\Lambda_{n,d}$ by narrowed versions because frequencies above $r_{0}/\lambda$ (instead of $D/\lambda$) are lost. All the remaining processing can be performed, leading to a fraction of $\lambda/r_{0}$ precision (instead of a fraction of $\lambda/D$). This accuracy should ensure that we enter the capture range of usual wavefront sensors over all the sub-apertures simultaneously, in order to operate an adaptive optics system.

5 Experimental validation

In the following, sub-PSFs of a multi-aperture telescope are brought to the parking position thanks to the ELASTIC algorithm, and then superimposed.

5.1 Implementation of ELASTIC

In order to test focal plane wavefront sensors, Onera built a dedicated bench called BRISE [33].

A fibered single-mode laser diode operating near $635$ nm is used as an unresolved source. It is collimated and provides an object at infinity to a plane segmented mirror with nineteen sub-apertures (Fig. 11 left). To introduce and correct piston/tip/tilt perturbations, each of these mirrors is supported by three piezoelectric actuators. These have no internal feedback control and suffer from hysteresis. In addition, the external corona is currently not fully functional, thus the validation is performed over the $6$ sub-apertures of the first corona. Downstream, a phase diversity module is used to simultaneously form a focused and a defocused image of the object on a $1300\leavevmode\nobreak\ \times\leavevmode\nobreak\ 1000$ pixels camera, from which we extract two $512\leavevmode\nobreak\ \times\leavevmode\nobreak\ 512$ images.

Because the active mounts are not perfectly deterministic, it is not possible to perform a single-step good correction from a randomly misaligned state. Therefore, a closed-loop sequence is performed to align the mirror. According to Eq. (A37), the minimal required diversity for ELASTIC is $0.94\lambda$. Measurements are then performed with a $1\lambda$ pseudo-defocus temporally applied on the 512$\times$512 pixel focal image. To have a full-field tip/tilt estimation, the FSCS is computed with a 2 pixel shift.

5.2 Loop closure on an unresolved source

The focal and pseudo-defocused images of an unresolved source are shown on Fig. 12 line 1. They include random tip/tilt errors on each sub-aperture. The aim is to obtain in the focal plane the same simulated parking position as Fig. 11 right. For the closed loop, we use an integrator with a gain control of $0.5$. Lines 2 and 3 of Fig. 12 show iterations 2 and 4 of the loop. Line 4 is the position obtained at iteration 20, close to the expected parking.

Fig. 13 presents the experimental evolution of the estimated tip/tilt errors during the closed-loop, for each iteration. The error is defined as the relative distance to the parking position. At first iteration, some tip/tilt errors are larger than $2\leavevmode\nobreak\ \lambda$. As from iteration 5, the estimated error is mainly less than $\lambda/2$. Moreover, the closed-loop control is stable as from around iteration 7 and keeps the sub-PSFs in the parking position until the end with an estimated error less than $\lambda/10$. It is to be noted that the estimation can be biased if two sub-PSFs come too close (iteration 2), but the error always converge to zero thanks to the closed-loop control.

It can be seen on iteration 20 that the sub-PSFs are not perfectly placed (they do not form the exact pre-defined parking). This is mainly due to the temporal modulation over each sub-aperture. Indeed, in open loop, if the diversity offset is consecutively switched on and off, the sub-PSFs do not go back to the exact same position they had before the offset. Hence, the instrumental limitations (plus a possible error in the telescope modeling and environment potential perturbations) can explain this visible positionning error on iteration 20. Nevertheless we can conclude that ELASTIC was succesfully implemented on BRISE bench, and was able to bring the telescope from a randomly misaligned configuration to a stable parking position. The error of the estimation, due to the above mentioned limitations, is now to be quantified.

5.3 Noise evaluation

To experimentaly evaluate the estimation standard deviation for several illumination conditions, the sub-PSFs are set in a fixed configuration (e.g. the previously mentioned parking position) and $\mathcal{N}_{o}$ open-loop diversity pairs of images are taken to estimate tip/tilt. Fig. 14 plots the estimated standard deviation for various illumination values (defined as the total flux in each image).

Results show that the estimation’s standard deviation is less than $\lambda/100$ as from a $10^{6}$ photo-electrons illumination, and decreases following a $1/\sqrt{\mathcal{N}_{ph}}$ slope. This highlights the stability of the estimation for reasonable flux conditions.

Despite the instrumental limitations (turbulence from the camera fan, hysteresis, possible calibration errors) that remain, the stability of the estimation is better than the capture range of the fine algorithms.

5.4 Final superimposition

After the closed-loop organisation of the sub-PSFs in parking position, an open-loop offset command is applied to finish the alignement (Fig. 15). In the focal plane, interference fringes are distiguishable. Eventhough the large piston errors were not corrected, they are smaller than the coherence length of the single-mode laser used as an unresolved source. Without hysteresis or calibrations errors, simulations showed in Section 4.4 that a $\lambda/100$ precision could be reached for high illumination, in the case of a 18 sub-aperture instrument. Experimentally, the result presented on Figure 15 shows that in the case of a 6 sub-aperture instrument, the open-loop offset command confines all sub-PSFs in an area which is less than $30$ pixels FWHM. Thus the sub-PSFs dispersion after superimposition is +/- 7 pixels peak-valley around the center, leading to an estimation of residual tip/tilt error smaller than $\lambda/8$ RMS, according to Appendix 3 since the FWHM of a single sub-PSF is $16$ pixels ($n_{s}=16$ on our bench).

This is confirmed by Fig. 16, which shows the result of a cophasing by phase diversity [17, 18] after superimposition by ELASTIC. Not only the sub-PSFs are all in the capture range of the fine phasing algorithm, but the width of the cophased pattern (2 pixel central core but around 30 pixel enveloppe) is comparable to the width of the aligned state after ELASTIC.

6 Conclusion

To estimate the tip/tilts (and pupil transmission amplitude) over each sub-aperture of a multi-aperture telescope, we introduced the ELASTIC method whose main features are to provide a polychromatic large amplitude estimation (up to the full camera field) with a simple hardware (only two images of a point source near the focal plane) and a closed-form unsupervised computation with small computing cost. Its typical use is to bring the sub-apertures from any distorted state up to a sufficiently aligned state (with residues much smaller than half the sub-PSF width, to enable their interference) to scan for pistons and ultimately enter the fine-phasing mode. Such a source (unresolved by each sub-aperture) can easily be found in a stellar field. For applications such as Earth observation from space, the ELASTIC algorithm may be extended to be insensitive to the object phase; however, it can be expected that large amplitude alignment is not frequent and can be operated by pointing the telescope at a star to operate ELASTIC with an unresolved source.

The main requirement of ELASTIC is to acquire two images, including known local tip/tilts offsets over the sub-apertures, so that each sub-PSF is uniquely identified by its position shift. These offsets need to be large enough to allow a clear separation between the sub-PSFs, hence are larger than those required by phase diversity in the fine phasing mode, which needs the sub-PSFs to remain superimposed in all diversity planes to interfere. This geometric diversity can be implemented by inserting a global aberration between the images (e. g. focal/defocused images) or by introducing a temporal modulation with the sub-aperture correction actuators themselves. This is also simpler since telecentricity is not required, and more efficient as only the tip-tilt part of the defocus is used by the algorithm, whereas higher-order modes only degrade the sub-PSFs.

ELASTIC is based on the computation of a frequency-shifted cross-spectrum, an operator detailed in appendix 6 introduced to simply invert the diversity process. We showed by means of numerical simulations with a JWST-like pupil how to optimize the free parameters (the diversity and the frequency shift) and estimated the limiting magnitude around $5\times 10^{3}$ photo-electons/sub-aperture/frame with $1024\times\leavevmode\nobreak\ 1024$ pixels images. Lastly, we performed an experimental validation that demonstrated the closed-loop alignment of a 6 aperture segmented mirror using a tip/tilt temporal diversity on the segments and $16\times\leavevmode\nobreak\ 10^{4}$ photo-electrons/sub-ap/frame with $512\times\leavevmode\nobreak\ 512$ pixels images.

Experience shows that the ELASTIC procedure is robust and is now regularly used to drive BRISE into the fine phasing mode. A short-term perspective is then to interface ELASTIC with a real-time fine phasing algorithm such as the one proposed in [34, 35].

Fundings

This research was partly funded by ONERA’s internal research project VASCO ; Thales Alenia Space co-funded S. Vievard’s PhD thesis; the BRISE bench was funded by French DGA

Acknowledgements

The authors would like to thank B. Denolle, A. Grabowski and C. Perrot for their early work on the algorithm and the experiments; J. Montri for the BRISE computer interface; and J.-P. Amans from GEPI (Galaxies, Etoiles, Physique et Instrumentation) of Observatoire de Paris-Meudon, for the design and manufacturing of the segmented mirror.

Appendix A Appendix

1 Definition of the Zernike modes

This appendix recalls the Zernike polynomials [32]: a basis of orthonormal modes $Z_{k}$ with normalized coordinates $\boldsymbol{u}_{n}=(u_{n},v_{n})$, where $|\boldsymbol{u}_{n}|\leq\leavevmode\nobreak\ 1$. The first 4 modes we will use are:

[TABLE]

2 From global to local aberrations

The global position $\boldsymbol{r}$ of any point in the pupil can be linked to a local coordinates $\boldsymbol{r}_{n}$ attached to each sub-aperture:

[TABLE]

with $\boldsymbol{c}_{n}=(c_{x,n},c_{y,n})$ the $n$th sub-aperture center normalized by $R$. We define $B$ as the global aperture diameter, and:

[TABLE]

with $\boldsymbol{u}_{n}$ the reduced coordinates over the sub-aperture and $\boldsymbol{u}^{\prime}$ the reduced coordinates over the global aperture. Therefore:

[TABLE]

Then, inserting in Eq. (A29), a defocus can be decomposed as:

[TABLE]

Therefore, the tip (or tilt) induced by a defocus on the $n^{th}$-subaperture is:

[TABLE]

3 Computation of the tip/tilt coefficients

A tilt disturbance can be expressed either by an angle $\theta$ in the object space, by a $k_{p}$ pixel shift over the detector (with a $p_{pix}$ pixel pitch) or by a $a_{m}$ coefficient, in waves, linked by:

[TABLE]

since from Eq. (A27) the peak-to-valley amplitude of $Z_{2}$ is 4. A consequence from Eq.(A36) is that the FWHM of a sub-PSF corresponds to as tilt variation of $\Delta\,a_{2}=1/4$ since the sub-PSF’s angular FWMH is $\Delta\theta=\lambda/2R$.

4 Typical diversity for ELASTIC

The goal of this appendix is to compute the amount of defocus to have the distance between the autospots just equal to their diameter. $P$ is defined as the distance between two aperture centers. It can be computed that the autocorrelation of an Airy pattern (which has an angular FWHM of $\lambda/(2R)$ rad) has an angular FWHM of $(2\lambda)/(3R)$ rad or 1/3 of wave from Eq. (A36). In addition, the separation between two adjacent sub-apertures with a defocus of amplitude $A_{4}$ (RMS value over the global full aperture) is from Eq. (A35) $\Delta a_{m}=2\sqrt{3}\frac{4R^{2}}{B^{2}}\,PA_{4}$. Equalizing this last expression with 1/3 gives the typical defocus value:

[TABLE]

5 Minimal distance between the sub-PSFs

The minimal distance between the sub-PSFs for the interspots not to pollute the signal is now defined, in the case of a single corona aperture. The separation between autospots from extremal apertures with a defocus of amplitude $A_{4}$ (RMS value over the full aperture) is $\Delta a_{m}=2\sqrt{3}\frac{4R^{2}}{B^{2}}(B-2R)A_{4}$ from Eq. (A35). Since the sub-PSF FWHM corresponds to a $1/4$ tilt from App. (3), the minimal separation is:

[TABLE]

6 The Frequency Shifted Cross-Spectrum (FSCS)

This appendix considers the various ways to correlate two signals $\boldsymbol{i}_{1}$ and $\boldsymbol{i}_{2}$, which here would stand for the two diversity PSFs.

Authors agree to define their cross-correlation or inter-correlation as $\boldsymbol{i}_{1,2}=\boldsymbol{i}_{1}\,\otimes\,\boldsymbol{i}_{2}$ (cf. Eq. 5) [36].

The cross-spectrum is sometimes defined as the FT ($\mathcal{F}$) of the cross-correlation: $\mathcal{F}\left[\boldsymbol{i}_{1,2}\right]=\mathcal{F}\left[\boldsymbol{i}_{1}\right]\times\mathcal{F}\left[\boldsymbol{i}_{2}\right]^{*}=\boldsymbol{s}_{1}\times\boldsymbol{s}_{2}^{*}$ and, in imaging through turbulence [37, 38, 39], as $\boldsymbol{s}_{1,1\boldsymbol{\delta}}=\left\langle\boldsymbol{s}_{1}\,\Delta_{m}\boldsymbol{s}_{1}^{*}\right\rangle$ where the average is performed over turbulent phase outcomes and $\Delta_{m}$ is the shift operator by a small increment $\boldsymbol{\delta}$ defined in Eq. (11), used by the Knox-Thompson method to estimate the object phase from turbulence-degraded images [40]. This cross-spectrum only holds for a stack of single images.

The product used by Eq. (12), defined as $\boldsymbol{s}_{1,2\boldsymbol{\delta}}=\boldsymbol{s}_{1}\,\Delta_{m}\boldsymbol{s}_{2}^{*}$, is a combination of these two products: it is performed in the spectral domain, between two different signals, with a small shift introduced to access some phase information as detailed by Fig. 5. To the best of our knowledge, it has not been introduced before so we name it the Frequency Shifted Cross-Spectrum (FSCS).

The FSCS can be equivalently understood in the direct domain from Fig. 4: the spectral shift of $\boldsymbol{s}_{2}$ is equivalent to the multiplication of $\boldsymbol{i}_{2}$ by a phase slope. In the real image of Fig. 4, line 1, right, each sub-PSF would inherit a phase, proportional to its $x$-coordinate. After what can be called a frequency-shifted cross-correlation (the FT of the FSCS), the positions of the autospots are still given by the differential positions of their related sub-PSFs in the two images as illustrated by Fig. 4, line 2, but their value would have a phase given by the associated phase (thus absolute position) in the sole $\boldsymbol{i}_{2}$ image.