## Abstract

Optical diffraction tomography (ODT) using Born or Rytov approximation suffers from severe distortions in reconstructed refractive index (RI) tomograms when multiple scattering occurs or the scattering signals are strong. These effects are usually seen as a significant impediment to the application of ODT because multiple scattering is directly linked to an unknown object itself rather than a surrounding medium, and a strong scatter invalidates the underlying assumptions of the Born and Rytov approximations. The focus of this article is to demonstrate for the first time that multiple scattering and high material contrast, if handled aptly, can significantly improve the image quality of the ODT thanks to multiple scattering inside a sample. Experimental verification using various phantom and biological cells substantiates that we not only revealed the structures that were not observable using the conventional approaches but also resolved the long-standing problem of *missing cones* in the ODT.

© 2017 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## 1. Introduction

Refractive index(RI), *n*, is the optical property that relates to the electromagnetic susceptibility (*χ _{e}* and

*χ*) representing whole light-matter interaction,

_{m}*D*) RI distribution of micro samples, especially biological cells, without any preprocessing such as staining, labeling, or protein tagging [1–25]. Recently, real-time observation [26] and compact optical setup [27] have also been achieved for routine monitoring in biological laboratories.

In optical diffraction tomography (ODT), there exist certain spatial frequency components that cannot be measured due to the limited projection angles imposed by objective lenses. This limitation, often coined as the *missing cone problem*, causes the under-estimation of RI values in tomograms and results in severe elongations of RI distributions along the optical axis. Towards this end, several iterative reconstruction algorithms have been introduced by exploiting *à priori* information such as positivity in RI differences or smoothness of samples. The interested readers are referred to [28] for a systematic comparison of various iterative reconstruction algorithms in ODT.

Aside from the missing cone problem, there exists a more fundamental, but relatively less studied, problem in ODT – *the multiple scattering problem*. Specifically, in ODT, the forward diffraction is usually modeled by the scalar Helmholtz equation,

*κ*= 2

*π/λ̄*is the wavenumber with

*λ̄*being the incident light wavelength in the free space,

*U*(

_{s}**r**) :=

*U*(

**r**) −

*U*

_{0}(

**r**),

**r**∈ ℝ

^{3}, is the scattered field from the scatterer with illumination dependent total electric field

*U*and the background electric field

*U*

_{0}in the surrounding medium without any scattering object. The density

*f*denotes the normalized scattering potential within a domain

*D*⊂ ℝ

^{3}representing the volume of the scattering object. It is defined by

*n*

_{0}(constant) is the RI of the surrounding homogeneous medium,

*χ*denotes the characteristic function of

_{D}*D*, and It is well-known that the solution of Eq. (2) can be represented by the nonlinear

*Lippmann-Schwinger*integral equation (see, e.g., [29, §13.1.1]),

^{3}is the entire imaging domain with smooth boundary

*∂*Ω,

*G*(

**r**,

**r′**) =

*e*

^{iκ|r−r′|}/4

*π*|

**r**−

**r′**| is the outgoing free-space Green’s function of 3

*D*Helmholtz equation and the operator

*K*:

_{D}*L*

^{2}(

*D*) →

*L*

^{∞}(ℝ

^{3}) is defined by

*f*(

**r**), for

**r**∈

*D*, from measured scattered field

*U*(

_{s}**r**) at

**r**∈ Γ ⊂

*∂*Ω. However, since Eq. (5) is nonlinear due to the coupling between the unknown potential

*f*and the unknown field

*U*(which itself depends on

*f*), the inverse problem becomes complicated. Accordingly, a common practice is to linearize the problem with respect to the size of the domain

*D*or the strength of its scatter. For example, the first-order Born approximation assumes that the total electric field in a sample is the same as the incident electric field inside

*D*, i.e.

*U*(

**r**) ≃

*U*

_{0}(

**r**), for

**r**∈

*D*(and therefore independent of

*f*thereby making Eq. (5) linear in

*f*) [29, §13.1.2]. The validity of these approximations is well studied in [4] under weak scatter assumption, i.e. when

*f*varies only slightly across the imaging domain. However, for high contrast scattering potential, or when multiple scattering occurs among individual weak scatterers, the Born approximation often causes severe distortions in reconstructed RI tomograms [30].

On the other hand, the resolution enhancement in highly scattering resonant media has been recently demonstrated in various far-field experiments [31–35]. The basic idea is that if the medium around the sources is engineered so that the point spread function (PSF) displays a much sharper peak than the homogeneous one, then one can focus on or resolve sub-wavelength details [31]. Therefore, the key research issue has been of designing the surrounding medium using various ingenious structures [31–37]. Moreover, it was mathematically shown that the associated resolution enhancement is due to the sub-wavelength resonance modes excited in the high contrast medium which can propagate into the far-field [38,39].

One of the key observations made in this paper is that these seemingly different topics are indeed closely related to each other. In fact, the *K _{D}* operator (cf. Eq. (6)) has been linked to the super-resolution effects in a resonant medium [38]. Moreover, while the

*K*operator in Ref. [38] is an engineered one by augmenting additional inhomogeneous layer with a known scattering potential, the similar operator in ODT problem is originated from the specimen itself (cf. Fig 1). The catch, though, is that the kernel of the operator (i.e.,

_{D}*f*) is unknown and needs to be estimated. Therefore, one of the main contributions of this paper is to provide a rigorous way to estimate the unknown

*K*operator in ODT and demonstrate for the first time that the unknown medium itself works as a self-resonator so that super-resolution imaging can be achieved from medium induced sub-wavelength resonance modes without any additional medium and/or PSF engineering. Specifically, in order to estimate

_{D}*K*operator, we show that Eq. (5) should be used without any Born or Rytov-type linearization, which is indeed possible thanks to the recent advances in joint sparse recovery in compressed sensing and its application to inverse scattering problems [40–49].

_{D}Based on extensive real experiments using various samples such as spherical beads, human red blood cells, hepatocyte cells, and multiple beads, we demonstrate clear resolution enhancement by the proposed method. Moreover, our approach is shown to simultaneously address the long-standing problem of missing cones in ODT.

## 2. Theory

#### 2.1. Sub-wavelength resonance modes

In order to facilitate ensuing discussion, let us briefly review the mathematical theory behind the resolution enhancement from the sub-wavelength resonance modes. We refer the interested readers to [38] for detailed arguments.

Consider an inverse source problem in a homogeneous medium

*s*∈

*L*

^{2}(Ω′) is the unknown density of the sources to be reconstructed which is compactly supported in a bounded smooth imaging domain Ω′ ⊂ ℝ

^{3}. Since the medium is homogeneous, the achievable diffraction limited resolution is basically

*d*≃

*λ̄*/2 =

*π/κ*, which is inversely proportional to the wavenumber. Thus, the key idea of resolution enhancement using resonant media is to insert tailored inhomogeneous material with known properties as surrounding medium to guide the wave to travel with respect to the altered wavenumber (see Fig. 1(a)) [31–35]. Precisely, the new wave propagates according to the modified Helmholtz equation,

*χ*denotes the characteristic function of a position

_{D′}*D′*⊂ Ω′ of the inserted material,

*f′*is a known positive function compactly supported in

*D′*, and

*τ′*≫ 1 is the contrast parameter. Fairly easy manipulations show that the perturbed field

*u*(

_{s}**r**) :=

*u*(

**r**) −

*u*

_{0}(

**r**), for

**r**∈ ℝ

^{3}, is given by

*K*is defined in (6) with

_{D′}*D*replaced by

*D′*and the scattering potential

*f*replaced by

*f′*.

The mathematical properties of the operator *K _{D′}* have been extensively studied in [38, Lemmas 3.1–3.5]. In Theorem 2.1, we summarize its key properties and refer the interested readers to [38] for further details. In the sequel,

*H*

^{2}(

*D′*) denote the space of functions

*φ*∈

*L*

^{2}(

*D′*) such that

*∂*,

_{i}φ*∂*∈

_{ij}φ*L*

^{2}(

*D′*) for all

*i*,

*j*= 1, 2, 3.

**Theorem 2.1**

*The operator K*_{D′}is compact from L^{2}(*D′*)*to L*^{2}(*D′*)*, bounded from L*^{2}(*D′*)*to H*^{2}(*D′*)*, and is Hilbert-Schmidt.**Let σ*(*K*)_{D′}*and σ*(_{p}*K*)_{D′}*denote the spectrum and the point spectrum of K*{_{D′}, respectively. The operator K_{D′}has a countable number of distinct eigenvalues*λ*}_{i}*which accumulate to zero, i.e.,*$$\sigma ({K}_{{D}^{\prime}})=\left\{0,{\lambda}_{1},{\lambda}_{2},\cdots ,{\lambda}_{n},\cdots \right\},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathit{where}\hspace{0.17em}\hspace{0.17em}\left|{\lambda}_{1}\right|\ge \left|{\lambda}_{2}\right|\ge \cdots ,\hspace{0.17em}\hspace{0.17em}\mathit{and}\hspace{0.17em}\hspace{0.17em}{\lambda}_{n}\to 0.$$*Moreover, σ*(*K*) \_{D′}*σ*(_{p}*K*) = {0}._{D′}*λ*∈*σ*(*K*)_{D′}*if and only if there is a non-trivial solution to*$$\{\begin{array}{ll}(\mathrm{\Delta}+{\kappa}^{2})v(\mathbf{r})=-\frac{{\kappa}^{2}}{\lambda}{f}^{\prime}(\mathbf{r})v(\mathbf{r}),\hfill & \mathbf{r}\in {D}^{\prime},\hfill \\ (\mathrm{\Delta}+{\kappa}^{2})v(\mathbf{r})=0,\hfill & {\mathbb{R}}^{3}\backslash {D}^{\prime},\hfill \\ v(\mathbf{r})\hspace{0.17em}\mathit{satisfies}\hspace{0.17em}\mathit{the}\hspace{0.17em}\mathit{Sommerfeld}\hspace{0.17em}\mathit{radiation}\hspace{0.17em}\mathit{condition}\hspace{0.17em}\mathit{as}\hfill & \left|\mathbf{r}\right|\to +\infty .\hfill \end{array}$$

The functions *v* satisfying (11) are termed as the *resonance modes*. They have sub-wavelength structure in *D′* for |*λ*| < 1 and can propagate into the far field. Furthermore, since the ensemble of generalized eigenspaces of the operator *K _{D′}* is dense in

*L*

^{2}(

*D′*),

*u*can be represented as

_{s}*a*denotes the expansion coefficient of

_{i}*u*

_{0}in {

*φ*} with

_{i}*φ*being the eigenvector associated with

_{i}*λ*. Therefore, if 1/

_{i}*τ′*is close to a specific

*λ*, the associated resonance mode

_{i}*v*:=

_{i}*K*[

_{D′}*φ*] =

_{i}*λφ*can be significantly amplified so that it can be measured in the far field, achieving sub-wavelength resolution.

_{i}#### 2.2. Super-resolution from intrinsically induced sub-wavelength resonance modes

Inspired by the theoretical finding of the sub-wavelength resonance modes, here we are interested in the understanding of the super-resolution in ODT. Note that the forward ODT scattering operator in Eq. (5) has exactly the same form as in Eq. (9). However, there exist fundamental differences between the two. As shown in Fig. 1(a), the *K _{D′}* operator in Eq. (9) is an engineered operator by adding inhomogeneous layer in surrounding medium at

*D′*whose optical property

*f′*is known. In contrast, the

*K*operator in ODT is generated from the specimen itself that is embedded in the

_{D}*homogeneous medium*at the location

*D*and with the scattering potential

*f*. Nevertheless, the general theoretical analysis summarized in Theorem 2.1 is still valid and the eigenmodes of the operator in ODT also satisfy Eq. (11), i.e., they are propagating sub-wavelength modes. As shown in Fig. 1(b), these sub-wavelength modes are intrinsically generated within the unknown medium from self-resonance phenomenon. Therefore, any additional engineering of medium or PSF is not required to achieve resolution enhancement. In addition, Eq. (12) reveals that more resolution enhancement can be achieved when the unknown object has a high contrast relative to the homogeneous background.

In spite of the aforementioned significant potential for resolution enhancement using self-resonance modes, an important technical issue is that neither the kernel of the *K _{D}* operator nor the domain

*D*is known

*à priori*in ODT. Therefore, different approximations have been sought. However, these conventional linearization approaches hitherto fail to exploit the self-resonance modes. For example, under the first order Born approximation, Eq. (5) is approximated by

*K*and the contribution of the sub-wavelength resonance modes is lost in the approximation process. Consequently, one cannot expect any resolution enhancement.

_{D}Over the last few decades, the problem of accurate estimation of the *K _{D}* operator has been regarded intractable due to the coupling of the two unknowns in Eq. (5), i.e., the scattering potential

*τκ*

^{2}

*f*(

**r**) and the support

*D*. Our group has demonstrated in [45–49] that this type of problems can be decoupled into a series of easier problems as follows.

- (P1) Estimate the induced current density
*I*(**r**) :=*f*(**r**)*U*(**r**), for**r**∈*D*, and its spatial support*D*from the measured scattered field*U*(_{s}**r**), for**r**∈ Γ, based on the forward relationship, - (P2) Estimate the total field
*U*(**r**), for**r**∈*D̂*, using the relationship,$$\widehat{U}(\mathbf{r})={U}_{0}(\mathbf{r})+\tau {\kappa}^{2}{\int}_{\widehat{D}}G(\mathbf{r},{\mathbf{r}}^{\prime})\widehat{I}({\mathbf{r}}^{\prime})d{\mathbf{r}}^{\prime},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathbf{r}\in D,$$where*Û*,*Î*, and*D̂*, respectively, denote the estimated total field, induced current, and the spatial support obtained by solving (P1) and (P2). - (P3) Estimate the scattering potential based on either$$\widehat{I}(\mathbf{r})=\widehat{U}(\mathbf{r})f(\mathbf{r}),\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathbf{r}\in \widehat{D},$$or$${U}_{s}(\mathbf{r})=\tau {\kappa}^{2}{\int}_{\widehat{D}}G(\mathbf{r},{\mathbf{r}}^{\prime})\widehat{U}({\mathbf{r}}^{\prime})f({\mathbf{r}}^{\prime})d{\mathbf{r}}^{\prime},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathbf{r}\in \mathrm{\Gamma}.$$

It is worthwhile remarking that the problems (P2) and (P3) are linear in their unknowns if (P1) is resolved which is nonlinear. Accordingly, an accurate solution of Eq. (14) is very crucial. Therefore, ODT measurements for all illumination angles should be used for the estimation of *D* and *I*(**r**) in Eq. (14). Towards this end, Eq. (5) for multiple illumination angles can be represented by

*M*≥ 1 is the number of illumination angles used. Due to the rapid scanning using galvanometer, the specimen remains still during the multiple illuminations; so

*f*(

**r**) is set identical in Eq. (18).

By concatenating the measurement vector ${U}_{s}^{(m)}(\mathbf{r})$, **r** ∈ Γ, for *m* = 1, · · ·, *M*, side-by-side into a matrix *𝒴*, we have

*ℐ*denotes the matrix composed of stacking current ${\left\{{\tilde{I}}^{(m)}(\mathbf{r})\right\}}_{m=1}^{M}$,

**r**∈ Ω with

*ℐ*has non-zero rows only corresponding to

**r**∈

*D*thanks to its construction. Since the support

*D*usually assumes a very small part of the entire volume,

*ℐ*has sparse non-zero rows. The estimation problem of

*ℐ*satisfying Eq. (19) under these conditions is often called the multiple measurement vector (MMV) problem in compressed sensing [50,51], which can be reformulated as the optimization problem,

*denotes the Frobenius matrix norm and ‖ · ‖*

_{F}_{0}denotes the number of non-zero rows. Thanks to the common sparsity pattern in

*ℐ*, it enables a further reduction in the number of required measurements [43,51,52], as the number of measurements required per sensor must account for the minimal features unique to that sensor. With the advances of high performance joint sparse recovery algorithms with the performance guarantee [40–44], the proposed

*K*operator estimation can be readily done in ODT by estimating the correct support

_{D}*D*: accordingly, we can realize super-resolution imaging from self-resonance sub-wavelength modes.

There exists another simplification. In ODT, the measurement plane, Γ, is in the far-field region. Therefore, one can consider the far-field expression of the Green’s function,

**r**/|

**r**| :=

**r̂**∈ 𝕊

^{2}for all

**r**∈ ℝ

*\ {*

^{d}**0**} with 𝕊

^{2}:= {

**r**∈ ℝ

^{3}: |

**r**| = 1}, and

*G*

_{∞}(

**r̂**,

**r′**) = exp(−

*iκ*

**r̂**·

**r′**)/4

*π*. Moreover, the incident field

*U*

_{0}(

**r**) is a plane wave, ${U}_{0}^{(m)}(\mathbf{r})={e}^{i{\mathbf{k}}_{m}\cdot \mathbf{r}}$, where

**k**

*=*

_{m}*κ*

**ŝ**

*is a real-valued spatial wavenumber determined by the directional unit vector*

_{m}**ŝ**

*∈ 𝕊*

_{m}^{2}for the

*m*-th illumination. Accordingly, Eq. (14) can be converted to a 3

*D*Fourier transform formula,

## 3. Methods

#### 3.1. Experimental setup

Fused silica microsphere (44054, Sigma-Aldrich Inc., USA), multiple fused silica microspheres (95581, Sigma-Aldrich Inc., USA), red blood cells (RBCs), and hepatocyte cells (Huh-7 cell line, Apath, Brooklyn, NY, USA) were used as experimental samples. All samples were prepared by following the standard protocols [53]. Fused silica microsphere was diluted in oil (*n* = 1.559 at 532 nm). RBCs, hepatocyte cells, and multiple fused silica microspheres were diluted in Dulbecco’s buffered saline (DPBS, *n* = 1.337 at 532 nm [54]) and sandwiched between two coverslips before loaded on a microscope (IX-73, Olympus Inc., Japan) which was modified for ODT.

For measuring the 3*D* RI tomograms of samples, we implemented an ODT system based on a Mach-Zehnder interferometer equipped with a two-axis galvanometric mirror using a laser source (532 nm, 50 mW, Cobolt, Solna, Sweden). We used a high-NA objective lens for the illumination (UPLFLN 60×, NA = 0.9 Olympus Inc., Japan) and for the detection (UPLSAPO 60×, NA = 1.42, oil immersion). The beam diffracted from a sample and a reference beam forms a spatially modulated hologram, which was recorded using a CMOS camera (Neo sCMOS, ANDOR Inc., Northern Ireland, UK). The illumination angle impinging onto a sample was systematically controlled using a two-axis galvanometric mirror (GVS012/M, Thorlabs, USA). Detalied information about the setup can be found elsewhere [55].

#### 3.2. Data preprocessing

In ODT, the specimens are not absorptive, so the main changes of the total field usually comes from phase variations. Therefore, when the optical path length is larger than the wavelength, the measurement phases are usually wrapped around, so we applied phase unwrapping algorithms. Moreover, the following preprocessing step was applied before the ODT reconstruction. Because the total field can be represented by

where*ϕ*(

_{s}**r**) denotes the phase perturbation, the scattered field becomes

*U*(

_{s}**r**) in Eqs. (14) and (17) after

*ϕ*(

_{s}**r**) is unwrapped.

#### 3.3. Joint support recovery

The first step of the reconstruction algorithm is to solve the MMV problem (21). Since the support of the scattering potential is invariant during various illuminations, the matrix *ℐ* commonly shares its support. Using this prior knowledge as well as far-field measurement model (23), we relax the non-convex optimization problem (21) to a more manageable convex problem using a mixed norm. More specifically, define the followings:

*w*≥ 0 denotes a weighting factor and the subscript

_{i}*g*denotes the index set that corresponds to the

_{i}*i*-th row of each vector

*I*

^{(m)},

*m*= 1, · · ·,

*M*. Now, using the primal and dual based alternating direction method (ADM) for group sparsity [56], we can solve the joint sparsity minimization problem.

#### 3.4. Scattering potential recovery

The scattering potential estimation was done using Eq. (16) for cases without phase wrapping, and using Eq. (17) for the ones with phase wrapping, respectively. Specifically, when the measurements do not require phase unwrapping, the following penalized regression was sufficiently accurate to obtain the scattering potential:

*𝒟*denotes the finite difference operator, and

*α*is the regularization parameter. We solve minimization problem (28) using the conjugate gradient (CG) method.

When Eq. (25) was used, the original forward model (17) was considered to account for the possible errors from phase unwrapping step. Precisely, we define

*D̂*denote the estimated total flux and the domain, respectively. Then, the inverse problem

*α*is the regularization parameter,

*𝒞*is the non-negativity real convex set, and

*ι*(

_{𝒞}*x*) is the indicator function that has value 0 when

*x*∈

*𝒞*, or ∞ otherwise. The minimization problem (30) was solved using alternating direction method of multiplier (ADMM) [57,58] after splitting

*f*in order to account for the positivity term.

#### 3.5. Computational resources and source codes

We performed the reconstruction using custom-made scripts in MatLab R2015a (MathWorks Inc., Natick, MA, USA) on a desktop computer (Intel Core i7-4770K CPU, 3.5 GHz, 32 GB RAM). To accelerate the ODT reconstruction speed, we utilized a graphics processing unit (GPU, GeForce GTX TITAN, nVidia Corp., Santa Clara, CA, USA); custom-made functions based on the Compute Unified Device Architecture (CUDA) were used. The source code scripts of the reconstruction algorithm are provided as Code 1 in supplementary material and are archived in [59].

## 4. Results

#### 4.1. Fused silica microsphere

Our reconstruction results of a 5-*μm*-diameter fused silica microsphere in oil were compared with those of the first order Rytov approximation. Figure 2 shows different slice images of RI tomograms of the sample. We observed that the proposed method successfully overcame the artifacts caused by the missing cone problems, such as the underestimation of RI values and the elongation of tomograms along the optical axis. Compared to the results using the Rytov approximation, the proposed method showed higher RI values which are closer to the manufacturer’s information. More importantly, it reduced the elongation in the measured RI tomograms along the *z*-axis, providing more reasonable shape that matches the ground-truth.

#### 4.2. Multiple microspheres

The proposed method is also applied for the case when multiple scattering occurs and the conventional Rytov approximation fails. The sample consists of multiple 0.2-*μm*-diameter fused silica microspheres which are aggregated. This is the case when the scatterers have weak scatters but are closely located to one another thereby causing multiple scattering. In this case, phase wrapping occurred due to the long optical path length from multiple scattering between 0.2-*μm*-diameter fused silica microspheres. Therefore, phase unwrapping was performed as a preprocessing. As expected, due to the multiple scattering, under the first order Rytov approximation, it was hard to see the structure of multiple beads (see Fig. 3). However, the proposed method revealed the structures of multiple microspheres improving the resolution of RI tomograms thanks to the resolution enhancement. Moreover, the proposed method shows high RI values and the elongation of RI tomograms along optical axis is reduced by mitigating the missing cone problem.

#### 4.3. Biological samples:RBCs and hepatocytes

We also applied the proposed method for the recovery of biological samples such as RBCs and hepatocytes. As shown in Fig. 4, the missing cone artifact was successfully reduced by the proposed method. Compared to the Rytov approximation, the proposed method results in higher RI values and reasonable shape, clearly revealing that the RBC have the narrow disk shaped morphology along the vertical axis.

In case of a hepatocyte, the proposed method not only mitigated the missing cone problem in terms of RI values and the shape of RI tomograms along the optical axis as in Fig. 5, but also achieved increased resolution. The 3*D* volume rendering images which correspond to two squares (one is dotted line and the other one is solid line) in Fig. 5 are shown in Fig. 6 (solid line). In case of 3*D* images, the volumes containing 3 slices in *z*-axis and *x*-axis are rendered using icy software [60]. Even though each scatterer has weak scatter, severe distortions could occur in RI tomograms when individual scatterers are closely located to one another [30]. Even, when the weak scatterers are overlapped along the optical axis as in Fig. 6, the proposed method reveals the structure that is hard to observe under the approximation as indicated by the white arrow.

## 5. Discussions

Similar nonlinear inverse problems have arisen in other types of inverse scattering problems such as diffuse optical tomography (DOT), electric impadence tomography (EIT), and magnetic resonance elastography. Based on the observation that the unknown constitutive parameters are usually sparsely distributed and do not vary during multiple illuminations in these problems, we have shown that the nonlinearity can be decoupled by exploiting the joint sparsity [46–49]. Specifically, it was shown that the unknown constitutive parameters can be reconstructed using a similar simple three-step approach as described in this paper. While superior reconstruction performance was demonstrated using numerical simulations, it was not clear why such a resolution enhancement was achieved and whether resolution enhancement can be achieved in real experimental setup.

In this regard, this work is the first theoretical explanation and empirical verification why such resolution enhancement occurs. In fact, another important contribution of this paper is to derive a general principle - resolution enhancement from sub-wavelength self-resonance modes - which can be universally applied in various inverse scattering problems. Specifically, for imaging unknown specimen with high contrast against the homogeneous background, the sub-wavelength self-resonance modes within the specimen can be amplified and propagated to the far-field, so without using any engineered medium or tailored point spread function (PSF), one can achieve super-resolution.

## 6. Conclusion

In this paper, we proposed a novel method to break through the Born and Rytov approximation to achieve resolution enhancement for ODT. Using the recent theory of super-resolution in resonant media, we showed that if the multiple scattering can be handled to estimate the *K _{D}* operator accurately, the resolution enhancement can be expected from the sub-wavelength self-resonance modes that propagate to the far-field. Moreover, we provided an algorithm to estimate

*K*operator using a joint sparse recovery algorithm from compressed sensing. We compared the proposed algorithm with the first order Rytov approximation, which is commonly used to linearize the nonlinear integral solution of system equation of ODT using various samples such as a 5-

_{D}*μm*-diameter fused silica microsphere, an RBC, a hepatocyte, and multiple fused silica microspheres. Results showed that the proposed method successfully overcame the missing cone problem by showing increased RI values and reduced elongation of the RI tomograms along the optical axis. More importantly, in case of the sample that consists of complex structures such as lipid oil droplets in hepatocytes and multiple fused silica microspheres, the proposed method revealed clear structures compared to the first order approximation, resulting in the improvement in the resolution.

## Funding

National Research Foundation of Korea (NRF) (NRF-2016R1A2B3008104); Korea Research Fellowship Program through the National Research Foundation (NRF) funded by the Ministry of Science and ICT (NRF-2015H1D3A1062400).

## Disclosures

Mr. Lee and Prof. Park have financial interests in TomoCube Inc., a company that commercializes optical diffraction tomography.

## References and links

**1. **F. Charrière, A. Marian, F. Montfort, J. Kuehn, T. Colomb, E. Cuche, P. Marquet, and C. Depeursinge, “Cell refractive index tomography by digital holographic microscopy,” Opt. Lett. **31**(2), 178–180 (2006). [CrossRef] [PubMed]

**2. **M. Fauver, E. Seibel, J. R. Rahn, M. Meyer, F. Patten, T. Neumann, and A. Nelson, “Three-dimensional imaging of single isolated cell nuclei using optical projection tomography,” Opt. Express **13**(11), 4210–4223 (2005). [CrossRef] [PubMed]

**3. **G. G. Levin, G. N. Vishnyakov, S. C. Zakarian, A. V. Likhachov, V. V. Pickalov, G. I. Kozinets, J. K. Novoderzhkina, and E. A. Streletskaya, “Three-dimensional limited-angle microtomography of blood cells: experimental results,” in *BiOS’98 International Biomedical Optics Symposium* (International Society for Optics and Photonics1998), pp. 159–164.

**4. **V. Lauer, “New approach to optical diffraction tomography yielding a vector equation of diffraction tomography and a novel tomographic microscope,” J. Microsc. **205**(2), 165–176 (2002). [CrossRef] [PubMed]

**5. **B. Simon, M. Debailleul, V. Georges, V. Lauer, and O. Haeberlé, “Tomographic diffractive microscopy of transparent samples,” Eur. Phys. J. Appl. Phys. **44** (01), 29–35 (2008). [CrossRef]

**6. **S. O. Isikman, W. Bishara, S. Mavandadi, W. Y. Frank, S. Feng, R. Lau, and A. Ozcan, “Lens-free optical tomographic microscope with a large imaging volume on a chip,” Proc. Natl. Acad. Sci. U.S.A. **108**(18), 7296–7301 (2011). [CrossRef] [PubMed]

**7. **L. Yu and M. K. Kim, “Wavelength-scanning digital interference holography for tomographic three-dimensional imaging by use of the angular spectrum method,” Opt. Lett. **30**(16), 2092 (2005). [CrossRef] [PubMed]

**8. **J. Kühn, F. Montfort, T. Colomb, B. Rappaz, C. Moratal, N. Pavillon, P. Marquet, and C. Depeursinge, “Submicrometer tomography of cells by multiple-wavelength digital holographic microscopy in reflection,” Opt. Lett. **34**(5), 653–655 (2009). [CrossRef] [PubMed]

**9. **M. Potcoava and M. Kim, “Optical tomography for biomedical applications by digital interference holography,” Meas. Sci. Technol. **19**(7), 074010 (2008). [CrossRef]

**10. **U. S. Kamilov, I. N. Papadopoulos, M. H. Shoreh, A. Goy, C. Vonesch, M. Unser, and D. Psaltis, “Learning approach to optical tomography,” Optica **2**(6), 517–522 (2015). [CrossRef]

**11. **J. A. Izatt, E. A. Wanson, J. G. Fujimoto, M. R. Hee, and G. N. Owen, “Optical coherence microscopy in scattering media,” Opt. Lett. **19**(8), 590–592 (1994). [CrossRef] [PubMed]

**12. **M. A. Choma, A. K. Ellerbee, C. Yang, T. L. Creazzo, and J. A. Izatt, “Spectral-domain phase microscopy,” Opt. Lett. **30** (10), 1162–1164 (2005). [CrossRef] [PubMed]

**13. **Z. Wang, D. L. Marks, P. S. Carney, L. J. Millet, M. U. Gillette, A. Mihi, P. V. Braun, Z. Shen, S. G. Prasanth, and G. Popescu, “Spatial light interference tomography (slit),” Opt. Express **19**(21), 19907–19918 (2011). [CrossRef] [PubMed]

**14. **O. Haeberlé, K. Belkebir, H. Giovaninni, and A. Sentenac, “Tomographic diffractive microscopy: basics, techniques and perspectives,” J. Mod. Opt. **57**(9), 686–699 (2010). [CrossRef]

**15. **R. Fiolka, K. Wicker, R. Heintzmann, and Stemmer, “Simplified approach to diffraction tomography in optical microscopy,” Opt. Express **17**(15), 12407–12417 (2009). [CrossRef] [PubMed]

**16. **S. S. Kou and C. J. Sheppard, “Image formation in holographic tomography: high-aperture imaging conditions,” Appl. Opt. **48**(34), H168–H175 (2009). [CrossRef] [PubMed]

**17. **Y. Ruan, P. Bon, E. Mudry, G. Maire, P. C. Chaumet, H. Giovannini, K. Belkebir, A. Talneau, B. Wattellier, S. Monneret, and A. Sentenac, “Tomographic diffractive microscopy with a wavefront sensor,” Opt. Lett. **37**(10), 1631–1633 (2012). [CrossRef] [PubMed]

**18. **S. Uttam, S. A. Alexandrov, R. K. Bista, and Y. Liu, “Tomographic imaging via spectral encoding of spatial frequency,” Opt. Express **21**(6), 7488–7504 (2013). [CrossRef] [PubMed]

**19. **A. Kuś, M. Dudek, B. Kemper, M. Kujawińska, and A. Vollmer, “Tomographic phase microscopy of living three-dimensional cell cultures,” J. Biomed. Opt. **19**(4), 046009 (2014). [CrossRef]

**20. **L. Su, L. Ma, and H. Wang, “Improved regularization reconstruction from sparse angle data in optical diffraction tomography,” Appl. Opt. **54**(4), 859–868 (2015). [CrossRef] [PubMed]

**21. **L. Ma, H. Wang, L. Su, Y. Li, and H. Jin, “Digital holographic microtomography with few angle data-sets,” J. Mod. Opt. **61**(14), 1140–1146 (2014). [CrossRef]

**22. **U. S. Kamilov, I. N. Papadopoulos, M. H. Shoreh, A. Goy, C. Vonesch, M. Unser, and D. Psaltis, “Optical tomographic image reconstruction based on beam propagation and sparse regularization,” IEEE Transactions on Computational Imaging **2**(1), 59–70 (2016). [CrossRef]

**23. **F. Merola, L. Miccio, P. Memmolo, G. Di Caprio, A. Galli, R. Puglisi, D. Balduzzi, G. Coppola, P. Netti, and P. Ferraro, “Digital holography as a method for 3d imaging and estimating the biovolume of motile cells,” Lab Chip **13**(23), 4512–4516 (2013). [CrossRef] [PubMed]

**24. **Y. Cotte, F. Toy, P. Jourdain, N. Pavillon, D. Boss, P. Magistretti, P. Marquet, and C. Depeursinge, “Marker-free phase nanoscopy,” Nature Photon. **7**(2), 113–117 (2013). [CrossRef]

**25. **Y. Sung and R. R. Dasari, “Deterministic regularization of three-dimensional optical diffraction tomography,” J. Opt. Soc. Am. A **28**(8), 1554–1561 (2011). [CrossRef]

**26. **K. Kim, K. Kim, H. Park, J. Ye, and Y. Park, “Real-time visualization of 3-D dynamic microscopic objects using optical diffraction tomography,” Opt. Express **21**(26), 32269–32278 (2013). [CrossRef]

**27. **K. Kim, Z. Yaqoob, K. Lee, J. Kang, Y. Choi, P. Hosseini, P. So, and Y. Park, “Diffraction optical tomography using a quantitative phase imaging unit,” Opt. Lett. **39**(24), 6935–6938 (2014). [CrossRef] [PubMed]

**28. **J. Lim, K. Lee, K. H. Jin, S. Shin, S. Lee, Y. Park, and J. Ye, “Comparative study of iterative reconstruction algorithms for missing cone problems in optical diffraction tomography,” Optics express **23**(13), 16933–16948 (2015). [CrossRef] [PubMed]

**29. **M. Born and E. Wolf, *Principles of Optics: Electromagnetic Theory of Propagation, Interference and Diffraction of Light.* (CUP Archive, 1999). [CrossRef]

**30. **M. Azimi and A. Kak, “Distortion in diffraction tomography caused by multiple scattering,” IEEE Trans. Med. Imaging **2**(4), 176–195 (1983). [CrossRef] [PubMed]

**31. **G. Lerosey, J. de Rosny, A. Tourin, and M. Fink, “Focusing beyond the diffraction limit with far-field time reversal,” Science **315**(5815), 1120–1122 (2007). [CrossRef] [PubMed]

**32. **F. Lemoult, G. Lerosey, J. de Rosny, and M. Fink, “Resonant metalenses for breaking the diffraction barrier,” Phys. Rev. Lett. **104**(20), 203901 (2010). [CrossRef] [PubMed]

**33. **F. Lemoult, A. Ourir, J. de Rosny, A. Tourin, M. Fink, and G. Lerosey, “Time reversal in subwavelength-scaled resonant media: beating the diffraction limit,” Int. J. Microw. Sci. Technol. **2011**, 425710 (2011). [CrossRef]

**34. **F. Lemoult, M. Fink, and G. Lerosey, “Acoustic resonators for far-field control of sound on a subwavelength scale,” Phys. Rev. Lett. **107**(6), 064301 (2011). [CrossRef] [PubMed]

**35. **S. Arhab, G. Soriano, Y. Ruan, G. Maire, A. Talneau, D. Sentenac, P. C. Chaumet, K. Belkebir, and H. Giovannini, “Nanometric resolution with far-field optical profilometry,” Phys. Rev. Lett. , **111**(5), 053902 (2013). [CrossRef] [PubMed]

**36. **J.-H. Park, C. Park, H. Yu, J. Park, S. Han, J. Shin, S. H. Ko, K. T. Nam, Y.-H. Cho, and Y. Park, “Subwavelength light focusing using random nanoparticles,” Nat. Photon. **7**, 454–458 (2013). [CrossRef]

**37. **C. Park, J.-H. Park, C. Rodriguez, H. Yu, M. Kim, K. Jin, S. Han, J. Shin, S. H. Ko, and K. T. Nam, “Full-field subwavelength imaging using a scattering superlens,” Phys. Rev. Lett. **113**, 113901 (2014). [CrossRef] [PubMed]

**38. **H. Ammari and H. Zhang, “Super-resolution in high-contrast media,” Proc. R. Soc. A. **471**(2178), 20140946 (2015). [CrossRef]

**39. **H. Ammari and H. Zhang, “A mathematical theory of super-resolution by using a system of sub-wavelength Helmholtz resonators,” Commun. Math. Phys. **337**(1), 379–428 (2015). [CrossRef]

**40. **S. F. Cotter, B. D. Rao, K. Engan, and K. Kreutz-Delgado, “Sparse solutions to linear inverse problems with multiple measurement vectors,” IEEE Trans. Signal Process. **53**(7), 2477–2488 (2005). [CrossRef]

**41. **D. Malioutov, M. Cetin, and A. S. Willsky, “A sparse signal reconstruction perspective for source localization with sensor arrays,” IEEE Trans. Signal Process. **53**(8), 3010–3022 (2005). [CrossRef]

**42. **D. P. Wipf and B. D. Rao, “An empirical Bayesian strategy for solving the simultaneous sparse approximation problem,” IEEE Trans. Signal Process. **55**(7), 3704–3716 (2007). [CrossRef]

**43. **J. M. Kim, O.K. Lee, and J. C. Ye, “Compressive MUSIC: revisiting the link between compressive sensing and array signal processing,” IEEE Trans. Inf. Theory **58**(1), 278–301 (2012). [CrossRef]

**44. **J. M. Kim, O. K. Lee, and J. C. Ye, “Improving noise robustness in subspace-based joint sparse recovery,” IEEE Trans. Signal Process. **60**(11), 5799–5809 (2012). [CrossRef]

**45. **J. C. Ye and S. Y. Lee, “Non-iterative exact inverse scattering using simultaneous orthogonal matching pursuit (S-OMP),” in *IEEE International Conference on Acoustics, Speech and Signal Processing* (IEEE2008), pp. 2457–2460.

**46. **O. Lee, J. M. Kim, Y. Bresler, and J. C. Ye, “Compressive diffuse optical tomography: noniterative exact reconstruction using joint sparsity,” IEEE Trans. on Medical Imaging **30**(5), 1129–1142 (2011). [CrossRef] [PubMed]

**47. **O. Lee and J. C. Ye, “Joint sparsity-driven non-iterative simultaneous reconstruction of absorption and scattering in diffuse optical tomography,” Opt. Express **21**(22), 26589–26604 (2013). [CrossRef] [PubMed]

**48. **O. K. Lee, H. Kang, J. C. Ye, and M. Lim, “A non-iterative method for the electrical impedance tomography based on joint sparse recovery,” Inverse Probl. **31**(7), 075002 (2015). [CrossRef]

**49. **J. Yoo, Y. Jung, M. Lim, J. C. Ye, and A. Wahab, “A joint sparse recovery framework for accurate reconstruction of inclusions in elastic media,” SIAM J. Imag. Sci. **10**(3), 1104–1138 (2017). [CrossRef]

**50. **O. K. Lee, J. M. Kim, Y. Bresler, and J. C. Ye, “Compressive diffuse optical tomography: non-iterative exact reconstruction using joint sparsity,” IEEE Trans. Med. Imag. **30**(5), 1129–1142 (2011). [CrossRef]

**51. **J. Chena and X. Huo, “Theoretical results on sparse representations of multiple measurement vectors,” IEEE Trans. Signal Process. **54**(12), 4634–4643 (2006). [CrossRef]

**52. **M. E. Davies and Y. C. Eldar, “Rank awareness in joint sparse recovery,” IEEE Trans. Inf. Theory **58**(2), 1135–1146 (2012). [CrossRef]

**53. **Y. Kim, H. Shim, K. Kim, H. Park, J. Heo, J. Yoon, C. Choi, S. Jang, and Y. Park, “Common-path diffraction optical tomography for investigation of three-dimensional structures and dynamics of biological cells,” Opt. Express **22**(9),10398–10407 (2014). [CrossRef] [PubMed]

**54. **O. Zhernovaya, O. Sydoruk, V. Tuchin, and A. Douplik, “The refractive index of human hemoglobin in the visible range,” Phys. Med. Biol. **56**(13), 4013–4021 (2011). [CrossRef] [PubMed]

**55. **K. Kim, H. Yoon, M. Diez-Silva, M. Dao, R. R. Dasari, and Y. Park, “High-resolution three-dimensional imaging of red blood cells parasitized by plasmodium falciparum and in situ hemozoin crystals using optical diffraction tomography,” J. Biomed. Opt. **19**(1), 011005 (2014). [CrossRef]

**56. **W. Deng, W. Yin, and Y. Zhang, “Group sparse optimization by alternating direction method,” Proc. SPIE **8858**, 88580R (2013). [CrossRef]

**57. **T. Goldstein and S. Osher, “The split Bregman method for l1-regularized problems,” SIAM. J. Imaging. Sci. **2**(2), 323–343 (2009). [CrossRef]

**58. **M. V. Afonso, J. M. Bioucas-Dias, and M. A. Figueiredo, “An augmented Lagrangian approach to the constrained optimization formulation of imaging inverse problems,” IEEE Trans. Image Process. **20**(3), 681–695 (2011). [CrossRef]

**59. **J. Lim, A. Wahab, G. Park, K. Lee, Y. Park, and J. C. Ye, “Source Code: Beyond Born-Rytov limit for super-resolution optical diffraction tomography,” figshare (2017) [Retrieved on Nov. 10, 2017], https://doi.org/10.6084/m9.figshare.5589787.

**60. **F. de Chaumont, S. Dallongeville, N. Chenouard, N. Hervé, S. Pop, T. Provoost, V. Meas-Yedid, P. Pankajakshan, T. Lecomte, Y. Le Montagner, T. Lagache, A. Dufour, and J.-C. Olivo-Marin, “Icy: an open bioimage informatics platform for extended reproducible research,” Nat. Methods **9**(7), 690–696 (2012). [CrossRef] [PubMed]