Characterization of three-dimensional field distribution of bowtie aperture using quasi-spherical waves and surface plasmon polaritons

We present an analytical formula to predict the three-dimensional field distribution of a nanoscale bowtie aperture using quasi-spherical waves (QSWs) and surface plasmon polaritons, which are excited by the fundamental waveguide mode and local plasmons of the aperture, respectively. Assuming two separate bowtie apertures in a metal film, we analysed the decay characteristics of QSWs using a finite difference time-domain method. To verify the formula, we recorded the spot patterns of the bowtie aperture on a photoresist film using various exposure times, and fit the patterns to the analytical formula in terms of the width and depth of the patterns. In addition, it was found that the formula successfully represented the dipole characteristics of the spot patterns, which were in agreement with the surface geometry, with a root-mean-square error of 9.4%. We expect that our theoretical formula will extend the potential applications of nanoscale bowtie apertures to plasmonic device fabrication, three-dimensional plasmonic lithography, and other technologies.

Scientific RepoRts | 7:45352 | DOI: 10.1038/srep45352 In this paper, we propose an analytical formula to predict the three-dimensional field intensity distribution around a bowtie aperture, which is one type of ridge aperture. In order to describe the surface waves by the bowtie aperture, we extended the QCW concept to quasi-spherical waves (QSWs) in a two-dimensional hole-type aperture, such as a circular hole or ridge aperture. In order to simplify the problem, we used the empirical solution to predict the QSW attenuation rate. We expected local plasmons to excite Hankel-type SPPs and QSWs, and for them to have a cosine dependence on the azimuthal angle. In a surface wave analysis using the finite difference time-domain (FDTD) method, we determined that the QSWs had an asymptotic attenuation rate affected by the distance from the aperture in a real metal. We quantitatively verified that our formula fit the surface geometry of spot patterns recorded on a photoresist (PR) thin film.

Results and Discussion
Intensity distribution by ridge aperture. When light is incident on a ridge aperture, as shown in Fig. 1, two emission sources are generated underneath the aperture: one is a horizontal dipole source, which is local plasmon located on the metal-dielectric interface, and the other is a point source created by the field transmitted through the aperture. The horizontal dipole excites surface waves composed of the evanescent mode of QSWs and SPPs.
For a perfect conductor, the evanescent mode of the QSW has a damping rate of ρ −1 , which is the same as the attenuation rate of a spherical wave. For a real metal, however, the decay rate of QSW is similar to the decay rate of a spherical wave in the vicinity of the aperture, and this wave has an asymptotic ρ −2 attenuation over a long distance from the aperture 7,20 . In the intermediate regime, we expect the evanescent mode of the QSW to fall off with increasing distance (as ρ −q ). Thus, with ϕ = 0, the magnetic field on the surface, H y = H SPP + H QSW , can be written as where δ is an arbitrary phase constant of the SPP, φ is an arbitrary phase constant of the QSW, A SPP is the SPP magnetic field amplitude, A QSW is the QSW magnetic field amplitude, and ρ, ϕ, z are the cylindrical coordinates. With Maxwell's curl equation applied to the evanescent mode of the QSW, the electric field on the surface is represented by where B SPP is the electric field amplitude of the SPP, and B QSW is the electric field amplitude of the QSW. We assumed that the electric field and magnetic field have a harmonic dependence on time, i.e., in the form of exp(−iwt). As shown in Fig. 1, the electromagnetic field of the surface wave has a cosine dependence of ϕ, and it can be achieved by multiplying both equations (1) and (2) with cos ϕ.
With equations (1) and (2), the intensity of the surface wave can be expressed as In Fig. 1, ψ represents the electromagnetic field vector defined by [E z , H y ]. The surface wave propagates on the metallic surface with cosine dependence, and the space wave propagates spherically because of the small ridge gap. The purple wave represents a field whose wavelength is equal to the incident wavelength; the red wave represents a field whose wavelength is shorter than the incident wavelength.
Scientific RepoRts | 7:45352 | DOI: 10.1038/srep45352 where E z is the z component of the electric field; H y * is the complex conjugate of the magnetic field on the surface, H y ; Φ is the phase delay between the evanescent modes of the QSW and SPP; and ρ, ϕ, z are the cylindrical coordinates. In order to derive equation (3), we used k SPP ≒ k 0 and disregarded the difference between the wave vectors of the SPP and QSW. The cosine term in equation (3) originates from the dipole source.
In addition to the surface wave, the space wave defined as the propagating mode of the QSW and transmitted field also contributes to the field distribution. Assuming that the ridge gap is so small that the source of the space wave can be regarded as a point source, the irradiance from the fundamental mode of the waveguide follows the inverse square law: where C is a constant related to the space wave 18 .
Evaluation of QSW decay rate. We calculated the QSW decay rate with FDTD (Lumerical 8.12.590) simulation. For experimentalists or some theoreticians in the field of nano-apertures, FDTD simulations are somewhat popularly used to analyse the electromagnetic field distribution of the bowtie aperture, slit aperture, and other waveguide because the FDTD solves Maxwell's equations numerically 16,[21][22][23] .
To simulate a situation in which both a QSW and SPP are generated underneath a bowtie-shaped aperture, we assumed a metal film with two bowtie apertures separated by a distance of x s . Figure 2 shows a schematic of the FDTD simulation geometry for extracting the QSW from the total field. In the FDTD analysis, we used an incident light wavelength of 405 nm, a metal layer of aluminium with a permittivity of − 23.9819 + 4.9508i, and a dielectric layer of silicon nitride with permittivity of 4.2845; x s was set at 4400 nm to meet the condition of a negligibly small QSW. We normalized the incident electric field by E = 1 and the magnetic field H = 1/Z 0 , where Z 0 is the vacuum impedance, denoted by (μ 0 /ε 0 ) 1/2 .
To estimate the QSW decay rate, we extracted the QSW from the total field. To separate the QSW from the surface wave, we used the interference of the surface wave by placing two ridge apertures separated by x s of 4400 nm, which is larger than 10 λ. In the vicinity of the central point between the two ridge apertures, we expect interference of SPP waves propagating in opposite directions along the surface where the SPP is much larger than the QSW, because of the latter's fast attenuation. The SPP field can then be obtained in the overall region. By subtracting the SPP from the total field, we can obtain the field distribution of the QSW in the FDTD simulation. Figure 3 shows the results of the FDTD simulation for the geometry shown in Fig. 2. As can be seen in Fig. 3a, light illumination on a subwavelength metallic ridge aperture results in horizontal electric dipole radiation, which depends on the azimuthal angle on the surface. On the contrary, the y component of the magnetic field generated by the aperture is dominated by the fundamental mode of its waveguide, and it is independent of the azimuthal angle. In order to analyse the surface wave along the metal surface, we focused on the field distribution along the x axis (y = 0) on the metal-air interface.
At the interface between the metal film and air, the magnetic field of the QSW falls off as a function of x −q , where q is a factor determining the decay rate of the field. In order to fit the FDTD simulation results with the fitting equation, we adopted the sinusoidal function instead of Euler's form. The electromagnetic field along the x axis on the surface can be written as where ψ is the electromagnetic field vector defined by [E z , H y ], ψ SPP is the electromagnetic field for the SPP, ψ QSW is the electromagnetic field for QSW, and the sign is + and − for increasing and decreasing, respectively, in the x direction. Far from the ridge aperture, the electromagnetic field is dominated by the SPP because of the fast attenuation of the QSW. We speculate that the QSW is much smaller than the SPP and can be disregarded at a distance above 1.6 μ m. For 1.6 μ m < x < 2.8 μ m, the electromagnetic field can be as approximated by

SPP S PP s
In order to predict the SPP field distribution, we used the FDTD simulation and fitted the numerical results with equation (6). The electric field of the SPP waves in the region of 1.6 μ m < x < 2.8 μ m can be described by   where δ is an arbitrary phase constant, and E SPP is the amplitude of the SPP electric field. Owing to the dipole characteristics, the fields from the apertures at x = 0 and x = 4400 nm have different signs. Figure 4 shows the electric field distribution of the QSW and SPP. It is found that the FDTD simulation data is well matched with the fitted curve with equation (7) as shown in inset of Fig. 4. From the fitting process, we obtained E SPP as 0.3544 and δ as − 0.1209. To obtain the decay rate of the QSW, we subtracted the SPP contribution from FDTD simulation results in the entire region. From equations (2) and (5), the QSW can be described by where δ − φ is a possible phase delay between the SPP and QSW, q is the decay rate, E QSW is the QSW electric field amplitude, and k 0 is the wave vector in free space. By fitting the subtracted data with equation (8) Assuming that the intensity distribution follows a Gaussian distribution, the exposure dose distribution by the surface wave on the PR can be described by where α is the attenuation coefficient, σ x is w /2 2 ln 10 x , w x is the full width at tenth of maximum (FWTM) in the x direction, and I 0,surf is the peak intensity due to the surface wave. For simplicity, we used the fact that α is the effective attenuation coefficient related with the SPP attenuation coefficient and the QSW evanescent mode. The dose distribution by the surface wave has a cosine dependence of the azimuthal angle because of the dipole radiation. Similarly, the exposure dose distribution by the propagating wave on the PR can be described by where σ y is w /2 2 ln 2 y , I 0,space is the peak intensity due to the space wave, and w y is the full width at half maximum (FWHM) in the y direction.
In order to demonstrate the decay rate of the space wave and surface wave obtained from the FDTD method, we used the calibration curve, such as w x and w y in equations (9) and (10). The calibration curve is defined as the geometrical dimension of the spot pattern with respect to exposure dose. By imposing the threshold dose of the PR, the theoretical calibration curves for the x-width (w x ) and y-width (w y ) can be obtained. Compared with the y-width, which can be described by w y = f(D ex ), the x-width is mathematically complex. Thus, by taking the inverse function D ex = g(w x ), we can obtain the calibration curve for the x-width. To calibrate the system, we used the method reported in ref. 24. By varying the exposure time in increments of 1 ms and using a fixed light intensity, we obtained a spot pattern array. By measuring the geometrical features of each spot pattern with an atomic force microscope (AFM), we obtained the calibration data for the x-width and y-width.
By using the decay rate of intensity with the threshold dose condition, the depth, x-width, and y-width can be obtained. In this study, the x-width is defined as the FWTM of the intensity distribution in the x direction, and the y-width is defined as the FWHM of the intensity distribution in the y direction. Because the x-width is determined near the metal surface, the x-width is mainly determined by E z , which is related to the surface wave due to the boundary condition of the conductor. On the contrary, because the surface wave excited by the horizontal electric dipole does not affect the field intensity in the y direction, the y-width is determined only by the space wave.
With the threshold dose and decaying rate of the spherical wave, the y-width can be calculated as where w y is the y-width, which is the FWHM in the y direction; D th is the threshold dose on the photoresist. By fitting the calibration data with equation (11), we obtained C/D th = 272.75 m 2 /s. In addition, an implicit function of the x-width can be obtained from the decay rate of the surface wave intensity (from equation (3)) and threshold condition as follows: where A SPP is the SPP magnetic field amplitude, B SPP is the SPP electric field amplitude, A QSW is the QSW magnetic field amplitude, B QSW is the QSW electric field amplitude, Φ is the phase difference between the SPP and QSW,  and w x is the x-width, which is the FWTM in the x direction. Before we fitted the results with equation (12), we used q = 1.18 and Φ = − 0.56, which were obtained from the FDTD simulation results with D th fixed at a constant value of 1. From the fitting process with equation (12), we obtained A SPP = 4.3, B spp = − 11.52, A QSW = − 1.64, and B QSW = 1.78 × 10 −4 . Figure 5 shows the calibration data and fitting curve for the x-width and y-width. With equations (9)-(12), we calculated the three-dimensional (3D) dose distribution by the bowtie aperture. By imposing the threshold dose condition on the total dose distribution, we obtained the three-dimensional spot shape in Fig. 6b. To assess the accuracy of the analytical formula, we used an ellipsoid described by an ellipsoidal Gaussian function as the reference model for comparison. Figure 6a shows an AFM image at a spot depth of 130 nm, and Fig. 6b shows the simulated result of our analytical formula, which considers the dipole radiation effects. Our analytical formula is close to reality, even though the effect of dipole radiation from the aperture creates an atypical shape of the spot. We defined the fitting error with the difference between the spot pattern and either the ellipsoidal Gaussian function or our analytical formula fitted to the spot pattern. Figure 6c and d show the error distribution of the ellipsoidal Gaussian function and our analytical formula, respectively. The analytical formula has a peak-to-valley (P-V) error value of 30.60 nm and a root-mean-square (RMS) error of 12.31 nm; the reference model has a P-V error of 67.65 nm and an RMS error of 18.75 nm. This means that the P-V and RMS errors in the analytical formula are 55% and 34% lower, respectively, than those of the reference model.
In summary, we presented an analytical formula of the 3D field intensity distribution by a bowtie aperture using the concept of QSWs and SPPs, which were generated by the fundamental waveguide mode and local surface plasmon polaritons of the aperture, respectively. Assuming two separate bowtie apertures on a metal film, we calculated the attenuation rate of the QSWs on the metallic surface using the FDTD method and obtained an asymptotic decay rate q of 1.18 for distance x < 11 λ. Using the spot patterns recorded on a photoresist film, we calibrated the formula in terms of the width and depth of the patterns and determined the coefficients for the analytical formula. Our analytical formula successfully reflected the dipole characteristics of the optical spots and was in close agreement with the surface geometry of the spot patterns, with an RMS error of 9.4%. We expect that our theoretical formula will extend the potential applications of nanoscale bowtie apertures to plasmonic device fabrication, three-dimensional plasmonic lithography, and other technologies.

Methods
Experimental system. The experimental setup of the plasmonic lithography system was shown in our previous work 24 . A continuous-wave diode laser with a wavelength of 405 nm was used as the light illumination source, and an objective lens (CFI LU Plan Epi El WE, Nikon; 100× , Na: 0.8) was used to focus the laser beam onto the aperture. The bowtie aperture was perforated in a silicon nitride layer with a thickness of 150 nm coated with a 150 nm aluminium layer. This configuration is well described in Fig. 2 for the case of one bowtie aperture.