Singleshot polychromatic coherent diffractive imaging with a high-order harmonic source

Singleshot polychromatic coherent diffractive imaging is performed with a highintensity high-order harmonic generation source. The coherence properties are analyzed and several reconstructions show the shot-to-shot fluctuations of the incident beam wavefront. The method is based on a multi-step approach. First, the spectrum is extracted from double-slit diffraction data. The spectrum is used as input to extract the monochromatic sample diffraction pattern, then phase retrieval is performed on the quasi-monochromatic data to obtain the sample’s exit surface wave. Reconstructions based on guided error reduction (ER) and alternating direction method of multipliers (ADMM) are compared. ADMM allows additional penalty terms to be included in the cost functional to promote sparsity within the reconstruction. © 2020 Optical Society of America under the terms of the OSA Open Access Publishing Agreement


Introduction
Extreme ultraviolet (XUV) and X-ray microscopy has shown the ability to image samples with nanometer spatial and femtosecond temporal resolutions. While high-resolution X-ray imaging made its breakthrough at synchrotron-radiation sources delivering the required high photon flux, newly developed laboratory sources have become increasingly attractive for such imaging applications. Generally, experiments using laboratory sources are motivated by the better accessibility and by more extensive possibilities to tailor the apparatus to the experiment compared to generic instruments at large-scale facilities. Even more importantly, short-pulse-laser-driven XUV and X-ray sources promise pump-probe time-resolved experiments with pulse durations down to few femtoseconds and typically very low arrival jitter between pump and probe. These sources particularly comprise plasma-based XUV lasers and high-order harmonic generation (HHG) sources. As both sources deliver inherently coherent radiation, coherent imaging methods including holography [1,2], coherent diffractive imaging (CDI) [3][4][5][6][7][8], and ptychography [9][10][11] are preferably applied. These methods have the additional advantage that they dispense of any X-ray imaging optics which typically have a limited photon efficiency and which can introduce aberrations.
Despite the excellent time-resolution delivered by the laboratory sources and despite continuous development of sources [12] and imaging methods, time-resolved imaging experiments have hardly been realized [13]. As an exception, nanosecond oscillations of sub-micron cantilevers have been imaged using an XUV laser source [14]. At HHG sources, single-shot femtosecond-pulse imaging has been demonstrated [8]. The development is mainly impeded by the low average photon flux of all laboratory sources which makes it difficult to realize spatially and time-resolved experiments on typically weakly scattering specimens. One strategy to overcome the fluence limitation is to exploit the photons delivered by the source more efficiently for the imaging process. In particular at HHG sources, the emitted spectrum, typically containing several, almost discrete harmonics of the driving laser, is typically spectrally filtered to a single wavelength for the imaging experiment. The monochromatization discards most of the photons generated in the HHG process. In order to fully utilize HHG sources for coherent imaging, it is important to develop approaches that can relax the temporal coherence and signal-to-noise ratio (SNR) requirements needed for successful reconstructions.
In this paper, we report results from single-shot CDI from ultra-short polychromatic pulses provided by a high-intensity HHG source. Specifically, we use HHG pulses containing several harmonics to image an aperture sample and shot-to-shot wavefront variations. The XUV spectrum at the sample location is characterized in double-slit measurements and additionally compared to values obtained from an XUV spectrometer. From the spectral intensities, we then recover a pseudo-monochromatic diffraction pattern which is fed into the CDI reconstruction algorithm [15]. There have been several experiments which have utilized the polychromatic spectra of the source [15][16][17] with recent progress being made that is based on a two-pulse approach [18].
Plane-wave CDI places the least restrictions on the sample and optics -it only needs the sample to be isolated to ensure a sufficiently high sampling of the measured intensity. To recover the phase lost in the detection process, numerical approaches have been devised which work by projecting the estimate between constraint sets in an iterative fashion [19]. The most common a priori knowledge used in these algorithms is the sample support (support constraint) and the measured diffraction intensities (Fourier or modulus constraint). The estimate is projected iteratively onto all the constraints in sequential order until a solution has been found corresponding to the intersection point of the constraint sets. Due to the nonlinear and non-convex nature of the problem coupled with noisy measurements, finding a solution can be challenging. The single-shot exposures used as input for the CDI algorithm in our experiment are characterized by low SNR. Therefore we here combine conventional CDI with the alternating directions method of multipliers (ADMM) [20]. We adapt ADMM for the CDI reconstruction procedure allowing for additional penalty terms to be used.

Experiment
The experiments were performed at the Intense XUV Beamline at the Lund High-Power Laser Facility, described in detail in [21] and with the relevant parts of the beamline depicted in Fig. 1. The XUV pulses were produced from HHG in argon, driven by an infrared (IR) Ti:sapphire femtosecond laser [22,23]. The vertically polarized driving laser pulses had a central wavelength of 815 nm, a pulse duration of 45 fs, a pulse energy of 25 mJ and a repetition rate of 10 Hz. After the Ar cell the IR and XUV pulses co-propagate until they reach a beam-separator (BS) consisting of two parallel Si plates aligned at the Brewster angle for the IR beam and used to separate the IR and XUV pulses in space. After reflecting off both Si plates both beams emerge parallel to, but slightly lower than, the original pulses with a theoretical transmission of 70 % and 0 % for the XUV and IR pulses, respectively. The XUV pulses are steered and focused onto the sample using two Sc/Si multilayer mirrors designed for 32.5 eV (38.1 nm). The first mirror (M1) is a flat mirror used at an incidence angle close to 45°and the second (M2) is a spherical mirror with 0.5 m focal length used close to normal incidence. Both mirrors are piezo-actuated to allow for precise steering of the beam, and the remaining IR pulses are absorbed by a 200 nm Al filter placed between M2 and the sample. The sample is placed on a three-axis stage controlled by three independent actuators allowing for the sample to be positioned within the focus which has a diameter of ∼100 µm (full width at half maximum). The diffraction data is measured with an in-vacuum Andor CCD camera cooled to −40°with 2048×2048 pixels and 13.5 µm×13.5 µm pixel size. In addition to the cooling we used 2×2 binning of the CCD to improve the SNR.
The beam-separator (BS) and the flat multilayer mirror (M1) reside on a breadboard which can be rotated in or out of the beam path. To obtain initial XUV spectra, the optics are rotated out and the beam passed to a flat-field XUV spectrometer [21]. The optics are rotated back in to acquire double-slit diffraction data. The double-slit is removed and replaced with a puzzle-piece sample. The double-slit and puzzle-piece samples were both fabricated with focused ion beam milling of Si 3 N 4 coated with a 270 nm thick layer of gold. The double-slit sample had 200 nm wide slits with 1, 2, 4, 8, and 16 µm separations. The sample-detector distance was 5 cm and each individual single-shot diffraction pattern was read out before the next pulse.
The ultrafast pulses and XUV wavelengths provided by the source are important for studying sample dynamics on the nanoscale level. The XUV wavelengths provide an inherent resolution improvement over visible light solely due to the short wavelength; however, material absorption in this regime plays a more significant role. Conventional microscope designs in the XUV typically rely on diffractive optics which can be difficult to align, restrict sample size and introduce aberrations into the system. The beam-separator and multilayer mirrors can be rotated out of the beam path allowing the source spectrum and intensity to be optimized prior to imaging. By avoiding image forming optics, this system increases the photon efficiency, provides easy alignment and allows for different sample environments while taking advantage of the source properties.

Spectrum and coherence
We briefly describe the methods used to extract the spectrum and spatial coherence properties. For more comprehensive reviews see [24][25][26]. Due to the detector's relatively long exposure time compared to the oscillations of the electromagnetic fields, any interference between harmonics is negligible and is therefore neglected in the following treatment. As a result, the beam will be treated as several independent fields each with partial spatial coherence. The cross-spectral density can be written as and after making the approximation that the field consists of several discrete wavelengths before the sample, where i, j indicate the harmonic order, U denotes the electric field and is considered as an ergodic random field and indicates an ensemble or time average. For a Young's double slit, the sample transmission is comprised of two rectangular openings t( After propagation through the sample the cross-spectral density becomes If the slits are sufficiently thin, the cross-spectral density can be approximated as being sampled only at the slit centers. In order to obtain an expression for the intensity at the detector, W(x 1 , x 2 , ω) needs to be propagated to the detector plane. Propagation of W(x 1 , x 2 , ω) to the far-field is calculated with a two-dimensional Fourier transform for each spatial coordinate (x 1 , x 2 ∈ R 2 ). The measured intensity is given by whereõ 1,2 (y) denotes the far-field propagation of o 1,2 (x), S i are the spectral intensities, µ 12 is the complex degree of coherence and y ∈ R 2 is a spatial coordinate in the detector plane. It is readily visible that this equation consists of an incoherent sum of several partially coherent far-field diffraction patterns. Upon taking the inverse Fourier transform of Eq. (4) we obtain several autocorrelation (Patterson) functions for each harmonic. In order to understand the autocorrelation function (or the inverse Fourier transform of the measured intensity) it is important to know how the reconstruction pixel size at the sample plane depends on the wavelength. This dependence will be used to extract the spectrum from the double-slit diffraction patterns. This relationship is derived by comparing the expression for propagating a continuous field with the discrete fast Fourier transform. The pixel size in the sample plane (∆ s i ) is related to the wavelength (λ i ), sample-detector distance (z), the detector pixel size (∆ d ) and the number of detector pixels along one direction (N) through the relation In discrete space with a fixed pixel size, each wavelength will give rise to a larger or smaller autocorrelation image depending on the wavelength. Using knowledge of the slit separation (d) it is then possible to determine the contributions that each harmonic makes to the total autocorrelation function. The distance between the i th peak and center of the autocorrelation is used to determine the harmonic spectrum at the sample using where n i is the number of pixels between the center and i th peak in the autocorrelation. The spectral intensities are then determined from the peak heights. From these relations, the absolute wavelengths and spectrum was calculated. This approach is similar to the one used in [27,28]. The multi-harmonic autocorrelation for a 16 µm slit separation is shown in Fig. 2(a) where each harmonic is separated vertically. Figure 2(b) shows a vertical lineout from the autocorrelation for several slit separations. It is apparent that the spectral resolution improves as the slit separation increases and the contributions from each harmonic begin to separate for distances above 4 µm. To confirm the accuracy of the method, the data shown in Fig. 2 is used to simulate the diffraction patterns. The comparison between data and extracted spectra are shown in Fig. 3. Diffraction from 2, 4, 8, and 16 µm slits are shown. The measured data (blue) from several double-slits are compared to the calculated diffraction (black) and appear to agree well.
Additionally, the spectrum is compared to that measured from an XUV flat-field spectrometer. The effect of the multilayer mirrors is readily apparent. An analytic expression was used to calculate the beam-separator transmission while reflectivity measurements were used for the multilayer mirrors resulting in a calculated spectrum shown by the red line in Fig. 4. Lastly, the visibility in the double-slit diffraction data is used to determine the magnitude of the complex degree of coherence and is shown in the inset of Fig. 4. As shown in prior experiments [27] the high-spatial coherence from the IR beam appears to translate to the XUV wavefronts.   4. The measured spectral intensities (data points) compared to the calculated (red line) and spectrum measured from the XUV spectrometer. The spectrometer measures data at a different position on the beamline. The difference in the spectrum can be accounted for through the multilayer mirror reflectivities which have been taken into account in the red curve. The energy uncertainties arise from uncertainty in the measured double-slit-todetector distance. A plot of the average complex degree of coherence magnitude is shown in the inset.
After the spectrum has been determined it is used as input for extracting the monochromatic diffraction pattern. This is accomplished by iteratively minimizing a cost function (J) defined as where I E is the experimental polychromatic intensity data, I is the monochromatic diffraction pattern variable, the residual is defined as R := i S i ψ i I − I E with the scaling operator given by The second term (βP(I)) can be additional penalty terms to increase the L 1 or gradient (total variation) sparsity or others from compressed sensing. The Fréchet derivative with respect to the unknown intensity given by ∇ I J = Σ i S i ψ † i is used to determine the proximal step for the ADMM algorithm [29]. ADMM is an iterative approach which solves the optimization problem by splitting it into smaller subproblems which are themselves easier to solve. The reference wavelength associated with the unknown intensity, λ 0 , can be chosen freely. Here, we chose the 37th harmonic. This does not affect the resolution, but merely samples the reconstruction at a finer scale. Using this approach each ψ i can be unique if the physical model requires. We only require that the adjoint (ψ † i ) be calculable in order to determine ∇ I J. This added flexibility may play an important role in treating each harmonic with its own propagation model [30]. The polychromatic and pseudo-monochromatic data are shown in Fig. 5. The data consists of diffraction from an open aperture shaped as a puzzle piece with a characteristic size of ∼5 µm×5 µm. A scanning electron micrograph of the puzzle piece is shown in the inset of Fig. 5. To combat the low SNR, the diffraction data is convolved with a Gaussian with standard deviation of 0.75 pixels. Typically, this would be a poor decision, but two factors made this a favorable approach. The first reason is that it keeps a simple signal processing pipeline. Penalty terms could have been introduced into the extraction algorithm to smooth the result and reduce the effect of the noise, but this would also introduce a constant which would need to be chosen appropriately. Second, due to several harmonics in the data, it is smoother than its monochromatic diffraction pattern. As a result, the method noise introduced due to Gaussian filtering is reduced [31].

Coherent diffractive imaging
Once the pseudo-monochromatic diffraction pattern has been calculated standard phase retrieval algorithms can be used to recover the phase. We implement a guided ER algorithm similar to [32] for these reconstructions combined with shrinkwrap to determine the exit surface wave and support simultaneously. The approach used six generations with 20 independent trials each with 1000 iterations. In addition, we implement ADMM using the ProxImaL framework [29]. This allows for greater flexibility and ability to incorporate additional penalty terms such as u 1 (L 1 ) or ∇u 1 (total variation (TV)). Due to the nonlinear nature of the CDI optimization problem it is possible for the iterate to get stuck in local minima. To prevent this, we employ a basin-hopping global optimization approach [33,34] for each proximal step of ADMM. We have found that this gives very repeatable reconstructions, but with the disadvantage that the algorithm takes much longer to compute. The time is reduced by optimizing only the pixels which are contained in a bounding box around the support. The reconstruction procedure takes about 1 min 20 sec which is shortened by having a fixed support. In comparison, the reconstructions based on the ER algorithm take slightly less than 5 min which could also be reduced if a fixed support was used. The computations were carried out using a standard computer with four processors running at about 3 GHz which allowed for four reconstructions to be calculated in parallel. Figure 6 shows a comparison between ER and ADMM with L 1 and TV penalty terms. The complex-valued exit surface wave is shown within a single image by mapping the amplitude and phase into value and hue respectively. All reconstructions show a phase ramp which is shown more clearly in Fig. 6(a). Introducing additional penalty terms can help reduce reconstruction artifacts if a beamstop is used or if, as in this case, the SNR is low. There is no missing data regions within our data which explains the similarity between reconstructions. The support is determined by using guided ER with the shrinkwrap method [35]. The L 1 and TV reconstruction have small coefficients of 10 −4 and 10 −3 respectively. The plots in Fig. 6(a) show the amplitude and phase for the TV reconstruction. It's tempting to analyze the sharp boundary for a resolution assessment, but this can be simply a result of the sharp support boundary. Instead we look at a knife edge cut across an edge of an inset where the support is constant. Unfortunately, the edges across these insets are not flat resulting in larger 10 %-90 % distances. Using this approach, the resolution is 3 pixels corresponding to 600 nm. The resolution here is constrained by the low number of detected photons. Even with the low SNR we are able to recover single-shot reconstructions. Figure 7 shows several of these reconstructions which reveal shot-to-shot wavefront fluctuations within the sample opening. Lowering the SNR and coherence requirements allows for a larger set of sources and datasets become viable for single-shot imaging. At the moment; however, we are still restricted to samples that have wavelength-independent transmission functions. Samples with wavelength-dependent absorption will introduce uncertainty into the spectral intensities. This uncertainty, in the simplest situation, will result in a partial recovery of the coherence in the data. To handle these situations, a generalized approach will need to be developed. Several reconstructions of the pseudo-monochromatic data for comparison between ER and ADMM with TV and L 1 penalty terms. The scale andbars apply to all three reconstruction images. The magnitude and phase are mapped to intensity (value) and hue respectively. In this way the complex-valued images can be displayed within a single image.

Conclusion
In conclusion, high-intensity HHG pulses were generated and used for single-shot polychromatic CDI. A direct method based on the cross-spectral density was described to recover the source spectrum from double-slit diffraction data. This method relies on the wavelength dependence of the sample-plane pixel size to separate the autocorrelation functions from each harmonic which is then used to determine the field's spectrum. The low SNR in the diffraction patterns can be compensated for by Gaussian filtering the image prior to recovering the pseudo-monochromatic diffraction pattern. The recovered spectrum was used to obtain a pseudo-monochromatic diffraction pattern and reconstructions based on guided ER and ADMM techniques were compared. The repeatability of ADMM was improved by using a global optimization step which helps to avoid local minima. Finally, single-shot diffraction patterns were used to image a puzzle-piece aperture and show the shot-to-shot wavefront fluctuations of the incident field.