Kernel formalism applied to Fourier based wave front sensing in presence of residual phases

TL;DR

This paper models Fourier-based Wave Front Sensors as linear integral operators with kernels, simplifying their analysis and optimization, especially in the presence of residual phases, and introduces a convolutional kernel approach for fast diagnostics.

Contribution

It introduces a convolutional kernel framework for Fourier-based WFSs, enabling simplified modeling and sensitivity analysis, even with residual phases.

Findings

Convolutional kernels allow concise WFS behavior representation.

Residual phases do not invalidate the convolutional model.

The approach simplifies sensitivity computation and optimization.

Abstract

In this paper, we describe Fourier-based Wave Front Sensors (WFS) as linear integral operators, characterized by their Kernel. In a first part, we derive the dependency of this quantity with respect to the WFS's optical parameters: pupil geometry, filtering mask, tip/tilt modulation. In a second part we focus the study on the special case of convolutional Kernels. The assumptions required to be in such a regime are described. We then show that these convolutional kernels allow to drastically simplify the WFS's model by summarizing its behavior in a concise and comprehensive quantity called the WFS's Impulse Response. We explain in particular how it allows to compute the sensor's sensitivity with respect to the spatial frequencies. Such an approach therefore provides a fast diagnostic tool to compare and optimize Fourier-based WFSs. In a third part, we develop the impact of the residual phases on the sensor's impulse response, and show that the convolutional model remains valid. Finally, a section dedicated to the Pyramid WFS concludes this work, and illustrates how the slopes maps are easily handled by the convolutional model.

Equations

$$\text{Output}{|}_X=\int \text{d}x\text{K}{|}_{X;x}\text{Input}{|}_x$$

$${\psi }_i={\mathbb{I}}_P{e}^{ı\varphi }$$

$$\text{WFS’s input:}{\mathbb{I}}_P\varphi$$

$${\varphi }_{\text{mod}}(\vec{\mathbf{\alpha }}){|}_{\vec{\text{r}}}=\frac{2\pi }{\lambda }\vec{\mathbf{\alpha }}.\vec{\text{r}}$$

$$\begin{array}[]{ccccc}w&:&\mathbb{R}^{2}&\to&\mathbb{R}_{+}\\ &&\vec{\boldsymbol{\alpha}}&\mapsto&w|_{\vec{\boldsymbol{\alpha}}}\\ \end{array}$$

$${\int }_{{\mathbb{R}}^2}w{|}_{\vec{\mathbf{\alpha }}}{\text{d}}^2\vec{\mathbf{\alpha }}=1$$

$$\psi_{d}(\vec{\boldsymbol{\alpha}})=\Big{[}\mathbb{I}_{P}\exp\big{(}\imath(\phi+\phi_{\text{mod}}(\vec{\boldsymbol{\alpha}})\big{)}\Big{]}\star\hat{m}$$

$$I(\phi)=\int_{\mathbb{R}^{2}}\text{d}^{2}\vec{\boldsymbol{\alpha}}~{}w|_{\vec{\boldsymbol{\alpha}}}~{}\\ \left|\Big{[}\mathbb{I}_{P}\exp\big{(}\imath(\phi+\phi_{\text{mod}}(\vec{\boldsymbol{\alpha}})\big{)}\Big{]}\star\hat{m}\right|^{2}$$

$$\exp (ı\varphi )=\sum _{q=0}^{\infty }\frac{{ı}^q{\varphi }^q}{q!}$$

$$I(\varphi )=I_{\text{constant}}+I_{\text{linear}}(\varphi )+I_{\text{quadratic}}(\varphi )+I_{\text{cubic}}(\varphi )+...$$

$$\Delta I(\varphi )$$

$$=$$

$$\Delta I(\varphi )\approx I_{\text{linear}}(\varphi )$$

$$\varphi <<1$$

$$I_{\text{linear}}({\color[rgb]{0,0,0}\phi})=2\textbf{Im}\bigg{[}\int_{\mathbb{R}^{2}}\text{d}^{2}\vec{\boldsymbol{\alpha}}~{}~{}w|_{\vec{\boldsymbol{\alpha}}}~{}\Big{(}\big{(}\mathbb{I}_{P}e^{\imath\phi_{\text{mod}}(\vec{\boldsymbol{\alpha}})}\big{)}\star\hat{m}\Big{)}\Big{(}\big{(}\mathbb{I}_{P}e^{-\imath\phi_{\text{mod}}(\vec{\boldsymbol{\alpha}})}{\color[rgb]{0,0,0}\phi}\big{)}\star\bar{\hat{m}}\Big{)}\bigg{]}$$

$$I_{\text{linear}}(\phi)|_{\vec{\textbf{R}}}=2\textbf{Im}\bigg{[}\int_{\mathbb{R}^{4}}\text{d}^{2}\vec{\textbf{r'}}\text{d}^{2}\vec{\textbf{r}}~{}~{}\mathbb{I}_{P}|_{\vec{\textbf{r'}}}(\mathbb{I}_{P}\phi)|_{\vec{\textbf{r}}}\hat{m}|_{\vec{\textbf{R}}-\vec{\textbf{r'}}}\bar{\hat{m}}|_{\vec{\textbf{R}}-\vec{\textbf{r}}}\int_{\mathbb{R}^{2}}\text{d}^{2}\vec{\boldsymbol{\alpha}}w|_{\vec{\boldsymbol{\alpha}}}e^{\frac{2\imath\pi}{\lambda}\vec{\boldsymbol{\alpha}}.(\vec{\textbf{r'}}-\vec{\textbf{r}})}\bigg{]}$$

$$I_{\text{linear}}(\phi)|_{\vec{\textbf{R}}}=2\textbf{Im}\bigg{[}\int_{\mathbb{R}^{4}}\text{d}^{2}\vec{\textbf{r'}}\text{d}^{2}\vec{\textbf{r}}~{}~{}\mathbb{I}_{P}|_{\vec{\textbf{r'}}}(\mathbb{I}_{P}\phi)|_{\vec{\textbf{r}}}\\ \hat{m}|_{\vec{\textbf{R}}-\vec{\textbf{r'}}}~{}\bar{\hat{m}}|_{\vec{\textbf{R}}-\vec{\textbf{r}}}\hat{w}|_{\vec{\textbf{r'}}-\vec{\textbf{r}}}\bigg{]}$$

$$I_{\text{linear}}(\varphi ){|}_{\vec{\text{R}}}={\int }_{{\mathbb{R}}^2}{\text{d}}^2\vec{\text{r}}\text{K}{|}_{\vec{\text{R}};\vec{\text{r}}}({\mathbb{I}}_P\varphi ){|}_{\vec{\text{r}}}$$

$$\text{K}{|}_{\vec{\text{R}};\vec{\text{r}}}=2\text{Im}\left[\bar{\hat{m}}{|}_{\vec{\text{R}}-\vec{\text{r}}}{\int }_{{\mathbb{R}}^2}{\text{d}}^2\overrightarrow{\text{r’}}{\mathbb{I}}_P{|}_{\overrightarrow{\text{r’}}}\hat{m}{|}_{\vec{\text{R}}-\overrightarrow{\text{r’}}}\hat{w}{|}_{\overrightarrow{\text{r’}}-\vec{\text{r}}}\right]$$

$$\text{Without modulation:}~{}~{}~{}~{}~{}\textbf{K}|_{\vec{\textbf{R}};\vec{\textbf{r}}}=2\textbf{Im}\Big{[}\bar{\hat{m}}|_{\vec{\textbf{R}}-\vec{\textbf{r}}}~{}(\mathbb{I}_{P}\star\hat{m})|_{\vec{\textbf{R}}}\Big{]}$$

$$\text{K}(\hat{m},\hat{w},{\mathbb{I}}_P){|}_{\vec{\text{R}};\vec{\text{r}}}$$

$$\text{K}{|}_{\vec{\text{R}};\vec{\text{r}}}=\text{IR}{|}_{\vec{\text{R}}-\vec{\text{r}}}$$

$$I_{\text{linear}}(\varphi )=({\mathbb{I}}_P\varphi )\star \text{IR}$$

$$\mathbb{I}_{P}|_{\vec{\textbf{r'}}}=\mathbb{I}_{P}|_{\vec{\textbf{r'}}-\vec{\textbf{r}}}+\vec{\textbf{r}}\cdot\left(\vec{\nabla}\mathbb{I}_{P}\big{|}_{\vec{\textbf{r'}}-\vec{\textbf{r}}}\right)+...$$

$$\textbf{K}|_{\vec{\textbf{R}};\vec{\textbf{r}}}=2\textbf{Im}\big{[}\bar{\hat{m}}(\hat{m}\star\big{(}\hat{w}\mathbb{I}_{P})\big{)}\big{]}\big{|}_{\vec{\textbf{R}}-\vec{\textbf{r}}}\\ +2\vec{\textbf{r}}\cdot\textbf{Im}\big{[}\bar{\hat{m}}(\hat{m}\star\big{(}\hat{w}\vec{\nabla}\mathbb{I}_{P})\big{)}\big{]}\big{|}_{\vec{\textbf{R}}-\vec{\textbf{r}}}+...$$

$${\mathbb{I}}_P{|}_{\overrightarrow{\text{r’}}}\approx {\mathbb{I}}_P{|}_{\overrightarrow{\text{r’}}-\vec{\text{r}}}$$

$$\hat{w}\partial {\mathbb{I}}_P=0$$

$$\textbf{IR}=2\textbf{Im}\big{[}\bar{\hat{m}}\big{(}\hat{m}\star(\hat{w}\mathbb{I}_{P})\big{)}\big{]}$$

$$\text{TF}\equiv \widehat{\text{IR}}$$

$$\text{K}{|}_{\vec{\text{R}};\vec{\text{r}}}=2\text{Im}\left[\bar{\hat{m}}{|}_{\vec{\text{R}}-\vec{\text{r}}}{\int }_{{\mathbb{R}}^2}{\text{d}}^2\overrightarrow{\text{r’}}\hat{m}{|}_{\vec{\text{R}}-\overrightarrow{\text{r’}}}(\hat{w}{\mathbb{I}}_P){|}_{\overrightarrow{\text{r’}}-\vec{\text{r}}}{\text{R}}_{\text{s}}{|}_{\overrightarrow{\text{r’}};\vec{\text{r}}}{\text{R}}_{\text{d}}{|}_{\overrightarrow{\text{r’}};\vec{\text{r}}}\right]$$

$${\text{R}}_{\text{s}}{|}_{\overrightarrow{\text{r’}};\vec{\text{r}}}$$

Full text

1 Introduction

Linear integral operators are the continuous version of matrices. They linearly transform an input into an output depending on a quantity called Kernel. Mathematically, such an operation may be written:

[TABLE]

$x$ (resp. $X$) is the variable of the input (resp. output) space, K making the link between those two spaces. Examples of integral operators are many: Fourier, Hilbert, Laplace transforms, etc. are the most famous. In a more general framework, it is absolutely relevant to try to describe a continuous linear system thanks to an integral operator. Advantages of such an approach are many. First of all, the Kernel allows to synthesize the system’s behavior into a compact and elegant mathematical quantity. Moreover, due to the fact linear integral operators may be seen as the continuous version of matrices (see Appendix C), the use of Kernel formalism may improve the understanding of the discrete description of systems that numerical approach often requires. We finally note that reconstructing the input from the output, i.e. inverting (1), depends on our capability to find the ”inverse Kernel” $\textbf{K}^{-1}$ of the system. Fortunately, depending on the nature of K, many methods already exist to calculate this inverse if it does exist. Knowing K is thus critical to build robust and relevant reconstruction algorithms.

The purpose of this article is to use this powerful formalism to study optical systems which probe the wavefront of the light. They are called Fourier-based Wave Front Sensors [1] and are essentially used in Adaptive Optics for Astronomy. The first part of this paper provides an optical description of Fourier-based WFSs. We make explicit the input and output of such sensors and introduce proper mathematical quantities to characterize the different optical elements. In a second part, we give the dependency of the Kernel with respect to these optical elements. The third part is dedicated to the particular case of convolutional Kernel. Phase reconstruction, in that situation, is greatly facilitated since the inverse Kernel has an explicit formulation. Moreover, it allows to define the WFS’s Impulse Response and Transfer Function which are compact and meaningful quantities allowing, for instance to rapidly visualize the WFS’s sensitivity (see Appendix A). Since the convolutional Kernel case is not systematic, we give the assumptions it relies on. Another part is dedicated to generalization of these results in presence of static and dynamic residual as it is often the case, in practice, in Adaptive Optics. These results are finally applied to the most used Fourier-based WFS which is the Pyramid WFSensor [2] in Appendix B.

2 Fourier-based Wave Front Sensing

2.1 Optical system

We consider the optical system shown in Fig. 1. The first plane corresponds to the pupil plane of the telescope. It contains a focusing device which is described as a perfect lens with a focal $f$, an aperture and another element that we call modulation. We describe it in detail in the next paragraph. The second plane is the focal plane which contains a filtering mask mathematically described by its transparency function $m$ and an imaging lens with focal $f/2$. The detector is finally placed in the next plane which is conjugated to the first pupil plane.

2.2 Optical propagation

The incoming field is called $\psi_{i}$, we assume that its light is monochromatic at wavelength $\lambda$. In the Adaptive Optics Wave Front Sensing context, the field is windowed by an entrance pupil. We call its phase $\phi$ for turbulent phase and describe the pupil geometry thanks to its indicator function $\mathbb{I}_{P}$. Besides, we assume that the total flux is unitary. Under these assumptions, we have:

[TABLE]

Let’s note that turbulent phase has a priori an infinite support. Nevertheless, due to windowing by the finite aperture, only the part which goes through the entrance pupil has a physical interest and may be coded by the WFS. As a consequence, the quantity to be measured, i.e. the sensor’s input is the turbulent phase windowed by the pupil:

[TABLE]

The incoming field is then shaped by the ”modulation”. Modulation consists in adding a dynamic aberration to the field via a controlled moving device placed in the entry pupil plane. The movement is regular and may be described as a cycle (although it can be generalised). The detector is then synchronized with it in order to have one image per cycle. This system was initially introduced by [2] in order to adjust the Pyramid WFS’s linearity range. It is specific to this sensor yet it should remain general. In this paper, we are interested in the particular case of ”tip/tilt modulation” which consists in adding a tip/tilt aberration in the field. This modulation phase is called $\phi_{\text{mod}}$ and defined as:

[TABLE]

The vector $\vec{\textbf{r}}$ codes the position in the pupil plane while $\vec{\boldsymbol{\alpha}}$ codes the amplitudes of the modulation tip and tilt. Brief clarifications about notations: $\phi_{\text{mod}}(\vec{\boldsymbol{\alpha}})|_{\vec{\textbf{r}}}$ means that the function $\phi_{\text{mod}}$ depends on one true spatial vector $\vec{\textbf{r}}$ but also on another non-spatial parameter $\vec{\boldsymbol{\alpha}}$. To describe the modulation, we also need to indicate the time spent for each modulation phase. To do so, we introduce the modulation weighting function $w$. Usually, this function directly depends on the time variable but we prefer here to use as variable the tip/tilt amplitudes $\vec{\boldsymbol{\alpha}}$:

[TABLE]

The weighting function is now a 2D function and $w|_{\vec{\boldsymbol{\alpha}}}$ codes the time spent for the modulation phase $\phi_{\text{mod}}(\vec{\boldsymbol{\alpha}})$. Such an approach presents several advantages. First, it allows to visualize the tip/tilt modulation in its natural focal plane. Indeed, in this plane $\vec{\boldsymbol{\alpha}}$ vector correspond to a spatial shift (that’s why, we do not use the notation $w(\vec{\boldsymbol{\alpha}})$ but $w|_{\vec{\boldsymbol{\alpha}}}$). Subsequently, the weighting function directly gives the profile of the modulation. Thus, it becomes possible to envision 2D modulation, as for instance a disk modulation (the classical description only allowed 1D modulation, as for instance ring modulation). To ensure the energy conservation, we finally enforce that the weighting function has a unitary 1-norm:

[TABLE]

Regarding the optical Fourier filtering stage, it is based on the fact that the focal plane corresponds to the reciprocal space (i.e. the spatial frequencies space) of the phase (or pupil) space. This fact is directly linked to Fraunhofer’s diffraction. Subsequently, a mask placed in this plane acts like a spatial frequencies filter. A relevant way to describe this mask consists in using its transparency function $m$. For a more detailed description of this quantity, the curious reader may refer to [1].

Describing mathematically the filtering process consists in getting the field in the detector’s plane, we call it $\psi_{d}(\vec{\boldsymbol{\alpha}})$ since it depends on the modulation phase. Optical propagation laws indicate that this field equals to the convolution product between the pupil plane field and the Fourier transform of the transparency function of the mask:

[TABLE]

where $\hat{.}$ is the Fourier transform and $\star$ the 2D convolution product. The last step consists in converting the detector’s field into photo-electrons. Physically, this detection consists in the integration during the modulation cycle of the squared modulus of $\psi_{d}$. The resulting intensity, which depends on the phase $\phi$, equals to:

[TABLE]

In the wave front sensing context, it is relevant to give the dependence of this intensity with respect to the phase. A way to efficiently describe this dependence has been shown in detail in [1]. It consists in doing a Taylor’s development on the phase term:

[TABLE]

which allows then to decompose the intensity into phase-constant, linear, quadratic, cubic, etc. terms:

[TABLE]

We note that the constant term does not contain any information about phase since it does not depend on $\phi$. Consequently, it is common practice to apply a return-to reference operation on the intensity which consists in numerically subtracting $I_{\text{constant}}$ to $I$. The resulting quantity which may be seen as the WFS’s output is called differential intensity:

[TABLE]

The first phase dependence of $\Delta I$ is thus the linear one. We note that this numerical operation is easy to do in practice since $I_{\text{constant}}$ equals to the intensity when there is no turbulent phase, i.e. $I(0)$. Finally, we specify that this return-to-reference operation is in practice associated to a flux normalization which is not necessary here since the flux has been assumed unitary.

3 Linear Model

A standard approximation consists in assuming that the WFS works in its linearity regime. Mathematically, it means that the differential intensity only contains the linear term:

[TABLE]

Calculating the linear regime of a WFS, i.e. the phase subspace where non-linear terms may be neglected, is a challenging task which requires the study of the WFS’s dynamic range but it is not the topic of this paper. Subsequently, we consider that the linear regime of a WFS essentially corresponds to the small phases domain:

[TABLE]

Such an assumption does make sense in Adaptive Optics context since most of the wave front sensing is done in closed loop, i.e. when phase-to-be-measured are residuals of the atmospheric turbulent phase.

3.1 Linear intensity

Within the small phases approximation framework, WFS’s output equals the linear intensity. We give an explicit expression of it initialy developed in [1]:

[TABLE]

where $\bar{.}$ means the complex conjugate. The linear intensity is a 2D map corresponding to the perfectly linear response of the WFS regarding to the phase $\phi$. By using the tip/tilt nature of the modulation (4), it becomes possible to reveal the spatial variations of $I_{\text{linear}}$:

[TABLE]

where $\vec{\textbf{R}}$ is the position vector in the intensity plane and $\vec{\textbf{r}}$ and $\vec{\textbf{r'}}$ the two integration variables coming from convolution products. We realize that, by integrating along tip/tilt amplitudes $\vec{\boldsymbol{\alpha}}$, the 2D Fourier transform of the weighting function appears:

[TABLE]

Such a result has already been observed for the 1D Pyramid WFS with linear tip/tilt modulation in [3]. Its 2D generalization for any kind of modulations and masks is here allowed thanks to the fact that weighting function is considered as depending on tip/tilt amplitude vector $\vec{\boldsymbol{\alpha}}$ and not on the time variable.

3.2 Fourier based WFS’s kernel

At this point, (17) is sufficiently developed to be interpreted as an integral transform performed on the input phase. Indeed, we can write

[TABLE]

where the Kernel K equals to:

[TABLE]

We observe that such a kernel depends on the WFS optical characteristics. More precisely, we see that it actually depends on ”pupil plane” functions: $\hat{m}$, $\hat{w}$ and $\mathbb{I}_{P}$. Moreover, it is worth noticing that this kernel may be also understood as the continuous interaction matrix with respect to the natural basis of the direct phase space, i.e. the ”Dirac phase basis”. Such a fact is obvious when looking at (18). This property of the Kernel is fundamental since the decomposition of an arbitrary phase on the Dirac phase basis is absolutely trivial. In other word, K is a very general descriptor of a WFS since it allows to compute naturally the WFS’s response to any set of phases.

Finally, since most of the Fourier-based WFSs work without modulation, we give the Kernel when the tip/tilt modulation mirror generates a motionless and centered focal spot. In this case, the weighting function is a Dirac function. Subsequently, its Fourier transform equals to the identity function. The resulting kernel is:

[TABLE]

In that specific case, we observe that the matrix associated to this kernel may be understood as an Hadamard’s product between a vertical matrix and a circulant one (see Appendix C).

As a conclusion we will remember that the Kernel of a Fourier based Wave Front Sensor may be written as a real function of three ”pupil plane” quantities and two spatial vectors, one in the input space the other in the output space.

[TABLE]

In other words, for a given optical configuration, K is a continuous 2D matrix. Such a fact explains why the Kernel requires in general, onerous simulations in terms of computational time. Using it to lead an exhaustive study of Fourier based WFSs seems thus unrealistic. Moreover, despite the elegance of its expression, the input/output relation (18) remains difficult to invert directly; K has to have additional properties to be exploited effectively. That is the goal of the next part where we will work with convolutional Kernels, i.e. Kernels which are not function of $\vec{\textbf{R}}$ and $\vec{\textbf{r}}$ but only of their difference $\vec{\textbf{R}}-\vec{\textbf{r}}$.

4 Convolutional model

From now, we will study the consequences of having a Kernel which would be convolutional. In that case, it exists a function, that we call IR for Impulse Response which allows to write:

[TABLE]

The input/output relation of the WFS becomes convolutional:

[TABLE]

It is worth noticing that convolutional systems are also called shift invariant. Indeed, we observe that if $\mathcal{T}$ is a translation operator, the convolution implies that $I_{\text{linear}}(\mathcal{T}[\phi])=\mathcal{T}[I_{\text{linear}}(\phi)]$. In other words, knowing the WFS’s response for a given phase is enough to know the response to any translation of this phase. Such a result explain why a convolutional Kernel may be more rapidly computed: the associated matrix is pure circulant, i.e. a strongly redundant matrix. On the other hand, (23) clarifies the notation ”Impulse Response”, indeed this quantity corresponds to the WFS response when the phase is a pure impulse, i.e. a centered Dirac.

4.1 Required assumptions

Since Kernel is not in general convolutional, we are now interested on the assumptions needed to have such a Kernel. Our starting equation is (19). By taking a closer look at it, we observe that the main difficulty comes from the combined presence of functions depending on $\vec{\textbf{r'}}$ and $\vec{\textbf{r'}}-\vec{\textbf{r}}$ variables; this fact is an obvious obstacle to shift invariance. In other words, converging toward the convolutional model consists in eliminating one of two previous dependency. To do so, we make a Taylor’s development of the pupil indicator function $\mathbb{I}_{P}|_{\vec{\textbf{r'}}}$ around $\vec{\textbf{r'}}-\vec{\textbf{r}}$:

[TABLE]

where $\vec{\nabla}$ is the 2D gradient operator. Such a development allows to write the Kernel as:

[TABLE]

The convolutional model consists in keeping the first term (which has a convolutional form) and neglecting the next ones (which do not). In most cases, such a simplification is not exact and constitute an approximation of the linear model. We call this assumption the sliding pupil approximation since it corresponds to the following assumption:

[TABLE]

Nevertheless, we note that (26) is not an assumption in two cases. First one is a bit unrealistic but worthy of mention: it is the infinite pupil case. Indeed we observe that $\mathbb{I}_{P}=\mathbb{I}$ implies a pure convolutional kernel. The second case is much more interesting since it shows how it is possible to choose the tip/tilt modulation in order to improve the accuracy of the convolutional model. It consists in using a weighting function $w$ which ensures:

[TABLE]

where $\partial\mathbb{I}_{P}$ corresponds to the area where the entrance pupil has discontinuities. In other words, a Fourier-based WFS is a convolutional sensor if the Fourier Transform of the weighting function $\hat{w}$ is barely null on the edge of the pupil. Such a condition is actually quite simple to reach. For a circular pupil and a ring modulation, for instance, $\hat{w}$ is the first Bessel function of the first kind: $J_{0}$. In order to ensure (27), we need adapt the modulation radius in such a way that the pupil radius corresponds to a Bessel function’s zero. An example is shown on Fig. 2. In other words, it would be possible to improve the accuracy of the convolutional model on the condition to properly set the modulation radius.

4.2 Impulse Response and Transfer function

Assuming the sliding pupil approximation, it is possible to identify the WFS’s Impulse Response IR in (25):

[TABLE]

The Impulse Response is a concise way to characterize a Fourier-based WFS while assuming it is shift invariant. We observe that the computational cost which is required to calculate the IR depending on the optical parameters (see Fig. 3) is absolutely negligible compared to calibration matrices’ approaches. Finally, we emphasize on the fact that the circulant matrix $\textbf{IR}|_{\vec{\textbf{R}}-\vec{\textbf{r}}}$ may be seen as the response of the sensor with respect to the Dirac phase basis, i.e. the natural basis of the phase space.

To get the response of the WFS in the phase spatial frequencies space, i.e. in a focal plane, we just have to consider the 2D Fourier transform of the impulse response that we call Transfer Function TF:

[TABLE]

For the sake of clarity, we do not give here the full expression of this important quantity. The reader is referred to the paragraph 5.3 for a comprehensive formulation of the TF depending on the optical parameters. Indeed, we will firstly generalize the convolutional model to more realistic wave front sensing contexts.

5 Wave front sensing in presence of residual phases

In this part, we complicate the sensing context by studying a typical Adaptive Optics case : the presence in the pupil plane of residual phases which are distinct from the phase-to-be-measured. They may be static (as for instance the vast majority of Non Common Path Aberrations (NCPA)) or dynamic (as for instance uncorrected turbulent phase). In the first case, we call these terms $\phi_{\text{s}}$ for static residuals. In the same way, we call the dynamic ones $\phi_{\text{d}}$ and we assume their statistical behavior to be known. In particular, the notation \big{\langle}.\big{\rangle} means ”average on the residual phases statistics”.

The residual phases are taken into account exactly as the tip/tilt modulation, i.e. by adding phases to the incoming field. We thus replicate the calculations of paragraph 3.1 while assuming the sliding pupil approximation (26) and using the linearity of the average operation. The WFS’s kernel associated to this new situation is:

[TABLE]

where functions $\textbf{R}_{s}$ and $\textbf{R}_{d}$ characterize the static and dynamic Residual phases:

[TABLE]

We observe unfortunately that residual phases break the convolutional model: it is impossible to write (30) as a convolutional relation without doing supplementary assumptions.

5.1 Dynamic residuals

We firstly consider dynamic residuals. In an Adaptive Optics context, we can assume that they are Gaussian-distributed and zero-mean phase, it allows us to simplify the function $\textbf{R}_{\text{d}}$ into:

[TABLE]

where D is the structure function of the residual phases [4]. This quantity is much more easy to compute and well-known concerning the full or residual atmospheric turbulence. In order to go even further, we can use the eventual stationarity of the dynamic residuals which says that the structure function only depends on the distance between the points where it is computed:

[TABLE]

Thus, the set of assumptions ”Gaussian-distributed + zero-mean + stationarity” allows to make $\textbf{R}_{d}$ function of $\vec{\textbf{r'}}-\vec{\textbf{r}}$ only. In other words, if these assumptions are valid, sensing in presence of dynamic residuals is compatible with the convolutional kernel. It also means that obtaining K does not require doing a large number of statistical realizations : knowing the structure function is enough. Such a fact implies a significant acceleration of numerical simulations. Finally, it is worth noticing that to take dynamic residuals into account, we just have to change the weighting function by using the structure function D:

[TABLE]

Such a fact is, afterwards, quite logical. The tip/tilt modulation itself may be understood as a ”dynamic residual phases”; its weighting function is a way to characterize its statistics.

5.2 Static residual

We are now interested in static residual. We observe in (31) that $\textbf{R}_{s}$ does not depend on $\vec{\textbf{r'}}-\vec{\textbf{r}}$ only. Thus the kernel is not convolutional due to the presence of static residual. To be able to use the convolutional model, we perform a Taylor’s development of the phases difference in (31) around $\vec{\textbf{r'}}-\vec{\textbf{r}}$:

[TABLE]

where $\mathcal{D}$ is a differential operator defined as:

[TABLE]

We observe that keeping the first term in (37) makes $\textbf{R}_{s}$ a function of $\vec{\textbf{r'}}-\vec{\textbf{r}}$ only. Generally this approximation is not valid but it is worth noticing that the only static residual phases for which it is not an approximation are pure tip/tilt aberrations. Physically, it means that the only static residuals which are rigorously compatible with the convolutional model are tip/tilt. This is not that surprising: adding such a static phase in the pupil plane is absolutely equivalent to using a shifted weighting modulation. Nevertheless, it can be shown that the first term of (37) are predominant compared to the next ones as soon as $\phi_{s}$ only contains low order aberrations. To summarize, if this ”low order approximation” is valid the convolutional model becomes usable again; taking into account static residual just consists in adapting the weighting function:

[TABLE]

5.3 Effective Modulation

We are now ready to give the Impulse Response in presence of dynamic and static residual phases. To do so, we examine the convolutional Kernel when the three following assumptions are valid:

• •

”Sliding pupil approximation” about the pupil geometry.

• •

”Gaussian-distribution + zero-mean + stationarity” about the dynamic residuals statistics.

• •

”Low order static aberration approximation” about static residual phase.

[TABLE]

At first sight, this Kernel depends on 5 distinct quantities: 3 related to the parameters of the optical Fourier filtering system ($m$,$w$,$\mathbb{I}_{P}$) and 2 which characterizes the residual phases ($\textbf{D},\mathcal{D}\phi_{s}$) but it is actually possible to reduce drastically this complexity by introducing a concise quantity that we call $\omega$ and which is defined via its Fourier transform:

[TABLE]

Physically, $\omega$ may be seen as an effective modulation weighting function when taking into account the finite size of the pupil and the residual phases. To be convinced of this, we will lead a small mental experiment by including just before the filtering mask, an imaginary camera. This one will allow us to visualize in the focal plane the effects of the tip/tilt modulation and residual phases. We start with the most basic sensing context, i.e. an infinite pupil without any residual phases. In that case, the fictive camera will record the weighting function $w$, exactly as if $w$ were the shape of an extended object. Thus the effective weighting function rigorously equals to the tip/tilt weighting function, cf. (42). Unfortunately, this case is not realistic: the pupil cannot have an infinite size, this one necessarily windows the incoming field which will result in a blurring of the tip/tilt weighting function. More precisely, $w$ is convolved with the Point Spread Function associated to the pupil. In the Fourier space, it exactly corresponds to (43). If we add now a static residual $\phi_{s}$, the Point Spread Function will be degraded. Thus, a new term has to be added in $\omega$; it is the (44). We can repeat this argumentation for dynamic residual phases only (45). The most general case corresponds to (46); it allows to understand the sensing in presence of dynamic and static residual phases as a tip/tilt modulation using the effective weighting function $\omega$. An example is given on Fig. 5 which shows how to compute the effective weighting function when sensing in open loop on a typical turbulent atmosphere. We may observe a blurring of the original weighting function due to the low-pass filtering induced by the $\mathbb{I}_{P}e^{-\frac{1}{2}\textbf{D}}$ term. Physically, the residual phases act like a supplementary modulation.

[TABLE]

In conclusion, to characterize a WFS in the convolutional model, a generic Impulse Response depending on the mask $m$ and the effective weighting function $w$ only is sufficient:

[TABLE]

To take into account the tip/tilt modulation, the pupil geometry and the residual phases, we just have to use (43) to (46). Such a result is also valid for the Transfer Function which has a concise expression when depending on $m$ and $w$.

[TABLE]

Indeed, we observe that the TF only requires basic mathematical operations: conjugation, convolution and double Fourier transform (i.e. symmetry operation). It is worth noticing that (48) is even more clear if the optical system is centro-symmetric. (It is for instance the case of the 4-faces Pyramid WFS with a ring modulation and isotropic residual phases.) As a matter of fact, we get:

[TABLE]

The computational simplicity of (48) or (49) is especially interesting because of its link with the WFS’s sensitivity regarding to the phase spatial frequencies (for more details, see appendix A). Indeed this quantity is the essential to know how noises are propagating in AO loops, see for instance [5] or to build efficient convolution based phase reconstructors [6].

6 Conclusion

This paper is dedicated to a natural extension of the theoretical work done in [1] about Fourier-based WFSs. Using the mathematical framework describing these sensors, we showed that they may be understood as 2D integral operators which are completely characterized by their Kernel.

In the first part, we provided an explicit expression of this quantity depending on the WFS’s optical parameters: pupil geometry, weighting function of the tip/tilt modulation and filtering mask and showed that it can be understood as the continuous calibration matrix with respect to the Dirac phase basis. Nevertheless, in spite of the generality of the kernel’s expression, this one remains difficult to calculate numerically because of the computational complexity introduced by the tip/tilt modulation.

The particular case of convolutional kernels, which correspond to highly redundant circulant matrices, solves this problem. In this case, the WFS’s output equals the convolution product between the phase-to-be-measured and a quantity called the Impulse Response of the sensor. It was found that this quantity directly depends on the optical parameters and then allows to calculate performance criterion like the sensor’s sensitivity with respect to the spatial frequencies in a much faster way than existing numerical simulations. The convolutional model then represents an ideal tool to compare and optimize in an exhaustive way the Fourier based WFSs, as for instance the very abundant Pyramid WFSs class. Indeed, if such a study has already been done for the filtering masks ($n$-faces Pyramids [7], Axicon [8, 9, 10, 11]; Flattened Pyramid [12]; masks with manufacturing errors, etc.) the tip/tilt modulation still remained to be explored due to the computing time it required. Thanks to the convolutional model, all the 1D modulations (ring, square, etc.) but also 2D modulations (disc, Gaussian, etc.) are easily taken into account. This also opens the way to random modulation or sensing with extended objects as Laser Guide Stars.

In Appendix B, we showed how the convolutional model could be applied to the classical 4-faces pyramid when using the slopes maps rather than differential intensity as output of the WFS. We also proved that the sensitivity derived from the convolutional model is in total agreement with experimental measurements [13].

Since the convolutional model is an approximation, we paid special attention to the required assumptions. In particular, we showed that the infinite pupil approximation [14, 15, 6] is not the only way to justify the use of that model: a weaker assumption called the sliding pupil assumption allows to be closer to the real model since it does take into account the pupil geometry’s. These results suggests that advances may be done in deconvolution-based phase reconstruction by adapting algorithms to sliding pupil approximation.

Finally, we extended the convolutional model to more realistic sensing contexts. The goal was to know how to adapt the convolutional Kernel in presence of residual phases which differ from the phase-to-be-measured. Two cases were addressed: static residual to handle for instance Non Common Path Aberrations and dynamic residuals to take into account non-corrected turbulent phases. We showed that under realistic assumptions about residuals statistics, one only has to modify the modulation weighting function, as for instance by introducing the structure function of the residual phases. These results are an important step in the theoretical understanding of wave front sensing in presence of residuals. Indeed, they demonstrate that these ”non pure tip/tilt phases” which vary in the pupil plane during the sensing may still be described via the tip/tilt modulation weighting function. Consequently, the convolutional model remains valid when sensing with residuals.

This result opens the path to detailed studies of the so-called ”optical gain”. These ”optical gains” correspond to variations of the Pyramid WFS sensitivity in presence of partially-corrected turbulent phases. A solution to this problem consists in adjusting modal gains in the AO loop to compensate the loss of sensitivity [16]. The determination of the modal gains can be derived from empirical methods [17], and the convolutional model now provides an alternative approach to estimate the optimal optical gain directly.

Funding

This project has been supported by the LABEX FOCUS, the WOLF ANR project (grant No 183707), the Action Spécifique Haute Résolution Angulaire (ASHRA) of CNRS/INSU co-funded by CNES and the European Union’s Horizon 2020 dearch and innovation programme under grant agreement (grant No 730890). It is also partly funded by the French Aerospace Lab (ONERA) in the frame of the VASCO Research Project and by the Laboratoire d’Astrophysique de Marseille (LAM). C. Correia received support from A*MIDEX (project no. ANR-11- IDEX-0001- 02) funded by the “Investissements d’Avenir” French Government program, managed by the French National Research Agency (ANR).

Acknowledgments

Victoria Hutterer (Johannes Kepler University, Linz) is acknowledged for her careful rereading and her judicious advice.

Appendix A Sensitivity with respect to phase spatial frequencies

This appendix aims to show the link between the Transfer Function and the sensitivity with respect to phase spatial frequencies. Such a basis describes phases in terms of Sines and Cosines:

[TABLE]

where $\vec{\textbf{k}}$ may be seen as a position vector in the spatial frequencies space. Usually, the sensitivity is defined thanks to a scalar which quantifies the ratio between WFS’s output and input. For instance, the sensitivity with respect to a phase mode $\phi$ is defined as:

[TABLE]

We showed in [1] that such a definition was consistent with [18] and [19] approaches. In our case, a more convenient way to represent the sensitivity consists in grouping the sine and cosine sensitivities into a unique number associated to $\vec{\textbf{k}}$.

[TABLE]

where we purposely did not write a scalar coefficient allowing to normalize these phases on the pupil support with respect to RMS norm. From now, we will assume that the convolutional model is valid to simplify this expression. We have for instance for the cos term:

[TABLE]

where we used the Plancherel theorem. Moreover, we know that

[TABLE]

As a consequence,

[TABLE]

[TABLE]

which gives

[TABLE]

With $|\widehat{\mathbb{I}_{P}}|^{2}$, we identify the Point Spread Function the diffraction limited imaging system. Since the indicator pupil function is real, the PSF is necessarily centro-symmetric. The sensitivity with respect to spatial frequencies finally equals to:

[TABLE]

This formula explains why the convolutional model is so efficient to explore the abundant world of Fourier-based WFSs. Indeed, it proves that it is possible to get very quickly an idea of the WFS’s the response with respect to spatial frequencies; there is no need to build an interaction matrix via time demanding end-to-end codes, using our knowledge of the system’s optical parameters ($m$, $w$ and $\mathbb{I}_{P}$) and the residuals nature (D and $\mathcal{D}\phi_{s}$) in (48) and (60) is enough.

Appendix B Application to the Pyramid WFS

In this appendix we apply the previous developments to the most famous and used [20, 21, 22, 23] Fourier-based WFS which is the Pyramid Wave Front Sensor. We study the following configuration: a 4-sided pyramid with an apex angle large enough to completely separate the pupil images that this optical device produces. In order to handle independently these 4 pupil images, we assume that they are generated by 4 Fourier based WFSs using 4 different filtering masks $m_{1}$, $m_{2}$, $m_{3}$ and $m_{4}$. Each of them corresponds to one quadrant of the cartesian tessellation performed by the Pyramid, see Fig. 6. Moreover, we assume that the effective modulation function $\omega$ is identical for these 4 WFSs.

[TABLE]

where

[TABLE]

The set of the differential intensities $\{\Delta I_{i}\}$ may constitute the Pyramid WFS’s output. Nevertheless historically we prefer to combine them to define two new signals called the slopes maps:

[TABLE]

The slopes maps have many advantages. They firstly improve the linearity range of the PWFS but they also condense in a smaller signal the phase information. Finally they allow to understand physically how the Pyramid performs the wave front sensing: as a matter of fact, $S_{x}$ (resp. $S_{y}$) may be seen as the phase derivative in the spatial frequencies space along the $x$-axis (resp. $a$-axis). Since (66) and (67) transform linearly the differential intensities, it is still possible to associate to them two impulses responses. Indeed, if we define

[TABLE]

[TABLE]

the WFS’s input/output relation remains convolutional:

[TABLE]

In order to check the validity of the convolutional model, we give on Fig. 7 the slopes maps when calculating for the linear (left) and convolutional (middle inserts) models. We choose an input phase which follows a typical atmospheric turbulence. We ensure that the WFS is working in its linearity range. We may observe on the residuals maps (right insert) that these models are in good agreement except at the edge of the pupil. Such a fact confirms the predictions of (27) which says that the sliding pupil approximation is not valid on the pupil discontinuities. These results are another proof in favor of phase reconstructors based on deconvolution.

We are now interested in the PWFS’s sensitivity with respect to the spatial frequencies. To do so, we just have to know the transfer functions associated to the slopes maps. These ones may be calculated directly from $\textbf{IR}_{x}$ and $\textbf{IR}_{y}$ but it is also possible to give the TFs’ dependency with respect to the masks and the effective weighting function thanks to (48):

[TABLE]

[TABLE]

It is then possible to use (60) to compute from $\textbf{TF}_{x}$ and $\textbf{TF}_{y}$ the sensitivities associated to $S_{x}$ and $S_{y}$ with respect to the spatial frequencies:

[TABLE]

We give on Fig. 8 these sensitivities for two modulation radii: 0.5 $\lambda/D$ on the left inserts and 2 $\lambda/D$ on the right ones. Color scale goes from 0-black which corresponds to non-seen frequencies to 1-white which codes the best seen frequencies.

For the slopes map $S_{x}$, we observe that the best sensitivity lies around the $x$-axis while the worst is around the $y$ one. In another word, $S_{x}$ does measure spatial frequencies in the $x$-direction but is blind to those in the $y$-direction. Fortunately, the other slopes map $S_{y}$ is complementary. The association of the two slopes maps is therefore able to all the spatial frequencies. Nevertheless, it is worth noticing that low spatial frequencies, i.e. those with are in the midst of the images, have degraded sensitivity in both slopes maps. The number of these badly seen frequencies is even increasing with the modulation radius. We also observe that the transition frequency between this degraded response and the flat-gray response exactly corresponds to the modulation radius. Such results are in perfect agreement with previous works [3] about the link between sensitivity and modulation.

It is also worth mentioning that the total sensitivity is maximum on the neighborhood of the edges of the pyramid mask. Such a fact is not surprising since we know from Foucault and his knife-edge test [24] that wave front sensing efficiently works where discontinuities are. We note that this area of maximum signal grows with the modulation radius. It comes from the fact that the modulated Point Spread Function enlarges with the modulation and thus communicates more with the edges of the pyramid mask.

The convolutional model is therefore able to reproduce well-known behaviors of the Pyramid WFS but it also constitutes a significiant improvement in the theoretical understanding of this widely used sensor. As a matter of fact, it is the only 2D model able to accurately mimic the cross structure which as been observed experimentally [13].

Appendix C From Integral transform to Matrix transform

We give in this appendix some elements allowing to understand the link between kernels and matrices. First, we use the lexicographic order for all the 2D images: phase, mask, weighting function, differential intensity, etc. handled in this article. As a consequence, 2D variables become 1D variable:

[TABLE]

With such notations, an integral transform may be written:

[TABLE]

In order to make appear a matrix transform, we discretize the spatial variables. Previous integral thus becomes a sum:

[TABLE]

We obviously identify a matrix operation:

[TABLE]

where K may be understood as:

[TABLE]

$x$ is thus the input space variable while $X$ is the output space variable.

Thanks to this parallel, we may interpret Kernels in the matrix formalism. A typical decomposition of Kernel is for instance:

[TABLE]

V, H and C are 1-variable functions and K and U have 2-variables. $\times$ is the usual scalar multiplication. In the matrix formalism, previous equation may be written:

[TABLE]

where $\circ$ is the Hadamard product which corresponds to a ”coefficient by coefficient” multiplication:

[TABLE]

Thanks to this interpretation we may detail the structure of V, H and C matrices (there is nothing to say about the Unspecified matrix U). H is a matrix only depending on the input variable $x$, it is thus an horizontal matrix:

[TABLE]

In the same way, V only depends on the output variable $X$. Is is a vertical matrix:

[TABLE]

C depends on $X-x$. In the matrix formalism, such a type of matrix is called circulant. It has the following form:

[TABLE]