Interferometric superlocalization of two incoherent optical point sources

Ranjith Nair; Mankei Tsang

doi:10.1364/OE.24.003684

1. Introduction

Rayleigh’s criterion for resolution of two incoherent point sources, which asserts that a minimum separation between the sources equal to the diffraction-limited spot size is necessary for them to be resolvable, has perhaps been the most influential resolution criterion in the history of optics despite its heuristic character [1,2]. Based as it is on the eye as a detection instrument, Rayleigh’s criterion neglects the possibility of better resolution using light of greater intensity or using a longer observation period. Pioneering studies on resolution limits as a function both of the source separation and the mean number of observed photons were made by Bettens et al. in [3] and more recently by Ram et al. in [4] (See also [5] in the context of electron microscopy). In [3], it was shown that for a given mean photon number, ideal continuum image-plane photon counting can locate the centroid, i.e., the midpoint, of the two radiating sources with a finite precision (depending on the size of the diffraction-limited spot in the image plane and on the number of photons detected) no matter how small their separation is. Based on the Cramér-Rao (CR) lower bound of estimation theory [6], it was also shown in [3,4] that any (unbiased) estimate of the separation between the sources based on image-plane photon counting must suffer a divergent mean squared error for a given mean photon number as the separation tends to zero. This phenomenon was dubbed Rayleigh’s curse in [7] as it suggests a fundamental limitation in resolving incoherent point sources even when the role of the average detected photon number is taken into account. In the past few decades, several far-field super-resolution techniques that circumvent the Rayleigh limit have revolutionized microscopy (see [8] for a review), but these rely on placing restrictions on the emissions from the sources and do not fundamentally challenge the criterion (or Rayleigh’s curse) for independently and incoherently radiating sources. Other developments, e.g., [9], use clever techniques to approach the performance of ideal continuum image-plane photon counting using imperfect detectors.

In [7], following preliminary work in [10], the problem of resolving two incoherent point sources was approached anew from the perspective of quantum estimation theory using the quantum Cramér-Rao (QCR) bound [11–13]. This bound provides a fundamental limit to the accuracy of estimating the source separation optimized over all possible measurement techniques allowed by quantum mechanics. Under a weak-source assumption, it was found in [7] that the QCR lower bound on the minimum mean squared error (MSE) of estimating the separation of two point sources is independent of that separation. Further, a linear optics-based measurement – spatial-mode demultiplexing (SPADE) – was proposed and shown in principle to approach the quantum bound for all values of the separation. These results are in stark contrast to the performance of image-plane photon counting mentioned above, as the divergent behavior of the minimum MSE with decreasing separation – Rayleigh’s curse – is conspicuously absent.

In this paper, we revisit the resolution problem from the point of view of the semiclassical theory of photodetection [14–16]. For sources emitting thermal radiation, for linear propagation of light through the imaging system, and for linear-optics processing followed by photon counting, the semiclassical theory of photodetection may be used to obtain the same results as a fully quantum analysis would give [14, 15]. This is because the thermal state has a positive P-representation, and its counting statistics may be obtained by statistical averaging over those of coherent states via Mandel’s rule [14]. The SLIVER scheme of this paper as well as the schemes of [7] satisfy the above necessary conditions for the applicability of the semiclassical theory. We emphasize that, for the state of radiation and measurements considered here and in [7], “semiclassical” does not imply “approximate” – the fully quantum and semiclassical treatments are in quantitative agreement.

Besides its greater familiarity to the optics community, a semiclassical analysis of the resolution problem offers several advantages. Firstly, the analysis can be carried out for arbitrary source strengths. Indeed, the work in [3, 4] assumes the sources are weak enough so that the counting statistics are nearly Poissonian. Similarly, the quantum analysis in [7] assumed, for mathematical tractability, that the sources are weak enough so that the state of light could be assumed to be confined to the zero and one-photon subspaces. Thus, a semiclassical treatment will incorporate the effects of multi-photon events neglected in [3, 4] and in the quantum analysis of [7]. Such an analysis would also reveal whether the superresolution effect predicted in [7] for weak sources persists as the source strength is increased. A second advantage of the semiclassical approach is that intuitions from semiclassical optics can be harnessed to both understand the reasons for superresolution as well as to design new measurement schemes which can also be analyzed semiclassically provided that they only involve linear optics and photodetection, possibly along with the use of other classical sources of light.

Accordingly, our contributions in this paper are three-fold:-

We set up the problem of resolving two incoherently radiating point sources in the framework of semiclassical photodetection theory. The sources can be of arbitrary strength and imaging systems with two-dimensional circularly-symmetric point-spread functions are studied. The thermal source model is applicable to a wide range of physical scenarios ranging from optical astronomy to fluorescence microscopy.
We propose a new interferometric scheme for estimating source separation that we call SLIVER (Super Localization by Image inVERsion interferometry) and that yields finite resolution for arbitrarily small values of the source separation and for arbitrary source strengths. Devices employing image inversion have been proposed, studied and demonstrated previously in the microscopy literature for various related applications [17–24], but the fundamental limits on the capability of such devices for improving lateral resolution seem not to have been explored before. Through an analysis of the photodetection statistics of SLIVER in the semiclassical framework, we show both using the Cramér-Rao bound and through explicit simulation that our method manifestly alleviates Rayleigh’s curse.
We offer an explanation for the superlocalization effect of SLIVER at small values of separation in the language of estimation theory applied to the photodetection statistics for thermal light.

2. Source and system model

Consider the depiction in Fig. 1 of the focused image of two incoherent point sources in the image plane of an imaging system coordinatized by ρ = (x,y). We assume that the imaging system is spatially invariant and coordinates have been rescaled so that the image is of unit magnification [25]. We also assume that the (possibly complex-valued) normalized amplitude point-spread function (PSF) ψ(ρ) satisfying

\int d ρ {| ψ (ρ) |}^{2} = 1

is circularly symmetric, so that ψ(ρ) depends only on the magnitude |ρ| for all ρ. It is straightforward make to ψ(ρ) spatially invariant and circularly symmetric in an imaging system using two lenses [25]. The prototypical example of such a PSF is of course the Airy disk pattern resulting from a circular pupil. The images of the sources 1 and 2 in the image plane are assumed to be centered at

\mp \frac{d}{2}

respectively, giving rise in the image plane to a combined field with complex amplitude

E (ρ) = A_{1} ψ (ρ + \frac{d}{2}) + A_{2} ψ (ρ - \frac{d}{2}) .

in normalized units of

\sqrt{photons \cdot m^{- 2}}

. We will assume throughout this paper, except in Section 5.3, that the centroid of the sources, i.e., their midpoint, has already been located, perhaps by image-plane photon counting [4]. It is taken to lie at the origin of the image plane. In Sec. 5.3, we consider the effect of a small error in centroid localization. Our initial results suggest that our measurement scheme is robust to such perturbations.

Fig. 1 A rendering of the focused image of two point sources 1 and 2 located at the tail and tip respectively of the vector d. The blue spots indicate the approximate extent of the point-spread function (PSF) centered at each source. The centroid of the sources is assumed to coincide with the intersection O of the optical axis with the image plane. The focus of this paper is the regime in which d = |d| is smaller than the width of the PSF, in which conventional image-plane photon counting performs poorly in estimating d.

Download Full Size | PDF

The pair of dimensionless complex numbers $A = (A_{1}, A_{2}) \in ℂ^{2}$ are the source amplitudes. The thermal and mutually incoherent nature of the sources dictates that the source amplitudes are circular-complex Gaussian random variables having the first and second moments [16]:-

E [A_{μ}] = 0

E [A_{μ} A_{ν}] = 0

E [A_{μ}^{*} A_{μ}] = ε_{μ}

E [A_{1}^{*} A_{2}] = 0

for µ,ν ∈ {1,2} ranging over the two sources. According to the above relations, the real and imaginary parts of the {A_µ} comprise four statistically independent zero-mean Gaussians of variance ε₁/2 (for the components of A₁) or ε₂/2 (for the components of A₂). Eqs. (3)–(6) are as appropriate for completely incoherent sources with respective strengths (i.e., mean photon number) ε₁ and ε₂, which can assume any values.

Given that the centroid has been located, the parameter of interest is the separation d = |d| between the two sources. The angle between the line joining the sources and the x-axis is another unknown parameter, but the measurement method proposed here works regardless of the value of that angle, which need not be known beforehand in practice. We will thus focus just on the separation d and apply the tools of single-parameter estimation theory in our analysis.

3. The SLIVER measurement and its statistics

The SLIVER measurement, illustrated schematically in Fig. 2 and described in detail below, essentially involves the separation of the image-plane field E(ρ) of Eq. (2) into its symmetric and antisymmetric parts with respect to inversion in the image plane about the optical axis. By definition, the symmetric and antisymmetric parts of E(ρ) are given respectively by

\begin{array}{l} E_{s} (ρ) : = \frac{E (ρ) + E (- ρ)}{2} = \frac{S}{2} [ψ (ρ + \frac{d}{2}) + ψ (ρ - \frac{d}{2})] \\ E_{a} (ρ) : = \frac{E (ρ) - E (- ρ)}{2} = \frac{D}{2} [ψ (ρ + \frac{d}{2}) - ψ (ρ - \frac{d}{2})] . \end{array}

Fig. 2 The SLIVER method for separation estimation:– The input field E(ρ) is separated into its symmetric (E_s(ρ)) and antisymmetric (E_a(ρ)) components with respect to inversion about the centroid of the sources, which is also the optical axis. This separation is accomplished by an image inversion interferometer as shown, where the image inversion box can be implemented by any of a variety of methods (see Sec. 3). The components impinge upon separate bucket detectors that collect spatially-unresolved photocounts N_s and N_a. The photocount record over a series of M detection windows of N_s and (or) N_a is processed to yield an estimate $\hat{d}$ of the separation d.

Download Full Size | PDF

In writing the right-hand sides of Eqs. (7), we have used the fact that ψ(ρ) is circularly symmetric and have defined the sum and difference

S = A_{1} + A_{2}

D = A_{1} - A_{2}

of the amplitudes of the two sources. The relations

\begin{array}{l} E [S] = E [D] = 0 \\ E [S^{2}] = E [D^{2}] = E [S D] = 0 \\ E [S^{*} S] = E [D^{*} D] = ε_{1} + ε_{2} \\ E [S^{*} D] = ε_{1} - ε_{2} \end{array}

hold for these new random variables, which are thus also identically distributed circular-complex Gaussians and statistically independent of each other when ε₁ = ε₂.

The field components of Eq. (7) may be obtained in spatially separated planes using an image inversion interferometer (III) that works by splitting the input field using a 50-50 beamsplitter, spatially inverting the field from one output about the optical axis and recombining the two beams at a second 50-50 beamsplitter. The two beams output from the second beamsplitter then have the field patterns of Eq. (7). The entire setup is thus essentially a balanced Mach-Zehnder interferometer with an inversion device in one arm (Fig. 2). The critical element in such setups is the inversion device – a device that implements the transformation E_out(ρ) = E_in(−ρ). Inversion devices have been proposed and demonstrated for various applications in the fluorescence microscopy community and a considerable literature on image inversion microscopy exists [17–24]. However, to the best of our knowledge, a fundamental statistical analysis of the kind we make here of the capabilities of an III for lateral resolution has not been done (See the discussion in Section 6). The spatial inversion itself can be accomplished in a variety of ways, e.g., using a 4f -arrangement with lenses [18], with a 3-D setup using plane mirrors [19], or using two concave mirrors [20]. We will be concerned here only with the outputs (7) of such a device rather than its detailed implementation, which may vary according to the application.

Suppose now that the output beams from the III are directed to different photodetectors. Conditioned on the values of S and D, the intensity patterns I_s(ρ_s) and I_a(ρ_a) of the fields on their respective detector planes are

I_{s} (ρ_{s}) = \frac{| S |^{2}}{4} [{| ψ (ρ_{s} + \frac{d}{2}) |}^{2} + {| ψ (ρ_{s} - \frac{d}{2}) |}^{2} + 2 I_{int} (ρ_{s}, d)],

I_{a} (ρ_{a}) = \frac{| D |^{2}}{4} [{| ψ (ρ_{a} + \frac{d}{2}) |}^{2} + {| ψ (ρ_{a} - \frac{d}{2}) |}^{2} - 2 I_{int} (ρ_{a}, d)],

where

I_{int} (ρ, d) = Re ψ^{*} (ρ + \frac{d}{2}) ψ (ρ - \frac{d}{2})

is an interference term. Conditioned on the values of S and D, semiclassical detection theory dictates that the photocounts on the two detector planes are independent inhomogeneous Poisson processes with the rate functions I_s(ρ_s) and I_a(ρ_a) respectively [14–16]. It follows that the spatially-unresolved integrated photocounts N_s and N_a at each detector are Poisson random variables with the means

E [N_{s} | S] = \frac{| S |^{2}}{2} \int d ρ_{s} I_{s} (ρ_{s}) = \frac{| S |^{2}}{2} [1 + δ (d)],

E [N_{a} | D] = \frac{| D |^{2}}{2} \int d ρ_{a} I_{a} (ρ_{a}) = \frac{| D |^{2}}{2} [1 - δ (d)],

where

\int d ρ I_{int} (ρ, d) = Re \int d ρ ψ^{*} (ρ) ψ (ρ - d)

= Re \int_{-} \int_{\infty}^{\infty} d x d y ψ^{*} (x, y) ψ (x - d, y)

\equiv δ (d) .

The circular symmetry of ψ(ρ) has been used to get Eq. (17), and we have retained only the separation d as argument. According to Eqs. (10), both |S|² and |D|² are exponentially distributed with mean (ε₁ + ε₂), so the photocounts N_s and N_a integrated over the photodetector surfaces are distributed according to Bose-Einstein distributions [14–16] with means

{\bar{N}}_{s} = \frac{ε_{1} + ε_{2}}{2} [1 + δ (d)] = ε_{ave} [1 + δ (d)],

{\bar{N}}_{a} = \frac{ε_{1} + ε_{2}}{2} [1 - δ (d)] = ε_{ave} [1 - δ (d)],

for ε_ave = (ε₁ + ε₂)/2. If ε₁ = ε₂, the photocounts are also statistically independent because S and D are independent in that case.

The SLIVER method consists of optimal statistical processing of the sequence of observed photocount measurements of N_s and (or) N_a over M independent shots using the same sample to obtain a good estimate $\hat{d}$ of the separation d – see Fig. 2 and Secs. 4 and 5. As we show below, the photocount from the antisymmetric component of the input field is much more informative in this regard.

4. Cramér-Rao bounds for the separation d

4.1. Number-resolved photon counting

Given the Bose-Einstein distribution

P_{s} (n) : = \Pr [N_{s} = n] = \frac{1}{{\bar{N}}_{s} + 1} {(\frac{{\bar{N}}_{s}}{{\bar{N}}_{s} + 1})}^{n},

of number-resolved photon counting of the inversion-symmetric field, we can calculate the fundamental Cramér-Rao (CR) lower bound of estimation theory [6] on the mean squared error (MSE)

E {[\hat{d} (N_{s}) - d]}^{2}

of any unbiased estimator

\hat{d} (N_{s})

of d that is a function of the photocount N_s. An unbiased estimator is one that satisfies

E [\hat{d} (N_{s})] = d

for all d in some interval of interest. The statistical expectations here are taken over the distribution of N_s. The CR bound is itself the reciprocal of the Fisher information

J_{d}^{(sym - pc)}

on the parameter d defined as

J_{d}^{(sym - pc)} = E_{N_{s}} {[\frac{\partial \log P_{s} (n)}{\partial d}]}^{2} .

The superscript indicates that the Fisher information pertains to photon counting of the inversion-symmetric field. This quantity turns out to equal

J_{d}^{(sym - pc)} = \sum_{n = 0}^{\infty} P_{s} (n) {[\frac{\partial \log P_{s} (n)}{\partial d}]}^{2}

= {[\frac{\partial {\bar{N}}_{s}}{\partial d}]}^{2} \sum_{n = 0}^{\infty} P_{s} (n) {[\frac{\partial \log P_{s} (n)}{\partial {\bar{N}}_{s}}]}^{2}

= {[\frac{\partial {\bar{N}}_{s}}{\partial d}]}^{2} \frac{1}{{\bar{N}}_{s} ({\bar{N}}_{s} + 1)}

= \frac{ε_{ave} γ^{2} (d)}{1 + δ (d)} {1 + ε_{ave} [1 + δ (d)]}^{- 1} .

Here we have defined

γ (d) = δ^{'} (d) = - Re \int_{-} \int_{\infty}^{\infty} d x d y ψ^{*} (x, y) \frac{\partial ψ (x - d, y)}{\partial x} .

In a similar fashion, we find the Fisher information of photon counting of the inversion-antisymmetric field to be

J_{d}^{(asym - pc)} = \frac{ε_{ave} γ^{2} (d)}{1 - δ (d)} {1 + ε_{ave} [1 - δ (d)]}^{- 1} .

The CR lower bounds for these individual measurements then respectively read

E {[\hat{d} (N_{s}) - d]}^{2} \geq \frac{1}{J_{d}^{(sym - pc)}} = \frac{1 + δ (d)}{ε_{ave} γ^{2} (d)} + {[\frac{1 + δ (d)}{γ (d)}]}^{2},

E {[\hat{d} (N_{a}) - d]}^{2} \geq \frac{1}{J_{d}^{(asym - pc)}} = \frac{1 - δ (d)}{ε_{ave} γ^{2} (d)} + {[\frac{1 - δ (d)}{γ (d)}]}^{2} .

For a fixed d, both bounds contain a term that displays the familiar shot-noise scaling with respect to the combined source strength ε_ave. Of greater interest here, however, is the behavior of the bounds as the separation reduces to zero for fixed source strengths. Since γ(d) →0 and δ (d) → 1 as d → 0, we have

\lim_{d \to 0} \frac{1 + δ (d)}{γ^{2} (d)} = \infty,

so that any unbiased estimator based on bulk direct detection of the inversion-symmetric field suffers a divergent mean squared error for separations d → 0. Thus, Rayleigh’s curse that plagues spatially-resolved image-plane photon counting remains in effect here [3,4]. However, note that

\lim_{d \to 0} \frac{1 - δ (d)}{γ^{2} (d)} = - \frac{δ^{'} (0)}{2 γ (0) γ^{'} (0)} = - \frac{1}{2 γ^{'} (0)} \equiv \frac{1}{2 {(Δ k)}^{2}} .

Here we have used

γ^{'} (d) = - Re \int_{-} \int_{\infty}^{\infty} d x d y \frac{\partial ψ^{*} (x + d, y)}{\partial x} \frac{\partial ψ (x, y)}{\partial x},

γ^{'} (0) = - \int_{-} \int_{\infty}^{\infty} d x d y {| \frac{\partial ψ (x, y)}{\partial x} |}^{2}

\equiv - {(Δ k)}^{2},

where (Δk)² is the squared spectral width of the PSF. Remarkably therefore, direct detection of the inversion-antisymmetric field does not suffer Rayleigh’s curse and a limiting MSE

\frac{1}{M J_{0}^{(asym - pc)}} = \frac{1}{2 M ε_{ave} {(Δ k)}^{2}} \equiv \frac{1}{N {(Δ k)}^{2}}

is potentially achievable for small d, where M is the number of independent measurements, and

N = 2 M ε_{ave}

is the the mean number of photons emitted over the M shots. The limiting Fisher information value 2Mε_ave(Δk)² is similar to the quantum limit obtained in [7] and achieved by the SPADE and binary SPADE schemes therein. An advantage of SLIVER over those schemes is that there is no need to couple the image-plane field into a waveguide(s) with mode profile(s) tailored to the PSF.

If ε₁ = ε₂, N_s and N_a are statistically independent so that the two Fisher information terms of eqs. (26) and (28) may be added to give the total Fisher information. The reciprocal of this quantity is then the CR bound for any unbiased estimate $\hat{d} (N_{s}, N_{a})$ of d based on both photocounts. In the remainder of the paper, however, we will only consider measurements in the antisymmetric output port.

4.2. On-off detection of the inversion-antisymmetric field

Consider bucket direct detection at the antisymmetric port with an on-off detector, i.e., a detector that only distinguishes between no photon and one or more photons. The no-click probability and click probability are respectively

P_{a} (0) = \frac{1}{1 + ε_{ave} [1 - δ (d)]}

P_{a} (> 0) = \frac{ε_{ave} [1 - δ (d)]}{1 + ε_{ave} [1 - δ (d)]} .

The Fisher information for this binary measurement can be shown to be

J_{d}^{(asym - oo)} = \frac{ε_{ave} γ^{2} (d)}{1 - δ (d)} {1 + ε_{ave} [1 - δ (d)]}^{- 2} < J_{d}^{(asym - pc)},

which is expected as on-off detection is a coarser measurement than photon counting. However, it is interesting that on-off detection also evades Rayleigh’s curse for small d and the Fisher information attains the optimal value at d = 0. Note also that the Fisher information of Eq. (40) is negligibly different from that of Eq. (28) for ε_ave ≪ 1. These behaviors are expected as multi-photon arrivals in the image-plane are rare in the regime ε_ave ≪ 1 while multi-photon detections at the antisymmetric output port are rare in the ε_ave [1 − δ (d)] ≪ 1 regime, i.e., the regime of d ≈ 0. Equations (28) and (40) indicate that the two antisymmetric detection methods should exhibit superresolution for arbitrary values of source strength beyond the ε_ave ≪ 1 regime.

We illustrate the above results for the circular Gaussian PSF

ψ_{G} (ρ) = \frac{1}{{(2 π σ^{2})}^{1 / 2}} \exp (- \frac{| ρ |^{2}}{4 σ^{2}}) .

The PSF-dependent quantities appearing in the CR bound then become:-

δ_{G} (d) = \exp (- \frac{d^{2}}{8 σ^{2}}),

γ_{G} (d) = - \frac{d}{4 σ^{2}} \exp (- \frac{d^{2}}{8 σ^{2}}),

{(Δ k)}_{G}^{2} = \frac{1}{4 σ^{2}} .

The Fisher information quantities of Eqs. (26), (28), and (40) pertaining to photon counting at the symmetric and antisymmetric ports, and on-off detection at the antisymmetric port of the interferometer are plotted as a function of the separation d in ^{Figures 3–4}. The source strength is ε_ave = 10⁻³ photons and 0.5 photons in Figures 3 and 4 respectively. The plots are normalized to the maximum value of $J_{d}^{asym - pc}$ attained at d = 0. We see that the information obtainable from the antisymmetric detection methods is maximum at d = 0 and decreases thereafter. The information from detection at the symmetric port is zero for d = 0 and reaches a peak around d = 2 – 2.5σ and decreases thereafter.

Fig. 3 Fisher information obtainable using the SLIVER method: - The Fisher information for photon counting $J_{d}^{asym - pc}$ at the symmetric output port of the image inversion interferometer (blue curve), photon counting at the antisymmetric output port of the image inversion interferometer ( $J_{d}^{asym - pc}$ – dashed green curve) and for on-off detection at the antisymmetric output port ( $J_{d}^{asym - oo}$ – dotted red curve) are plotted against the source separation d. The plots are normalized with respect to the maximum value 2ε_ave (Δk)² of $J_{d}^{asym - pc}$ and $J_{d}^{asym - oo}$ , attained at d = 0. The superlocalization effect, consistent with the quantum bound of [7], is evident for d ~ 0. The circular Gaussian PSF of Eq. (41) has been assumed and the plots are independent of the value of the half-width σ. The average source strength ε_ave = 0.5 photons.

Download Full Size | PDF

Fig. 4 Fisher information obtainable using the SLIVER method: - The Fisher information for photon counting $J_{d}^{sym - pc}$ at the symmetric output port of the image inversion interferometer (blue curve), photon counting at the antisymmetric output port of the image inversion interferometer ( $J_{d}^{asym - pc}$ – dashed green curve) and for on-off detection at the antisymmetric output port ( $J_{d}^{asym - oo}$ – dotted red curve) are plotted against the source separation d. The plots are normalized with respect to the maximum value 2ε_ave (∆k)² of $J_{d}^{asym - pc}$ and $J_{d}^{asym - oo}$ , attained at d =0. The superlocalization effect, consistent with the quantum bound of [7], is evident for d ∼ 0. The circular Gaussian PSF of Eq. (41) has been assumed and the plots are independent of the value of the half-width σ. The average source strength ε_ave = 0.5 photons. so that superresolution persists outside the regime ε_ave ≪ 1 and is more marked for number-resolved measurements.

Download Full Size | PDF

5. Maximum-Likelihood (ML) estimates of d: Monte-Carlo analysis

To further support the predictions from the CR bounds, we present Monte-Carlo simulations of the MSE for on-off detection and photon counting at the antisymmetric output port of the interferometer. The circular Gaussian PSF of Eq. (41) is assumed and the simulation results are independent of the half-width σ.

5.1. On-off detection in the antisymmetric output port

Consider direct detection of the output E_a(ρ) in M detection windows using an on-off detector. The ensuing measurement record is a bit string k = (k₁,…,k_M), with k_m = 0 if the detector did not fire in the m-th detection window and k_m = 1 if it did. The maximum likelihood (ML) estimator [6] for d is then:-

{\hat{d}}_{ML} = {\begin{array}{l} 2 σ \sqrt{- 2 \ln (1 - \frac{K}{(M - K) ε_{ave}})} & if \frac{K}{M - K} < ε_{ave} . \\ 2 σ & otherwise . \end{array}

Here, $K = \sum_{m = 1}^{M} k_{m}$ is the total number of clicks observed in the detector during the measurement, which is a sufficient statistic for generating the estimate. The second case above is necessary because the equation for the maximum likelihood estimate of d has no solution if K/(M − K) ≥ ε_ave, in which case the estimate is set to an arbitrary value. As M increases, the probability of such “large deviations” goes to zero.

Figure 5 shows the simulated MSE of the ML estimate for a source strength of ε_ave = 0.2 photons. The MSE is shown scaled relative to the minimum value ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ of the CR bound for M-values of 100, 200 and 400. As the separation increases beyond the half-width, the curves for different M approach the Cramér-Rao bound for on-off detection in the antisymmetric output port. We see that the ML estimate actually beats the CR bound for small d ≲ 0.3σ. This is because the ML estimate is biased for finite M, as may be expected from the highly nonlinear nature of the estimator of Eq. (45). We recall here that bias and MSE are independent characteristics of an estimator. Many minimum MSE estimators are biased [6], and techniques for reducing MSE by introducing bias are well-known in the signal processing literature [26]. Since MSE is the relevant metric for most applications, the behavior of the MSE observed here is a welcome – though intriguing – feature. For all values of the separation, the performance – when not superior to that predicted by the CR bound – is within a factor of 2 of it. The behavior of the simulated MSE of SLIVER with respect to the CR bound is very similar to that reported for the SPADE and binary SPADE schemes in [7].

Fig. 5 Simulated mean squared errors for the maximum likelihood estimator ${\hat{d}}_{ML}$ of Eq. (45) as a function of the source separation for M = 100,200,400 measurements. Each measurement consists of a binary value indicating whether or not an on-off detector in the antisymmetric output port of the interferometer registered a click. The MSE curves are each scaled by ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ , the value of the Cramér-Rao bound for that M at d = 0. Each data point was obtained as an average of 10⁵ Monte-Carlo runs. The Cramér-Rao bound for on-off detection in the antisymmetric output port (solid blue curve) normalized to unity at its minimum value of ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ is also shown. A circular Gaussian PSF (Eq. (41)) was assumed and the source strength ε_ave = 0.2 photons. The simulated results are independent of the half-width σ.

Download Full Size | PDF

For PSFs other than the circular Gaussian, the ML estimator takes a different form than that in Eq. (45), but depends on the measurement result only through the quantity $(\frac{K}{M - K})$ .

5.2. Photon counting in the antisymmetric output port

For number-resolved photon counting in the antisymmetric output port, the measurement record p = (p₁,…,p_M) is now a vector of non-negative integers. The total number of photons counted $- P = \sum_{m = 1}^{M} p_{m}$ – is a sufficient statistic and the ML estimator for d is:-

{\hat{d}}_{ML} = {\begin{array}{l} 2 σ \sqrt{- 2 \ln (1 - \frac{P}{M ε_{ave}})} & if \frac{P}{M} < ε_{ave} . \\ 2 σ & otherwise, \end{array}

where the last value is again an arbitrary assignment needed in the event that a solution to the equation for the ML estimate for d does not exist.

Figure 6 shows the results for the MSE of the ML estimate for a source strength ε_ave = 5 photons, well outside the ε_ave ≪ 1 regime. For such a source, an on-off detector would need a very large M to collect useful statistics because of detector saturation. The simulated MSE is shown scaled relative to the minimum value ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ of the CR bound for M-values of 50, 100 and 150. The ML estimate again beats the CR bound for small d ≲ 0.1σ owing to the biasedness of the estimator of Eq. (46). However, the curves for each M closely approach the CR bound for photon counting in the antisymmetric output port for separations d ≳ 0.5σ. For all values of the separation, the observed performance – when not superior to that predicted by the CR bound – is within a factor of 2 of it.

Fig. 6 Simulated mean squared errors for the maximum likelihood estimator ${\hat{d}}_{ML}$ of Eq. (46) as a function of the source separation for M = 50,100,150 measurements. Each measurement counts the number of photons in the antisymmetric port, i.e., is a measurement of N_a of Sec. 2. The MSE curves are each scaled by ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ , the value of the Cramér-Rao bound for that M at d = 0. Each data point was obtained as an average of 10⁵ Monte-Carlo runs. The Cramér-Rao bound for photon counting in the antisymmetric output port (solid blue curve) normalized to unity at its minimum value of ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ is also shown. A circular Gaussian PSF (Eq. (41)) was assumed, the source strength ε_ave = 5 photons, and the simulated results are independent of the half-width σ.

Download Full Size | PDF

For PSFs other than the circular Gaussian, the ML estimator takes a different form than that in Eq. (46), but depends on the measurement result only through the quantity $(\frac{P}{M})$ .

5.3. Effect of centroid mismatch: On-off detection in the antisymmetric output

We have so far assumed that the centroid has been located at the optical axis, so that the vector c from the origin of the image plane to the midpoint of d is zero (Fig. 1). A detailed study of the operation and performance of SLIVER when this assumption is removed will be given elsewhere. Here, we briefly evaluate the performance of on-off detection in the antisymmetric port for nonzero c. It can be shown that the no-click and click probability are given by

\begin{array}{l} P_{a} (0) = {[1 + α ε_{1} + β ε_{2} + (α β - γ^{2}) ε_{1} ε_{2}]}^{- 1}, \\ P_{a} (> 0) = 1 - P_{a} (0), \end{array}

where

\begin{array}{l} α = \frac{1 - δ (| 2 c - d |)}{2}, \\ β = \frac{1 - δ (| 2 c + d |)}{2}, \\ γ = \frac{δ (| d |) - δ (| 2 c |)}{2} . \end{array}

Consider the case of c parallel to the separation vector d (Fig. 7). For the Gaussian PSF of Eq. (41), let |c| = ξσ, so that ξ is a normalized misalignment factor. For such a configuration, and for fixed c with ξ = 0.1, simulation results for the MSE of the estimator of Eq. (45) are shown in Fig. 8. Note that the estimator used cannot depend on the centroid position, as this is unknown in practice. We have assumed ε₁ = ε₂ = ε_ave = 0.2 photons. Comparing to Fig. 5, we see that the misalignment has little effect on the MSE for separations d ≳ 0.5σ. For smaller separations, the MSEs in Fig. 8 are greater than those in Fig. 5 and the MSE does not go to zero for zero separation. However, the MSEs still beat the CR bound for the aligned SLIVER system, suggesting they also beat the CR bound – yet to be determined – for the misaligned system.

Fig. 7 The case of centroid mismatch considered in Section 5.3: The vector OM = c from the optical axis to the centroid is nonzero but parallel to the separation vector 12 = d. The ordering of the points (1−O−M−2) shown is illustrative – the simulation results are valid even if O does not lie between 1 and 2.

Download Full Size | PDF

Fig. 8 Simulated mean squared errors for the estimator Eq. (45) for a misaligned SLIVER system using on-off detection in the antisymmetric port shown as a function of the source separation for M = 100,200,400 measurements. We have assumed that the centroid is a distance |c| ≡ξσ away from the inversion axis, that c and d are parallel, and that ξ = 0.1. The MSE curves are each scaled by ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ , the value of the Cramér-Rao bound for the ideal aligned measurement (ξ = 0) and that M at d = 0. Each data point was obtained as an average of 10⁵ Monte-Carlo runs. The Cramér-Rao bound for on-off detection in the antisymmetric output port of an aligned SLIVER system (solid blue curve) normalized to unity at its minimum value of ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ is also shown. A circular Gaussian PSF (Eq. (41)) was assumed and the source strengths ε₁ = ε₂ = ε_ave = 0.2 photons. The simulated results are independent of the half-width σ.

Download Full Size | PDF

Figure 9 displays the simulated MSEs for the same source and system parameters, with M = 100 and for the misalignment factor values ξ = 0.1,0.2,0.3 and 0.4. We see that the performance degrades gracefully with increasing misalignment, remains within an order of magnitude of the CR bound of the aligned system, and does not display the divergent behavior of image-plane photon counting [3, 4, 7].

Fig. 9 Simulated mean squared errors for the estimator Eq. (45) for misaligned SLIVER as a function of the source separation for M = 100 measurements and for misalignment factors ξ = 0.1,0.2,0.3 and 0.4. The MSE curves are each scaled by ${[2 M ε_{ave} {(Δ k)}^{2}]}^{- 1}$ , the value of the Cramér-Rao bound for the ideal aligned measurement (ξ = 0) at d = 0. Each data point was obtained as an average of 10⁵ Monte-Carlo runs. The Cramér-Rao bound for on-off detection in the antisymmetric output port of an aligned SLIVER system normalized to unity at its minimum value is also shown. A circular Gaussian PSF (Eq. (41)) was assumed and ε₁ = ε₂ = ε_ave = 0.2 photons. The simulated results are independent of the half-width σ.

Download Full Size | PDF

A detailed study of SLIVER under misalignment – including derivations of the formulas in this section, Cramér-Rao bounds, and performance analyses – will be given elsewhere. Here, we simply note the important point that the MSE of a misaligned SLIVER system is still finite for arbitrarily small separations, and Rayleigh’s curse has not returned.

6. Discussion and conclusions

Motivated by the greater familiarity and facility of analysis afforded by the semiclassical theory of photodetection over the full quantum theory, we have set up the problem of resolving two incoherently radiating point sources in this framework. This theory is quantitatively exact for such sources and for measurements involving only linear-optics processing prior to photodetection. The sources are assumed to be of arbitrary strength and circularly-symmetric PSF’s directly related to real imaging systems are considered. Since the analysis in [7] was largely confined to a one-dimensional model of the point-spread function with a weak-source assumption, the model used here is of interest both theoretically and in practice (See also the remarks in Sec. V therein). We have not made full semiclassical analyses of the binary SPADE and SPADE protocols of [7] here, as these are rather involved and will be given elsewhere. Preliminary calculations indicate that the results from [7] for ε_ave ≪ 1 agree with those from the semiclassical analysis to $O (ε_{ave}^{2})$ .

We have also proposed and analyzed a novel interferometric scheme – the SLIVER method of Sec. 3 – for estimating the separation between the sources once their centroid has been located. The scheme alleviates Rayleigh’s curse and provides superresolution at all values of source strengths. The essential ingredient of the scheme is an image inversion interferometer (III) – a device that separates the symmetric and antisymmetric parts of the image-plane field relative to inversion about the centroid of the sources. The SLIVER method, unlike the binary SPADE and SPADE schemes, does not require coupling the image-plane field into one or more waveguides with mode profiles tailored to the PSF of the imaging system. As such, it appears to be more readily implementable in practice.

An explanation of the superlocalization effect of SLIVER for estimating d ≈ 0 may be given in terms of photon statistics as follows. If both sources were superimposed at the centroid (so that d = 0), the field E_a(ρ) in the antisymmetric output port of the III is identically zero despite their incoherence. If d is small but nonzero, the mean photon number ${\bar{N}}_{a}$ is finite but still small as per Eq. (20). The Fisher information on the mean of a Bose-Einstein distribution equals ${[{\bar{N}}_{a} ({\bar{N}}_{a} + 1)]}^{- 1}$ and hence is very large in this region of d. This extreme sensitivity to changes in d is partially offset by the weak sensitivity of ${\bar{N}}_{a}$ on d via γ(d) which is close to zero for d ≈ 0 (cf. the term in square brackets in Eq. (25)). These two effects compensate each other to give a finite value of the Fisher information for separations close to zero. In the symmetric port, on the other hand, the mean output is large at d ≈ 0 and the Fisher information is correspondingly small, so that the second effect dominates and superlocalization is not obtained. Similar remarks apply to on-off detection in the two output ports. It is rather remarkable that appropriate linear-optics processing followed by spatially-unresolved on-off detection can alleviate the divergent MSE behavior plaguing ideal spatially-resolved image-plane photon counting for closely separated sources. Somewhat intriguingly, receivers using the same toolkit of operations (linear-optics processing with bulk on-off detection) also approach the quantum limits of binary [27, 28] and M-ary [29] communication with laser (coherent-state) light.

As mentioned in Sec. 3, various imaging modalities employing image inversion have been proposed and demonstrated previously in the microscopy literature [17–24]. The importance of the antisymmetric output of the III for improving lateral resolution has also been noticed, e.g., in [18]. However, the claimed resolution enhancements in these works are in terms of the reduction of effective widths of the PSF by fixed factors, which in turn are related to the intensity patterns at the two outputs of the III (or their difference). Thus, from a statistical viewpoint, it appears that these works focus on the spatial structure of the first moment of the photocount random process (i.e, the intensity) at the two output channels of the interferometer and on the information about the separation d extractable from it. In contrast, we have focused on spatially unresolved detection at the antisymmetric output port of the III and on the Cramér-Rao bound as the fundamental limit to the mean squared error of any (unbiased) estimate of the separation. By considering statistically optimal estimators, we can exploit the information available in the full probability distribution of the measurement rather than just its first moment. Such subtler estimators are necessitated by the limited number of photons that can be collected in most imaging applications. The maximum likelihood estimator, which approaches the Cramér-Rao bound in the limit of large M [6], was shown to yield performance in close agreement to the bound and sometimes exceeding it.

In order to focus on the spatial aspects of the resolution problem, we have, as in [7], suppressed the explicit temporal dependence of the field in our analysis. In effect, we have assumed that the light emitted by the sources is in a single quasimonochromatic temporal mode in each detection window. Such a situation can be easily realized using a pseudothermal light source [16], but is in general unrealistic for many astronomical and biological imaging situations in which the coherence time of the emitted radiation is typically much shorter than the duration of the detection window. An approximate treatment of such cases can be given using the concept of the number of coherence cells [16]. This approach effectively redefines the number of shots M according to the coherence properties of the light, with exact results obtainable in certain cases. These issues will be explored in detail elsewhere.

The effects of loss and non-unity quantum efficiency of the detectors can be incorporated in our analysis by appropriate scaling of the source strengths. Since spatially homogeneous losses simply scale the semiclassical field amplitudes, thermal states remain thermal states under such attenuation. Thus, a system resolving sources of strength ε_ave using a detector with quantum efficiency η can be modeled as a system with source strength ηε_ave and an ideal detector.

We have assumed through most of this paper that the centroid of the sources has been located at the optical axis. In optical astronomy applications, the centroid position may be available from previous observational records, or else a long observation time should be available to determine the centroid by, say, image-plane photon counting. In microscopy applications, the time available for measurements may be more limited, e.g., because of photobleaching in fluorescent samples [30], and lack of knowledge of the centroid position can be a significant impediment to implementing SLIVER. To address this problem, one can first perform image-plane photon counting using a portion of the available light to determine the centroid position to within a small fraction of the PSF width before implementing SLIVER. The simulations in Sec. 5.3 indicate that the MSE of separation estimation should still be well below that achievable using image-plane photon counting. A more general study of the effects of centroid mismatch is best carried out in the framework of multi-parameter estimation theory and will be given elsewhere.

Apart from exploring these issues, various extensions and refinements of the schemes of this paper suggest themselves. For example, spatially-resolved detection in one or both output ports of the interferometer can only improve the performance of SLIVER. In any case, the results here should hopefully spur investigations into these issues in order to reap the superlocalization gain promised by the SLIVER method.

Acknowledgments

This work is supported by the Singapore National Research Foundation under NRF Grant No. NRF-NRFF2011-07.

References and links

1. Lord Rayleigh, “XXXI. Investigations in optics, with special reference to the spectroscope,” The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science 8, 261–274 (1879). [CrossRef]

2. M. Born and E. Wolf, Principles of Optics: Electromagnetic Theory of Propagation, Interference and Diffraction of Light (Cambridge University, 1999). [CrossRef]

3. E. Bettens, D. Van Dyck, A. DenDekker, J. Sijbers, and A. Van den Bos, “Model-based two-object resolution from observations having counting statistics,” Ultramicroscopy 77, 37–48 (1999). [CrossRef]

4. S. Ram, E. S. Ward, and R. J. Ober, “Beyond Rayleigh’s criterion: a resolution measure with application to single-molecule microscopy,” Proceedings of the National Academy of Sciences of the United States of America 103, 4457–4462 (2006). [CrossRef]

5. S. V. Aert, A. den Dekker, D. V. Dyck, and A. van den Bos, “High-resolution electron microscopy and electron tomography: resolution versus precision,” Journal of Structural Biology 138, 21–33 (2002). [CrossRef] [PubMed]

6. H. L. Van Trees, Detection, Estimation, and Modulation Theory : Part I, 1st ed. (Wiley-Interscience1st Ed, 2001).

7. M. Tsang, R. Nair, and X.-M. Lu, “Quantum theory of superresolution for two incoherent optical point sources,” (2015). http://arxiv.org/abs/1511.00552.

8. S. Weisenburger and V. Sandoghdar, “Light microscopy: an ongoing contemporary revolution,” Contemporary Physics 56, 123–143 (2015). [CrossRef]

9. J. Chao, S. Ram, E. S. Ward, and R. J. Ober, “Ultrahigh accuracy imaging modality for super-localization microscopy,” Nat. Meth. 10, 335–338 (2013). [CrossRef]

10. M. Tsang, “Quantum limits to optical point-source localization,” Optica 2, 646–653 (2015). [CrossRef]

11. C. W. Helstrom, Quantum Detection and Estimation Theory (Academic Press, 1976).

12. A. S. Holevo, Probabilistic and Statistical Aspects of Quantum Theory (Edizioni Della Normale, 2011). [CrossRef]

13. M. G. Paris, “Quantum estimation for quantum technology,” International Journal of Quantum Information 7, 125–137 (2009). [CrossRef]

14. L. Mandel and E. Wolf, Optical Coherence and Quantum Optics (Cambridge University, 1995). [CrossRef]

15. J. H. Shapiro, “The quantum theory of optical communications,” IEEE Journal of Selected Topics in Quantum Electronics 15, 1547–1569 (2009). [CrossRef]

16. J. W. Goodman, Statistical Optics (John Wiley & Sons, 1985).

17. N. Sandeau and H. Giovannini, “Increasing the lateral resolution of 4pi fluorescence microscopes,” J. Opt. Soc. Am. A 23, 1089–1095 (2006). [CrossRef]

18. K. Wicker and R. Heintzmann, “Interferometric resolution improvement for confocal microscopes,” Opt. Express 15, 12206–12216 (2007). [CrossRef] [PubMed]

19. K. Wicker, S. Sindbert, and R. Heintzmann, “Characterisation of a resolution enhancing image inversion interferometer,” Opt. Express 17, 15491–15501 (2009). [CrossRef] [PubMed]

20. D. Weigel, R. Foerster, H. Babovsky, A. Kiessling, and R. Kowarschik, “Enhanced resolution of microscopic objects by image inversion interferometry,” Opt. Express 19, 26451–26462 (2011). [CrossRef]

21. D. Weigel, H. Babovsky, A. Kiessling, and R. Kowarschik, “Investigation of the resolution ability of an image inversion interferometer,” Opt. Commun. 284, 2273–2277 (2011). [CrossRef]

22. D. Weigel, T. Elsmann, H. Babovsky, A. Kiessling, and R. Kowarschik, “Combination of the resolution enhancing image inversion microscopy with digital holography,” Opt. Commun. 291, 110–115 (2013). [CrossRef]

23. D. Weigel, H. Babovsky, A. Kiessling, and R. Kowarschik, “Imaging properties of different types of microscopes in combination with an image inversion interferometer,” Opt. Commun. 332, 301–310 (2014). [CrossRef]

24. D. Weigel, H. Babovsky, A. Kiessling, and R. Kowarschik, “Widefield microscopy with infinite depth of field and enhanced lateral resolution based on an image inverting interferometer,” Opt. Commun. 342, 102–108 (2015). [CrossRef]

25. J. W. Goodman, Introduction to Fourier Optics (Roberts and Company Publishers, 2005).

26. S. Kay and Y. C. Eldar, “Rethinking biased estimation,” IEEE Signal Processing Magazine 25, 133–136 (2008). [CrossRef]

27. R. S. Kennedy, “A near-optimum receiver for the binary coherent state channel,” Tech. Rep. 108, MIT RLE Quarterly Progress Report (1973).

28. S. J. Dolinar, “An optimum receiver for the binary coherent state channel,” Tech. Rep. 111, MIT RLE Quarterly Progress Report (1973).

29. R. Nair, S. Guha, and S.-H. Tan, “Realizable receivers for discriminating coherent and multicopy quantum states near the quantum limit,” Physical Review A 89, 032318 (2014). [CrossRef]

30. U. Kubitscheck, (Ed.), Fluorescence Microscopy: from Principles to Biological Applications (John Wiley & Sons, 2013). [CrossRef]

Interferometric superlocalization of two incoherent optical point sources

Abstract

1. Introduction

2. Source and system model

3. The SLIVER measurement and its statistics

4. Cramér-Rao bounds for the separation d

4.1. Number-resolved photon counting

4.2. On-off detection of the inversion-antisymmetric field

5. Maximum-Likelihood (ML) estimates of d: Monte-Carlo analysis

5.1. On-off detection in the antisymmetric output port

5.2. Photon counting in the antisymmetric output port

5.3. Effect of centroid mismatch: On-off detection in the antisymmetric output

6. Discussion and conclusions

Acknowledgments

References and links

Cited By

Figures (9)

Equations (48)

Optics Express