Far-field linear optical superresolution via heterodyne detection in a higher-order local oscillator mode

The Rayleigh limit has so far applied to all microscopy techniques that rely on linear optical interaction and detection in the far field. Here we demonstrate that detecting the light emitted by an object in higher-order transverse electromagnetic modes (TEMs) can help achieving sub-Rayleigh precision for a variety of microscopy-related tasks. Using optical heterodyne detection in TEM01, we measure the position of coherently and incoherently emitting objects to within 0.0015 and 0.012 of the Rayleigh limit, respectively, and determine the distance between two incoherently emitting slits positioned within 0.28 of the Rayleigh limit with a precision of 0.019 of the Rayleigh limit. Extending our technique to higher-order TEMs enables full imaging with resolution significantly below the Rayleigh limit in a way that is reminiscent of quantum tomography of optical states.

Introduction. Since the invention of the optical microscope, there has been a quest for enhancing its resolution. The Rayleigh criterion [1] establishes the minimum resolvable distance in a direct image of a pair of sources to be limited by diffraction according to d R = 1.22λ/2NA, where λ is the wavelength and NA is the numerical aperture of the objective lens. In the past century, a number of techniques for circumventing the Rayleigh limit have emerged [2]. These methods rely, for example, on using shorter-wavelength radiation [3], near-field probing [4], nonclassical [5] or nonlinear [6] optical properties of the object or switching the object's emission on and off [7][8][9]. However, these approaches are often expensive and not universally applicable. Therefore finding a linear-optical microscopy technique that is operational in the far-field regime remains an important outstanding problem.
A promising approach to addressing this problem is by detecting the light emitted by the object in higherorder transverse electromagnetic modes (TEMs). A point source emits primarily into the fundamental TEM 00 mode. However, when the emitter is displaced from the center of the fundamental mode, higher-order TEMs are illuminated. In this way, null measurements of small displacements are possible. For example, homodyne detection in TEM 01 has been used to detect small displacements of a laser beam [10] and to tracking particles with a nanometer resolution [11]. These experiments were performed for single emitters only.
More recently, Tsang et al. showed theoretically that sub-Rayleigh distances between two identical point sources can be estimated by measuring the photon count rate in TEM 01 [12,13]. Remarkably, for distances below the Rayleigh limit, the per-photon uncertainty of this measurement is much less than that associated with direct imaging. This feature is particularly valuable when the number of photons the object can radiate is limited, such as in the case of photobleaching.
Inspired by these developments, we overcome the Rayleigh limit by means of heterodyne detection, tak-ing advantage of the fact that the heterodyne detector is only sensitive to the electromagnetic field in the mode of the local oscillator (LO) [14]. We apply the technique in a variety of settings. First, we determine the positions of single objects emitting coherent and incoherent light. Second, we determine the distance between two identical incoherent objects separated by a distance below the Rayleigh limit. In both cases, we demonstrate a measurement precision significantly below that limit.
In addition, we utilize a mathematical analogy between decomposing an image into TEMs and representing the quantum state of a harmonic oscillator in the Fock basis to propose a new microscopy technique. By measuring intensity and phase of the object's emission into all TEMs one can completely reconstruct its image, in principle with an arbitrarily high precision. In an experiment, the imaging resolution depends on the number of TEMs in which the detection can be realized.
Concept. Consider an optical microscope objective lens that is used to image a plane object with transverse spatial field distribution E(x). The field distribution in the image plane is then given by the convolution where T (·) is the transfer function of the objective lens and we have assumed the magnification to be unity for convenience (Supplementary). The transfer function can be approximated by a Gaussian with the width σ ≈ 0.21λ/NA [15]. The narrower the aperture in the lens, the wider the transfer function, and the stronger the distortion of the image. The heterodyne detector generates a current that is proportional to the overlap between the LO and the sig- where E LO (x ) is the spatial profile of the LO. The LO is prepared in the TEM 01 mode such that the corresponding fundamental TEM 00 mode is matched to the image of a point source located at x = 0: For this conceptual discussion, we assume that the source is a point located at position x p , so that E(x) = δ(x − x p ), in which case we have which for a Gaussian transfer function reduces to J(x p ) = 1 2σ x p e −x 2 p /8σ 2 and the corresponding electronic power The signal vanishes at x p = 0, enabling sensitive null measurement of the source position with respect to the central point [10]. Especially useful is the enhancement associated with determining the distance d between two incoherent point sources. Suppose the sources are located at x = ±d/2, where d is the distance between them. The signal in TEM 01 is given by P (d/2) + P (−d/2), which is proportional to d 2 in the leading order. The intensity in the direct image of this object, on the other hand, is given by I(x ) = S(x +d/2)+S(x −d/2), where S(·) = |T (·)| 2 is the point spread function of the microscope. Approximating . The effect on the image due to the separation of the slits is also of the second order in d, but with a macroscopic zeroth-order background. Any noise in this background will have a deleterious effect on the measurement precision of the term of interest.
Experiment. The experimental setup is shown in Fig. 1. Both the signal and LO are obtained from a home-made external cavity diode laser, producing ∼ 45 mW at a wavelength of λ = 780 nm. We prepare the LO in the desired TEM by transmitting it through a temperature stabilized monolithic cavity with a finesse of about 275 [16]. The laser frequency is locked to the cavity resonance by means of the Pound-Drever-Hall technique [17]. The signal beam passes through an acousto-optic modulator operating at 40 MHz in order to reduce the effect of flicker noise in balanced detection. The modulated beam is collimated to an about 5 mm diameter and sent to a diaphragm with four pairs of slits of 0.15 mm width whose centers are separated by d = 0.25, 0.50, 0.75 and 1.00 mm (3B Scientific U14101). The power of the beam transmitted through the slits is ∼ 200 µW.For the measurements with incoherent light, we place a white paper library card before the slits. During the data acquisition, the card is moved in the transverse plane by a motorized translation stage to achieve averaging over the incoherent light statistics. The incoherent optical power behind the slits is ∼ 10 µW. After the slits, the light propagates in free space for L = 84 cm and passes through an iris diaphragm. The diameter of the diaphragm is measured using an optical microscope as 0.8 ± 0.1 mm and independently estimated with a higher precision from the fit to the data as 0.87 ± 0.01 mm. This latter value corresponds to a numerical aperture of NA= 0.52 × 10 −3 and the Rayleigh distance of d R = 0.912 mm.
The field transmitted through the diaphragm is refocused by an objective lens and subjected to heterodyne detection by means of a balanced detector (ThorLabs PDB150A-SP). To align the detector, we first match the signal's mode to the LO prepared in TEM 00 . Subsequently, the monolithic cavity temperature is changed to transmit the TEM 01 mode and the LO is adjusted slightly to minimize the interference with the signal. The output photocurrent from that detector is observed with a spectrum analyzer set to a zero span mode at 40 MHz with a resolution bandwidth of 1 kHz and a video bandwidth of 10 kHz. An average of 100 traces is acquired for each measurement.
Results and discussion. We first measure the position of a single light source. In this case, the beam only passes through a single 0.15-mm slit, whose position x p is controlled by a translation stage. We conduct the measurement for both coherent and incoherent light.
In the coherent case, our setting resembles that of Hsu et al. [10]. We acquire the data for each point in Fig. 2(a) once, except for a set of points around x p = 0 shown in the inset. For these latter points, the signal value is acquired ten times to estimate the stochastic experimental error. Two types of experimental imperfections con- tribute to that error. First, the orthogonality between the LO and signal modes is imperfect and fluctuates due to air movements and vibration of the optics. As a result, the magnitude of the signal power fluctuation at x p = 0 is comparable to the signal itself. Second, the power of the LO fluctuates due to the instability of the lock of the laser to the cavity resonance. These fluctuations translate directly into those of the detector's output, so the errors at high signals are proportional to the signal.
In view of this analysis, we model the error of our measurement to behave as σ 2 = c 2 + (gP ) 2 , where c is the uncertainty due to the mode mismatch fluctuations and bP is due to the fluctuations of the LO. The actual experimental behavior of the error is consistent with this model. We use this model to find the uncertainty of estimating the slit position x p from the signal power. To this end, we calculate the Fisher information as a function of and find that this function is maximized at x p = 0.012 mm. The value of the Fisher information at this point corresponds to the Cramér-Rao uncertainty bound [18] of δx p = 1.4 µm, almost three orders of magnitude below the Rayleigh limit (Supplementary).
In the experiment with incoherent light [ Fig. 2(b)], the signal's spatial structure is a speckle pattern, changing in time as we move the library card. The intensity of the signal field in the mode being detected is then governed by thermal statistics, so large fluctuations, whose magnitude is on the scale of the signal, are present even for large signals. While the effect of these fluctuations is reduced due to averaging, it is still the dominant source of error for high signal powers. For low signal powers, similarly to the coherent case, the contribution to the uncertainty associated with non-constant mode matching between the signal and LO becomes significant. The general behavior of the experimental error is therefore still consistent with σ 2 = c 2 + (gP ) 2 , albeit with different values of c and g compared to the coherent case. Experimentally, we evaluate these coefficients by acquiring the heterodyne signal at each x p three times. We find that the Fisher information is maximized and the uncertainty of x p is minimized for x p = 0.14 mm, corresponding to δx p = 11 µm.
Next, we measure the distance between two incoherent light sources. Each slit pair is centered on the laser beam and the output signal is sampled twelve times to estimate the error. The resulting data are shown in Fig. 3(b). The error analysis for this setting is similar to that for a single incoherent slit and yields a minimum of δd = 17 µm for d = 0.18 mm. It is interesting to compare the precision of our measurement with what can be achieved by conventional imaging. The images of pairs of slits acquired with a conventional CCD camera are shown in Fig. 3(a), bottom row. For separations significantly below the Rayleigh limit, these images do not resolve the slits and cannot be used to determine the separations. Our technique, on the other hand, permits this determination with the precision comparable to the camera pixel size (7 µm).
At a more quantitative level, the above error analysis implies that the Fisher information decreases to zero, and the estimation uncertainty tends to infinity, for d → 0. This is a common feature of both our technique and conventional imaging, arising because of a nonvanishing uncertainty in evaluating the signal at d = 0 (Supplementary). The fundamental reason for this feature in both cases is the shot noise. In practice, the limitation in our experiment is the fluctuation of the mode matching between the LO and the signal, which can be reduced by using more stable optics. Furthermore, the shot noise level can be reduced by using squeezing [10,11]. The conventional imaging technique, on the other hand, is often limited by technical fluctuations of the electronic signal produced by individual pixels in the camera, which complicates precise measurement of minute image width variations associated with varying d. Note that the techniques of Refs. [12,13] utilize photon counting in higher order modes and do not suffer from the fundamental precision limitation associated with the shot noise. These methods are however still vulnerable to practical noise sources such as detector dark counts.
Application to imaging. Performing heterodyne detection in higher-order Hermite-Gaussian modes TEM 0n permits reconstruction of the full image of the object with sub-Rayleigh resolution. To see this, we write the heterodyne detector output photocurrent (3) as where J(x) is the photocurrent in response to a point source at x given by Eq. (4). We assume, as previously, that the transfer function is given by Eq. (2) while the LO is in TEM 0n of width σ: where H n (·) are the Hermite polynomials. Integral (4) then corresponds to a Weierstrass transform of that poly-nomial and is given by a remarkably simple expression We see that, for objects of size σ, photocurrent J 0n gives the nth moment of the field in the object plane.
The set of photocurrents acquired for multiple modes can be further utilized to find the decomposition of E(x) into the Hermite-Gaussian basis and thereby reconstruct the full image of the object with a sub-Rayleigh resolution. Let α kn be the coefficients of the Hermite polynomial of degree k, so that H k (x/2σ) = n α kn (x/2σ) n . Then, according to Eqs. (6) and (8), we have (9) Because Hermite-Gaussian functions form an orthonormal basis in the Hilbert space of one-dimensional functions, it follows that Knowing all β k , we can calculate E(x). This approach is reminiscent of representing a quantum state of a harmonic oscillator as a superposition of Fock states, whose wavefunctions in the position basis are given by Hermite-Gaussian functions. We name the above-described imaging technique Hermite-Gaussian microscopy, or HGM. Acquiring heterodyne photocurrents for a sufficiently high number of TEMs in principle allows HGM to reconstruct the image with arbitrarily high resolution. The first few tens of TEMs, which are attainable in experimental optics, permit significant improvement of imaging with respect to the Rayleigh limit, as evidenced by Fig. 4.
The above conceptual description can be readily extended to practically relevant cases. Two-dimensional imaging is possible by scanning over both indices of TEM mn and measuring the photocurrent J mn for each pair (m, n) up to a desired maximum. Acquiring both the sine and cosine quadratures of the heterodyne photocurrent permits phase-sensitive reconstruction of the object's field.
A somewhat less trivial extension is to incoherent images. In this case, the output power of the heterodyne detector is given by (see Supplementary for the derivation). Similarly to the coherent case, we obtain moments of the field distribution in the object plane. However, now we obtain only the even coefficients of the decomposition of I(x) into the Hermite-Gaussian basis. Therefore the information about only the even component of function I(x) is retained. An image reconstructed with these data will be a sum of the original intensity profile I(x) with the "ghost" image I(−x). For two-dimensional microscopy, three ghost images, I(x, −y), I(−x, y) and I(−x, −y), will be added to the original image I(x, y). Their effect can be eliminated by placing the entire object into a single quadrant of the x-y plane. Summary. We have used heterodyne detection in the TEM 01 Hermite-Gaussian mode to obtain sub-Rayleigh precision in the measurement of the position of single coherent and incoherent sources, as well as that of the separation between two incoherent sources. With numerical apertures below 10 −3 , our measurement precision is on a scale of a few microns. If our technique is used with state-of-the art microscopes with NA ∼ 1, precision on nanometer scales can be expected. By utilizing higher Hermite-Gaussian modes, the technique can be extended to full imaging with sub-Rayleigh resolution.
Note added. While working on this manuscript, we became aware of similar research being pursued by Sheng at al. [19] and Tham et al. [20].

SUPPLEMENTARY INFORMATION
Theoretical prediction for the signal. Here we calculate the heterodyne detector signal in response to the electromagnetic field generated by an object of a specific shape. We used this method to obtain the theoretical curves in Figs. 2 and 3 of the main text. The calculation is for one-dimensional objects, but is readily extended to two dimensions.
We start by briefly reviewing Abbe's microscope resolution theory.Ẽ where k ⊥ is the orthogonal component of the wavevector and constant normalization factors are neglected throughout the calculation. The objective lens is located in the far field at distance L from the object plane. The position X in the lens plane is then related to k ⊥ according to The lens, in turn, generates the inverse Fourier image in its focal plane: whereT (k ⊥ ) is the transmissivity of the lens as a function of the transverse position in its plane. If this transmissivity is a constant, corresponding to an infinitely wide lens, the image is identical to the object: E (x ) = E(x). If the lens is of finite width, the image is distorted according to where image of the lens. In other words, the image is a convolution of the object with T (·). Heterodyne detection of the image field yields the electronic signal given by Eq. (3) in the main text. If the field E(x) is coherent, that equation is sufficient to calculate the output signal.
If the object is incoherent, we must take the average of the signal over all realizations of E(x) to find the output power of the heterodyne detector photocurrent: For an incoherent image, and hence Our goal is to determine the heterodyne detector signal as a function of the object configuration, E(x) in the coherent case and I(x) in the incoherent case. These calculations are simplified if we take into account a specific shape of the transmissivity function associated with a round diaphragm: T (k ⊥ ) = θ(k ⊥max − |k|) with k ⊥max = 2πR/(Lλ) according to Eq. (S2), R = 0.4 ± 0.05 mm is the radius of the diaphragm and θ(·) is the Heaviside step function. In the Fourier domain this translates into where J 1 (·) is a Bessel function. This is further approximated as T (x −x) ≈ e −(x−x ) 2 /4σ 2 with σ ≈ 0.21λL/R = 0.31 ± 0.03 mm. The LO field in the TEM 00 mode is optimized to match the mode E(x ) in the coherent case for E(x) = δ(x), so E LO (x ) = e −x 2 /4σ 2 . Subsequently it is switched to the TEM 01 mode with E LO (x ) = x e −x 2 /4σ 2 . Substituting this mode into the expression (3) in the main text for the coherent case, we find (S11) For the incoherent case, the mean signal power (S8) equals While the above results are valid for coherent and incoherent objects of any shape, a simple analysis leading to Eq. (5) in the main text is sufficient for the purposes of our experiment, as evidenced by Fig. S1. The only visible difference between the two models is that the curve calculated for the incoherent case using the complete model does not reach zero at the slit position x p = 0. This is because an incoherent slit of a finite width, which can be seen as a combination of multiple mutually incoherent point sources positioned around x p = 0, makes a nonzero contribution to TEM 01 . Error analysis. Our experimental setup yields the electronic signal power P (x p ) with root mean square (rms) uncertainty σ(P (x p )) as a function of the slit position x p . Suppose the task is to estimate x p from the observed power. Below, we present a method for determining and minimizing the error of this estimation.
We can treat the observed power as a random variable whose probability distribution f (P, x p ) = 1 √ 2πσ e −[P −P (xp)] 2 /2σ 2 (S13) is Gaussian with rms width σ centered on the mean power P (x p ). The problem then reduces to that of estimating parameter x p from this random variable. The uncertainty of this estimation can be found using the Cramér-Rao error bound,δx p ≥ 1/F , where is the Fisher information. In the neighborhood of x p = 0, the power P (x p ) can be approximated as ax 2 p , while the uncertainty σ 2 = c 2 + [gP (x p )] 2 as discussed in the main text, with constants a, c and g evaluated from the experimental data. Accordingly, we find for the Fisher information in the limit g 1 F = 4a 2 θ 2 a 2 g 2 θ 4 + c 2 . Next, we determine the value of x p where the Fisher information is maximized so the measurement is the most sensitive. Taking the derivative of F , we find that the maximum value F max = 2a/gc is reached at x p = c/ag. It is these values that we use to evaluate the estimation uncertainties of x p in the main text.
Derivation of Eq. (10) in the main text. We use Eq. (S7) to write the output power of the heterodyne detector as the average