Surface Plasmon Resonance (SPR) Sensing Simulator

Theory & Methods: Surface Plasmon Resonance Sensing

1. SPP Dispersion and Momentum Matching

Surface plasmon-polaritons (SPPs) are hybrid electromagnetic waves—coupled transverse-magnetic (TM) photons and collective electron-density oscillations—that propagate along the interface between a metal (with a complex relative permittivity, Re[ε_m] < 0) and a dielectric (Re[ε_d] > 0). (The operator Re denotes the real part.) Because the electromagnetic field is coupled to the massive charge carriers, the SPP wave travels slower than light in the bulk dielectric, possessing a larger in-plane momentum (wavevector k_x). For a single metal-dielectric interface, the SPP dispersion relation is defined as:

k_x^SPP = (ω / c) [ (ε_d ε_m) / (ε_d + ε_m) ]^(1/2)

where ω is the angular frequency of the wave and c is the speed of light in vacuum. Consequently, freely propagating light from the dielectric medium cannot excite SPPs directly because it lacks sufficient longitudinal momentum (k_x < Re[k_x^SPP]). To bridge this momentum gap, the Kretschmann-Raether configuration employs attenuated total reflection (ATR).

2. The Kretschmann-Raether Configuration

In this optical geometry, laser light enters a high-index coupling prism (n₁) and strikes the base at an internal angle θ_int greater than the critical angle for total internal reflection (TIR). Under TIR conditions, the incident wave does not propagate into the lower-index dielectric medium (n₄); instead, it generates an evanescent field that decays exponentially into the thin active metallic film (Gold, n₃).

If the metallic layer is sufficiently thin (e.g., 50 nm), this evanescent field tunnels through the metal and reaches the other metal-dielectric interface (e.g., Gold-Air). Resonance occurs when the in-plane wavevector of the light in the prism exactly matches the propagation wavevector of the SPP at the outer interface:

(ω / c) n₁ sin(θ_int^res) ≈ Re[ k_x^SPP ]

At this specific angle, the incident photon energy is resonantly transferred into the SPP mode, which subsequently dissipates via Joule heating in the metal (intrinsic damping, γ_in) and leakage radiation back into the prism (radiative damping, γ_rad). This extraction of energy manifests as a sharp dip in the macroscopic reflectance curve.

3. The 4-Layer Recursive Fresnel Model

Predicting the macroscopic reflectance profile R(θ_int), including its minimum intensity and Full-Width at Half-Maximum (FWHM), requires a solution to Maxwell's equations across the stratified boundary. This SPR sensor utilizes a 4-layer architecture: Prism (1), Adhesion Layer (2), Active Metal (3), and Dielectric Environment (4).

The total amplitude reflection coefficient is computed by recursively expanding the generalized Fresnel equations. For an interface between an incident medium i and a transmitted medium j, the single-interface p-polarization reflection coefficient r_ij^p and s-polarization reflection coefficient r_ij^s are:

r_ij^p = (n_j²β_i − n_i²β_j) / (n_j²β_i + n_i²β_j) r_ij^s = (β_i − β_j) / (β_i + β_j)

Here, n is the refractive index and r_ij is the single-interface reflection coefficient. The variables β_i and β_j represent the complex z-components of the wavevector in their respective layers. Governed by the conservation of in-plane momentum across all boundaries, the general form for any layer j is:

β_j = [ (ω / c)² n_j² − k_x² ]^(1/2)

where k_x = (ω / c) n₁ sin(θ_int) is the incident in-plane wavevector. Because s-polarized light lacks a normal electric-field component, it cannot induce the surface charge density necessary to excite SPPs. The net reflection coefficient (r_1,4) and total reflectivity (R_SPR = |r_1,4|²) for the 4-layer stack incorporates the multiple internal reflections and accumulated phase shifts through the intermediate films:

r_1,4 = (r₁₂ + r₂₃ exp[2iβ₂d₂] + r₃₄ exp[2i(β₂d₂ + β₃d₃)] + r₁₂r₂₃r₃₄ exp[2iβ₃d₃]) / (1 + r₁₂r₂₃ exp[2iβ₂d₂] + r₁₂r₃₄ exp[2i(β₂d₂ + β₃d₃)] + r₂₃r₃₄ exp[2iβ₃d₃]) R_SPR = |r_1,4|²

where i² ≡ −1. The variables d₂ and d₃ represent the physical thicknesses of the adhesion layer and the active metallic layer, respectively. This simulation computes the recursive Fresnel expansion.

4. Discrete Convolution for Beam Divergence

Because real lasers emit beams with a finite angular divergence, the measured macroscopic reflectance is a continuous convolution of the ideal plane-wave reflectance and the laser's Gaussian angular intensity distribution. Since an analytical integration of this convolution is intractable, the simulation approximates it using a finite Riemann sum:

R_measured(θ₀) ≈ Σ [ w_j R_ideal(θ₀ + Δθ_j) ] / Σ w_j

To balance real-time computational performance with numerical accuracy, this uniformly sampled discrete convolution is evaluated over an angular window of ±3σ (where σ is the standard deviation of the beam profile) using 21 evenly spaced evaluation nodes (index j). The parameter w_j represents the statistical weight of the j-th node, determined by the underlying Gaussian intensity distribution.

Note: To maintain computational efficiency, the simulation treats the angular beam divergence as invariant across the refracting boundaries.

5. Macroscopic Prism Geometry and Facet Losses

In a practical laboratory setup, the macroscopic geometry of the coupling prism dictates the transformation between the internal angle of incidence (θ_int) and the external angle of incidence (α) relative to the entrance facet normal. For an isosceles triangular prism with base angle A, this transformation is governed by Snell's law:

θ_int = A − sin⁻¹[ (n_ambient / n₁) sin α ]

Within this geometric framework, an increasingly positive external angle (α > 0), measured relative to the facet normal, indicates that the incident beam is on the prism apex side of the normal and rotated toward the prism apex, yielding a smaller internal angle (θ_int < A). Conversely, a negative external angle (α < 0) designates that the incident beam is rotated toward the prism base and resides on the base side of the facet normal, producing a larger internal angle (θ_int > A).

The incident beam undergoes reflection losses at both the entrance and exit facets. For an isosceles geometry, the internal angle of refraction at the entrance facet and the internal angle of incidence at the exit facet are symmetric: θ_t = A − θ_int. The corresponding amplitude reflection coefficients for p-polarized and s-polarized light are:

r_facet^p = (n₁ cos α − n_ambient cos θ_t) / (n₁ cos α + n_ambient cos θ_t) r_facet^s = (n_ambient cos α − n₁ cos θ_t) / (n_ambient cos α + n₁ cos θ_t)

The total macroscopic reflectance integrates these transmittance losses via T_in = T_out = 1 − |r_facet|², such that the final signal evaluates to:

R_macro = (1 − |r_facet|²)² R_SPR

For the hemicylindrical coupling geometry, the beam enters and exits strictly normal to the curved facets (α = 0, θ_t = 0), yielding constant, angle-independent transmission losses. For a mixed-polarization beam containing an s-polarization intensity fraction f_s, the total reflectance is the weighted superposition of the orthogonal states:

R_mix = (1 − f_s) R_macro^p + f_s R_macro^s

6. Bulk Sensitivity of the SPR Sensor

The overall sensitivity S of an SPR sensor can be expressed as the product of the sensitivity of the output signal (reflectance R) to the effective index of the SPP mode (n_eff^SPP) and the sensitivity of the effective index to the refractive index of the dielectric medium (n_d). Applying the chain rule yields:

S = (∂R / ∂n_eff^SPP) (∂n_eff^SPP / ∂n_d)

The first factor depends on the method of SPP excitation (e.g., prism coupling) and modulation, while the second factor represents the intrinsic properties of the SPP mode. To quantify the intensity sensitivity S_I directly from the macroscopic reflectance, a finite-difference approach is utilized:

S_I = ΔR / Δn = [ R(θ_int; n_d + Δn) − R(θ_int; n_d) ] / Δn

Evaluating S_I across the angular domain reveals that the maximum sensitivity occurs at the inflection points of the SPR dip, pinpointing the optimal operating angles θ_opt for differential sensing. The sign of S_I provides direct physical intuition regarding the sensor's operating point relative to the resonance minimum. An increase in the dielectric index (n_d) increases the SPP propagation wavevector; therefore, the momentum-matching condition requires a correspondingly larger incident angle, shifting the entire reflectance curve to the right. Thus, if the sensor operates on the descending left flank of the dip (θ < θ_res), a rightward shift of the curve moves the minimum further from the observation point, causing the measured reflectance to increase (S_I > 0). Conversely, operating on the ascending right flank (θ > θ_res) pulls the resonance minimum closer to the observation point, decreasing the measured reflectance (S_I < 0).

7. Phase Interrogation and Ellipsometry

While amplitude interrogation tracks the dip in reflected intensity, phase interrogation evaluates the differential phase shift between the orthogonal polarization states. The complex reflection coefficient is r = |r| exp(iδ). The differential phase shift Δ is defined as:

Δ = δ_p − δ_s

where the individual phase delays are extracted from the complex Fresnel equations via δ_j = tan⁻¹[ Im(r_j) / Re(r_j) ]. At the SPP resonance angle, the phase of the p-polarized component undergoes a steep topological transition. This sharp phase gradient generally offers a lower limit of detection (LOD) than intensity mapping.

If the incident beam is strictly p-polarized, a direct physical measurement of Δ is impossible without the s-polarized internal reference wave; experimental extraction of δ_p then necessitates an external split-path interferometer. Computationally, however, r_p and r_s operate as intrinsic transfer functions of the stratified media. Simulating Δ under pure p-polarization visualizes the theoretical polarization-dependent phase retardation inherent to the boundary conditions, irrespective of the incident optical field amplitudes.

8. Sensor Figure of Merit (FOM)

The resolving power of the SPR sensor is quantified by its Figure of Merit (FOM), which normalizes the maximum sensitivity against the spectral width of the resonance:

FOM = |S_I|_max / Δθ_FWHM

where Δθ_FWHM is the Full-Width at Half-Maximum of the SPR dip. Narrower resonance curves minimize the overlap of adjacent signals in differential measurements, yielding a superior FOM.

9. Synthetic Measurement Noise

To replicate the stochastic fluctuations of physical photodetectors (shot noise and thermal noise), a pseudo-random Gaussian distribution is applied to the ideal reflectance profile utilizing the Box-Muller transform:

R_measured = R_ideal + σ [ -2 ln(U₁) ]^(1/2) cos(2π U₂)

where σ is the standard deviation (noise amplitude), and U₁ and U₂ are independent uniformly distributed random variables on the interval (0, 1].