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+i1.8, in reasonable accord with the interpolated empirical value [16] of -4.9+i2.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 H, three additional field components associated with a current density vector, J. The complete system of time-domain equations for E, H, and J is designed to be 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)}$ ⁽¹⁾ (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^SPPz=ẑ ·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 ∇·E^SPP=0 is not exactly satisfied by Eq. (4), it is satisfied in the large 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 |² 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.

 

Fig. 1. Time-averaged |E_x |² for an isolated d=200 nm diameter nanohole in a 100 nm thick gold film on top of a glass layer: (a) z=4 nm, x-y profile; (b) x-z cut for y=0 of result in (a); (c) same as (b) but for the case of no hole present; (d) result of subtracting the amplitudes in (b) and (c) and taking the absolute square of the result. The field shown is the Fourier transform at incident wavelength λ=532 nm of our FDTD result.

Download Full Size | PDF

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 ₁ exp(-ρ′/L), and Im[E_z ]∝Y₁ 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₁ exp(-ρ′/L) and Im[E_x ]∝- J₁ exp(-ρ′/L). Figure 2 is indeed consistent with these relations. For x<100 nm, the simple point source SPP model is not applicable.

 

Fig. 2. Profiles of scaled field components, E_x′=E_x k_spp /(Aα) and E_z ′=E_z/A, along the incident polarization direction (x with y=0) just above (z=4 nm) the metal film. Symbols are the FDTD results and curves correspond to the Bessel function fits described in the text.

Download Full Size | PDF

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 |²=|E_x |²+|E_y |² 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 |² 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 |², rather than the total intensity which is also shown. The probe-free |E_parallel |² 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.

 

Fig. 3. Comparison of theoretically estimated NSOM signal (solid red symbols) with the total intensity (dashed black line) and the intensity of parallel electric field component (solid blue line) on the metal surface of a single nanohole.

Download Full Size | PDF

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/λ⁴ 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.

 

Fig. 4. Transmission spectrum of a single d=200 nm nanohole in a 100 nm thick metal film on glass, treating the metal as gold (solid curve) and as a perfect electrical conductor (dashed curve). The metal film is sandwiched between glass and air, with the light incident from the glass side.

Download Full Size | PDF

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

 

Fig. 5. Transmission spectra for 2-D periodic square arrays of d=200 nm nanoholes with square lattice spacing D=600 nm in a 100 nm thickmetal film. (a) Treating the metal film as gold. (b) Treating the metal film as a PEC.

Download Full Size | PDF

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.

 

Fig. 6. | E_z |² for the 2-D square hole array in a gold film of Fig. 5(a). (a)–(f) correspond to λ=610 nm, 620 nm, …, 660 nm.

Download Full Size | PDF

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 |² in Fig. 6), with the hole excitation being most intense near the transmission maximum which occurs near the Fig. 6(d) wavelength.

 

Fig. 7. | E_x |² for the 2-D hole array in the gold film of Fig. 5(a). (a)–(f) correspond to λ=610 nm, 620 nm, …, 660 nm.

Download Full Size | PDF

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 |² 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.

 

Fig. 8. |E_z |² for the 2D hole array in (a) gold and (b) PEC films at wavelength 600 nm. The contrast has been increased by a factor of 3 comparing to other figures in order to show the weak Wood’s anomaly features.

Download Full Size | PDF

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² ) 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.

 

Fig. 9. Multiple resonance Fano model fit (solid curve) to the FDTD transmission data (symbols) for the 2D array of holes in a gold film.

Download Full Size | PDF

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.

 

Fig. 10. Charge distribution of the (1,0)_air mode at (a) the transmission minimum, and (b) the transmission maximum.

Download Full Size | PDF

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.