(Top)
See also
- Faulkner, H. M. L. & Rodenburg, J. M. Movable Aperture Lensless Transmission Microscopy: A Novel Phase Retrieval Algorithm. Phys. Rev. Lett. 93, 023903 (2004).
- Rodenburg, J. M. & Faulkner, H. M. L. A phase retrieval algorithm for shifting illumination. Applied Physics Letters 85, 4795–4797 (2004).
- Rodenburg, J. M. et al. Hard-X-Ray Lensless Imaging of Extended Objects. Phys. Rev. Lett. 98, 034801 (2007).
The Ptychographic Iterative Engine (PIE) is the first iterative algorithm used for ptychography by Rodenburg and Faulkner (1,2). PIE extends the traditional Gerchberg-Saxton algorithm (3,4) to use the multiple probe positions of ptychography. Like the Gerchberg-Saxton algorithm, PIE works by repeatedly moving between the real and reciprocal spaces and applying the respective constraints in each space, solving piecewise for the whole object.
Real space constraint:
\[ \psi_j(\mathbf{x}) = O(\mathbf{x}) \cdot P(\mathbf{x} - \mathbf{X}_j) \]
Reciprocal space constraint:
\[ I_j(\mathbf{k}) = |{\Psi}_j(\mathbf{k})|^2 = | \mathcal{F} [\psi_j(\mathbf{x})] |^2 \]
Notation
- Coordinates
- Real-space coordinates: \(\mathbf{x}\)
- Reciprocal-space coordinates: \(\mathbf{k}\)
- Reconstruction
- The “object” (sample): \(O\)
- The “probe”: \(P\)
- The reconstructions are denoted with a carat: \(\hat{O}\), \(\hat{P}\)
- The \(n\)th iteration is in parentheses: \(\hat{O}^{(n)}\), \(\hat{P}^{(n)}\)
- The scan positions are indexed with \(j\): \(\mathbf{X}_j\)
- The exit waves (“views”): \(\psi_j(\mathbf{x})\)
- The measured intensities: \(I_j(\mathbf{k})\)
Algorithm
Let us outline the basics of the original PIE algorithm given by Rodenburg and Faulkner in (2).
Algorithm summary
- Start with a guess: \(\hat{O}^{(n)}\)
- Create the exit waves: \(\hat{\psi}_j^{(n)}(\mathbf{x}) = P(\mathbf{x} - \mathbf{X}_j) \cdot \hat{O}_j^{(n)}(\mathbf{x})\)
- Propagate: \(\hat{\Psi}_j = \mathcal{F}[\ \hat{\psi}_j\ ]\)
- Correct amplitudes: \(\hat{\Psi}^{\prime}_j = (\sqrt{I_j} / |\hat{\Psi}^{(n)}_j|) \times \hat{\Psi}^{(n)}_j\)
- Backpropagate: \(\hat{\psi}^{(n+1)}_j = \mathcal{F}^{-1}[\ \hat{\Psi}^{\prime} \ ]\)
- Update object according to: \[\hat{O}_j^{(n+1)} = \hat{O}_j^{(n)} + \left[ \frac{|P|}{P_{\mathrm{max}}} \right] \times \frac{P^{*}}{|P|^2 + \alpha} \times \beta \times (\hat{\psi}^{(n+1)}_j - \hat{\psi}^{(n)}_j)\]
- Select different scan position and repeat steps 2 to 6 until error converges \[\epsilon = \sum_j{\left|\left|\ I_j - |\hat{\Psi}_j|^2 \ \right|\right|^2}\]
PIE assumes a known probe \(P(\mathbf{x})\) at the object plane, given probe positions \(\mathbf{X}_j\), and their downstream intensity measurements \(I_j\). Usually these measurements will be at the diffraction plane connected via a Fourier transform.
The algorithm starts with a guess of the object \(\hat{O}^{(n)}\). The exit wave of the \(j_\mathrm{th}\) scan position is calculated with the usual model. (The lowercase \(\hat{o}\) denotes a ‘piece’ of the object illuminated by the probe)
\[ \hat{\psi}_j^{(n)}(\mathbf{x}) = P(\mathbf{x} - \mathbf{X}_j) \cdot \hat{o}_j^{(n)}(\mathbf{x}) \]
The exit wave is propagated to the diffraction plane \(\hat{\Psi}(\mathbf{k}) = \mathcal{F} [\ \hat{\psi}(\mathbf{x})\ ]\), then the reciprocal space constraint is applied.
\[ {\hat{\Psi}^{\prime}_j} = \frac{\sqrt{I_j}}{|\hat{\Psi}_j^{(n)}|}\ \hat{\Psi}_j^{(n)} \]
This “corrects” the amplitudes while keeping the phases of the guessed wave. The corrected wave is then backpropagated to become the new guess for the exit wave.
\[ \hat{\psi}_j^{(n+1)} = \mathcal{F}^{-1}[\ \hat{\Psi}_j^\prime\ ] = P \cdot \hat{o}_j^{(n+1)} \]
Solving for \(\hat{o}\) gives an update equation:
\[ \hat{o}_j^{(n+1)} = \hat{o}_j^{(n)} + \frac{1}{P} (\hat{\psi}^{(n+1)}_j - \hat{\psi}^{(n)}_j) \]
Importantly, PIE tweaks this update equation with a weighting function related to \(P\), and the hyperparameters \(\alpha\) and \(\beta\) :
\[ \hat{o}_j^{(n+1)}(\mathbf{x}) = \hat{o}_j^{(n)}(\mathbf{x}) + \left[ \frac{|P(\mathbf{x} - \mathbf{X}_j)|}{P_{\mathrm{max}}} \right] \times \frac{P^{*}(\mathbf{x} - \mathbf{X}_j)}{|P(\mathbf{x} - \mathbf{X}_j)|^2 + \alpha} \times \beta \times (\hat{\psi}^{(n+1)}_j(\mathbf{x}) - \hat{\psi}^{(n)}_j(\mathbf{x})) \]
- The weighting function \(\left[\ {|P|}\ /\ {P_{\mathrm{max}}}\ \right]\) gives credence to where the illumination is brightest to update those areas the most.
- \(\alpha\) ensures that a division-by-zero does not occur.
- \(\beta\) controls the feedback; i.e. the size of the update at each iteration.
These steps are repeated for the next piece (next index \(j\)) until the sum squared error \(\epsilon\) converges. In practice, it is best to randomize the order of \(j\) for each loop across the scan positions.
\[ \epsilon = \sum_j{\left\|\ I_j - |\hat{\Psi}_j|^2 \ \right\|^2} \]
Simultaneously solving the probe: ePIE
See also
- Maiden, A. M. & Rodenburg, J. M. An improved ptychographical phase retrieval algorithm for diffractive imaging. Ultramicroscopy 109, 1256–1262 (2009).
The simultaneous reconstruction of the probe has been a key part in the robustness of ptychography. Maiden and Rodenburg introduced the extended PIE (ePIE) algorithm in 2009 (5), alongside the development of other probe-solving algorithms such as the Difference Map or Conjugate Gradient Descent.
Steps 1 to 5 are identical to PIE, but a new update equation is used in ePIE. The weighting function uses the intensity instead of the amplitude, which cancels the \(|P(\mathbf{x} - \mathbf{X}_j)|^2\) term in the denominator:
\[ \hat{o}_j^{(n+1)}(\mathbf{x}) = \hat{o}_j^{(n)}(\mathbf{x}) + \frac{{P^{(n)}}^*(\mathbf{x} - \mathbf{X}_j)}{{P^{(n)}}^2_{\mathrm{max}}} \times \beta_O \times (\hat{\psi}^{(n+1)}_j(\mathbf{x}) - \hat{\psi}^{(n)}_j(\mathbf{x})) \]
Importantly, ePIE introduces an update equation for the probe, which is symmetric to the object update:
\[ \hat{P}^{(n+1)}(\mathbf{x}) = \hat{P}^{(n)}(\mathbf{x}) + \frac{{o^{(n)}}_j^{*}(\mathbf{x} + \mathbf{X}_j)}{{o^{(n)}}^2_{j,\ \mathrm{max}}} \times \beta_P \times (\hat{\psi}^{(n+1)}_j(\mathbf{x}) - \hat{\psi}^{(n)}_j(\mathbf{x})) \]
As the problematic denominator has been canceled, \(\alpha\) is no longer needed. The hyperparameters \(\beta_O\) and \(\beta_P\) can be independently chosen for the update size; the original ePIE paper uses \(\beta_O = \beta_P = 1\).
Comment:
Note that the real space constraint \(\psi = O(\mathbf{x}) \cdot P(\mathbf{x} - \mathbf{X}_j)\), implies \(O\) and \(P\) are on equal footing. The only requirement is that \(\psi\) can be separated into a product, and the factors can shift their positions independently. An algorithm cannot identify which of the factors is truly the ‘probe’ and which is the ‘object’. To the algorithm, reconstructing the probe and reconstructing the object is identical, which is reflected in the symmetry of the update equations of ePIE.
In practice, because we want to image a large field of view from a localized probe, we set the dimensions of the probe as the array detector, and calculate the object size from the shift positions. However, from the algorithm’s point of view, it is indistinguishable to having a large ‘probe’ to image a small ‘object’.