Phase retrieval is the process ofalgorithmically finding solutions to thephase problem. Given a complex spectrum${\displaystyle F(k)}$, of amplitude${\displaystyle |F(k)|}$, and phase${\displaystyle \psi (k)}$:

${\displaystyle F(k)=|F(k)|e^{i\psi (k)}={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}f(x)e^{-2\pi ik\cdot x}dx}$

wherex is anM-dimensional spatial coordinate andk is anM-dimensionalspatial frequency coordinate. Phase retrieval consists of finding the phase that satisfies a set of constraints for a measured amplitude. Important applications of phase retrieval includeX-ray crystallography,transmission electron microscopy andcoherent diffractive imaging, for which${\displaystyle M=2}$.^[1] Uniqueness theorems for both 1-D and 2-D cases of the phase retrieval problem, including the phaseless 1-Dinverse scattering problem, were proven by Klibanov and his collaborators (see References).

This sectiondoes notcite anysources. Please helpimprove this section byadding citations to reliable sources. Unsourced material may be challenged andremoved.(September 2024) (Learn how and when to remove this message)

Here we consider 1-Ddiscrete Fourier transform (DFT) phase retrieval problem. The DFT of a complex signal${\displaystyle f[n]}$ is given by

${\displaystyle F[k]=\sum _{n=0}^{N-1}f[n]e^{-j2\pi \frac{kn}{N},}=|F[k]|\cdot e^{j\psi [k]}k=0,1,\dots ,N-1}$,

and the oversampled DFT of${\displaystyle x}$ is given by

${\displaystyle F[k]=\sum _{n=0}^{N-1}f[n]e^{-j2\pi \frac{kn}{M},}k=0,1,\dots ,M-1}$,

where${\displaystyle M>N}$.

Since the DFT operator is bijective, this is equivalent to recovering the phase${\displaystyle \psi [k]}$. It is common recovering a signal from its autocorrelation sequence instead of its Fourier magnitude. That is, denote by${\displaystyle \hat{f}}$ the vector${\displaystyle f}$ after padding with${\displaystyle N-1}$ zeros. The autocorrelation sequence of${\displaystyle \hat{f}}$ is then defined as

${\displaystyle g[m]=\sum _{i=max\{1,m+1\}}^N{\hat{f}}_i\overline{{\hat{f}}_{i-m}},m=-(N-1),\dots ,N-1}$,

and the DFT of${\displaystyle g[m]}$, denoted by${\displaystyle G[k]}$, satisfies${\displaystyle G[k]=|F[k]{|}^2}$.

The error reduction is a generalization of theGerchberg–Saxton algorithm. It solves for${\displaystyle f(x)}$ from measurements of${\displaystyle |F(u)|}$ by iterating a four-step process. For the${\displaystyle k}$th iteration the steps are as follows:

Step (1):${\displaystyle G_k(u)}$,${\displaystyle {\varphi }_k}$, and${\displaystyle g_k(x)}$ are estimates of, respectively,${\displaystyle F(u)}$,${\displaystyle \psi }$ and${\displaystyle f(x)}$. In the first step we calculate theFourier transform of${\displaystyle g_k(x)}$:

${\displaystyle G_k(u)=|G_k(u)|e^{i{\varphi }_k(u)}=\mathcal{F}(g_k(x)).}$

Step (2): The experimental value of${\displaystyle |F(u)|}$, calculated from the diffraction pattern via the signal equation^{[clarification needed]}, is then substituted for${\displaystyle |G_k(u)|}$, giving an estimate of the Fourier transform:

${\displaystyle G_k^{\prime }(u)=|F(u)|e^{i{\varphi }_k(u)},}$

where the ' denotes an intermediate result that will be discarded later on.

Step (3): the estimate of the Fourier transform${\displaystyle G_k^{\prime }(u)}$ is then inverse Fourier transformed:

${\displaystyle g_k^{\prime }(x)=|g_k^{\prime }(x)|e^{i{\theta }_k^{\prime }(x)}={\mathcal{F}}^{-1}(G_k^{\prime }(u)).}$

Step (4):${\displaystyle g_k^{\prime }(x)}$ then must be changed so that the new estimate of the object,${\displaystyle g_{k+1}(x)}$, satisfies the object constraints^{[clarification needed]}.${\displaystyle g_{k+1}(x)}$ is therefore defined piecewise as:

where${\displaystyle \gamma }$ is the domain in which${\displaystyle g_k^{\prime }(x)}$ does not satisfy the object constraints. A new estimate${\displaystyle g_{k+1}(x)}$ is obtained and the four step process is repeated.

This process is continued until both the Fourier constraint and object constraint are satisfied. Theoretically, the process will always lead to aconvergence,^[1] but the large number of iterations needed to produce a satisfactory image (generally >2000) results in the error-reduction algorithm by itself being unsuitable for practical applications.

The hybrid input-output algorithm is a modification of the error-reduction algorithm - the first three stages are identical. However,${\displaystyle g_k(x)}$ no longer acts as an estimate of${\displaystyle f(x)}$, but the input function corresponding to the output function${\displaystyle g_k^{\prime }(x)}$, which is an estimate of${\displaystyle f(x)}$.^[1] In the fourth step, when the function${\displaystyle g_k^{\prime }(x)}$ violates the object constraints, the value of${\displaystyle g_{k+1}(x)}$ is forced towards zero, but optimally not to zero. The chief advantage of the hybrid input-output algorithm is that the function${\displaystyle g_k(x)}$ containsfeedback information concerning previous iterations, reducing the probability of stagnation. It has been shown that the hybrid input-output algorithm converges to a solution significantly faster than the error reduction algorithm. Its convergence rate can be further improved through step size optimization algorithms.^[2]

Here${\displaystyle \beta }$ is a feedback parameter which can take a value between 0 and 1. For most applications,${\displaystyle \beta \approx 0.9}$ gives optimal results.{Scientific Reports volume 8, Article number: 6436 (2018)}

For a two dimensional phase retrieval problem, there is adegeneracy of solutions as${\displaystyle f(x)}$ and its conjugate${\displaystyle f^{\ast }(-x)}$ have the same Fourier modulus. This leads to "image twinning" in which the phase retrieval algorithm stagnates producing an image with features of both the object and itsconjugate.^[3] The shrinkwrap technique periodically updates the estimate of the support by low-pass filtering the current estimate of the object amplitude (by convolution with aGaussian) and applying a threshold, leading to a reduction in the image ambiguity.^[4]

The phase retrieval is an ill-posed problem. To uniquely identify the underlying signal, in addition to the methods that adds additional prior information likeGerchberg–Saxton algorithm, the other way is to add magnitude-only measurements like short time Fourier transform (STFT).

The method introduced below mainly based on the work of Jaganathanet al.^[5]

Given a discrete signal${\displaystyle \mathbf{x}=(f[0],f[1],...,f[N-1]{)}^T}$ which is sampled from${\displaystyle f(x)}$. We use a window of lengthW:${\displaystyle \mathbf{w}=(w[0],w[1],...,w[W-1]{)}^T}$ to compute the STFT of${\displaystyle \mathrm{f}}$, denoted by${\displaystyle \mathbf{Y}}$:

${\displaystyle Y[m,r]=\sum _{n=0}^{N-1}f[n]w[rL-n]e^{-i2\pi \frac{mn}{N}}}$

for${\displaystyle 0\le m\le N-1}$ and${\displaystyle 0\le r\le R-1}$, where the parameter${\displaystyle L}$ denotes the separation in time between adjacent short-time sections and the parameter${\displaystyle R=\lceil \frac{N+W-1}{L}\rceil }$ denotes the number of short-time sections considered.

The other interpretation (called sliding window interpretation) of STFT can be used with the help of discrete Fourier transform (DFT). Let${\displaystyle w_r[n]=w[rL-n]}$ denotes the window element obtained from shifted and flipped window${\displaystyle \mathbf{w}}$. Then we have

${\displaystyle \mathbf{Y}=[{\mathbf{Y}}_0,{\mathbf{Y}}_1,...,{\mathbf{Y}}_{R-1}]}$, where${\displaystyle {\mathbf{Y}}_r=\mathbf{x}\circ {\mathbf{w}}_r}$.

Let${\displaystyle Z_w[m,r]=|Y[m,r]{|}^2}$ be the${\displaystyle N\times R}$ measurements corresponding to the magnitude-square of the STFT of${\displaystyle \mathbf{x}}$,${\displaystyle {\mathbf{W}}_r}$ be the${\displaystyle N\times N}$diagonal matrix with diagonal elements${\displaystyle \left(w_r[0],w_r[1],\dots ,w_r[N-1]\right).}$ STFT phase retrieval can be stated as:

Find${\displaystyle \mathbf{x}}$ such that${\displaystyle Z_w[m,r]={\left|⟨{\mathbf{f}}_m,{\mathbf{W}}_r\mathbf{x}⟩\right|}^2}$ for${\displaystyle 0\le m\le N-1}$ and${\displaystyle 0\le r\le R-1}$, where${\displaystyle {\mathbf{f}}_m}$ is the${\displaystyle m}$-th column of the${\displaystyle N}$-point inverseDFT matrix.

Intuitively, the computational complexity growing with${\displaystyle N}$ makes the method impractical. In fact, however, for the most cases in practical we only need to consider the measurements corresponding to${\displaystyle 0\le m\le M}$, for any parameter${\displaystyle M}$ satisfying${\displaystyle 2W\le M\le N}$.

To be more specifically, if both the signal and the window are notvanishing, that is,${\displaystyle x[n]\ne 0}$ for all${\displaystyle 0\le n\le N-1}$ and${\displaystyle w[n]\ne 0}$ for all${\displaystyle 0\le }$${\displaystyle n\le W-1}$, signal${\displaystyle \mathbf{x}}$ can be uniquely identified from its STFT magnitude if the following requirements are satisfied:

1. ,
2. ${\displaystyle 2W\le M\le N}$.

The proof can be found in Jaganathan' s work,^[5] which reformulates STFT phase retrieval as the following least-squares problem:

${\displaystyle \underset{\mathbf{x}}{min}\sum _{r=0}^{R-1}\sum _{m=0}^{N-1}{\left(Z_w[m,r]-{\left|⟨{\mathbf{f}}_m,{\mathbf{W}}_r\mathbf{x}⟩\right|}^2\right)}^2}$.

The algorithm, although without theoretical recovery guarantees, empirically able to converge to the global minimum when there is substantial overlap between adjacent short-time sections.

To establish recovery guarantees, one way is to formulate the problems as a semidefinite program (SDP), by embedding the problem in a higher dimensional space using the transformation${\displaystyle \mathbf{X}=\mathbf{x}{\mathbf{x}}^{\ast }}$ and relax the rank-one constraint to obtain a convex program. The problem reformulated is stated below:

Obtain${\displaystyle \hat{\mathbf{X}}}$ by solving:for${\displaystyle 1\le m\le M}$ and${\displaystyle 0\le r\le R-1}$

Once${\displaystyle \hat{\mathbf{X}}}$ is found, we can recover signal${\displaystyle \mathbf{x}}$ by best rank-one approximation.

Phase retrieval is a key component ofcoherent diffraction imaging (CDI). In CDI, the intensity of thediffraction pattern scattered from a target is measured. The phase of the diffraction pattern is then obtained using phase retrieval algorithms and an image of the target is constructed. In this way, phase retrieval allows for the conversion of a diffraction pattern into an image without anoptical lens.

Using phase retrieval algorithms, it is possible to characterize complex optical systems and their aberrations.^[6] For example, phase retrieval was used to diagnose and repair the flawed optics of theHubble Space Telescope.^[7][8]

Other applications of phase retrieval includeX-ray crystallography^[9] andtransmission electron microscopy.

• Phase problem
• Crystallography
• X-ray crystallography
• Coherent diffraction imaging
• Transport-of-Intensity Equation
• Phase correlation

• Klibanov, M. V. (1985). "On uniqueness of the determination of a compactly supported function from the modulus of its Fourier transform".Soviet Mathematics - Doklady.32:668–670.
• Klibanov, M.V. (1987). "Determination of a function with compact support from the absolute value of its Fourier transform and an inverse scattering problem".Differential Equations.22:1232–1240.
• Klibanov, M.V. (1987). "Inverse scattering problems and restoration of a function from the modulus of its Fourier transform".Siberian Math. J.27 (5):708–719.doi:10.1007/bf00969199.S2CID 120840929.
• Klibanov, M. V. (1989). "Uniqueness of the determination of distortions of a crystal lattice by the X-ray diffraction in a continuous dynamical model".Differential Equations.25:520–527.
• Klibanov, M.V. & Sacks, P.E. (1992). "Phaseless inverse scattering and the phase problem in optics".J. Math. Phys.33 (11):2813–3821.Bibcode:1992JMP....33.3813K.doi:10.1063/1.529990.
• Klibanov, M. V.; Sacks, P.E. (1994). "Use of partial knowledge of the potential in the phase problem of inverse scattering".J. Comput. Phys.112 (2):273–281.Bibcode:1994JCoPh.112..273K.doi:10.1006/jcph.1994.1099.
• Klibanov, M. V.; Sacks, P.E.; Tikhonravov, A.V. (1995). "The phase retrieval problem".Inverse Problems.11 (1):1–28.Bibcode:1995InvPr..11....1K.doi:10.1088/0266-5611/11/1/001.S2CID 250916850.
• Klibanov, M. V. (2006)."On the recovery of a 2-D function from the modulus of its Fourier transform".J. Math. Anal. Appl.323 (2):818–843.doi:10.1016/j.jmaa.2005.10.079.

• Crystallography
• Mathematical physics
• Mathematical chemistry
• Inverse problems

• Articles with short description
• Short description matches Wikidata
• Articles needing additional references from September 2024
• All articles needing additional references
• Wikipedia articles needing clarification from January 2022