## Abstract

Extensive 3-D finite-difference time-domain simulations are carried out to elucidate the nature of surface plasmon polaritons (SPPs) and localized surface plasmon polaritons (LSPs) generated by nanoscale holes in thin metallic films interacting with light. Both isolated nanoholes and square arrays of nanoholes in gold films are considered. For isolated nanoholes, we expand on an earlier discussion of Yin et al. [Appl. Phys. Lett. 85, 467–469 (2004)] on the origins of fringe patterns in the film and the role of near-field scanning optical microscope probe interactions. The associated light transmission of a single nanohole is enhanced when a LSP excitation of the nanohole itself is excited. Periodic arrays of nanoholes exhibit more complex behavior, with light transmission peaks exhibiting distinct minima and maxima that can be very well described with Fano lineshape models. This behavior is correlated with the coupling of SPP Bloch waves and more directly transmitted waves through the holes.

©2005 Optical Society of America

## 1. Introduction

The pioneering experimental work of Ebbesen et al. [1] on optical transmission through subwavelength hole arrays has stimulated much work, with the state-of-the-art being well described by several articles in a recent focus issue of Optics Express [2] and references therein. Despite so much work, there continues to be uncertainty in the mechanistic details. For example, the precise role played by surface plasmon polaritons (SPPs) in achieving light transmission that is more than expected based on independent holes is still not completely resolved. We should point out that not only light transmission, but the general behavior of SPPs in the plane of the metal film are of interest from the point of view of potential device applications. Recently, for example, experimental near-field scanning microscopy (NSOM) coupled with realistic 3-D time-domain simulations showed that isolated nanoholes could act as point sources of SPPs in thin Au films [3].

The present contribution extends the isolated nanohole simulations [3] to include discussion of transmission and gives simulation results for square arrays of nanoholes in gold films. Analogous calculations corresponding to a perfect metal where no plasmons can exist are also carried out. These calculations allow us to develop a more complete understanding of both SPP generation and the light transmission. Among the issues we address are the roles played by localized surface plasmon (LSP) excitations around the nanoholes [4–6], Wood’s anomalies [7–10], SPP Bloch waves [11,12], and Fano resonances [8,13]. In the case of an isolated nanohole we also present a more detailed analysis of the effect of the NSOM probe [3]. Section 2 below outlines computational details, Section 3 discusses isolated nanoholes, Section 4 discusses arrays of nanoholes, and Section 5 summarizes our conclusions and places them in the context of other analyses.

## 2. Computational details

With Cartesian coordinates *x, y*, and *z*, a system consists of an underlying glass layer for *z*<-*h*, a gold (Au) film for *z*=-*h* to 0, and air for *z*>0. One or more cylindrical nanoholes of diameter *d* from *z*=-*h* to 0, containing air, are present. We examine surface plasmon generation and propagation on the film surfaces, and transmission of light into the air. While we investigated numerous variations on the system parameters, all of the numerical results presented correspond to either gold (Au) or perfect electrical conductor (PEC) films with *h*=100 nm, *d*=200 nm and, for square arrays, lattice spacing *D*=600 nm.

The finite-difference time-domain (FDTD) method [14] is used to solve Maxwell’s equations. It represents the electric, magnetic and current density fields, * E*(

*x,y,z,t*),

*(*

**H***x,y,z,t*), and

*(*

**J***x,y,z,t*), on evenly spaced, staggered spatial grids and involves a simple time-stepping algorithm. Each system is defined by specifying the relative dielectric constant,

*ε*

_{R}(

*x, y, z, ω*). For air we set

*ε*

_{R}=

*ε*

_{air}=1, for glass we set

*ε*

_{R}=

*ε*

_{glass}=2.25, and for Au we set

*ε*

_{R}=

*ε*

_{Au}(

*ω*), with a Drude model,

for the frequency dependence. We choose the parameter set (*ε*_{∞}*, ω*_{D}*, γ*_{D}
) to *fit* empirical data for the real and imaginary parts of the metal’s dielectric constant over a given frequency or wavelength range of interest [15]. Specifically, to fit the Au empirical dielectric constant information [16] for wavelengths close to 532 nm, we find set D-Au-1, (*ε*_{∞}*, ω*_{D}*, γ*_{D}
)=(12.9965, 9.8528 eV, 0.2401eV), and for a more global description of wavelengths in the 500–1000 nm we find set D-Au-2, (*ε∞, ω*_{D}*, γ*_{D}
)=(11.4577, 9.4027 eV, 0.08314 eV). Set D-Au-1 is the same set as used in our preliminary study [3]; near 532 nm it gives *ε*_{metal}
=-4.9+*i*1.8, in reasonable accord with the interpolated empirical value [16] of -4.9+*i*2.3. Set D-Au-2 is less accurate than D-Au-1 in the 500–550 nm range (mostly due to underestimating the imaginary part), but is significantly better than D-Au-1 in the 550–1000 nm range, yielding real parts within 10 % and imaginary parts within 35% of the empirical values. (When treating the metal film as a PEC we set the electric field to be zero within the film.)

In order to describe the frequency-dependence in Eq. (1), an auxiliary differential equations scheme is employed [14,15]. This involves not explicitly including Eq. (1) in the FDTD update equations but propagating, in addition to the three components each of * E* and

*, three additional field components associated with a current density vector,*

**H***. The complete system of time-domain equations for*

**J***, and*

**E, H***is designed to be*

**J***equivalent*to the frequency-domain problem with metallic regions described by Eq. (1) [15]. The total field/scattered field formalism is used to generate appropriate incident plane waves from below the metal film, and the grids are truncated with uniaxial perfect matching layers to accomplish absorption of field components approaching grid edges where appropriate [14].

For the single hole case, the incident field is launched from the glass substrate with a flat circular profile of diameter 6 µm which gives a plane wave like source for the majority of the film. The edges are truncated with a Gaussian profile. This incident field is not intended to mimic an experimental laser beam, but to mimic normally incident infinite plane waves close to the hole. Far from the hole, there is some wave vector dispersion, but we verified that it does not affect any of the results displayed here. Oblique incidence is discussed in [3].

In the case of periodic structures, a perfect plane wave is launched from the glass substrate. The FDTD method can be used to obtain spectra of various types over a range of wavelengths from a single calculation. A calculation is then performed with the incident wave being multiplied by a Gaussian in time in order to obtain a pulse that contains information over a desired range of wavelengths. We often display wavelength or frequency resolved results, which are obtained by complex Fourier transformation of our time-domain field components at one or more frequencies of interest to yield complex, frequency-resolved fields, the absolute magnitude of which can be associated with time-averaged intensities.

Total transmission spectra *T*_{tot}*(λ)* are obtained by Fourier transforming the electric and magnetic fields on a surface above the holes and constructing the surface integral of the outward Poynting vector. The transmission spectra presented in our figures are obtained by subtracting out the transmission spectrum *T*_{film}
(λ) for the case of a pure film with no holes from the total transmission spectrum, and normalizing by the incident power over the hole,

where *I*_{inc}
is the incident intensity and *a*=*d*/2 is the hole radius.

All results presented are based on a 3-D, parallel FDTD implementation that uses the Message Passing Interface library similar in manner to Ref. [17]. (Note that for the isolated nanohole case, with no NSOM probe present, cylindrical coordinates could be used to simplify the problem and reduce the numerical effort [18].) The grid spacing in each spatial dimension is 4 nm. For the isolated nanohole calculations, the computational grid is 7.5 µm×7.5 µm×0.4 µm.

The calculations that include the NSOM probe are similar to those in Ref. [3]. An experimental aperture NSOM probe in collection mode is explicitly modeled (see Fig. 3 of Ref. [3]). It involves a glass fiber with 100 nm thick Al side coating. The tip aperture diameter is 80 nm, and opens up to a cylinder of diameter 240 nm. The total height of the tip structure is 0.6 µm and it is truncated by the UPML. The surface integral of the Poynting vector over the cylinder top is a measure of the NSOM signal. We verify that the flux calculation is at a sufficient height to be out of the evanescent tip region. The experimental tip is estimated about 5~10nm above the metal surface [3]; we take it to be 8 nm above the surface. In this case the total computational grid is 7.5 µm×7.5 µm×1 µm and a separate calculation is run for each position of the probe.

For the periodic square arrays with period of 600 nm, the computational grid is 0.6 µm×0.6 µm×4 µm. The extra length in the vertical direction allows us to resolve Wood’s anomalies, if present. Periodic boundary conditions are applied to simulate an infinite square array. Outgoing field components are absorbed towards the edge of the relevant grids with Uniaxial Perfect Matching Layers (UPML) [14] adapted to include Drude metals.

## 3. Isolated nanoholes

Conventional SPPs on a metal film/dielectric interface propagate in the *x-y* plane with wave vector magnitude [19]

With incident light propagating from below along *z* as we assume, there are no incident field components in the *x-y* plane to couple to SPPs. However, SPPs are still generated due to diffraction by the hole edges. They are particular evanescent diffraction modes [10,20] with real wave vector magnitude greater than *ω/c* and decay factor *exp*(-*α z*), *α*≥0, associated with an imaginary component for the upward or *z*-component of the wave vector. Evanescent diffraction modes can also occur if the metal is described as a perfect electrical conductor (PEC) where no plasmons can exist and, indeed, in hole systems in dielectric films. As emphasized recently by Lezec and Thio [10], one must not immediately jump to the conclusion that SPPs are at the origin of all enhanced transmission observations. However, strong evidence for their role in the particular systems we study here will be presented.

#### 3.1 Point source SPP waves and fringe patterns

SPP waves created along the rim of a nanohole, if the nanohole is sufficiently small, act as if the center of the hole is a point source for waves in the plane of the surface [3]. This makes them behave somewhat differently from conventional SPPs. We give an approximate complex phasor form for them with cylindrical coordinates, (*ρ, φ, z*), and unit vectors, (*$\widehat{\rho}$
, $\widehat{\phi}$
, ẑ*). Assuming a small hole centered at the origin, and incident *x*-polarized light, the frequency-resolved complex form is, for *z*≥0,

where *A* is a constant. ${H}_{1}^{\left(1\right)}$
^{(1)} (*k*_{SPPρ}
) is the *m*=1 Hankel function [21],

*J*_{m}
and *Y*_{m}
are integer-order Bessel functions of the first and second kind, respectively. *Y*_{m}
diverges at the origin representing the point source. It is easy to see that *E*^{SPP}*z=ẑ* ·* E^{SPP}* exactly satisfies Helmholtz’s equation for

*z*>0; the assumed evanescent

*z*-dependence leads to

*k*

_{SPP}>ω/c, consistent with Eq. (3). While the relation ∇·

*=0 is not exactly satisfied by Eq. (4), it is satisfied in the large*

**E**^{SPP}*k*

_{SPP}

*ρ*limit where [21]

is valid. The Cartesian *E*_{x}
and *E*_{y}
components inferred from the cylindrical ones in Eq. (4) do not exactly satisfy the vector Helmholtz equation except in this limit. In practice, Eq. (6) is quite good except very close to the origin. More exact, but less transparent forms can also be developed by inserting ${E}_{\phi}^{\mathit{\text{SPP}}}$
=0 and ${E}_{z}^{\mathit{\text{SPP}}}$
from Eq. (4) into ∇·
^{SPP}
=0, and integrating to obtain ${E}_{\rho}^{\mathit{\text{SPP}}}$
.

Expressions similar to ${E}_{z}^{\mathit{\text{SPP}}}$
from Eq. (4), but with any other ${H}_{m}^{\mathit{\left(}\mathit{1}\mathit{\right)}}$
function multiplied by either *cos(mφ)* or *sin(mφ)* [21], and linear combinations thereof, also exactly satisfy Helmholtz’s equation. It is the specific form for the incident wave that leads to just a *cos(φ)* dependence. Similarly, ${E}_{x}^{\mathit{\text{SPP}}}$
can have only isotropic and *cos(2φ)* angular dependencies, consistent with *m*=0 and *m*=2 Hankel outgoing contributions, and ${E}_{y}^{\mathit{\text{SPP}}}$
can have only a *sin(2φ)* dependence, consistent with an *m*=1 Hankel contribution. (This is shown by expanding the field components in series involving *sin(mφ)* and *cos(mφ)*, insertion into Maxwell’s equations as in Ref. [18], and deducing which terms couple into the incident wave.) Again, the ${E}_{x}^{\mathit{\text{SPP}}}$
and ${E}_{y}^{\mathit{\text{SPP}}}$
components inferred from Eq. (5) can be shown to be consistent with these angular dependencies if one is in the limit where Eq. (7) applies.

Figure 1(a) displays the time-averaged |*E*_{x}
|^{2} in the *x-y* plane near *z*=0 inferred from Fourier transformation of our FDTD calculations at *λ*=532 nm. Radial fringes that correlate reasonably well with periodicity *λ*_{SPP}
=2*π/k*_{SPP}
≈480 nm, particularly for large *ρ*, are evident. Figure 1(b) corresponds to an *x-z* cut with *y*=0 of the same intensity as in Fig. 1(a). The fringes are seen as structured intensity patterns near top metal surface. However, time-averaged field intensities, if only SPPs as above contribute, would not display significant fringes along *ρ* because, the absolute square of either Eq. (6) or (7) exhibits a simple *exp*[- *ρ/L*]/*ρ* decay for increasing *ρ*. Here we have also introduced an SPP radial decay factor [1/*L*≈2 Im(*k*_{SPP}
)]. As suggested by Wannemacher [4], interference between SPP waves and directly penetrating waves through the metal film leads to such fringes. Ref. [3] presented a simple analysis of this. We give here a more explicit pictorial demonstration. To estimate the directly penetrating waves, we carried out a comparable calculation to that in the previous sections but for an Au film with no hole in it, Fig. 1(c). Interestingly, in addition to directly penetrating waves there are still some small fringes evident in Fig. 1(c). They arise because our incident light is a circular spot, not a perfect plane wave (Sec. 2), leading to a small amount of SPP excitation even *without a* hole. The creation of surface plasmons in structureless metal films due to focused laser beams has been observed experimentally [22]. Figure 1(d) is the result of subtracting the amplitude of the field in Fig. 1(c) from that of Fig. 1(b), and taking the absolute square of the result. The large fringes evident in Fig. 1(b) are absent in Fig. 1(d). The SPP intensity now appears as more smoothly varying intensity on the metal/air interface. Figure 1(d) also shows the presence of relatively larger SPP intensity on the bottom metal/glass interface. These bottom SPPs were not visible in Fig. 1(b), owing to saturation by the incident intensity.

Equations (3) and (4) should provide a good description of the SPP field contributions when the hole diameter *d*≪*λ*. This condition is only weakly satisfied with *d*=200 nm and the *λ*=532 nm wavelength of interest. In this case, along a given radial line, it is better to view the edge of the hole as source, i.e. to replace *ρ* with *ρ′=ρ-d/2* in Eq. (4). Figure 2 displays, as symbols, the *y*=0, *z*=4 nm cuts along *x* of the real and imaginary parts of the complex *E*_{x}
and *E*_{z}
fields at *λ*=532 nm inferred from our FDTD calculations that result after subtracting the directly penetrating amplitudes discussed above. The result of fits based on Eqs. (3) and (4) with *ρ* replaced by *ρ*′, and also including an extra exp[-*ρ*′/L] loss term, are given as curves in Fig. 2. We see that for *x*>100 nm, it is possible to fit the FDTD results very well. Thus, Re[*E*_{z}
]∝*J*
_{1}
*exp*(-*ρ*′/*L*), and Im[*E*_{z}
]∝*Y*_{1}
*exp*(-*ρ*′/*L*). Since *φ*=0 we have *E*_{x}
=*E*ρ, and *E*_{x}
should be exp[*iπ*]=-*i* away from *E*_{z}
according to Eq. (4). This type of phase relation is also seen in ordinary SPPs emanating from lines [19]. In our case it implies that Re[*E*_{x}
]∝*Y*_{1}
*exp*(-*ρ*′/*L*) and Im[*E*_{x}
]∝- *J*_{1}* exp*(-*ρ*′/*L*). Figure 2 is indeed consistent with these relations. For x<100 nm, the simple point source SPP model is not applicable.

The SPP behavior of Eq. (4) is consistent with zero field in the y-z plane (*φ=π*/2). However, Fig. 1(a) indicates this is not so. The reasons for this are non-SPP features such as dipolar and LSP radiation.

#### 3.2 NSOM probe interactions

Reference [3] presented preliminary results of our FDTD calculations that included an NSOM probe, demonstrating a correlation of the upward probe-free FDTD Poynting vector with the both the experimental results and the FDTD results based on including the probe. Here we extend these results and demonstrate a remarkable correlation of the intensity contribution from the field components parallel to the surface, |*E*_{parallel}
|^{2}=|*E*_{x}
|^{2}+|*E*_{y}
|^{2} with the results allowing for the probe. The nature of the calculations is outlined in Sec. 2 above and also in Ref. [3]. Figure 3 shows| *E*_{parallel}
|^{2} calculated by FDTD without the NSOM probe (solid curve) and the collected NSOM signal inferred from our FDTD calculation including the NSOM probe (symbols). The good agreement shows two important facts. First, the fringing is still preserved even with the presence of the NSOM probe. Second, the signal picked up by the NSOM probe is well described |*E*_{parallel}
|^{2}, rather than the total intensity which is also shown. The probe-free |*E*_{parallel}
|^{2} is a better indicator of the NSOM measurement than the upward Poynting vector component of Ref. [3] which deviates from the probe result at the hole edge due to energy absorption by the hole. We also carried out the same calculation for the case of a nanohole in a PEC film. No fringes were observed.

#### 3.3 Transmission from isolated nanoholes

Figure 4 shows transmission spectra calculated treating the metal film as (i) Au with the appropriate Drude model above and as (ii) a PEC. See Sec. 2 for further system details. In the PEC case, the transmission monotonically decreases as the wavelength becomes larger. This is qualitatively consistent with Bethe-Bouwkamp theory [23, 24] which predicts that transmission should scale as *1/λ*^{4}
in the small hole limit. (Note that we are treating a finite metal thickness, not the infinitely thin PEC layer of Bethe and Bouwkamp. The magnitude of the transmission from this latter case is significantly larger than the finite height case.) Fig. 4 shows that the Au film case is significantly different from the PEC case, with noticeably larger transmission over the whole 500–1000 nm wavelength range. The maximum in the transmission is near 580 nm and can be associated with the excitation of a localized surface plasmon (LSP). The LSP resonance structure near 580 nm corresponds to a dipolar charge oscillation along the polarization direction of the incident field in the hole walls. Inspection of the field components (not shown) indicates the near fields are most intense near the hole edges near the bottom metal/glass interface, indicating the character of the nanohole LSP resonance. Figure 4 also shows a smaller peak/shoulder structure near 500 nm which is associated with an LSP nanohole resonance that is more concentrated at the hole edges of the top air/metal interface. Calculations with different system parameters (not shown) indicate that the position of the transmission maximum red shifts with increasing nanohole diameter, and blue shifts with increasing film thickness. Such behavior is also seen in the LSP resonances of metal nanoparticle disks [5, 25]. The electrostatic approximation found in [26] for ellipsoid nanoparticles gives a similar physical interpretation of such resonant shifts. The correlation of the hole and particle excitations is consistent with Babinet’s principle [27].

In contrast to a metal nanoparticle, however, the nanohole excitation decays not only into scattered light, but into surface plasmon polaritons (SPPs) emanating from the hole in the plane of the film, as shown in Ref. [3]. We have also monitored the flux of the SPPs on the top (air/metal) and bottom (metal/glass) interfaces and find, near 580 nm that bottom SPPs are most intense whereas near 500 nm top SPPs are most intense, i.e. SPP production on the film surfaces does correlate with transmission.

## 4. Arrays of nanoholes

A square, 2-D array of nanoholes with lattice spacing D can act as a grating, providing wave vector components 2*πn*_{x}*/D* along x and 2*πn*_{y}*/D* along y, with *n*_{x}
and *n*_{y}
integers, sufficient to excite SPP wave vectors of magnitude *k*_{SPP}
, Eq. (3). This implies *k*_{SPP}
=(*2π/D*) (${{n}_{x}}^{\mathit{2}}$ + ${{n}_{y}}^{\mathit{2}}$
)^{1/2}, which rearranges to [9]

Eq. (8) is an implicit equation for *λ*_{SPP}
because the metallic dielectric constant depends on wavelength. SPPs consistent with Eq. (8) are different in nature from the isolated hole SPPs of Sec. 3: (i) they occur at only certain, discrete wavelengths and (ii) they are best thought of as SPP-Bloch waves (SPP-BWs) [11,12], i.e. standing waves consistent with Bloch’s theorem. If one repeats this grating analysis but replaces *k*_{SPP}
with the wave vector for ordinary light, *k*=${\epsilon}_{d}^{\mathit{1}\mathit{/}\mathit{2}}$* ω/*_{c}
, one arrives at *λ*_{W}
=*Dε*^{1/2}
/_{d}(${{n}_{x}}^{\mathit{2}}$ + ${{n}_{y}}^{\mathit{2}}$
)^{1/2}, the condition for Wood’s anomalies, i.e., light waves diffracted to move in the plane of the surface. Both SPP-BWs and Wood’s anomalies can exist in the same problem and, depending on system details, lie close to one another.

While SPP-BWs were originally associated with light transmission maxima and Wood’s anomalies with the nearby minima [9], different interpretations have been suggested including those in Refs. [7–13, 28]. The interpretation [8, 13] we have found most appropriate of our systems draws parallels with Fano’s discrete + continuum state model of resonance scattering in quantum mechanics [29]. It associates the SPP-BW with the discrete state and direct hole transmission with the continuum state. Depending on the system studied, the Wood’s anomalies might also act as discrete states in the same fashion. Our FDTD studies below largely confirm and solidify this interpretation. Our results also point to an extension of this picture to include contributions from LSP excitations as discussed in Sec. 3.

#### 4.1 Transmission from nanohole arrays

The solid curve in Fig. 5(a) is the calculated transmission spectrum for a periodic 2-D array of *d*=200 nm holes in a gold film with lattice spacing *D*=600 nm and other details as in Sec. 2.

The transmission is calculated over the cell area for each individual hole and normalized by the incident energy flux over one hole area. This allows direct comparison with the isolated hole transmission. The corresponding results, but treating the metal film as a PEC with no possibility of any SPPs, are given in Fig. 5(b).

The 2-D array transmission in Fig. 5(a) has maxima that are higher than the isolated hole results. However, only the peak near *λ*=945 nm represents a significant relative enhancement of ≈10. This modest enhancement is also consistent with recent results of Lezec and Thio [10]. Note that the isolated hole result in Fig. 5(a) still has nonperiodic plasmon effects present -- LSPs and propagating SPPs as discussed in Sec. 3. The overall shape of the 2-D array transmission is in fact similar to that of the isolated hole case (apart from the 945 nm region), indicating that such independent hole plasmon excitation is also important in the 2-D case. If one compares the 2-D gold film results to the PEC film results of Fig. 5(b) where no plasmons of any sort can exist, the relative enhancement inferred can be significantly larger than 10. In this sense it is still true for the 2-D array studied here that there can be large transmission enhancements when a metal that supports plasmon excitation is used. Regarding the 2-D PEC array results in Fig. 5(b), while the magnitude of the transmission is considerably less than the gold case of Fig. 5(a), there are also discrete maxima and the PEC array lead to increased transmission compared to the isolated hole result.

The lower wavelength minima adjacent to each maximum in the 2-D array transmission peak structures in Fig. 5 can be correlated with either SPP-BWs or Wood’s anomalies. For example, the minimum to the left of the 945 nm maximum occurs at 922 nm. From Eq. (8), the *n*_{x}
=±1, *n*_{y}
=0 (or vice versa) glass/metal interface SPP-BW, which we will denote as (1,0)
_{glass}
is predicted to occur at 928 nm, just 6 nm of the observed minimum. Similarly, the transmission minima at 676 nm, 628 nm, and 550 nm can be correlated with the (1,1)
_{glass}
, (1,0)
_{air}
, and (2,0)
_{glass}
SPP-BWs at 688 nm, 629 nm, and 558 nm, respectively. As wavelength decreases, and consistent with Eq. (8), the density of SPP-BWs increases such that overlapping resonance features lead to more complicated patterns. This is already the case by the 550 nm minimum in Fig. 5(a). It is likely that the (2,1)
_{glass}
SPP-BW feature, which is predicted to have a minimum near 537 nm and a maximum likely in the 550 nm region, is overlapping with it. For the PEC case, Fig. 5(b), there are just three Wood’s anomalies predicted in the spectral range shown: 900 nm, 636 nm and 600 nm, corresponding to (1,0)
_{glass}
and (1,1)
_{glass}
and (1,0)
_{air}
, respectively. Fig. 5(b) shows PEC resonance features within 1–2 nm these predictions, with transmission minima 898 nm, 636 nm, and 600 nm, in good accord with the anomaly positions. However in this case the resonance features are extremely narrow, on the order of just a few nm. We therefore carried out our FDTD integrations to a much longer time compared to the gold case to obtain the result shown. Some artificial ringing is evident in the FDTD results of Fig. 5(b) due to the finite propagation length and narrow spectral features. The gold case might also have Wood’s anomaly contributions, but it is dominated by the stronger LSP resonance and thus stronger SPP-BW features. Choosing system parameters (thicker films or different hole sizes) to lessen the importance of the LSP resonance could increase the relative importance of the Wood’s anomaly features.

The above correlation of the zero order SPP-BW positions with transmission minima may seem surprising. However, it has been previously noticed in slits [31, 32].

#### 4.2 Field patterns

The correlation of zero-order SPP-BW positions, Eq. (8), with transmission minima in Fig. 5(a) is not sufficient proof of their role in the transmission resonance structures. Equation (8) is an approximate relation; higher level considerations will not only shift such positions but lead to band structures with band gaps [12]. Furthermore, coupling of such states with external light leads to broadening and the existence of states with SPP-BW character over continuous frequency ranges. We can, however, simply inspect how the actual electric fields near the metal surfaces change as one scans through the wavelengths of a resonance feature and deduce explicitly the SPP-BWs and other features that are present.

Figure 6 has *x-z* cuts of the |*E*_{z}
|, the surface normal component of the electric field as we scan through wavelengths encompassing the resonance structure in the 610–660 nm range of Fig. 5(a). We associated this structure with a zero-order (1,0)
_{air}
SPP-BW from Eq. (8) at 628 nm. Indeed, one sees from this figure the growth and disappearance of evanescent components that can be assigned as *n*_{x}
=±1, *n*_{y}
=0 on the metal/air side. Inspection of the field in the *x-y* plane for a fixed *z* near the top of the metal (not shown), further confirms that there is an *n*_{x}
=1 nodal structure along *x* but not *y*. We find waves that can be characterized as having a significant (1,0)
_{glass}
SPP-BW component actually exist over the entire, continuous 620 nm–655 nm range. Notice that Figs. 6e and 6f are showing the start of a similar range of states but characterized by (1,1)
_{glass}
on the glass side, consistent with the zero-order SPP-BW wave predicted to be near *λ*=688 nm. |*E*_{z}
| for wavelengths near and including the minimum in the transmission at 628 nm all look almost identical to Fig. 6(c), which is for 630 nm, i.e. there is nothing qualitatively different about this wavelength region insofar as |*E*_{z}
| is concerned. |*E*_{z}
| for wavelengths near the transmission maximum at 646 nm all look very similar to Fig. 6(d) which is for *λ*=640 nm. While qualitatively similar to all the other |*E*_{z}
| panels exhibiting SPP-BW character, the maximum region can be correlated with the most intense and extended such waves. One expects the system to display absorption maxima where the SPP-BWs are most intense, and this is consistent with the experimental results of Ref. [30]. We will see next that in the maximum region these waves are also strongly coupled to hole excitations.

Figure 7 is similar to Fig. 6 but corresponds to the electric field component, *E*_{x}
, along the incident polarization axis. As one scans through the 610 nm–660 nm range we see a localized nanohole excitation, also coupled to SPP-BWs (see |*E*_{z}
|^{2} in Fig. 6), with the hole excitation being most intense near the transmission maximum which occurs near the Fig. 6(d) wavelength.

Along with SPP-BW features Wood’s anomaly features can also occur. As noted in the discussion of Fig. 5, Wood’s anomalies are very narrow spectrally and SPP-BWs, when present, tend to have a greater impact. Nonetheless the frequency-resolved fields near Wood’s anomalies do indeed show their presence for both metal and PEC cases. Figure 8 shows |*E*_{z}
|^{2} at *λ*=600 nm, where a (1,0)
_{air}
Wood’s anomaly feature at the top air/metal interface is predicted. The 2-D hole array in gold case, Fig. 8(a), indeed shows a relatively weak, extended wave along the top surface, characteristic with this Wood’s anomaly [31]. Notice that the bottom surface, however, also shows SPP-BW character on the bottom metal/glass interface. The corresponding PEC result is given in Fig. 8(b); in this case the Wood’s anomaly dominates since no plasmon features are possible.

#### 4.3 Fano Lineshape Analysis

Many of the nearly isolated resonance features of Fig. 5 are not Lorentzian in form but more typical of Fano resonances [29] which can be asymmetric with distinct minima and maxima. Refs. [13] and [8] introduced use of the Fano model in the context of 2-D hole array problems. The model arises by assuming a zero-order discrete state, *φ*_{a}
, embedded in a zero-order continuum of states. For an isolated resonance, such as the *λ*=922-945 nm structure in Fig. 5(a), the transmission is then expected to be of the form

*T*_{b}
is associated with transmission that is uncoupled to the discrete state. We define our transmission, Eq. (9), to be relative to the direct penetration, which is essentially *T*_{b}*. T*_{a}
is associated with transmission due to a zero-order continuum state that, upon allowing for coupling, will mix with the zero-order discrete state (see [29] for details). If *T*_{a}
is a constant, the maximum relative transmission is *T*_{a}
(*1+q*^{2}
) which occurs for *ε*=1/
_{q}
, and the minimum is 0, which occurs for *ε*=-*q*. Notice that the nearly isolated resonances in Fig. 5(a) have transmission minima very close to zero, which is a strong indicator of the validity of the Fano model. The parameters of the model of relevance to us are thus *T*_{a}
, *ω*_{r}
and *γ*_{r}
. The parameter *ω*_{r}
lies between the minimum and maximum and characterizes the resonance position. It can be thought of as the position of the zero-order discrete state *plus* a shift due to interaction with the continuum: *ω*_{r}
=*ω*_{a}
+Δ [29]. The parameter *γ*_{r}
characterizes the lifetime of the discrete state, 1/*γ*_{r}
, and the linewidth. Actually since the resonances are not perfectly isolated, a better approach is to adopt a multiple resonance model, which can be inferred from Eq. (65) of Ref. [29], and adapted to our case to be:

We fit Eq. (10) to the *λ*=550 nm - 1000 region of Fig. 5(a), allowing for four resonances. The result, Fig. 9, agrees remarkably well with our FDTD transmission results. Recently, Eq. (10) was also shown to fit reasonably well experimental hole array data for infrared transmission, although the experimental data differed from Eq. (10) by not having zero or near zero transmission minima [33]. The present demonstration, comparing Eq. (10) with rigorous FDTD simulation results that do not have experimental imperfections shows conclusively that the Fano resonanc e form is indeed correct.

The Fano model is consistent with the interaction of a zero-order discrete state and a zero-order continuum state. Constructive and destructive interference of these states leads to the Fano lineshapes. The precise nature of these states is, however, somewhat elusive. However, the fact that Eq. (8) can be used to predict the approximate locations of the resonance features, as well as the field patterns of Sec. 4.2 which display SPP-BW features, indicate that SPP-BWs, when they can exist are good approximations to the discrete states. Wood’s anomalies are the likely candidates for the discrete states in the case of PEC metals and they can also contribute in the spectral regions for which they are allowed in the metal case. The continuum state in this Fano picture is an isolated hole excitation state which, in addition to having Bethe-Bouwkamp type contributions, could also include contributions from LSP excitation when plasmon excitations are possible.

The interference effect that produces the minima and maxima can be explicitly visualized from our FDTD calculations. Figure 10 shows the charge distribution inferred from ∇·**D** at the wavelength of the transmission minimum and maximum for the (1,0)
_{air}
mode. At the transmission minimum, the charge distribution on the side wall of each individual nanohole shows an anti-symmetric pattern between the top and bottom portion of the hole. We can imagine this pattern arising from the interference of SPP-BWs on the top surface with LSP or isolated hole excitation. In contrast, a symmetric charge distribution is observed at the transmission maximum, consistent with constructive interference of the SPP-BW and isolated hole excitation states.

The bandgap interpretation [12] of the transmission spectra is similar in spirit to the Fano model, involving mixing of discrete and continuum states. In essence the bandgap picture imagines perturbing the zero-order SPP-BWs consistent with Eq. (8) by considering some sort of more realistic model of the system. This splits some of the degeneracies associated with the zero-order states. Associated with each transmission resonance feature are upper and lower frequency modes that have resulted from the splitting. The upper frequency modes are long-lived and couple weakly to outside light. The lower frequency modes have short lifetimes and couple better to outside light. The upper state position is thus associated with the transmission minimum and the lower state position with the transmission maximum, i.e. the maximum is red-shifted relative to the minimum. We have performed some FDTD calculations to verify some aspects of the bandgap interpretation. Rather than light incident from the glass side, we have generated SPPs using a total field/scattered field technique [14] on the top of the gold surface and analyzed their interaction with the holes. The steady state pattern shows that the SPP wave stays in the metal film region at the wavelength of the transmission minimum and is most intense near the hole edge at the wavelength of transmission maximum. This qualitatively agrees with the bandgap model expectation where at the transmission minimum (the top bandgap edge) the lifetime of SPP-BW is infinitely long and thus less likely to couple to holes or the outside. While for the lower bandgap edge, the lifetime is finite and thus likely to couple. Because of the lossy characteristics of the SPP waves, there is no open bandgap for the SPP-BWs in our case.

## 5. Concluding remarks

We carried out FDTD simulations of nanoscale holes in thin gold films on glass substrate with incident light from the glass side. We focused on both SPP generation in the film and light transmission. Calculations treating the metal film as gold with a dielectric constant that supports surface plasmon generation and as a PEC that cannot have any surface plasmons were presented.

Regarding isolated nanoholes, we correlated transmission maxima with LSP excitations of the holes. Traveling SPPs on the top and bottom surfaces are also generated by the holes, and their most intense production was correlated with the transmission maxima. The corresponding PEC transmission results were less structured and significantly less than the gold results, indicating that for a single nanohole plasmon excitation clearly enhances transmission. These results complement results obtained by Wannemacher for silver with different methods. We also augmented our previous discussion [3] of the origin of SPP fringe patterns and the role of NSOM probes by carrying out calculations that explicitly include the probe. We found that the flux measured by the probe is more reflective of the square of the electric field component parallel to the incident polarization axis and not the total intensity.

Periodic arrays of nanoholes display significantly richer structure. For the particular system considered, we find only relatively modest transmission enhancement relative to isolated nanoholes. This result is consistent with results for different systems discussed by Lezec and Thio [10]. Despite the occurrence of numerous maxima and minima, the transmission spectrum retains overall features of the isolated nanohole spectrum, indicating the importance of LSP excitations. We investigated how SPP-BWs and Wood’s anomalies, i.e. special modes due to the periodic grating structure, lead to both transmission minima and transmission maxima through coupling with more directly transmitting continuum states. This was accomplished with analysis of the spectra based on a multiple-resonance Fano model and detailed inspection of the field components.

The relative contributions of SPP-BW, Wood’s anomaly, LSP as well as other features (e.g., waveguide modes) of course depend on the system parameters, which accounts for the sometimes conflicting theoretical models for the enhanced transmission. For example, in gold film cases where the hole size is far away from an LSP resonance, e.g. thicker films or very different hole diameters, the Wood’s anomaly features could be comparable to the LSP/SPP-BW Fano profiles. We have carried out some calculations to see how changing the system parameters influence the results. For example, a thicker 250 nm metal film leads to a blue shifting of the transmission peaks by about 10~20 nm and lower transmission. If the holes are filled with high dielectric material, we find that the transmission is enhanced. Finally, we note that actual experimental condition will involve smoothed over edges. Preliminary calculations indicate that smoothed edges lead to somewhat weaker SPP-BW and Wood’s anomaly features and slightly higher transmission. A more thorough investigation of such features, as well as possible deviations of material properties from their bulk values is planned.

## Acknowledgments

This work was supported by the U.S. Department of Energy, Basic Energy Sciences, under contract W-31-109-ENG-38. GCS was supported by the US Department of Energy under grant No. DEFG02-03-ER15487. The calculations were performed at the High Performance Computing Research Facility, Mathematics and Computer Science Division, Argonne National Laboratory.

## References and Links

**1. **T.W. Ebbesen, H.J. Lezec, H.F. Ghaemi, T. Thio, and P.A. Wolff, “Extraordinary optical transmission through sub-wavelength hole arrays,” Nature **391**, 667–669 (1998). [CrossRef]

**2. **Focus Issue: Extraordinary light transmission through sub-wavelength structured surfaces, Opt. Express12, pp. 3618–3706 (2004), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-12-16-3618 [PubMed]

**3. **L. Yin, V.K. Vlasko-Vlasov, A. Rydh, J. Pearson, U. Welp, S.-H. Chang, S.K. Gray, G.C. Schatz, D.E. Brown, and C.W. Kimball, “Surface plasmons at single nanoholes in Au films,” Appl. Phys. Lett. **85**, 467–469 (2004). [CrossRef]

**4. **R. Wannemacher, “Plasmon-supported transmission of light through nanometric holes in metallic thin films,” Opt. Commun. **195**, 107–118 (2001). [CrossRef]

**5. **J. Prikulis, P. Hanarp, L. Olofsson, D. Sutherland, and M. Kall, “Optical spectroscopy of nanometric holes in thin gold films,” Nano Lett. **4**, 1003–2007 (2004). [CrossRef]

**6. **A. Degiron, H.J. Lezec, N. Yamamoto, and T.W. Ebbesen “Optical transmission properties of a single subwavelength aperture in a real metal,” Opt. Commun. **239**, 61–66 (2004). [CrossRef]

**7. **M.M.J. Treacy, “Dynamical diffraction in metallic optical gratings,” Appl. Phys. Lett. **75**, 606–608 (1999). [CrossRef]

**8. **M. Sarrazin, J.-P. Vigneron, and J.-M. Vigoureux, “Role of Wood anomalies in optical properties of thin metallic films with a bidimensional array of subwavelength holes,” Phys. Rev. B **67**, 085415 (1–8) (2003). [CrossRef]

**9. **H. F. Ghaemi, T. Thio, D. E. Grupp, T. W. Ebbesen, and H. J. Lezec, “Surface plasmons enhance optical transmission through subwavelength holes,” Phys. Rev. B **58**, 6779–6782 (1998). [CrossRef]

**10. **H.J. Lezec and T. Thio, “Diffracted evanescent wave model for enhanced and suppressed optical transmission through subwavelength hole arrays,” Opt. Express **12**, 3629–3651 (2004), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-12-16-3629 [CrossRef] [PubMed]

**11. **L. Salomon, F. Grillot, A.V. Zayats, and F. de Fornel, “Near-field distribution of optical transmission of periodic subwavelength holes in a metal film,” Phys. Rev. Lett. **86**, 1110–1117 (2001). [CrossRef] [PubMed]

**12. **S.A. Darmanyan and A.V. Zayats, “Light tunneling via resonant surface plasmon polariton states and the enhanced transmission of periodically nanostructured metal films: An analytical study,” Phys. Rev. B **67**, 035424(1–7) (2003). [CrossRef]

**13. **C. Genet, M.P. van Exter, and J.P. Woerdman, “Fano-type interpretation of red shifts and red tails in hole array transmission spectra,” Opt. Commun. **225**, 331–336 (2003). [CrossRef]

**14. **A. Taflove and S.C. Hagness, *Computational Electrodynamics: The Finite-difference Time-Domain Method*, Second Edition, (Artech House, Boston, 2000).

**15. **S.K. Gray and T. Kupka, “Propagation of light in metallic nanowire arrays: Finite-difference time-domain results for silver cylinders,” Phys. Rev. B **68**, 045415(1–11) (2003). [CrossRef]

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

**17. **G. Guiffaut and K. Mahdjoubi, “A parallel FDTD Algorithm using the MPI library,” IEEE Ant. and Prop. Mag. **43**, 94–103 (2001). [CrossRef]

**18. **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]

**19. **H. Raether, *Surface Plasmons on Smooth and Rough Surfaces and on Gratings*, (Springer-Verlag, New York, 1988).

**20. **R.D. Grober, T. Rutherford, and T.D. Harris, “Model approximation for the electromagnetic field of a near-field optical probe,” Appl. Optics **19**, 3488–3495 (1996). [CrossRef]

**21. **G.B. Arfken and H.J. Weber, *Mathematical Methods for Physicists*, (Academic Press, New York, 1995).

**22. **H. Kano, S. Mizuguchi, and Satoshi Kawata, “Excitation of surface-plasmon polaritons by a focused laser beam,” J. Opt. Soc. Am. B , **15**, 1381 (1998) [CrossRef]

**23. **H.A. Bethe, “Theory of diffraction by small holes,” Phys. Rev. **66**, 163–182 (1944). [CrossRef]

**24. **C.J. Bouwkamp, “On Bethe’s theory of diffraction by small holes,” Philips Res. Rep. **5**, 321–332 (1950).

**25. **K.L. Kelly, E. Coronado, L.L. Zhao, and G.C. Schatz, “The optical properties of metal nanoparticles: The influence of size, shape, and dielectric environment,” J. Phys. Chem. B **107**, 668–677 (2003). [CrossRef]

**26. **C. F. Bohren and D. R. Huffman, *Absorption and Scattering of Light by Small Particles*, Second Edition (Wiley, New York, 1983) p. 344.

**27. **J.D. Jackson, *Classical Electrodynamics*, Second Edition (Wiley, New York, 1975) pp. 438–441.

**28. **A. Krishnan, T. Thio, T.J. Kim, H.J. Lezec, T.W. Ebbesen, P.A. Wolff, J. Pendry, L. Martin-Moreno, and F.J. Garcia-Vidal, “Evanescently coupled resonance in surface plasmon enhanced transmission,” Opt. Commun. **200**, 1–7 (2001). [CrossRef]

**29. **U. Fano, “Effects of Configuration Interaction on Intensities and Phase Shifts,” Phys. Rev. **124**, 1866–1878 (1961). [CrossRef]

**30. **W. L. Barnes, W. A. Murray, J. Dintinger, E. Devaux, and T.W. Ebbesen, “Surface plasmon polaritons and their role in the enhanced transmission of light through periodic arrays of sub-wavelength holes in a metal film,” Phys. Rev. Lett. **92**, 107401(1–4) (2004). [CrossRef] [PubMed]

**31. **J. M. Steele, C.E. Moran, A. Lee, C.M. Aguirre, and N.J. Halas, “Metallodielectric gratings with subwavelength slots: Optical properties,” Phys Rev. B **68**, 205103(1–7) (2003). [CrossRef]

**32. **Q. Cao and P. Lalanne, “Negative Role of Surface Plasmons in the Transmission of Metallic Gratings with Very Narrow Slits,” Phys. Rev. Lett. **88**, 057403(1–4) (2002). [CrossRef] [PubMed]

**33. **W. Fan, S. Zhang, B. Minhas, K.J. Malloy, and R.J. Brueck, “Enhanced infrared transmission through subwavelength coaxial metallic arrays,” Phys. Rev. Lett. **94**, 033902(1–4) (2005). [CrossRef] [PubMed]