1. INTRODUCTION

Superresolution optical far-field microscopy has undergone a tremendous evolution since roughly two decades ago, when it was shown that the classical resolution limit [1,2] posed by diffraction can be overcome [3,4], resulting in the development of a large variety of methods achieving superresolution. One group of methods relies on stimulated ground- or excited-state depletion and the non-linear response of fluorescence markers to given excitation intensities to deterministically engineer the effective excitation point-spread function (PSF) [3–6]. Other methods stochastically localize single photoswitchable molecules with an accuracy of a few tens of nanometers via centroid fitting of the PSF [7–10]. Another branch of methods makes use of higher-order intensity cross correlations in the Fourier plane [11–13] or autocorrelations in the image plane of a microscope [14–17]. For the latter group of correlation microscopy (CM) techniques, either super-Poissonian bunched-light emission due to statistical fluctuations [14] or sub-Poissonian antibunched-light emission of fluorescence markers can be used to enhance the resolution, both in wide-field [16] and confocal microscopy [17]. Finally, structured illumination microscopy (SIM) leads to a doubled resolution by use of the principle of moiré fringes and linear wave optics [18], and its non-linear derivative saturated SIM (SSIM) leads to further improvements and, in principle, unlimited resolution, though at the cost of necessitating high intensities [19]. Other derivatives combine SIM with the third-order process of CARS or with graphene plasmons to access more higher spatial frequency information than ordinary SIM [20–22]. Note that sub-wavelength phenomena can also be found in other fields of physics, for example, in sub-wavelength atom localization due to the non-linear behavior of coherent population trapping (CPT) [23–25] and sub-wavelength lithography via Rabi oscillations [26,27]. CPT was also proposed to highly increase the resolution in a microscopy-themed derivative [28].

Here, we report on a novel superresolution method that relies on intensity correlation measurements in the image plane of a microscope in combination with structured illumination to image fluorophores that exhibit antibunching. We therefore term it structured illumination quantum correlation microscopy (SIQCM). Linear low-intensity excitation and linear detection suffice such that the technique holds promise to highly enhance the resolution in biological imaging. The detrimental effects due to high intensities that occur with many superresolution techniques, leading to phototoxicity and photobleaching in fluorophores, do not arise. We demonstrate that already very low correlation orders $m$ provide highly enhanced superresolution that scales favorably as $m+\sqrt{m}$. The present manuscript focuses on the highly enhanced lateral resolution using a simple wide-field microscopic setup; however, one can easily extend the scheme, as CM as well as SIM each on their own already provide optical sectioning capability for 3D imaging [16,29,30].

The technique makes use of $m$th-order correlations and antibunched photon emission, inherently present in most common fluorophores, even at room temperature [31–34]. Hence, the required quantum emitters are already broadly in use in fluorescence microscopy. Furthermore, SIM is a well-established technique in biological imaging with commercial microscopes, attaining the theoretically predicted resolution enhancement, being widely spread.

2. THEORY

Let $h(\mathbf{r})$ be the PSF of a given microscope, where $\mathbf{r}$ denotes the position in the image plane, and $H(\mathbf{k})\equiv \text{FT}\{h(\mathbf{r})\}$ is the corresponding optical transfer function (OTF) obtained by a Fourier transform (FT), where $\mathbf{k}$ denotes the spatial frequency in reciprocal space. The full width at half-maximum (FWHM) of the PSF determines the resolution power a microscope provides to discern individual close-by emitters. Narrowing the FWHM of the PSF results in an increased resolution power. Equally, in reciprocal space, the OTF represents the resolution power via the size of the observable region. A larger observable region equals a higher resolution power. Later [cf. Eq. (3)], in the $m$ th-order correlation microscopy signals ${\mathrm{CM}}_m$, the effective PSF reads $h_m(\mathbf{r})\equiv {[h(\mathbf{r})]}^m$, and its corresponding OTF shall be defined as $H_m(\mathbf{k})$. In general, $h_m(\mathbf{r})$ gets narrower for increasing correlation order $m$ and its FWHM approximately scales as $m^{-1/2}$. Vice versa, the radius of the observable region in reciprocal space $H_m(\mathbf{k})$ is increased by the factor $\sqrt{m}$. Microscopes usually possess the circularly symmetric Airy disk ${(2J_1(r)/r)}^2$, with $r=|\mathbf{r}|$ as the PSF [35], which allows one to resolve individual incoherent emitters as individual sources of radiation as long as their separation $d$ is at least on the order of $d\ge \lambda /2$ or, more precisely, $d\ge 0.61\lambda /\mathcal{A}$ [1,2]. We denote the Rayleigh limit as $d_{\mathrm{R}}\equiv 0.61\lambda /\mathcal{A}$, with $\mathcal{A}$ as the numerical aperture of the microscope objective and $\lambda$ as the wavelength of the emitted fluorescent light. Without loss of generality, we assume a magnification of one (or rather minus one) throughout our theoretical treatment such that the coordinates in the object and image planes, $\mathbf{R}$ and $\mathbf{r}$, respectively, can be regarded as equal, i.e., $\mathbf{R}\equiv \mathbf{r}$. Further, the image plane coordinate $\mathbf{r}$ will be given in units of the Rayleigh limit $d_{\mathrm{R}}$, resulting in a dimensionless quantity. This intrinsically connects the image plane coordinate to the resolution power of the imaging system (determined via the wavelength $\lambda$ and the numerical aperture $\mathcal{A}$), which then is given by $\mathbf{r}\equiv 1$.

To measure fluorescence photons in the image plane, the fluorophores in the object plane need to be driven by an excitation light field. In classical linear optics, the fluorophores respond linearly to a given excitation intensity $I_0$. We treat the fluorophores quantum mechanically as a two-level system with ground $|g⟩$ and excited states $|e⟩$; however, this is only the case for intensities $I_0\ll I_{\mathrm{sat}}$, where the saturation intensity $I_{\mathrm{sat}}\propto 1/{\tau }_l^2$ depends on the excited states’ lifetime ${\tau }_l$. The general expression for the intensity emitted by a two-level system driven by a given excitation intensity ($I_0$ in units of $I_{\mathrm{sat}}$) reads [36]

Let us first assume a continuous and spatially uniform excitation illumination in the object plane with intensity $I_{\mathrm{str}}(\mathbf{r},t)=I_0$ and the fluorophore density distribution $n(\mathbf{r})\propto \sum _{i=1}^N\delta (\mathbf{r}-{\mathbf{r}}_i)$ to be comprised of individual point-like sources at positions ${\mathbf{r}}_i$ that emit statistically independent, i.e., incoherent radiation. Note that we can also assign (relative) weights to the independent emitters in case their photon emission rates differ. Differences in (relative) emission rates would be enhanced in the (higher-order) intensity autocorrelations. However, usually, fluorophores emit sufficiently uniformly, and our technique does not require very high correlation orders to achieve highly enhanced superresolution, in contrast to superresolution optical fluctuation imaging (SOFI) [14]. Further, in SOFI, this problem is resolved by using balanced cumulates [37], and our higher-order correlation signals can be adapted accordingly. Therefore, and to keep the analysis illustrative, we consider uniform emission rates here.

Considering a linear response, i.e., $I_0\ll I_{\mathrm{sat}}$, the intensity in the image plane reads

Taking the square of the intensity ${[G^{(1)}(\mathbf{r})]}^2=\sum _{i=1}^N{[h(\mathbf{r}-{\mathbf{r}}_i)]}^2+\sum _{i\ne j}^{N}h(\mathbf{r}-{\mathbf{r}}_i)h(\mathbf{r}-{\mathbf{r}}_j)$, we obtain an incoherent sum of narrowed PSFs $h_2(\mathbf{r}-{\mathbf{r}}_i)$; however, in addition, we also obtain detrimental cross terms between emitters. To understand the effect of the cross terms, consider the simple example of three close-by emitters, equidistantly separated by $0.7d_{\mathrm{R}}$, and the resulting intensity distribution in the image plane schematically depicted by the black line in Fig. 1(a). Due to the closeness, the intensity shows no dips in between the emitter positions, such that the emitters cannot be discerned as individual sources of radiation. The same is true for the squared intensity, shown by the blue line. However, when subtracting the cross terms (see the red line) from the squared intensity, the resulting distribution (green line) displays clear dips in between the emitter positions, since the cross terms’ main contribution arises in between the emitter positions. A simple a posteriori subtraction of the cross terms, though, is not possible, since one would need to know the emitter positions in advance. That is, a priori knowledge would be required.

 

Fig. 1. (a) Three equidistant emitters, separated by $d=0.7d_{\mathrm{R}}$, were chosen, where the black, blue, green, and red curves, respectively, represent the intensity $G^{(1)}(\mathbf{r})$, the intensity squared ${[G^{(1)}(\mathbf{r})]}^2$, the ${\mathrm{CM}}_2(\mathbf{r})$ signal of Eq. (3), and the scaled second-order correlation function $\frac{1}{2}{G}^{(2)}(\mathbf{r})$. (b) The same setup with added structured illumination $I_{\mathrm{str}}(\mathbf{r})=[\frac{1}{2}+\frac{1}{2}\cos ({\mathbf{k}}_0\mathbf{r}+\phi )]$ for phase $\phi =0$, such that the black and green curves now represent Eqs. (4) and (7). For the green curve, the effectively squared structured illumination $I_{\mathrm{str}}{(\mathbf{r}=0.7)}^2=0.65$ addresses the outer emitters, which are then weaker than the central emitter. The latter one can therefore be better resolved. Varying the phase $\phi$ addresses different emitters sequentially. Note that the image plane coordinate $\mathbf{r}$ is given in units $d_{\mathrm{R}}$.

Download Full Size | PDF

However, without any a priori knowledge, the cross terms can be removed by subtracting the second-order intensity autocorrelation function $G^{(2)}(\mathbf{r})\equiv G^{(2)}(\mathbf{r},\mathbf{r})=⟨{\widehat{E}}^{(-)}(\mathbf{r}){\widehat{E}}^{(-)}(\mathbf{r}){\widehat{E}}^{(+)}(\mathbf{r}){\widehat{E}}^{(+)}(\mathbf{r})⟩\propto 2\sum _{i\ne j}^{N}h(\mathbf{r}-{\mathbf{r}}_i)h(\mathbf{r}-{\mathbf{r}}_j)$. Here, the squared terms $h_2(\mathbf{r})$ vanish, as each two-level system can emit only one photon simultaneously, that is, $⟨{\widehat{\sigma }}_i^{+}{\widehat{\sigma }}_i^{+}{\widehat{\sigma }}_i^{-}{\widehat{\sigma }}_i^{-}⟩=0$. Antibunching of the individual emitters is thus the crucial ingredient here. Conducting the subtraction of the two signals, we obtain

 

Fig. 2. Schematic setup to combine SIM and CM to obtain the new SIQCM technique. The independent quantum emitters radiate fluorescent light (red) after being excited by a structured illumination standing-wave pattern (green). The SIQCM signal is obtained by post-processing a series of images that are captured by a CCD camera in the image plane. For details, see the text in Section 4.

Download Full Size | PDF

 

Fig. 3. (a)–(c) Comparison of the observable regions in reciprocal space provided by ordinary SIM, ${\mathrm{SIQCM}}_2$, and ${\mathrm{SIQCM}}_3$ (left to right), where $\alpha =3$, $\alpha =4$, and $\alpha =6$ orientations of the structured illumination pattern were chosen, respectively. The axes are normalized by the modulus of the highest spatial frequency $k_{\max }\equiv 1$ transmitted through the imaging system. The green disks stem from the densities $\tilde{n}(\mathbf{k}\pm {\mathbf{k}}_0)$ [cf. Eq. (5)], while the red and black disks originate from the higher harmonics arising in the ${\mathrm{SIQCM}}_m$ signals [cf. Eq. (8)]. Note that the size of individual disks is enlarged from (a) to (c).

Download Full Size | PDF

Now, considering a two-dimensional structured illumination $I_{\mathrm{str}}(\mathbf{r},t)=I_0[\frac{1}{2}+\frac{1}{2}\cos ({\mathbf{k}}_0\mathbf{r}+\phi )]$ with spatial frequency ${\mathbf{k}}_0$ and adjustable phase $\phi$, a linear response of the fluorophores, and ordinary intensity measurements, one obtains a doubled resolution by the principle of moiré fringes. The illumination pattern and the investigated sample produce beat patterns in the object and the image plane such that initially unobservable spatial frequencies in reciprocal space are shifted by the amount $k_0=|{\mathbf{k}}_0|={(k_x^2+k_y^2)}^{1/2}$ into the observable region and thus can be accessed [cf. Fig. 3(a)]. In general, it is useful to define $I_{\mathrm{str}}=I_{\mathrm{str}}(\mathbf{r},\alpha ,\phi )$, where $\alpha =\tan (k_y/k_x)$ is the orientation of the adjustable illumination pattern. Note that a larger $k_0$ effectively enlarges the observable region in reciprocal space and thus the resolution by a higher amount; however, $k_0$ is limited by diffraction and the given numerical aperture $\mathcal{A}$ of the microscope objective. Hence, by use of far-field wave optics, infinitesimally dense fringe spacings in the source plane cannot be produced. According to the image plane coordinate, we assume dimensionless spatial frequencies (normalized by the highest spatial frequency transmitted by the imaging system) and the illumination pattern to take the maximum value $k_0=k_{\max }=k_{\max }(\lambda ,\mathcal{A})\equiv 1$. Following the derivation for Eq. (2) with the new modulated $I_{\mathrm{str}}(\mathbf{r})$, the resulting signal reads (see also Ref. [18])

One image does not allow one to separate the three independent components, such that three images with three different phases $\phi =0,\frac{2\pi }{3},\frac{4\pi }{3}$ are required, creating the linear system $A\mathbf{n}=\mathbf{G}$, where the matrix $A$ describes the resulting system, and $\mathbf{n}$ denotes a vector with entries $\tilde{n}(\mathbf{k})$, $\tilde{n}(\mathbf{k}-{\mathbf{k}}_0)$, and $\tilde{n}(\mathbf{k}+{\mathbf{k}}_0)$. The vector $\mathbf{G}$ on the right-hand side of the system possesses the entries $\tilde{G}(\mathbf{k},\alpha ,0)$, $\tilde{G}(\mathbf{k},\alpha ,\frac{2\pi }{3})$, and $\tilde{G}(\mathbf{k},\alpha ,\frac{4\pi }{3})$, which represent the Fourier transforms of the (experimentally) measured data. The system is solved by applying the inverse matrix $\mathbf{n}=A^{-1}\mathbf{G}$. To sufficiently cover the enlarged area in reciprocal space, it is necessary to chose at least three orientations $\alpha =0,\frac{1\pi }{3},\frac{2\pi }{3}$ [cf. Fig. 3(a)], resulting in a total of 9 measurements.

Taking a non-linear fluorophore response into account, higher harmonics of $\cos ({\mathbf{k}}_0\mathbf{r}+\phi )$ arise, enabling access to higher spatial frequencies in reciprocal space via SSIM. The arising higher harmonics can be read out easily when plugging $I_0\cos ({\mathbf{k}}_0\mathbf{r})$ into Eq. (1) as the excitation illumination and compiling the Fourier cosine series, which reads

Our new approach combines the strength of both methods (CM and linear SIM) to enhance the already superresolving signals tremendously within the linear low-intensity regime. A schematic sketch of the setup is displayed in Fig. 2. The method utilizes the structured illumination $I_{\mathrm{str}}(\mathbf{r},\alpha ,\phi )$ and considers a linear response of the fluorophores. Then, the intensity as well as intensity autocorrelations $G^{(m)}(\mathbf{r})$ in the image plane of the microscope need to be evaluated. Combining the signals—equal to ${\mathrm{CM}}_2$ in Eq. (3)—we obtain the (second-order) SIQCM signal (illustrated by the green line in Fig. 1(b) for the previously chosen example),

The result of Eq. (8) is depicted in Fig. 3(b), where, due to the additional higher harmonic, images for five different phases $\phi =0,\frac{2\pi }{5},\frac{4\pi }{5},\frac{6\pi }{5},\frac{8\pi }{5}$ per orientation are required. In general, for the $m$ th-order signal, $2m+1$ components arise, such that $2m+1$ images are required per orientation. The need for a large number of orientations $\alpha$ is, however, relaxed due to the enlarged disks. Considering the maximum speed, we chose four orientations $\alpha =0,\frac{1\pi }{4},\frac{2\pi }{4},\frac{3\pi }{4}$, resulting in a total of 20 images to sufficiently cover the highly enlarged observable area. When speed is not the major goal, one can chose more orientations $\alpha$ to obtain a higher quality, which is also considered in regular SIM.

After obtaining the individual Fourier components $\tilde{n}(\mathbf{k})$, $\tilde{n}(\mathbf{k}\pm {\mathbf{k}}_0)$, and $\tilde{n}(\mathbf{k}\pm 2{\mathbf{k}}_0)$, they need to be assembled properly in reciprocal space, that is, by applying the same procedure which is conducted in SIM. The extracted raw components are so far scaled by the circularly symmetrical OTF $H_2(\mathbf{k})$ or by its shifted versions $H_2(\mathbf{k}\pm {\mathbf{k}}_0)$ and $H_2(\mathbf{k}\pm 2{\mathbf{k}}_0)$. To obtain an approximately homogeneous disk, we divide the enlarged observable area [cf. Fig. 3(b)] into subregions and rescale the components by use of a Wiener filter ${\tilde{n}}_{new}(\mathbf{k})=\tilde{n}(\mathbf{k})/(H_2(\mathbf{k})+\gamma )$, where the constant $\gamma >0$ prevents division by zero. In general, the modulus of $\gamma$ depends on the signal-to-noise ratio (SNR) a given measurement provides. After the assembly, we apply a triangular apodization, resembling the Fourier transform of an Airy disk, to the homogeneous disk to reduce ringing in the final image [19]. The final image is obtained by taking the modulus of the inverse Fourier transform of the assembled (and post-processed) disk in reciprocal space. Using the more advanced deconvolution methods proposed and applied in SIM would result in an even further enhanced resolution [39,40].

3. SIMULATIONS

For the simulations, we chose masks with point-like emitters and calculated data as it would be detected by a CCD with discrete and finite pixels. Note that here we are assuming perfect data, i.e., discrete intensity values matching the theoretical calculations without noise. The experimental requirements to obtain the sought-after SIQCM signals with preferably high SNRs, i.e., sufficient statistics for the second- and higher-order correlations, will be discussed later. For rescaling by use of the Wiener filter, we used $\gamma =0.05$. The pixel size and post-processed area were chosen in such a way that the offsets in reciprocal space approximately match the integer numbers, as we considered a real-valued sinusoidal modulation. To remove the necessity for integer numbers in reciprocal space (which is challenging to realize in a real experiment), one can use a complex wave vector in real space [19].

Simulations illustrating the resolution power of the ordinary intensity measurements $G^{(1)}(\mathbf{r})$, ${\mathrm{CM}}_2(\mathbf{r})$, SIM, and ${\mathrm{SIQCM}}_2(\mathbf{r})$ are presented in Fig. 4. For this, three $3\times 3$ arrays with grid constants $d=1.0d_{\mathrm{R}}$, $d=0.5d_{\mathrm{R}}$, and $d=0.29d_{\mathrm{R}}$ and six independent emitters each are imaged by use of the enlisted techniques. Resolving all details of the arrays corresponds to resolving the individual grid constants. The first array is resolved by every method, as the chosen distance corresponds to the classical resolution limit; however, $G^{(1)}(\mathbf{r})$ barely resolves the individual emitters. ${\mathrm{CM}}_2$ provides a moderately increased resolution and SIM the second-best resolution power. We want to point out that even though ${\mathrm{CM}}_2(\mathbf{r})$ and SIM already provide superresolved images, ${\mathrm{SIQCM}}_2$ outperforms both methods by far and provides the highest resolution power. Reducing the source separation to $d=0.5d_{\mathrm{R}}$, only SIM and ${\mathrm{SIQCM}}_2$ can resolve the individual emitters, and, finally, for $d=0.29d_{\mathrm{R}}$, only our new method resolves the array. Resolving the last array corresponds to a resolution improvement of 3.45, approximately matching the theoretical prediction, and indeed, even slightly overcoming it. This can be explained by (i) the phenomenological nature of Rayleigh’s criterion, (ii) the perfect-quality data used in the simulation, and (iii) the fact that the size of the observable region defined by $H_m(\mathbf{k})$ is only approximately adhering to the $\sqrt{m}$ scaling. Additional smaller improvements could be obtained by advanced deconvolution techniques, as pointed out at the end of Section 2.

 

Fig. 4. Comparison of the resulting final images utilizing ordinary intensity measurements $G^{(1)}(\mathbf{r})$, ${\mathrm{CM}}_2(\mathbf{r})$, SIM, and ${\mathrm{SIQCM}}_2(\mathbf{r})$, imaging three $3\times 3$ arrays with grid constants $d=1.0d_{\mathrm{R}}$, $d=0.5d_{\mathrm{R}}$, and $d=0.29d_{\mathrm{R}}$, and six independent emitters each; see the masks at the top. The bar within each mask represents the Rayleigh limit $d_{\mathrm{R}}$. The depicted areas in the final images differ from top to bottom, as the emitters are distributed over a smaller area. However, the areas are not shrunk according to relative distances, as the Airy disk’s size in the intensity measurements $G^{(1)}(\mathbf{r})$ remains the same for each run.

Download Full Size | PDF

To show that the resolution power of our technique scales very favorably as $m+\sqrt{m}$ compared to ordinary CM, which merely scales as $\sqrt{m}$, we also present simulations for third-order SIQCM [see the illustration in Fig. 3(c)]. We chose six orientations $\alpha =0,\frac{1\pi }{6},\frac{2\pi }{6},\frac{3\pi }{6},\frac{4\pi }{6},\frac{5\pi }{6}$, resulting in a total of 42 images as seven phases $\phi$ are required per orientation. In Figs. 5(a)–5(e), the resulting final images for the same $3\times 3$ array with $d=0.29d_{\mathrm{R}}$, imaged by the use of $G^{(1)}(\mathbf{r})$, ${\mathrm{CM}}_3(\mathbf{r})$, and with three different reconstruction approaches for ${\mathrm{SIQCM}}_3(\mathbf{r})$, are presented, together with the corresponding observable region in reciprocal space in Figs. 5(f)–5(i). The emitter geometry that was previously just resolved by ${\mathrm{SIQCM}}_2$ is expectedly not resolved by ${\mathrm{CM}}_3$, but clearly resolved by the ${\mathrm{SIQCM}}_3$ signal. For all three reconstruction approaches, we first subdivided the observable regions into sections and rescaled the individual Fourier components by use of the Wiener filter, resulting in a mostly homogenous disk in reciprocal space. Note that, obviously, the homogeneity is meant in terms of the scaling of the Fourier components, as the specific structure originating from the emitter geometry remains [cf. Figs. 5(h) and 5(i)]. For the approach depicted in Fig. 5(c), the resulting image is simply the modulus of the inverse Fourier transform of the homogenous disk. To remove ringing, in Fig. 5(d), after the inverse Fourier transform, we omitted the small, imaginary parts acquired throughout the numerical evaluation (which should be zero in theory) and cropped negative values, instead of simply taking the modulus. Note, though, that this approach might not readily be used on real data with noise. For the reconstruction method illustrated by Fig. 5(e), we applied triangular apodization to the homogenous disk before calculating the modulus of the inverse Fourier transform to obtain the final image. This method has also been used to produce the images in the last column of Fig. 4. In the given simulation, the approach that crops negative values to remove ringing performs best and provides the smallest FWHM of the effective PSF.

 

Fig. 5. Imaging the same $3\times 3$ array as used in Fig. 4, with grid constant $d=0.29d_{\mathrm{R}}$. The upper illustrations (a)–(e) show the resulting final images in real space, and the lower illustrations (f)–(i) show the corresponding reciprocal space for the images in (a)–(c) and (e). For image (d), the reciprocal space distribution of image (h) was used, however, with a reconstruction strategy different from the one used in (c). The ordinary intensity measurement $G^{(1)}(\mathbf{r})$ is shown in (a), ${\mathrm{CM}}_3(\mathbf{r})$ in (b), and ${\mathrm{SIQCM}}_3(\mathbf{r})$ with three different reconstruction approaches in (c)–(e). In (c), a homogenous disk was used for the reconstruction. In (d) the same disk was used, though negative values were cropped in the real-space distribution before taking the modulus. In (e), a triangular apodization was applied to the homogenous disk, the same procedure that was utilized for ${\mathrm{SIQCM}}_2(\mathbf{r})$ in Fig. 4.

Download Full Size | PDF

4. EXPERIMENTAL PROCEDURE, DATA EVALUATION, AND SIGNAL-TO-NOISE ANALYSIS

To implement the results of the theory and simulation, one would need to consider the following experimental geometry and data evaluation procedure in a wide-field microscopy setup. The first diffraction orders of a coherently illuminated grating are focused into the specimen by a microscope objective, where, due to interference, a standing wave pattern with an adjustable phase and orientation is created (both have to be known for the reconstruction). The subsequent fluorescence emission is captured by the same objective and imaged onto a pixelated detector, e.g., an electron multiplying charged-coupled decive (EM-CCD), intensified CCD (ICCD), or intensified sCMOS camera (IsCMOS) that records intensity frames.

More specifically, these need to be temporally and spatially resolved photon count maps. The time resolution of individual frames can easily reach subnanosecond regimes via gating, such that photons measured within one frame possess full temporal correlation. Another possibility would be to introduce time gating via ultrashort pulses, as previously implemented in CM [16,17]. Noise from the detection process can be mainly circumvented via thresholding the detector, which decreases the quantum efficiency and count rate but increases the SNR. In the work by Schwartz et al. [16], 4 million frames have been acquired at a pulse repetition rate of 1 kHz, which resulted in a 1.1 h measurement time to obtain high-quality third-order correlation maps. However, as the authors already pointed out, the rapid progress in detector technology with respect to low noise, high quantum efficiency, and frame rate will significantly improve practical aspects of quantum correlation microscopy techniques. Our SIQCM technique equally relies on the evaluation of correlations in the image plane, where structured illumination has simply been added to the signal, thus posing the same requirements toward an experimental implementation.

Once the SIQCM maps are obtained with a given SNR, the post-processing in Fourier space—that is, combining images with different phases and orientations into one superresolved image—is exactly the same as in SIM and SSIM. Hence, we simply have to take into account the SNR of the $m$th-order SIQCM signal and the strength of the arising higher harmonics to determine whether all harmonics can be used in the Fourier space assembly or only a limited number. Further, for the given experimental data, it is necessary to determine which correlation order can reasonably be determined from the data with a still sufficiently high SNR. Here, one should find the balance between correlation order $m$ to be evaluated and the SNR, since increasing the correlation order will increase the noise, though the strength of the same higher harmonic is increased as well. This can be seen from Table 1, where, for correlation orders $m=1,\dots ,5$, the theoretically expected coefficients of the higher harmonics are given in relation to the constant term $0{\mathbf{k}}_0\mathbf{r}$, represented by the central blue circles in Fig. 3.

Table 1. Coefficients of Higher Harmonics Relative to the Constant Term $0{\mathbf{k}}_0\mathbf{r}$ for Correlation Orders $m=1,\dots ,5$

View Table

To illustrate the influence of noise, we further conducted simulations where random noise was added. The results are shown in Fig. 6, where the $3\times 3$ array of Fig. 4 with grid constant $d=0.29d_{\mathrm{R}}$ was imaged by use of $G^{(1)}(\mathbf{r})$, ${\mathrm{CM}}_2(\mathbf{r})$, SIM, and ${\mathrm{SIQCM}}_2(\mathbf{r})$. As before, each emitter (represented by a delta peak) is convolved by the Airy disk $(2J_1(|\mathbf{r}|)/|\mathbf{r}|{)}^2\in [\mathrm{0,1}]$ with $(2J_1(|{\mathbf{r}}_i-{\mathbf{r}}_i|)/|{\mathbf{r}}_i-{\mathbf{r}}_i|{)}^2=1$ for emitter $i$. For each pixel, the noise was added by summing a random positive number $\mathrm{\Delta }(\mathbf{r})$ to the theoretically calculated values, which are given by the right-hand sides of Eqs. (2)–(4) and (7). For the simulations, we chose $\mathrm{\Delta }(\mathbf{r})\in [0,0.1]$, $\mathrm{\Delta }(\mathbf{r})\in [0,0.3]$, and $\mathrm{\Delta }(\mathbf{r})\in [0,1.0]$. Note that for the given simulations, we did not alter the reconstruction process, and in particular for the rescaling via the Wiener filter, we still utilized $\gamma =0.05$. For the lowest noise value, the resulting final image is hardly affected by the noise. For the second value, noise is apparent, though the emitter geometry is still clearly reconstructed. For the largest value $\mathrm{\Delta }(\mathbf{r})=1.0$, which constitutes a large amount of the given signal, the emitter positions are barely resolved due the high noise in the final image. This can be explained by (i) the underestimated rescaling constant $\gamma =0.05$, where by rescaling, the noise contributions are enhanced too strongly and (ii) by the use of the outer Fourier disks whose contribution is then on the order of the noise. However, taking care of noise in the reconstruction process, i.e., in the assembly of the individual Fourier components, would enhance the final image quality.

 

Fig. 6. Imaging the same $3\times 3$ array as used in Fig. 4, with grid constant $d=0.29d_{\mathrm{R}}$ and added noise $\mathrm{\Delta }$. The illustrations show $G^{(1)}(\mathbf{r})$, ${\mathrm{CM}}_2(\mathbf{r})$, SIM, and ${\mathrm{SIQCM}}_2(\mathbf{r})$ (from left to right) and noise levels $\mathrm{\Delta }=0.1$, $\mathrm{\Delta }=0.3$, and $\mathrm{\Delta }=1.0$ (from top to bottom), respectively.

Download Full Size | PDF

5. CONCLUSION AND OUTLOOK

We introduced a new quantum imaging technique we call SIQCM, which is based on the profitable merger of linear SIM with antibunching CM. For a linear low-intensity standing wave illumination pattern and linear detection of photon autocorrelations in the image plane of a microscope, our technique provides, in theory, unlimited superresolution, with the improvement scaling favorably as $m+\sqrt{m}$ with the correlation order $m$. Hence, it has the potential to increase the spatial resolution in imaging a variety of samples and, in particular, biological ones. Further, we anticipate the SIQCM concept to be applicable to super-Poissonian bunched-light emission (e.g., used in SOFI, due to the on–off blinking of fluorophores), where autocorrelations in the image plane can be combined into cumulants that equally lead to a signal with a narrowed PSF. Adding structured illumination will not only introduce offsets by $\pm {\mathbf{k}}_{\mathbf{0}}$ but also higher harmonics with offsets up to $\pm m{\mathbf{k}}_{\mathbf{0}}$. The optical sectioning capability provided by CM as well as SIM can also be implemented, enabling three-dimensional imaging with increased axial resolution compared to 3D-SIM [30], due to the higher harmonics along the $z$-axis.

Our new SIQCM approach would bring similar benefits to two-photon microscopy [41], and vice versa. Considering a standing-wave excitation pattern with a wavelength within the red or near-infrared part of the spectrum, short wavelength photons from the UV or blue part of the spectrum are emitted by fluorophores due to two-photon absorption. Since the absorption cross section is inherently dependent on the squared excitation intensity, the resulting effective illumination structure is of the form $I_{\mathrm{str}}(\mathbf{r},t)=I_0{[\frac{1}{2}+\frac{1}{2}\cos (\frac{1}{2}{\mathbf{k}}_0\mathbf{r}+\phi )]}^2$, where $\frac{1}{2}{\mathbf{k}}_0$ corresponds to the Fourier-limited standing-wave pattern for the near-infrared illumination wavelength, and shifts by $\pm {\mathbf{k}}_0$ (corresponding to the fluorescence wavelength), as used in regular SIM, already appear in the fluorescence intensity signal [42]. Evaluating correlations additionally would result in taking the $2m$-th power of the structured illumination, resulting in higher harmonics with offsets up to $\pm m{\mathbf{k}}_0$. The well-known advantages of two-photon microscopy, high penetration depth, energy deposition (and thus photobleaching) only within the vicinity of the focal plane, and inherent optical sectioning capability would be added to our highly improved superresolution.

Using bunched-light emission, our approach should be applicable with the current technology and at a reasonable speed, as SOFI already provides acquisition times of a few seconds. To obtain the sought-after CM and SIQCM signals, autocorrelations in the image plane can also be determined by evaluating the cross correlations of neighboring pixels, which reduces the experimental requirements and introduces an effectively denser sampling in the image plane [15,16]. The latter fact is of practical importance, as resolutions achievable with SIQCM would often exceed the sampling density of CCD cameras in use, and thus, interpolation can be circumvented. Further, autocorrelations and cross correlations can be evaluated in combination to enhance the available statistics [15].