## Abstract

Abstract

In this paper, a simple confocal laser scanning microscopic (CLSM) image mapping technique based on the finite-difference time domain (FDTD) calculation has been proposed and evaluated for characterization of a subwavelength-scale three-dimensional (3D) void structure fabricated inside polymer matrix. The FDTD simulation method adopts a focused Gaussian beam incident wave, Berenger’s perfectly matched layer absorbing boundary condition, and the angular spectrum analysis method. Through the well matched simulation and experimental results of the *xz*-scanned 3D void structure, we first characterize the exact position and the topological shape factor of the subwavelength-scale void structure, which was fabricated by a tightly focused ultrashort pulse laser. The proposed CLSM image mapping technique based on the FDTD can be widely applied from the 3D near-field microscopic imaging, optical trapping, and evanescent wave phenomenon to the state-of-the-art bio- and nano-photonics area.

© 2007 Optical Society of America

## 1. Introduction

Ultrashort pulsed lasers have opened up a new possibility for laser nanofabrication in three-dimensions (3D) that would impact hugely on applications such as high-density 3D optical data storage, micro-fluidic devices and laser nanosurgery of cells and tissues [1–11]. The physical processes of the pulsed laser nanofabrication are rich in nature and often not limited to a single process, but in the fabrication in a dielectric transparent material it involves melting, evaporation, ionizations and/or optical dielectric breakdown, which cause a refractive index change, void formation, explosion and multiples of those effects.

In particular, void formation has been under intense investigation for its importance in optical device applications, and recent efforts to reduce the void size has pushed its limit from the micron-size regime to the sub-micron or nano-size regime (i.e., nanovoids) [11]. As the size regime gets smaller, the characterization of the voids becomes more and more challenge. Currently, there are two major techniques used for the void characterization, scanning electron microscopy (SEM) and confocal laser scanning microscopy (CLSM) [12–14]. While SEM can provide much better spatial resolution, it requires a complicated sample preparation such as milling and cutting to reveal the voids inside the sample. The CLSM, on the other hand, provides fast and easy three-dimensional optical sectioning of the voids without much effort in sample preparation. However, because the size regime of the nanovoids are in the same order as the laser focal spot, the observed images are complicated convolution of the laser focal spot point spread function and the void structure, which often leaves arguable interpretations.

In order to extract information about the void shape and size, a detailed analysis of the observed image is required. Until now, only rough estimates of void size and shape are extracted, and no detailed analysis has been conducted. In this paper, we numerically simulate the image formation of the nanovoids (size regime 100~500 nm) in a 2D CLSM configuration, and provide an insight into the interpretation and correlation of the image to the actual void size and shape.

In order to simulate the beam focusing, propagation and detection, we employed the finite difference time domain (FDTD) technique [15]. The FDTD method has been widely used for analyzing the EM phenomena in all the ranges from microwave to x-ray optics, including near-field optics [16, 17], photonic crystal [18], lasing [19], optical tweezer [20], scattering from cells [21], optical storage [22], and diffractive optics [23]. While there are several other numerical methods for electromagnetic (EM) fields in an arbitrary geometry [24, 25], the FDTD method could provide an exact time-domain vectorial solution to the Maxwell’s equations. Here, we provide a method to apply the FDTD technique to the transmission- and reflection-type confocal detection geometry with a pinhole.

The paper is outlined as follows. In Section 2, the FDTD simulation methods for CLSM image formation of a nanovoid shall be given. From Sections 3 to 5, we shall present the CLSM images with variation in detection numerical aperture (NA), size, shape and index change. In Section 6, an experimental validation of the simulation will be presented, followed by a conclusion in Section 7.

## 2. FDTD calculation for confocal laser scanned microscopic image formation

The FDTD method, as one of the most powerful computational methods in EM problems, has been widely used to model the wave propagation, scattering, and radiation since it was first introduced by Yee [26] in 1966. The EM problems can be solved by the following Maxwell’s differential equations for a discretized arbitrary 3D geometrical space.

$$+\left(\frac{\Delta t}{\Delta u\epsilon (i,j,k+1\u20442)}\right)\left[\begin{array}{c}{H}_{y}^{n+1\u20442}(i+1\u20442,j,k+1\u20442)-{H}_{y}^{n+1\u20442}\left(i-1\u20442,j,k+1\u20442\right)\\ +{H}_{x}^{n+1\u20442}(i,j-1\u20442,k+1\u20442)-{H}_{x}^{n+1\u20442}(i,j+1\u20442,k+1\u20442)\end{array}\right]$$

$${H}_{z}^{n+1\u20442}\left(i+1\u20442,j+1\u20442,k\right)=\left(1-\frac{\Delta t\sigma (i+1\u20442,j+1\u20442,k)}{\epsilon (i+1\u20442,j+1\u20442,k)}\right){H}_{z}^{n-1\u20442}(i+1\u20442,j+1\u20442,k)$$

$$+\left(\frac{\Delta t}{\Delta u\epsilon (i+1\u20442,j+1\u20442,k)}\right)\left[\begin{array}{c}{E}_{x}^{n}(i+1\u20442,j+1,k)-{E}_{x}^{n}(i+1\u20442,j,k)\\ +{E}_{y}^{n}(i,j+1\u20442,k)-{E}_{y}^{n}(i+1,j+1\u20442,k)\end{array}\right]$$

where *µ* is the magnetic permeability, *σ* is the electric conductivity, *ε* is the dielectric permittivity, *Δt* is the temporal step size, and *Δu* is the spatial step size. Figure 1 shows the schematic diagram for a CLSM system and our simplified FDTD calculation structure. As shown in Fig. 1, we analyze two different modes, transmission and reflection-type CLSM.

Due to a limitation in computation power available, a 3D nanovoid structure is estimated with a 2D structure. While it is true that the linear polarization in 3D case is more complicated than mere TE and TM modes in the 2D case, it nevertheless can provide us with a qualitative image forming behavior in a confocal geometry. The rotational symmetry of the void and focal spot also helps the reduction in dimensions.

The calculation space is divided into small units called Yee cells. Complex dielectric constants are assigned to each cell for given materials. The FDTD method calculates the EM fields in each cell by integrating the Maxwell’s equations in a “leap-frog” procedure until the steady state is reached. In the case of a sinusoidal plane wave excitation source, the steady state is reached when all scattered fields vary sinusoidally in time. In our case, the tightly focused laser excitation source is applied by a Gaussian beam with a ramped delay as given in the following Eq. (3),

$${\tau}_{d}=\frac{1}{2}\times \left[1+\mathrm{erf}\left\{\frac{\left(m-20\right)}{5\sqrt{2}}\right\}\right]$$

where *f* is the focal length and *D _{x}* is the aperture diameter of the objective lens,

*w*is the full width half maximum (FWHM) value of the incident Gaussian beam,

_{0}*k*(=2

*π/λ*) is the wave number,

*m*is the iteration number, and

*erf*is the error function. In order to obtain an accurate field distribution, the cell size must be much less than the size of relevant features because the FDTD algorithm is stable only if

*Δt*≤

*Δu/(c√2)*, where c is the speed of light. In our case we used the cell size

*Δu*of 20 nm and the time step

*Δt*of 3.3349×10

^{-17}s, which sufficiently satisfies the Courant stability condition.

The 2D Maxwell equations decouple into two polarization sets: the TM set, *E _{z}H_{x}H_{y}*, and the TE set

*H*. Each field set can be modeled independently of the other. For absorbing boundary condition (ABC) in order to truncate the FDTD computational domain on homogeneous media we employed Berenger’s highly effective ABC, the perfectly matched layer (PML), which gives zero reflection at the absorbing boundary for all frequencies and all angles of incidence [27]. For the TM case, if

_{z}E_{x}E_{y}*E*is known, then

_{z}*H*and

_{x}*H*can be calculated at the cross-sectional plane of the FDTD grid using Eq. (1). Using wave propagation techniques, all the field components can be determined at a plane of interest well past the FDTD computational grid.

_{y}For obtaining the transmission- and reflection-type CLSM image we first performed each FDTD simulation with respect to all of the predefined inside positions of the *xz*-scanning window with the minimum distance variation of 0.1 µm. To analyze the frequency (amplitude and phase) response of the FDTD simulation results we can simply take the Fourier transform of the complex electric field for the time-domain data at every point of interest with satisfying the optical sectioning (axial) condition of a confocal imaging system. This can be done by the equation

where *f _{m}* and

*x*are the frequency and detector positions satisfying the lateral condition of the confocal imaging system and

*t*and

_{T}*t*are the final time step number and the period of the optical sectioning distance determined by the NA of the collection lens in the FDTD iteration, respectively. This Eq. (4) can be divided into its real and imaginary parts. The next step is to determine the scattered field distribution at the image plane after going through the pinhole position. In other words, we need to forward or backward propagate of the complex electric field at the detection plane to the confocal imaging plane using the exact solution of the Helmholtz equation for a single scalar field component. For obtaining the full-wave vector description of the subwavelength-scale light propagation in transmission- and reflection-type CLSM, we used an angular spectrum analysis method [28, 29]. The exact solution of the Helmholtz equation in homogeneous media may be expressed in the form

_{PH}$$\beta =\{\begin{array}{cc}\sqrt{{({n}_{m}\u2044\lambda )}^{2}-\left({f}_{x}^{2}\right)}& \mathrm{when}\left({f}_{x}^{2}\right)\le {({n}_{m}\u2044\lambda )}^{2}\\ i\sqrt{\left({f}_{x}^{2}\right)-{({n}_{m}\u2044\lambda )}^{2}}& \mathrm{otherwise}\end{array}$$

where *n _{m}* is the refractive index of the homogeneous media and

*β*can be either real or imaginary. If we assume

*β*is real and the reflection

*R(f*is zero, the solution (5) represents the full vector waves that propagate in forward direction

_{x})*k*=2

*π(fx, β)*. When

*β*is imaginary, the waves either decay or increase exponentially. It follows from the radiation condition that the field must assume the form of an outgoing spherical wave in the far-zone. If

*U(x, z*is known,

_{0})*U(x, z)*can be determined by the following procedure: Firstly, take the Fourier transform of

*U(x, z*to determine

_{0})*T(f*. Then, multiply

_{x}, z_{0})*T(f*with the propagation kernel

_{x}, z_{0})*exp*[

*i2π{f*}] to get

_{x}x+β(z-z_{0})*T(f*. Finally, the observation field

_{x}, z)*U(x, z)*can be determined by taking the inverse Fourier transform of

*T(fx, z)*. To get the power at the photomultiplier tube (PMT) receiver, we take summing up the absolute square value of each point inside the pinhole area. The reflection CLSM image can be acquired by the same manner placing the transmission

*T (f*equals zero.

_{x})FDTD simulation is performed using a personal computer with a Pentium4 3.6 GHz CPU machine. Our MatLab^{™}-based FDTD program has the FDTD cell grid sizes of 960×400 ≅ 19.2 µm×8 µm for NA1.4 (θ_{inc}=67.84°) and 300×400≅6 µm×8 µm for NA0.9 (θ_{inc}=36.87°), the PML depth of 12 layers, and the total time slots of 2640 Δt. The amplitude of the focused Gaussian incident wave *E _{z}* (or H

_{z}) is assumed as 1 V/m (A/m).

The angular spectrum method is a very simple, fast, and elegant rigorous model since it uses the plane wave decomposition method. Using the proposed FDTD-based CLSM image mapping technique we could get the exact matched transmission- and reflection-type CLSM image between the FDTD simulation results and experimentally captured results. Also, we could extract the exact position and topological shape factors of the subwavelength-scale void structure, which was fabricated by femtosecond laser microexplosion under the tightly focused Gaussian beam illumination.

## 3. Effect of numerical aperture on CLSM void imaging

Figure 2 shows FDTD simulated snapshots of 2D TE/TM focused wave propagation in a void structure with different NAs of a single subwavelength-scale void structure inside the predefined xz-scanning window under the focused He-Ne laser (λ=632.8 nm) illumination. For the case of wave propagation in various other positions of the void, animations were prepared where the position of the void structure is scanned relative to the focal spot position - See Supplimentary Animation A (NA 0.9, TM) and B (NA 1.4, TM) attached. In these simulations, we assumed that the refractive indices of the void and polymer matrix are 1.0 and 1.5, respectively. The electric and magnetic field distributions after going through a single void structure with the size of 0.56 µm×1.38 µm (aspect ratio=2.4643) under the different NA values of 0.9 (air) and 1.4 (oil-immersion).

In order to build up a confocal scanned image, the void structure was scanned in the focal field as in the animations, and then the transmitted or reflected wave was collected by the same NA lens and a pinhole according to the method specified in Section 2 above. Figure 3 shows the transmission- and reflection-type *x-z* CLSM images with respect to different NAs using our proposed FDTD-based CLSM image mapping techniques. In each simulation the minimum pixel resolution is 0.1 µm and the scanning window size is 3.1 µm×5.1 µm. As shown in Figs. 3(a) and 3(b) one cannot locate the exact void structure and position in the case of low NA objective lens. This is exacerbated with the aberration caused by index mismatch between the objective lens immersion medium (air) and the polymer matrix.

However, one can extract the position and topological shape factors of the voids from the calculated reflection- and transmission-type CLSM images using the high NA of 1.4 (oil-immersion). In particular, in transmission images [Fig. 3(d)] the axial length of the void structure can be obtained by measuring the length between the maximum bright point and minimum dark point of the *z*-scanned cross sectional profile. Such bright and dark points arise due to the reflection and focusing characteristics of the elliptical void shape at two different focal points. For example, the bright point in the transmission image arises due to the lensing effect of the bottom interface of the void (more transmission), and the dark point arises due to the reflection at the top of the void interface (less transmission). The conjugated reflection CLSM images show the opposite information — a bright peak at the tip of the void, demonstrating that the reflection signal is maxima when the focus is at the interface between the void and the matrix. This maximum point matches well the dark point of the corresponding transmission image [red lines in the Fig. 3(d)].

From this result, one can conclude that the axial size of the void can be accurately extracted by observing and measuring the bright and dark points of the transmission images, or bright points in the transmission-reflection conjugated images. In reality, the bright point of the reflection image is much sharper than the dark point of the transmission image and therefore the latter method would yield a better measure of the axial size.

Furthermore, the transversal (lateral) size can be acquired by measuring full-width half-maximum points of the dark region at the transmission image. From Fig. 3(d) we can obtain the shape factor of ~2.5 with the lateral and axial size of ~0.56 µm and 1.4 µm, respectively. These values match well with the pre-designed simulation input parameters of 2.4643 µm, 0.56 µm, and 1.38 µm, respectively.

## 4. Effect of refractive index difference on CLSM imaging

The laser induced subwavelength-scale 3D structure fabrication involves many physical changes of matter inside the focused volume, starting from a small change of the refractive index induced by phase change (melting, solidification), to a gas or plasma production leading to a void formation. In order to see the difference in the various index states of fabricated structure in CLSM, the FDTD simulation was conducted for index variations and the resulting images are shown in Fig. 4(a)–4(f). In Figs. 4(g) and 4(h) show the contrast ratio (defined as (*I _{max} - I_{min}*)/(

*I*+

_{max}*I*), where

_{min}*I*is the maximum value of the detected power and

_{max}*I*is the background level) between the bright and dark points (defined in the preceding section) of the transmission- and reflection-type CLSM with respect to the refractive index variations.

_{min}From the Fig. 4(g) one can notice that the contrast of the transmission-type confocal images in the TM mode do not show much differences along the refractive index variations, while that in the TE mode increases as enlarging the refractive index differences. However, the shape factor and exact position information can be extracted through TM transmission mode CLSM imaging. In addition, we could measure out the real-depth scale of the microexploded nanoscale void structure from each brightest point on the z-axis’ cross-sectional profile of the combined transmission- and reflection-type CLSM image up to the 80×200 nm size limitation. As shown in the Fig. 4(h) the reflection-type confocal mode, the image contrasts are mostly below 0.1, which may not be sufficient enough to provide information the shape factors of the fabricated structure. The detailed contrast values are shown in Table 1.

## 5. Effect of size and shape factors on CLSM void imaging

Figure 5 shows the simulated 2D TM FDTD-based transmission- and reflection-type CLSM images of the subwavelength-scale 3D structure with different size ratios, shapes, and multilayered (void shelled structure with 20% and 40% high density solidification) refractive index variations in the PVA polymer matrix. The 20–40% densification shells were considered to accommodate the possibility of variation in densification of polymer, although the recent work of Gamaly et al. [30, 31] reports the typical densification observed was up to 20%. Figures 5(a) and 5(b) show the simulated confocal transmission- and reflection-type CLSM images of a small sphere with the radius of 0.28 µm and a large sphere with the radius of 0.56 µm.

Also the surrounding region of the microexploded void structure has increased density (refractive index) than that of the original matrix following the mass (energy) conservation law [32, 33]. For simulating the confocal images including with more exact physical phenomena of our microexploded void structure, we imposed the high densitity refractive index of 1.55 with the expanded ratio of 1.2 and 1.4. Figures 5(c)–5(f) show the simulated 2D TM FDTD-based transmission- and reflection-type CLSM image with respect to different size (0.28 µm×0.79 µm, 0.4 µm×0.96 µm, 0.56 µm×1.38 µm, and 0.72 µm×1.725 µm) variations of the void structure. From the shadow region of the simulated transmission confocal images in the Figs. 5(c)–5(f) we can notice that the exact void position and topological shape factors can be extracted.

## 6. Experimental validation of the CLSM void imaging

Figure 6 shows our laser microexploded structure inside of the PVA polymer and ablated structure on the PVA polymer with the thickness of ~30 µm under the tightly focused (NA=1.4) femtosecond (780 nm, 100 fs Ti:sappire pulse, 80 MHz) laser illumination without any amplifications.

The fabricated pattern has the subwavelength-scale 3D void structure with the average dimension of ~0.56 µm×1.38 µm (aspect ratio of ~ 2.46). Figure 6(a) shows the input pattern with the total dimension of 40 µm×40 µm and the period of each void of 1.29 µm. The SEM image of the ablated subwavelength-scale patterns, which was recorded at the surface region between the PVA polymer and air interfaces, is also shown in Fig. 6(b). From the ablated subwavelength-scale 3D void patterns we can estimate the lateral dimension changes from 0.4110 µm to 0.5955 µm.

The most disadvantages of SEM and TEM is that they cannot provide any 3D information although they can resolve a structure up to nanometer-scale. However, the CLSM system provides the 3D information of the structure. Figures 6(c)–6(e) are showing the transmission and reflection CLSM images at the recording plane of a 3D character pattern inside the PVA polymer volume and theirs cross-sectional depth-scanned images using commercial transmission-type (Fluoview BX50, Olympus) and home made reflection-type CLSM system. From the black shadow region in Fig. 6(e) we estimate that the lateral size and the axial depth of the microexploded 3D void structure have 0.56 µm×0.56 µm and 1.38 µm, respectively. From Fig. 6(f) we can also figure out the position of the fabricated character pattern and axial dimension inside the PVA polymer matrix. The transmission- and reflection-type CLSM images are acquired with the high NA objective lens (UPlanApo 60X/1.4 oil-immersion).

For detailed comparison between the experimentally observed and calculated CLSM images we measure the *xz*-directional cross-section profiles at the center position of a fabricated single void structure. Figures 7(a) and 7(b) show the *xy*-plane view and the *xz*-scanned transmission-type CLSM images of a microexploded single void structure inside the PVA polmer matrix. The xz-direction cross-sectional profiles at the center position of the simulation and experiment result of a single void structure are shown in Fig. 7(c).

Figure 8 shows the reflection-type CLSM simulation and experimental images of a microexploded single void structure inside the PVA polmer matrix. From the Fig. 7(c) and 8(c) we see that the simulation results are similarly well matched with the experimentally obtained result. The slight mismatch between the simulation and experimental profiles may be resulted from the 2D simulation of 3D problem, absence of spherical aberration of the high NA objective lens and nonlinear laser-plasma interactions at the microexploded subwavelength-scale void structure. From the comparison results, we verified that our proposed FDTD-based transmission- and reflection-type CLSM image mapping technique matches well with the experimental result. And also, we can easily extract the exact position and topological shape factor of the microexploded subwavelength-scale single void structure inside the PVA polymer matrix.

## 7. Conclusions

We have first, for authors’ knowledge, proposed and evaluated a FDTD-based transmission- and reflection-type CLSM image mapping technique for analyzing the subwavelength-scale 3D void structure microexploded inside the PVA polymer matrix. In the FDTD simulation a focused Gaussian beam incident wave condition and Berenger’s PML ABC were imposed with respect to a homogenous refractive index medium and different void structurs. For optical sectioning and confocal imaging, we used a sectionized time integration method and an angular spectrum analysis method. From the simulation and experimental results of the transmission- and reflection-type CLSM images of the microexploded subwavelength-scale single void structure, we characterized the exact position and topological shape factor of the void, which was fabricated by a tightly focused ultrashort pulse laser with the high NA 1.4 objective lens (oil-immersion). We believe that our proposed FDTD-based CLSM image mapping technique can be assisted for analyzing the 3D near-field microscopic imaging, optical trapping, and evanescent wave phenomenon to the state-of-the-art bio- and nano-photonics area.

## Acknowledgments

We are grateful thanks to the Centre for Astrophysics and Supercomputing (CAS) at Swinburne University of Technology for supporting the computing resources. The authors acknowledge the support for the Australian Research Council.

## References and links

**1. **H. Misawa and S. Juodkazis, *3D Laser Microfabrication: Principles and Applications*, (Wiley-VCH Verlag, Weinheim, 2006). [CrossRef]

**2. **D. A. Parthenopoulos and P. M. Rentzepis, “Three-dimensional optical data storage memory,” Science **245**, 843–845 (1989). [CrossRef] [PubMed]

**3. **S. Kawata, H. -B. Sun, T. Tanaka, and K. Takada, “Finer features for functional microdevices,” Nature **412**, 697–698 (2001). [CrossRef] [PubMed]

**4. **B. H. Cumpston, S. P. Ananthavel, S. Barlow, D. L. Dyer, J. E. Ehrlich, L. L. Erskine, A. A. Heikal, S. M. Kuebler, I. -Y. S. Lee, S. M. Maughon, J. Qin, H. Röckel, M. Rumi, X. L. Wu, S. R. Marder, and J. W. Perry, “Two-photon polymerization initiators for three-dimensional optical data storage and microfabrication,” Nature **398**, 51 (1999). [CrossRef]

**5. **E. N. Glezer, M. Milosavljevic, L. Huang, R. J. Finlay, T. -H. Her, J. P. Callan, and E. Mazur, “Three-dimensional optical storage inside transparent materials,” Opt. Lett. **21**, 2023–2025 (1996). [CrossRef] [PubMed]

**6. **W. Watanabe, T. Toma, K. Yamada, J. Nishii, K. -I. Hayashi, and K. Itoh, “Optical seizing and merging of voids in silica glass with infrared femtosecond laser pulses,” Opt. Lett. **25**, 1669–1671 (2000). [CrossRef]

**7. **K. Minoshima, A. M. Kowalevicz, I. Hartl, E. P. Ippen, and J. G. Fujimoto, “Photonic device fabrication in glass by use of nonlinear materials processing with a femtosecond laser oscillator,” Opt. Lett. **26**, 1516–1518 (2001). [CrossRef]

**8. **D. Day and M. Gu, “Formation of voids in a doped polymethylmethacrylate polymer,” Appl. Phys. Lett. **80**, 2404–2406 (2002). [CrossRef]

**9. **K. Furusawa, K. Takahashi, H. Kumagai, K. Midorikawa, and M. Obara. “Ablation characteristics of Au, Ag, and Cu metals using a femtosecond Ti:sapphire laser,” Appl. Phys. **A 69**S359-S (1999).

**10. **Y. Cheng, K. Sugioka, and K. Midorikawa, “Freestanding optical fibers fabricated in a glass chip using femtosecond laser micromachining for lab-on-a-chip application,” Opt. Express **13**, 7225–7232 (2005). [CrossRef] [PubMed]

**11. **A. Vogel, J. Noack, G. Huttmann, and G. Paultauf, “Mechanisms of femtosecond laser nanosurgery of cells and tissues,” Appl. Phys. **B 81**, 1015–1047 (2005).

**12. **M. Gu, *Principles of Three-Dimensional Imaging in Confocal Microscopy* (World Scientific, Singapore, 1996). [CrossRef]

**13. **T. R. Corle and G. S. Kino, *Confocal Scanning Optical Microscopy and Related Imaging Systems*, (Academic Press, San Diego, 1996).

**14. **
Olympus Fluoview Resource Center, http://www.olympusfluoview.com/.

**15. **A. Taflove and S. C. Hagness, *Computational Electrodynamics: The Finite-Difference Time-Domain Method*, 2nd ed., (Artech House, Boston, 2000).

**16. **R. X. Bian, R. C. Dunn, and X. S. Xie, “Single molecule emission characteristics in near-field microscopy,” Phys. Rev. Lett. **75**, 4772–4775 (1995). [CrossRef] [PubMed]

**17. **T. C. Chu, W.-C. Liu, and D. P. Tsai, “Enhanced resolution induced by random silver nanoparticles in near-field optical disks,” Opt. Commun. **246**, 561–567 (2005). [CrossRef]

**18. **J. M. López-Alonso, J. M. Rico-García, and J. Alda, “Photonic crystal characterization by FDTD and principal component analysis,” Opt. Express **12**, 2176–2186 (2004). [CrossRef] [PubMed]

**19. **S. Chang and A. Taflove, “Finite-difference time-domain model of lasing action in a four-level two-electron atomic system,” Opt. Express **12**, 3827–3833 (2004). [CrossRef] [PubMed]

**20. **R. C. Gauthier, “Computation of the optical trapping force using an FDTD based technique,” Opt. Express **13**, 3707–3718 (2005). [CrossRef] [PubMed]

**21. **R. Drezek, A. Dunn, and R. Richards-Kortum, “Light scattering from cells: finite-difference time-domain simulations and goniometric measurements,” Appl. Opt. **38**, 3651–3661 (1999). [CrossRef]

**22. **J. B. Judkins, R. W. Ziolkowski, and C. W. Haggans, “Two-dimensional finite-difference time-domain simulation for rewritable optical disk surface structure design,” Appl. Opt. **35**, 2477–2487 (1996). [CrossRef] [PubMed]

**23. **D. W. Prather and S. Shi, “Formulation and application of the finite-difference time-domain method for the analysis of axially symmetric diffractive optical elements,” J. Opt. Soc. Am. **A 16**, 1131–1142 (1999). [CrossRef]

**24. **O. J. F. Martin, C. Girard, and A. Dereux, “Generalized field propagator for electromagnetic scattering and light confinement,” Phys. Rev. Lett. **74**, 526–529 (1995). [CrossRef] [PubMed]

**25. **C. Hafner, *The Generalized Multiple Multipole Technique for Computational Electromagnetics*, (Artech House, Boston, 1990).

**26. **K. S. Yee, “Numerical solution of initial boundary value problems involving Maxwell’s equations,” IEEE Trans. Antennas Propag. **14**, 302–307 (1966). [CrossRef]

**27. **J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic waves,” J. Comput. Phys. **114**, 185–200 (1994). [CrossRef]

**28. **M. Born and E. Wolf, *Principles of Optics*, 7th ed. (Pergamon, 1999).

**29. **S. Mellin and G. Nordin, “Limits of scalar diffraction theory and an iterative angular spectrum algorithm for finite aperture diffractive optical element design,” Opt. Express **8**, 705–722 (2001). [CrossRef] [PubMed]

**30. **E. G. Gamaly, S. Juodkazis, K. Nishimura, H. Misawa, and B. Luther-Davies, “Laser-matter interaction in the bulk of a transparent solid: Confined microexplosion and void formation,” Phys. Rev. **B 73**, 214101 (2006).

**31. **S. Juodkazis, K. Nishimura, S. Tanaka, H. Misawa, E. G. Gamaly, B. Luther-Davies, L. Hallo, P. Nicolai, and V. T. Tikhonchuk, “Laser-induced microexplosion confined in the bulk of a sapphire crystal: Evidence of multimegabar pressures,” Phys. Rev. Lett. **96**, 166101 (2006). [CrossRef] [PubMed]

**32. **M. Straub, M. Ventura, and M. Gu, “Multiple higher-order stop gaps in infrared polymer photonic crystals,” Phys. Rev. Lett. **91**, 043901 (2003). [CrossRef] [PubMed]

**33. **M. Straub, M. Ventura, and M. Gu, “Microvoid channel polymer photonic crystals with large infrared stop gaps and a multitude of higher-order bandgaps fabricated by femtosecond laser drilling in solid resin,” Thin Solid Films **453**–**454**, 522–526 (2004). [CrossRef]