A complete and computationally efficient numerical model of aplanatic solid immersion lens scanning microscope

Rui Chen; Krishna Agarwal; Colin J. R. Sheppard; Jacob C. H. Phang; Xudong Chen

doi:10.1364/OE.21.014316

1. Introduction

The need for high resolution imaging necessitates the use of high numerical aperture (NA) optical systems. Solid immersion lens (SIL) microscopy provides higher spatial resolution, improved light collection efficiency, and the capability of subsurface imaging [1–3]. Experimental implementations of this technique as a wide-field microscope [2, 4] and a scanning microscope [5–9] have been reported in the field of imaging solid state devices. For the SIL microscope, only the focusing of incident light through SIL has been theoretically investigated in detail [3, 8, 10–12]. However, the remaining optical path of the microscope also has a significant influence on the resolution and image quality, as studied for aplanatic solid immersion lens (ASIL) wide-field microscope in [13]. It indicates the importance of such a complete study and modelling of ASIL scanning microscope (ASIL-SM) as well. A computational model allows a controlled study of the impact of each individual or a set of parameters, which is often difficult to achieve in experiments.

A paraxial scanning microscope without SIL was investigated in [14–16]. An immersion type confocal scanning light microscope (CSLM) with numerical aperture NA = 1.3 was investigated experimentally in [17]. Furthermore, [18, 19] considered image formation of dielectric point scatterers in conventional and confocal microscopes with high NA. In addition, a numerical method, viz., finite-difference time-domain (FDTD) was used for modeling the optical microscope [20, 21]. A four-step process, i.e., illumination, scattering, resampling, and image formation, was designed to treat the numerical model of the optical microscope. Due to the huge computational load of FDTD, examples of only small scatterers were given for a scanning microscope [20] and some simplifications were assumed in order to deal with larger object structures [21].

The problem of developing a computation efficient model for ASIL-SM is very challenging. The case of ASIL microscope is very different from other free space or oil-immersion microscopes [22]. It is not only a high NA system (higher than oil immersion lens for silicon substrate), it also involves refraction along a spherical interface which is offset from the focal point of the objective. Thus, the integrals are quite complicated and a minor approximation may lead to significant error in computation. Further, the computation of the vector diffraction integrals in ASIL microscope [11, 23] suffers from a big computational challenge due to the highly oscillatory behavior of the integrands. In the meanwhile, the computational requirements of a scanning system is several times more than its wide-field counterpart for any microscopy system, depending upon the scanning resolution. Thus, the development of a numerical model that is fast as well as accurate is of prime importance and significant challenge for ASIL-SM. Addressing this challenge is the main motivation of this paper.

In addition to the above, we also seek to address the properties and engineering parameters of ASIL-SM. We report the theoretical resolution capability for small scatterers. We also study the effect of polarization, NA, and detector pinhole size on the resolution and imaging characteristics.

While using the basic theoretical premise developed in [13], we propose a complete and computationally efficient model of ASIL-SM with finite-sized detector. The highlights of our model are as follows:

This model is a complete model in the theoretical sense, since it considers the entire process from focusing of the incident light to the image generation with minimal and practical assumptions.
Our model has been made computationally efficient by the use of the conjugate gradient fast Fourier transform (CG-FFT) [24] in the second subsystem, chirp z transform (CZT) in the first subsystem and CZT with fast convolution (CZT-FC) in the third subsystem.
The modular nature of the model (Fig. 1) makes it easy to simulate and engineer several system parameters, such as the polarization of the illumination, the object structures in the focal plane, or the pinhole size.
In our knowledge, this paper is the first paper demonstrating the imaging properties of the ASIL-SM microscope as a complete system. Besides studying the theoretical resolution using small (point-like) scatterers, we also demonstrate the capability of simulating images of large and complicated object structures. Interesting and practical object patterns are used for investigation. Especially, practically important examples of USAF target and an object structure with various materials used in silicon fabrication industry are also illustrated.
The effect of polarization, NA, and detector pinhole sizes are studied as important system design parameters. Such a study is of practical importance for predicting and understanding the microscopy results better.

Fig. 1 Block diagram of the computational model.

Download Full Size | PDF

The structure of this paper is as follows. The proposed model of ASIL-SM is presented in Section 2. It is followed in Section 3 by numerical simulations and analysis for small scatterers and large objects. Finally, we conclude this paper in Section 4.

2. Computational model of ASIL scanning microscope with finite sized detector

2.1. Background

The diagram of ASIL-SM is shown in Fig. 2. We assume that the radius R of the ASIL is far larger than the wavelength λ in air, and that the refractive index of ASIL n_sil is same as that of the substrate. The NA of the ASIL is NA_sil = (n_sil/n_obj)²NA_obj, where NA_obj is the NA of the objective and n_obj is the refractive index of material between ASIL and objective. Further, f_obj and f_ccd are the focal lengths of the objective lens and the detector lens respectively, and R_PH denotes the radius of pinhole. With reference to Fig. 1, the first subsystem comprises of the computation of the electric fields in the focal region formed due to the focusing of the incident beam (coherent, collimated, and of a certain polarization). In terms of microscope components, it involves the objective lens and ASIL. The second subsystem comprises of the electromagnetic interaction of the focal fields with the object structure, which results in the induction of secondary sources. It should include multiple scattering effect between different portions of the discretized object structure. More details about the second subsystem can be found in [13, 24]. In this paper, we refrain from going into detail about the subsystem 2, and concentrate more on the subsystems 1 and 3. The third subsystem comprises of the radiation of electric fields from the induced sources, which pass through the ASIL, objective, detector lens and the pinhole, in order to compute a net intensity value corresponding to the focal point in the ASIL region.

Fig. 2 Diagrammatic description of ASIL scanning microscope. The refractive index of ASIL is the same as that of the substrate where the object structures are present.

Download Full Size | PDF

The basic computational models of subsystems 1 and 3 are the same as in [13]. It is notable that our model of subsystem 1 is in line with the extensive work done on modeling the focusing of incident light. Thus, like other works on focusing, our model of subsystem 1 also includes the evaluation of the diffraction integrals [10, 11]. Another approach is to use a vectorial spherical harmonics theory for focusing [12]. Both direct evaluation of the diffraction integral as well as the use of vectorial spherical harmonics theory suffer from a big computational challenge due to highly oscillatory behavior of the integrands. Further, even for the subsystem 3, the computation of the dyadic Green’s function (DGF) also involves the computation of diffraction integrals with highly oscillatory integrands [23, 25]. The computational load for subsystem 3 is significantly larger than subsystem 1 since the DGF should be computed for each pair of points in the object and image regions.

Our model uses an alternative and relatively less-used approach of the Fourier transform (FT) to evaluate the diffraction integrals of subsystems 1 and 3. Some relevant works on using FT in diffraction integrals for optical setups different from ASIL-SM appear in [22, 26–29]. However, none of these works can be applied directly for the case ASIL-SM. This is because these works do not deal with the presence of a spherical interface such as the ASIL behind the objective lens. The problem is further complicated by the fact that the focal point of the objective is not coincident with the center of the ASIL, the refractive index mismatch is higher, and the numerical aperture of ASIL-SM is large (non-paraxial).

While the use of FT helps in dealing with the oscillatory integrand issue of the diffraction integral and reducing the computational complexity, we augment the computation speed further by using the CZT [22, 30]. Finally, we use the convolution theorem [31] to further speed up the computation of the DGF. The details are provided in sections 2.2 to 2.5.

2.2. Focusing of incident light through ASIL

The optical diagram of the first subsystem is given in Fig. 3. A non-paraxial monochromatic wave $E_{inc} (θ) = {[E_{inc}^{x}, E_{inc}^{y}, E_{inc}^{z}]}^{T}$ is incident on the objective lens with ${N A}_{obj} = n_{obj} sin θ_{obj}^{\max}$ , where $sin θ_{obj}^{\max}$ is the maximum value of θ_obj. After refraction at the objective lens and the spherical interface of the ASIL, the local polarization components immediately inside the ASIL surface $E_{sil} = {[E_{sil}^{x}, E_{sil}^{y}, E_{sil}^{z}]}^{T}$ for the ray along the angular direction (θ_obj, ϕ_obj) can be described by the following equation [13]:

E_{sil} = {(\frac{n_{0}}{n_{obj}} cos θ_{obj})}^{1 / 2} \frac{f_{obj}}{r_{OA}} exp (i k_{obj} f_{obj} - i k_{obj} r_{O A}) T \cdot E_{inc} (θ)

where the matrix, T, is written as

T = \frac{1}{2} [\begin{matrix} (t_{s} + t_{p} cos θ_{sil}) - (t_{s} - t_{p} cos θ_{sil}) cos 2 ϕ_{obj} & - (t_{s} - t_{p} cos θ_{sil}) sin 2 ϕ_{obj} & - 2 t_{p} sin θ_{sil} cos ϕ_{obj} \\ - (t_{s} - t_{p} cos θ_{sil}) sin 2 ϕ_{obj} & (t_{s} + t_{p} cos θ_{sil}) + (t_{s} - t_{p} cos θ_{sil}) cos 2 ϕ_{obj} & - 2 t_{p} sin θ_{sil} sin ϕ_{obj} \\ 2 t_{p} sin θ_{sil} cos ϕ_{obj} & 2 t_{p} sin θ_{sil} sin ϕ_{obj} & 2 t_{p} cos θ_{sil} \end{matrix}] .

The factor

{(\frac{n_{0}}{n_{obj}} cos θ_{obj})}^{1 / 2}

is the apodization function to account for the sine condition and energy conservation. t_s and t_p are transmission coefficients at the interface of the ASIL when the wave is traveling from the objective to the ASIL. The term of

\frac{f_{obj}}{r_{O A}} exp (i k_{obj} f_{obj} - i k_{obj} r_{O A})

denotes the change in the magnitude and phase due to the wave propagating from the Gaussian reference sphere (GRS) to the outer surface of the ASIL. It is worth mentioning that Eq. (1) is more accurate for ASIL case [13], whereas [11] modelled the focusing effect for an SIL of arbitrary thickness.

Fig. 3 Equivalent representation of subsystem I. The details can be found in [13].

Download Full Size | PDF

Taking the point F in Fig. 3 as the reference point, we can express the electric field at a point P(x_sil, y_sil, z_sil) in the focal region of ASIL by the following expression [13]:

E (x_{sil}, y_{sil}, z_{sil}) = \iint_{Ω} \frac{i r_{F A} exp (i k_{sil} r_{F A})}{2 π k_{z_{sil}}} E_{sil} exp (i k_{x_{sil}} x_{sil} + i k_{y_{sil}} y_{sil} + i k_{z_{sil}} z_{sil}) d k_{x_{sil}} d k_{y_{sil}},

where Ω denotes the integration domain

{(k_{x_{sil}}^{2} + k_{y_{sil}}^{2})}^{1 / 2} \leq k_{sil} sin θ_{sil}^{\max}

. Considering the relationship n_objr_OA = n_silr_FA, we can rewrite Eq. (3) as

E (x_{sil}, y_{sil}, z_{sil}) = α \iint_{Ω} E_{f f t} (k_{x_{sil}}, k_{y_{sil}}) exp (i k_{x_{sil}} x_{sil} + i k_{y_{sil}} y_{sil}) d k_{x_{sil}} d k_{y_{sil}}

where

E_{f f t} = {(cos θ_{obj})}^{1 / 2} \frac{exp (i k_{z_{sil}} z_{sil})}{cos θ_{sil}} T \cdot E_{inc} (θ)

and

α = - \frac{i f_{obj}}{2 π k_{sil}} \frac{n_{obj}}{n_{sil}} {(\frac{n_{0}}{n_{obj}})}^{1 / 2} exp (i k_{obj} f_{obj})

. E(x_sil, y_sil, z_sil) can be understood as the two dimensional inverse Fourier transform (2D-IFT) of E_fft(k_{x_sil}, k_{y_sil}), except for the constant factor α. After discretization of Eq. (4), we have

\begin{matrix} E (p, q, z_{sil}) = α \sum_{m = - M / 2}^{M / 2 - 1} \sum_{n = - N / 2}^{N / 2 - 1} E_{f f t} (m, n) exp [2 π i (\frac{p m}{M} + \frac{q n}{N})] Δ k_{x_{sil}} Δ k_{y_{sil}}; \\ (- M / 2 \leq p \leq M / 2 - 1), (- N / 2 \leq q \leq N / 2 - 1), \end{matrix}

where

M = 2 sin θ_{sil}^{\max} k_{sil} / Δ k_{x_{sil}}

and

N = 2 sin θ_{sil}^{\max} k_{sil} / Δ k_{y_{sil}}

are the numbers of samples in the x and y directions in both spatial and spatial frequency domains, respectively. Therefore, the electric field in the focal region of the ASIL using FFT is given by

E (p, q, z_{sil}) = α IFFT [E_{f f t} (m, n)] Δ k_{x_{sil}} Δ k_{y_{sil}} .

For the FFT algorithm, the spatial domain (x_sil, y_sil) has the same sample size (M and N) as the spatial frequency domain (k_{x_sil}, k_{y_sil}) respectively. However, in the spatial domain, the region of interest is near and around the focal point only, which is only a small portion compared with the required sample size of the Fourier frequency domain for capturing the highly oscillatory behavior of the integrand. In order to sufficiently resolve the details of the region of interest, zero-padding is often used to increase the output sampling size, and thus the resolution of the spatial domain is determined by

Δ x_{sil} = 2 π / (M_{f f t} Δ k_{x_{sil}}), Δ y_{sil} = 2 π / (N_{f f t} Δ k_{y_{sil}}),

where M_fft and N_fft is the sample size after zero-padding. It subsequently increases the size of the input matrix and the computation time.

Compared to FFT, CZT allows us to avoid the use of the zero-padding, and select only the region of interest as the output space. Due to this, the number of computations is reduced significantly and the CZT scheme is much faster than FFT. For CZT, we choose the same spatial frequency domain as FFT and the spatial domain as x_sil ∈ [x₁, x₂] and y_sil ∈ [y₁, y₂], then Eq. (5) is given by [30, 31]

\begin{matrix} E (p, q) = α \sum_{m = - M / 2}^{M / 2 - 1} \sum_{n = - N / 2}^{N / 2 - 1} E_{f f t} (m, n) A_{x}^{- m} A_{y}^{- n} W_{x}^{p m} W_{x}^{q n} Δ k_{x_{sil}} Δ k_{y_{sil}}; \\ (- P / 2 \leq p \leq P / 2), (- Q / 2 \leq q \leq Q / 2), \end{matrix}

where A_x = exp(−iΔk_{x_sil}x₁), A_y = exp(−iΔk_{y_sil}y₁), W_x = exp[iΔk_{x_sil}(x₂ −x₁)/P] and W_y = exp[iΔk_{y_sil} (y₂ − y₁)/Q]. The flexibility of the CZT allows the choice of an arbitrarily shaped and sized output region, and expresses the CZT as a convolution, permitting the use of the well-known convolution algorithm. For the details, please refer to [30, 31].

2.3. Fast computation of DGF using CZT

The DGF for ASIL microscope has been presented recently in [23]. By means of the DGF, we can obtain the intensity distribution in the CCD/detector region corresponding to an induced current distribution in the focal plane of ASIL. In this section, the fast calculation of DGF for the imaging of scattered light using CZT is presented. The optical diagram is shown in Fig. 4. For a current dipole Il(r_sil) located at r_sil(x_sil, y_sil, z_sil) close to the aplanatic point of ASIL, O_sil, the electric field immediately after the interface of the Gaussian reference sphere representing the CCD lens is given by [23]

E_{c c d} = i ω μ \frac{exp (i k_{obj} f_{obj})}{8 π f_{obj}} exp (- i k_{sil} \cdot r_{sil}) {(\frac{n_{obj} cos θ_{c c d}}{n_{c c d} cos θ_{obj}})}^{\frac{1}{2}} [G_{x}, G_{y}, G_{z}] \cdot Il (r_{sil})

where k_sil = [k_{x_sil}, k_{y_sil}, k_{z_sil}]^T denotes the propagation vector of light inside the ASIL and

G_{x} = [\begin{matrix} (t_{s} + t_{p} cos θ_{c c d} cos θ_{sil}) - (t_{s} - t_{p} cos θ_{c c d} cos θ_{s i l}) cos 2 ϕ_{obj} \\ - (t_{s} - t_{p} cos θ_{c c d} cos θ_{sil}) sin 2 ϕ_{obj} \\ 2 t_{p} sin θ_{c c d} cos θ_{sil} cos ϕ_{obj} \end{matrix}],

G_{y} = [\begin{matrix} - (t_{s} - t_{p} cos θ_{c c d} cos θ_{sil}) sin 2 ϕ_{obj} \\ (t_{s} + t_{p} cos θ_{c c d} cos θ_{sil}) + (t_{s} - t_{p} cos θ_{c c d} cos θ_{sil}) cos 2 ϕ_{obj} \\ 2 t_{p} sin θ_{c c d} cos θ_{sil} sin ϕ_{obj} \end{matrix}],

G_{z} = - 2 t_{p} sin θ_{sil} [\begin{matrix} cos θ_{c c d} cos ϕ_{o b j} \\ cos θ_{c c d} sin ϕ_{o b j} \\ sin θ_{c c d} \end{matrix}] .

t_s and t_p are redefined transmission coefficients at the ASIL interface for a wave travelling from the ASIL to the objective [23], which are different from the transmission coefficients in Eq. (2). The far-field E_ccd can now be introduced into the angular spectrum representation and thus the electric field at a point r_ccd(x_ccd, y_ccd, z_ccd) in the CCD region is given by [23]

E (r_{c c d}) = i ω μ {\bar{\bar{G}}}_{P S F} \cdot Il (r_{sil}) .

In the present paper, the DGF is written differently from that in [23], such that it is amenable to a 2D-IFT format (similar to Eq. (4))

{\bar{\bar{G}}}_{P S F} = β \iint_{Ω_{c c d}} \bar{\bar{g}} exp [i k_{x_{c c d}} (x_{c c d} + M_{l a t} x_{sil}) + k_{y_{c c d}} (y_{c c d} + M_{lat} y_{sil})] d k_{x_{c c d}} d k_{y_{c c d}},

\bar{\bar{g}} = \frac{[G_{x}, G_{y}, G_{z}]}{cos θ_{c c d}} {(\frac{cos θ_{c c d}}{cos θ_{obj}})}^{1 / 2} exp (i k_{z_{c c d}} z_{c c d} - i k_{z_{sil}} z_{sil}),

where

β = \frac{- i}{16 π^{2} k_{c c d}} \frac{f_{c c d}}{f_{obj}} {(\frac{n_{obj}}{n_{c c d}})}^{1 / 2} exp (i k_{c c d} f_{c c d} + i k_{obj} f_{obj})

is a constant. Ω_ccd is the integration domain

{(k_{x_{c c d}}^{2} + k_{y_{c c d}}^{2})}^{1 / 2} \leq k_{c c d} sin θ_{c c d}^{\max}

, and

θ_{c c d}^{\max}

is the maximum value of θ_ccd. We also note that k_{x_sil}= −k_{x_ccd} M_lat and k_{y_sil} = −k_{y_ccd} M_lat, where

M_{l a t} = {(\frac{n_{sil}}{n_{obj}})}^{2} \frac{n_{obj} f_{c c d}}{n_{c c d} f_{obj}}

is transverse magnification of ASIL. Based on the similarity of Eq. (12a) with Eq. (4), we can follow the processes in Section 2.2 to compute the DGF efficiently and accurately using CZT.

Fig. 4 Equivalent representation of subsystem III. The details can be found in [13].

Download Full Size | PDF

While the adaptation of the fast techniques for focusing and propagation may appear to be similar to each other, we highlight that the adaptation of the propagation integral in subsystem 3 is more involved and requires care in computing the Fourier terms of each component of the DGF. We note that while the remaining procedure of applying CZT is similar to section 2.2, the modification of the DGF into the form of Eq. (12a) is crucial for being able to apply CZT.

2.4. Image formation using convolution theorem

Based on the analysis of Sections 2.2 and 2.3, the first and third subsystem can be evaluated using CZT instead of direct integration. We assume that an object structure is present in the focal plane of ASIL, as shown in Fig. 5(b). After we obtain the electric field distribution in the focal plane, shown in Fig. 5(a), the induced current distribution Il(x_sil, y_sil) in the focal plane of ASIL can be obtained by MoM [13]. The electric field in the focal plane of CCD region is then given by

E (x_{c c d}, y_{c c d}) = i ω μ {\bar{\bar{G}}}_{P S F} (x_{sil} + x_{c c d} / M_{lat}, y_{sil} + y_{ccd} / M_{lat}) \cdot I l (x_{sil}, y_{sil}),

where

{\bar{\bar{G}}}_{P S F} = \frac{β}{M_{lat}^{2}} \iint \bar{\bar{g}} exp [- i k_{x_{sil}} (x_{sil} + \frac{x_{c c d}}{M_{lat}}) - i k_{y_{sil}} (y_{sil} + \frac{y_{c c d}}{M_{lat}})] d k_{x_{sil}} d k_{y_{sil}}

which is just a modified expression of Eq. (12a). For Eq. (13), the numerical implementation is straightforward. In the domain where the object structure is located, an equidistant sampling x_sil = mΔx_sil and y_sil = nΔy_sil with m = −M/2,...,M/2 and n = −N/2,...,N/2, is made (or simply reused from the discretization in subsystem 1). The sampling in the focal plane of detector region is x_ccd = pΔx_ccd and y_ccd = qΔy_ccd with p = −M/2,...,M/2 and q = −N/2,...,N/2, where we choose Δx_ccd = M_latΔx_sil, Δy_ccd = M_latΔy_sil, such that the grid points in the detector region are image conjugates of the grid points in the object region.

Fig. 5 The description of the location of object structure when it is scanned relative to the optical microscope. (a) The focal spot and (b) the object domain in the focal plane of ASIL. The location of the object domain when (c) the first pixel (1,1) and (d) the (m,n) pixel are scanned, i.e,. the corresponding pixel is at the center of focal spot. The red arrows in (c) and (d) denote the scanning directions.

Download Full Size | PDF

Taking the element G_xx of DGF and Il_x(x_sil, y_sil) for example, we can express the electric field produced by these two elements in the discrete form

E_{x x} (p, q) = i ω μ \sum_{m} \sum_{n} {I l}_{x} (m, n) G_{x x} (p + m, q + n),

which is the two-dimensional discrete convolution of matrices Il_x(−p, −q) and G_xx(p, q), which can be expressed by

E_{x x} (p, q) = i ω μ ({I l}_{x} (- p, - q) * G_{x x} (p, q)) .

where the symbol ‘*’ denotes the convolution operator. According to the convolution theorem [31], FFT and IFFT can be used to evaluate this equation:

E_{x x} (p, q) = i ω μ IFFT {FFT [{I l}_{x} (- p, - q)] FFT [G_{x x} (p, q)]},

Similarly, other elements of DGF can also be evaluated using the fast convolution theorem and then we can obtain the image of object structure by considering all of the elements of DGF.

2.5. Scanning system implementation

For the scanning system, the object structure is assumed to be scanned relative to the optical system, and the total intensity of the light passing through a finite sized detector pinhole is collected at the detector. Noting the effects of employing a finite sized detector, the detected signal can be written as the integration of the intensity (the modulus square of the electric field vector) over the detector area S (i.e., the pinhole size),

I = \int_{S} {| E (x_{c c d}, y_{c c d}) |}^{2} d S = \int_{0}^{R_{P H}} \int_{0}^{2 π} {| E (x_{c c d}, y_{c c d}) |}^{2} d r_{c c d} d ϕ .

After discretization, the equidistant sampling of the object in the focal plane of ASIL are x_sil = mΔx_sil and y_sil = nΔy_sil with m ∈ [−M/2, M/2] and n ∈ [−N/2, N/2], as shown in Fig. 5(b). When the object is scanned laterally as shown in Figs. 5(c) and 5(d), the detector pinhole is fixed at the focal point of the detector lens and the object is shifted in the focal plane of ASIL. The scanning action makes each pixel of the object locate at the center of focal spot in turn. When the (m, n)th pixel is moved to the focal point of ASIL, the induced current distribution inside the focal spot in Fig. 5(d) produces an electric field in the detector plane. Therefore, the signal detected on the (m, n)th pixel in the detector region, I(m, n), determined by Eq. (18), corresponds to the (m, n)th pixel in the object region. After scanning is done, I(m, n), a matrix with the same dimension as discretization M × N, is the image of the original object structure using ASIL-SM with the pinhole of radius, R_PH.

2.6. Example of enhancement in computational efficiency

In order to understand the accuracy and computational efficiency of the proposed model, we take an example in which we compare the direct integral approach (DI) used in [13]. The setup of the numerical example is similar to the setup used in the numerical examples in section 3. In Table 1, we show the relative error and the time taken by the direct integral approach (DI), the FFT-approach (FFT), and the proposed CZT approach (CZT) using a personal computer with an Intel(R) Core2 Duo 3.16GHz processor and 3GB of RAM. The output spatial domain in the focal plane is [−1.0λ, 1.0λ] with 101 × 101 pixels. We also compare the effect of the value of M, N. It is seen that even with very large value of M, N, M = N = 800, the computation using CZT is about 36 times lesser than DI. Further, for an error of less than 10⁻², M = N = 100 is sufficient and the speed up is by a factor of 2000. We have used M = N = 200 for both subsytem 1 and subsystem 3 in the remaining part of the paper. We also summarize the computation efficiency for the whole system in Table 2 using the same computer. It is seen that the computation of a few hours in DI approach get reduced to a few seconds using the proposed approach.

Table 1. Comparison of calculation times and relative error for DI, FFT and CZT methods

View Table | View all tables in this article

Table 2. Comparison between computational time using direct integral (DI) in [13] and CZT methods presented in this paper for each scan point of the ASIL-SM.

View Table | View all tables in this article

3. Numerical simulations and analysis

In all the numerical simulations below, we consider the following settings. The SIL is made of silicon with a refractive index of n_sil = 3.5 and the operating wavelength is 1340 nm in free space, both of which are appropriate for subsurface imaging of silicon chips. The radius of the SIL is R = 1.5 mm. The incident light has unity magnitude. The focal length of the objective and detector lens are f_obj = 0.01 m and f_ccd = 0.1 m, respectively. The objective and detector media are free space, i.e., n_obj = 1 and n_ccd = 1, respectively. Thus, the lateral magnification of the microscope is M_lat = 122.5. The discretization (scanning pixel size) is Δx_sil = Δy_sil = 0.02λ for both small scatterer and large object examples. The semiconductor industry typically uses planar resolution target chips (like Metrochip [32] and MRS-5 [33]) for characterizing the resolution of the microscopes used in silicon failure analysis. Thus, we consider planar object structures only and thickness of the plane is taken as 0.02λ.

In our numerical examples, we have considered two polarizations, x-polarization and circular polarization. Further, for small scatterers, we consider two values of NA, NA_sil = 2.4 and 3.3. Throughout the numerical examples, we have considered various pinhole sizes like 1 μm, 25 μm, and 100 μm. We highlight that pinhole radius of 1 μm and 100μm correspond to the confocal and wide-field cases respectively. In our opinion, polarization and detector pinhole radius are two important system parameters that can be used ASIL-SM engineering. We make a note that pinhole radius of 1 μm reduces the signal strength considerably, but this effect can be practically dealt with by using high sensitivity and low noise InGaAs detector.

3.1. Imaging of small scatterers

For the example of small scatterers [13], we show the normalized intensity distribution of the image of the small scatterer for NA_sil = 3.3 and 2.4 in Figs. 6 and 7 respectively. The comparison for pinhole radius of R_PH = 1μm and 100μm and different polarizations are included. The paraxial approximation [34], referred to as PA, is also provided.

Fig. 6 Normalized intensity distribution of the image of a small scatterer along a lateral axis for different pinhole radii for NA_sil = 3.3, (a) R_PH = 1μm and (b) R_PH = 100μm, using x-polarized illumination (two curves: y_sil = 0 and x_sil = 0) and circularly polarized illumination, respectively. PA denotes paraxial approximation.

Download Full Size | PDF

Fig. 7 The caption is the same as that of Fig. 6 except for NA_sil = 2.4.

Download Full Size | PDF

For x-polarized illumination, it is seen that when high NA (NA_sil = 3.3, Fig. 6) is employed, the image asymmetry of a small scatterer is much larger than that using low NA (NA_sil = 2.4, Fig. 7). It is also seen that the curve corresponding to symmetric paraxial theory is quite close to the curve corresponding to x_sil = 0 (see Fig. 6 for clear illustration of this fact). The curve corresponding to circularly polarized illumination, as expected, lies in between the curves corresponding to x_sil = 0 and y_sil = 0 of x- polarized illumination, as shown in Figs. 6 and 7, since it is a superposition of x- and y- polarized illuminations. Further, it is observed that the bigger pinhole radius results in more pronounced sidelobes as seen in Figs. 6(b) and 7(b).

Fig. 8 shows the resolution (FWHM) of imaging a small scatterer as a function of radius of pinhole, R_PH, for different values of NA_sil using x-polarized and circularly polarized illumination, respectively. It is seen that the resolution of ASIL-SM gradually deteriorates with the increase of pinhole size and then remains at almost a constant value when the pinhole is big enough to capture most of the signal produced by the small scatterers. Although the results shown here correspond to substantially large NA, such observations were reported for the paraxial microscope also (non-SIL type) in [19, 34].

Fig. 8 Resolution of imaging a small scatterer as a function of radius of pinhole, R_PH, for different NA_sil using x-polarized illumination ((a) y_sil = 0 and (b) x_sil = 0) and (c) circularly polarized illumination, respectively.

Download Full Size | PDF

A quantitative comparison of the resolution of ASIL confocal microscope (RCM - FWHM for R_PH = 1μm), ASIL wide field microscope (RWFM - FWHM for R_PH = 100μm), and the resolution predicted by subsystem 1 (FWHM of focal spot) is provided in Table 3. For x-polarized illumination, it is found that the ratio of RCM and RWFM along y-axis is close to 1.4, as shown in Table 3 (red color numbers). This is similar to the observation made in [14,34] for a paraxial non-SIL microscope that the theoretically expected improved resolution for confocal microscope is about 1.4 times better than wide-field microscope for a point object. Further, the RWFM agrees with the prediction obtained using the FWHM of the focal spot for x- polarization along the y-axis (x_sil = 0). However, these observations do not hold for the x-polarized wave along the x-axis (y_sil = 0), as well as for the circularly polarized illumination.

Table 3. Comparison of RWMF, RCM and FWHM of subsystem 1.(Unit:λ)

View Table | View all tables in this article

3.2. Imaging of large object structures

The proposed model is capable of simulating and analyzing the images of large object structures. Several examples are presented in this section to illustrate this: an annular ring pattern (Fig. 9(a)), a ‘08’ digital pattern (Fig. 9(b)), a USAF resolution target pattern (Fig. 9(c)) and a pattern with materials similar to integrated circuits (called IC pattern for simplicity) (Fig. 9(d)). In the first three examples, the refractive index of the structures is 1.5, which corresponds to silica (silicon dioxide) in silicon substrate. The last example closely mimics the patterns and materials used in integrated circuits. It considers features with different materials, viz., silicon nitride (n = 2.0), silica (n = 1.5), cobalt silicide (n = 1.3), and gold (n = 0.41 + i9.11), all refractive indices at 1340 nm [35].

Fig. 9 Refractive index (magnitude) distribution of the cross section of the object structures, (a) an annular ring, (b) a ‘08’ digital pattern, (c) a designed three-bar resolving power test target and (d) the IC pattern.

Download Full Size | PDF

The annular ring pattern and the ‘08’ digital pattern are discussed first. The inner and outer radii of the annular ring are R_i = 0.12λ and R_o = 0.24λ, respectively. The ‘08’ digital pattern, is sized such that each line has a width of t = 0.1λ and the distance between any two adjacent parallel lines is d = 0.2λ, as shown in Fig. 9(b). Fig. 10 shows the simulated images for the annular ring and the ‘08’ digital pattern in the top and bottom rows respectively.

Fig. 10 Image of the annular ring (the top row) and ‘08’ digital pattern (the bottom row) using proposed model with NA_sil = 3.3 and different pinhole radius, R_PH = 1μm (a,d,g and j), R_PH = 25μm (b,e,h and k), and R_PH = 100μm (c,f,i and l). The first three columns are for x- polarized incidence while the last three columns for the circular polarized incidence. The horizontal and vertical coordinates are [−0.7, 0.7] and [−0.55, 0.55] for x_sil(λ) and y_sil(λ) in the bottom row, respectively. Both of coordinates are [−0.3, 0.3] in the top row.

Download Full Size | PDF

For x-polarized illumination and R_PH = 1μm, 25μm, we can see only two symmetric spots away from the optical axis in the x-direction in Figs. 10(a) and 10(b) for the annular ring pattern and only the vertical lines in Figs. 10(g) and 10(h) for the ‘08’ digital pattern. For x-polarized illumination and R_PH = 100μm in Figs. 10(c) and 10(i), the image quality of both patterns deteriorates even further. Thus, it can be concluded that x- polarized illumination (and likewise other linear polarizations) is not suitable for imaging such structures, irrespective of the size of pinhole. Comparatively, circular polarized illumination results in better images. It is observed that small pinhole radius results in reasonable imaging (see Figs. 10(d),10(e), 10(j) and 10(k)). However, there are some artifacts in the form of a small dot in the centre (annular pattern) or in the form of the lines in the middle (’08’ pattern). For pinhole radius R_PH = 25μm, we see a bleaching effect for both the patterns (Figs. 10(e) and 10(k)). The images are of poor quality for high values of pinhole radius (Figs. 10(f) and 10(l)).

The next example corresponds to resolution targets similar to the USAF resolution test chart. The target consists of nine groups of pattern decreasing in size from group 1 to 9, as shown in Fig. 9(c). Each group consists of two similar patterns (two sets of lines) at right angle to each other. Each pattern consists of three lines separated by spaces of the same size as the line width. The line length is five times the line width. The line width is uniform for a group and changes from t₁ = 0.2λ to t₉ = 0.04λ at a step 0.02λ. In addition to the aforementioned nine groups, four squares corresponding to the pattern size of the smallest four groups are also included.

Fig. 11 shows the simulated images for the x-polarized and circular polarized illuminations. It is found that we can resolve the target element of line width 0.10 λ for circularly polarized illumination (Fig. 11(e)), and 0.12 λ for x-polarized illumination (Fig. 11(b)), using a pinhole radius of 25μm. On the other hand, for the pinhole radius of 1μm, the images are not good for either circular or x-polarized illumination. For the circular polarized illumination, both the horizontal and vertical bar patterns have artifacts. For the x-polarized illumination, the vertical bar patterns have artifacts while the horizontal bar patterns have very low intensity.

Fig. 11 The caption is the same as Fig. 10 except that the three-bars pattern is imaged. The horizontal and vertical coordinates are [−2.55, 2.55] and [−2.10, 2.10] for x_sil(λ) and y_sil(λ), respectively.

Download Full Size | PDF

In the last example, the six white stubs in Fig. 9(d) are of gold. The dimension of each stub is 0.06λ × 0.02λ. The other structures are silica, silicon nitride, and cobalt silicide from outer to inner as shown in Fig. 9(d). The thickness of each line is 0.06λ and the distance between the lines is 0.20λ. Fig. 12 shows the simulated images of the IC pattern for the x-polarized and circular polarized illuminations. As expected, irrespective of the size of pinhole, the gold stubs have very strong intensity (Figs. 12(a)–12(c) and Figs. 12(g)–12(i)) due to the large refractive index (magnitude) of gold. However, the horizontal gold stubs are most clearly imaged structures for x- polarized incidence as shown in Figs. 12(a)–12(c). This is different from dielectric structures where vertical features were more prominent in Figs. 10 and 11. We note that while other features are difficult to visualize in the images in (Figs. 12(a)–12(c) and Figs. 12(g)–12(i)), they can be easily visualized by using a log scale plot of the intensity, as shown in (Figs. 12(d)–12(f) and Figs. 12(j)–12(l)). In the log scale, we note that all the features are visible. The artifacts similar to those noted in the previous examples are also clearly visible. Further, we see that circular polarization with pinhole radius of 15μm gives the best image.

Fig. 12 The caption is the same as Fig. 10 except that the IC pattern is imaged with the pinhole size R_PH = 1μm, 15μm and 25μm. Both the horizontal and vertical coordinates are [−1.60, 1.60] for x_sil(λ) and y_sil(λ). The last three columns create images using a base 10 logarithmic scale for intensity.

Download Full Size | PDF

This example raises a concern about the visibility and fidelity of the features in the microscope image in presence of noise, especially for the case of small pinhole like 1μm, where the measured intensity is very low and the signal-to-noise ratio (SNR) is very poor. In order to study the effect of noise on the image quality (when viewed in logarithmic scale), we add successively increasing noise to the detected intensity (circular polarization, 1μm pinhole radius) and plot the logarithmic image in Fig. 13. We see that even in the presence of Gaussian noise with SNR = 20dB, as shown in Fig. 13, all the features are clearly visible.

Fig. 13 The logarithmic image of the IC pattern in the presence of Gaussian noise with (a) no noise, (b) SNR = 40dB, (c) SNR = 30dB and (d) SNR = 20dB, using circular polarization and 1μm pinhole radius.

Download Full Size | PDF

4. Conclusion

This paper presents a complete and computationally efficient model of ASIL-SM and provides the ability to analyze image formation in a more systematic way than can be performed experimentally. For example, imaging different object structures would each require their own fabrication, whereas the presented model can predict results by simply changing a few parameters and then if needed, some structures can be checked by experiment [36].

The computational efficiency is derived by using fast algorithms like CG-FFT and CZT-FC for each subsystem in the computational model. This is an important contribution to the research and development in ASIL based microscopy technology since modeling and simulating ASIL microscopy system is very challenging due to spherical refraction interface, highly oscillatory diffraction integrals, and pressing requirement of reducing the scanning resolution. It is found that the presented model significantly reduces the computational time, especially for the large object structures from few hours to a few seconds. We highlight that the model has been experimentally tested and simulated results match the experimental results very well. Comparison of the simulation model with experimental results was reported in [36].

Further, engineering parameters of ASIL-SM, like the polarization, NA, and the detector pinhole size have been identified. Their effects on resolution and imaging quality have been studied. We hope that this study provides a preliminary analysis of suitable imaging setup for ASIL-SM and incites ASIL-SM system engineering using these parameters. We have already begun working on designing more complicated polarizations and pupil filter designs for achieving better resolution.

Acknowledgments

This work was supported by the Singapore Ministry of Education (MOE) grant under Project No. MOE2009-T2-2-086.

References and links

1. S. M. Mansfield and G. S. Kino, “Solid immersion microscope,” Appl. Phys. Lett. 57, 2615–2616 (1990) [CrossRef] .

2. Q. Wu, L. P. Ghislain, and V. B. Elings, “Imaging with solid immersion lenses, spatial resolution, and applications,” Proc. IEEE 88, 1491–1498 (2000) [CrossRef] .

3. S. B. Ippolito, B. B. Goldberg, and M. S. Ünlü, “High spatial resolution subsurface microscopy,” Appl. Phys. Lett. 78, 4071–4073 (2001) [CrossRef] .

4. F. H. Köklü, J. I. Quesnel, A. N. Vamivakas, S. B. Ippolito, B. B. Goldberg, and M. S. Ünlü, “Widefield subsurface microscopy of integrated circuits,” Opt. Express 16, 9501–9506 (2008) [CrossRef] [PubMed] .

5. Semicaps, “Optical fault localization system,” http://www.semicaps.com/innovations.htm (2011).

6. L. P. Ghislain and V. B. Elings, “Near-field scanning solid immersion microscope,” Appl. Phys. Lett. 72, 2779–2781 (1998) [CrossRef] .

7. D. A. Fletcher, K. B. Crozier, C. F. Quate, G. S. Kino, K. E. Goodson, D. Simanovskii, and D. V. Palanker, “Near-field infrared imaging with a microfabricated solid immersion lens,” Appl. Phys. Lett. 77, 2109–2111 (2000) [CrossRef] .

8. A. N. Vamivakas, R. D. Younger, B. B. Goldberg, A. K. Swan, M. S. Ünlü, E. R. Behringer, and S. B. Ippolito, “A case study for optics: The solid immersion microscope,” Am. J. Phys. 76, 758–768 (2008) [CrossRef] .

9. F. H. Köklü and M. S. Ünlü, “Subsurface microscopy of interconnect layers of an integrated circuit,” Opt. Lett. 35, 184–186 (2010) [CrossRef] [PubMed] .

10. S. H. Goh and C. J. R. Sheppard, “High aperture focusing through a spherical interface: Application to refractive solid immersion lens (RSIL) for subsurface imaging,” Opt. Commun. 282, 1036–1041 (2009) [CrossRef] .

11. K. M. Lim, G. C. F. Lee, C. J. R. Sheppard, J. C. H. Phang, C. L. Wong, and X. Chen, “Effect of polarization on a solid immersion lens of arbitrary thickness,” J. Opt. Soc. Am. A 28, 903–911 (2011) [CrossRef] .

12. T. X. Hoang, X. Chen, and C. J. R. Sheppard, “Multipole theory for tight focusing of polarized light, including radially polarized and other special cases,” J. Opt. Soc. Am. A 29, 32–43 (2012) [CrossRef] .

13. R. Chen, K. Agarwal, Y. Zhong, C. J. R. Sheppard, J. C. H. Phang, and X. Chen, “Complete modeling of subsurface microscopy system based on aplanatic solid immersion lens,” J. Opt. Soc. Am. A 29, 2350–2359 (2012) [CrossRef] .

14. C. J. R. Sheppard and A. Choudhury, “Image formation in scanning microscope,” Opt. Acta 24, 1051–1073 (1977) [CrossRef] .

15. C. J. R. Sheppard and T. Wilson, “Image formation in scanning microscopes with partially coherent source and detector,” Opt. Acta 25, 315–325 (1978) [CrossRef] .

16. C. J. R. Sheppard and T. Wilson, “The theory of the direct-view confocal microscope,” J. Microsc.-Oxf. 124, 107–117 (1981) [CrossRef] .

17. G. J. Brakenhoff, P. Blom, and P. Barends, “Confocal scanning light-microscopy with high aperture immersion lenses,” J. Microsc.-Oxf. 117, 219–232 (1979) [CrossRef] .

18. T. Wilson, R. Juskaitis, and P. Higdon, “The imaging of dielectric point scatterers in conventional and confocal polarisation microscopes,” Opt. Commun. 141, 298–313 (1997) [CrossRef] .

19. P. Török, P. D. Higdon, and T. Wilson, “Theory for confocal and conventional microscopes imaging small dielectric scatterers,” J. Mod. Opt. 45, 1681–1698 (1998) [CrossRef] .

20. P. Török, P. R. T. Munro, and E. E. Kriezis, “High numerical aperture vectorial imaging in coherent optical microscopes,” Opt. Express 16, 507–523 (2008) [CrossRef] [PubMed] .

21. R. L. Coe and E. J. Seibel, “Computational modeling of optical projection tomographic microscopy using the finite difference time domain method,” J. Opt. Soc. Am. A 29, 2696–2707 (2012) [CrossRef] .

22. M. Leutenegger, R. Rao, R. A. Leitgeb, and T. Lasser, “Fast focus field calculations,” Opt. Express 14, 11277–11291 (2006) [CrossRef] [PubMed] .

23. L. Hu, R. Chen, K. Agarwal, C. J. R. Sheppard, J. C. H. Phang, and X. Chen, “Dyadic Green’s function for aplanatic solid immersion lens based sub-surface microscopy,” Opt. Express 19, 19280–19295 (2011) [CrossRef] [PubMed] .

24. Y. Zhong and X. Chen, “An FFT twofold subspace-based optimization method for solving electromagnetic inverse scattering problems,” IEEE Trans. Antennas Propag. 59, 914–927 (2011) [CrossRef] .

25. R. Chen, K. Agarwal, C. J. R. Sheppard, J. C. H. Phang, and X. Chen, “Resolution of aplanatic solid immersion lens based microscopy,” J. Opt. Soc. Am. A 29, 1059–1070 (2012) [CrossRef] .

26. C. W. McCutchen, “Generalized aperture and the three-dimensional diffraction image,”J. Opt. Soc. Am. 54, 240–242 (1964) [CrossRef] .

27. J. Lin, X. C. Yuan, S. S. Kou, C. J. R. Sheppard, O. G. Rodriguez-Herrera, and J. C. Dainty, “Direct calculation of a three-dimensional diffracted field,” Opt. Lett. 36, 1341–1343 (2011) [CrossRef] [PubMed] .

28. C. J. R. Sheppard and K. G. Larkin, “Vectorial pupil functions and vectorial transfer functions,” Optik 107, 79–87 (1997).

29. J. Lin, O. G. Rodriguez-Herrera, F. Kenny, D. Lara, and J. C. Dainty, “Fast vectorial calculation of the volumetric focused field distribution by using a three-dimensional fourier transform,” Opt. Express 20, 1060–1069 (2012) [CrossRef] [PubMed] .

30. L. Rabiner, R. Schafer, and C. Rader, “The chirp z-transform algorithm,” IEEE Trans. Acoust. Speech 17, 86–92 (1969).

31. A. Oppenheim, R. Schafer, and J. Buck, Discrete-Time Signal Processing, 2nd ed (Prentice Hall, 1999).

32. http://www.tedpella.com/metrochip_html/metrochip-calibration-target.htm.

33. http://www.2spi.com/catalog/magnifiers/magnification-standard-geller-MRS-5.php.

34. T. Wilson and C. Sheppard, Theory and Practice of Scanning Optical Microscopy, vol. 1 (London: Academic Press, 1984).

35. P. B. Johnson and R. W. Christy, “Optical constants of the noble metals,” Phys. Rev. B 6, 4370–4379 (1972) [CrossRef] .

36. K. Agarwal, R. Chen, L. S. Koh, C. J. R. Sheppard, J. C. H. Phang, and X. Chen, “Experimental validation of the computational model of aplanatic solid immersion lens scanning microscope,” presented at Focus on microscopy 2013, Maastricht, The Netherlands, 24–27 Mar. 2013

	M = N = 50	M = N = 100	M = N = 200	M = N = 400	M = N = 800
DI - time (s)	36	36	36	36	36
FFT- time (s) Relative error	0.40 1.8 × 10⁻²	1.62 6.9 × 10⁻³	6.38 2.2 × 10⁻³	26.72 8.9 × 10⁻⁴	107.7 2.7 × 10⁻⁴
CZT- time (s) Relative error	0.008 1.4 × 10⁻²	0.018 6.0 × 10⁻³	0.068 2.2 × 10⁻³	0.226 8.3 × 10⁻⁴	0.906 2.6 × 10⁻⁴

System	Small structure (One small sphere)		Large structure (Annular ring)
System	DI	CZT	DI	CZT
Subsystem I Subsystem III	8.36 s 20.0 s	0.04 s 0.09 s	8.36 s around 2 h	0.05 s 0.10 s
Subsystem II (CG-FFT)	0.28 s		0.91 s
Total	28.64 s	0.41 s	around 2 h	1.06 s

	x-polarization (y_sil = 0)				x-polarization (x_sil = 0)				Circular polarization

	RWFM	RCM	FWHM	$\frac{R W F M}{R C M}$	RWFM	RCM	FWHM	$\frac{R W F M}{R C M}$	RWFM	RCM	FWHM	$\frac{R W F M}{R C M}$
NA_sil = 3.3	0.210	0.093	0.219	2.26	0.152	0.109	0.152	1.39	0.177	0.100	0.181	1.77
NA_sil = 3.0	0.196	0.108	0.214	1.81	0.168	0.120	0.168	1.40	0.181	0.114	0.189	1.59
NA_sil = 2.7	0.204	0.125	0.222	1.63	0.188	0.135	0.188	1.39	0.196	0.130	0.204	1.51
NA_sil = 2.4	0.222	0.144	0.239	1.59	0.212	0.152	0.212	1.39	0.216	0.148	0.225	1.46

	M = N = 50	M = N = 100	M = N = 200	M = N = 400	M = N = 800
DI - time (s)	36	36	36	36	36
FFT- time (s) Relative error	0.40 1.8 × 10⁻²	1.62 6.9 × 10⁻³	6.38 2.2 × 10⁻³	26.72 8.9 × 10⁻⁴	107.7 2.7 × 10⁻⁴
CZT- time (s) Relative error	0.008 1.4 × 10⁻²	0.018 6.0 × 10⁻³	0.068 2.2 × 10⁻³	0.226 8.3 × 10⁻⁴	0.906 2.6 × 10⁻⁴

System	Small structure (One small sphere)		Large structure (Annular ring)
System	DI	CZT	DI	CZT
Subsystem I Subsystem III	8.36 s 20.0 s	0.04 s 0.09 s	8.36 s around 2 h	0.05 s 0.10 s
Subsystem II (CG-FFT)	0.28 s		0.91 s
Total	28.64 s	0.41 s	around 2 h	1.06 s

A complete and computationally efficient numerical model of aplanatic solid immersion lens scanning microscope

Abstract

1. Introduction

2. Computational model of ASIL scanning microscope with finite sized detector

2.1. Background

2.2. Focusing of incident light through ASIL

2.3. Fast computation of DGF using CZT

2.4. Image formation using convolution theorem

2.5. Scanning system implementation

2.6. Example of enhancement in computational efficiency

3. Numerical simulations and analysis

3.1. Imaging of small scatterers

3.2. Imaging of large object structures

4. Conclusion

Acknowledgments

References and links

Cited By

Figures (13)

Tables (3)

Equations (21)

Optics Express